[XLS-SVY-01]: Prediksi Pasut dengan VBA Excel

Referensi : Harmonic Analysis and Prediction of Tides
Platform : Excel 2003
Lokasi File : <<available upon request>>

 

RUMUS YANG DIGUNAKAN

Prediksi pasang surut (pasut) menggunakan rumus Harmonic Analysis dengan metode least square untuk mencari 9 constituents  utama pembangkit pasut dan 9 phase.

Dalam referensi tertulis,

“The ideal tide curve for any given port is represented as an average height Z0 plus a sum of terms (“constituents”) each of which is of the form f(t) = H cos(at + \phi). The time t is measured in hours, and f comes out in feet. The numbers H,a, \phi are the amplitude, the speed and the phase of the constituent.”

Karena data hasil bacaan pasut yang tersedia adalah dalam meter dan dicatat tiap jam, maka nilai f tidak lagi dinyatakan dalam feet, melainkan dalam meter.

Rumus kemudian dikembangkan menjadi:

image

dalam format matrix dapat ditulis:

image

Jika MatrixX telah didapat, maka dengan menggunakan rumus(4) dan rumus(5),amplitudo dan phase tiap constituent dapat dihitung.

INPUT (DATA HASIL PENGUKURAN MUKA AIR)

Hasil pengukuran muka air didapat dari pengamatan selama 29 hari dengan pencatatan tiap jam, dimulai jam 0:00 sampai jam 23:00 tanggal 1 Maret 2008 sampai jam 23:00 tanggal 29 Maret 2008 yang ditabelkan sebagai berikut:

image

PROSES PEMBACAAN DATA DAN PERHITUNGAN

Proses pembacaan dan perhitungan dengan menggunakan program Visual Basic Application (macro) yang ada di microsoft excel 2003.

Khusus untuk prosedure atau program perkalian matirx dan inverse didapat dari http://www.alglib.net/

Pada file excel spread sheet saya, data pengukuran muka air pertama (1.850m) ada di cell $C$10 sedangkan data terakhir (2.110m) ada di $Z$38, sehingga range data pengukuran adalah di $C$10:$Z$38.

Pastikan bahwa di bagian Declaration di awal program diset :

Option Explicit ‘berguna untuk mendeteksi definisi variable
Option Base 1 ‘hitungan matrix dimulai dari 1, kalo tidak diset, maka default index matrix adalah 0

Sub PrediksiPasutDenganLeastSquare()

‘1. membaca data pengukuran muka air
dim cr as Range, rgData as Range

set rgData=Activesheet.Range(“$C$10:$Z$38”)

For Each cr In rgData.Columns(1).Cells
i = i + 1
j = 1 + (i – 1) * 24
ReDim Preserve MatrixL(1 To j + 23)
For jm = 0 To 23
MatrixL(j + jm) = cr.Offset(0, jm).Value
Next jm
Next

‘2. membaca atau set variable periode tiap Constituents

Dim w(1 To 9) As Double ‘periode dari 9 Constituents pasut
Dim pi As Double
pi = 4 * Atn(1)
w(1) = 2# * pi / 12.4206 ‘M2
w(2) = 2# * pi / 12# ‘S2
w(3) = 2# * pi / 12.6582 ‘N2
w(4) = 2# * pi / 11.9673 ‘K2
w(5) = 2# * pi / 23.9346 ‘K1
w(6) = 2# * pi / 25.8194 ‘O1
w(7) = 2# * pi / 24.0658 ‘P1
w(8) = 2# * pi / 6.2103 ‘M4
w(9) = 2# * pi / 6.1033 ‘MS4

‘3. membuat atau membetuk MatrixA() atau Matrix coeffisien

Dim it As Integer,MatrixA() as Double
ReDim Preserve matrixA(1 To UBound(MatrixL), 1 To 19)

For i = 1 To UBound(MatrixL)
matrixA(i, 1) = 1
it = i
For j = 1 To 9
matrixA(i, 2 * j) = Cos(w(j) * it)
matrixA(i, 2 * j + 1) = -Sin(w(j) * it)
Next j
Next i

‘4. proses least square

Dim weight() As Double: ReDim Preserve weight(UBound(MatrixL, 1))
For i = LBound(weight) To UBound(weight): weight(i) = 1: Next i ‘matrix bobot=matrix identitas
MatrixX() = clsLSQ.LSPAR(matrixA, MatrixL, weight)

‘clsLSQ adalah class yang saya buat untuk proses perhitungan least square, tidak dibahas di session ini

‘5. menghitung amplitudo dan phase 9 konstituent dan sekaligus menampilkan hasilnya di excel
Dim A As Double, B As Double, H(1 To 9) As Double, Phase(1 To 9) As Double
Dim ph As Double, Zo As Double

Const addPrint As String = “H66”

With Range(addPrint)
.Offset(0, 3) = MatrixX(1)’mencetak Zo atau mean sea level
For i = 2 To 19 Step 2
j = i / 2
A = MatrixX(i): .Offset(j, 0) = A
B = MatrixX(i + 1): .Offset(j, 1) = B
ph = Atn(B / A) ‘kwadran I
If A < 0 Then
ph = ph + pi ‘kwadran II dan III
Else
If B < 0 Then ph = ph + 2 * pi ‘kwadran IV
End If
‘phase dikonversi ke derajat
Phase(j) = ph * 180 / pi: .Offset(j, 2) = Phase(j)
H(j) = Sqr(A * A + B * B): .Offset(j, 3) = H(j)

Next i
End With

‘6. Membandingkan muka air hasil pengukuran  dengan muka air hasil hitungan
Dim Ht As Double, SumHCos As Double
Dim Cetak() As Double

ReDim Preserve Cetak(1 To UBound(MatrixL), 1 To 5)

Zo = MatrixX(1)
For i = LBound(MatrixL) To UBound(MatrixL)
Cetak(i, 5) = CDbl(i)’mencetak nomer urut
Cetak(i, 1) = Range(addFirstDate) + (i – 1) / 24’mencetak hari dan jam
Cetak(i, 2) = MatrixL(i)’mencetak muka air pengukuran
SumHCos = 0
For j = 1 To 9
SumHCos = SumHCos + H(j) * Cos(w(j) * i + (Phase(j) * pi / 180))’rumus harmonic
Next j
Ht = Zo + SumHCos ‘rumus(1)
Cetak(i, 3) = Ht’cetak muka air hasil hitungan
Cetak(i, 4) = Ht – MatrixL(i)’muka air hitungan – pengukuran
Next i
Range(“A90:E785”) = Cetak

End Sub

Dengan menggunakan data perbandingan antara muka air hasil pengukuran dengan perhitungan, maka dapat dibuat graphic sebagai berikut:

image

dari hasil hitungan didapat standard deviasi sebesar 0.11meter.

Advertisements

18 thoughts on “[XLS-SVY-01]: Prediksi Pasut dengan VBA Excel

    1. mengacu ke rumus.1, jika akan digunakan 15 komponen, maka diperlukan informasi tambahan periode untuk sebanyak 15buah. Pada tulisan ini saya menggunakan informasi periode sebanyak 9 buah untuk 9 konstituent

      1. terkait hasil akhir dari 9 komponen dengan 15 komponen, apakah amplitudo setiap komponen memiliki perbedaan? karena setelah saya menggunakan invers 15 komponen (dengan matriks 31×31) hasil amplitudonya sama semua :(

      2. Jika menggunakan 15 komponen, maka diperlukan informasi tambahan kecepatan sudut untuk masing-masing komponen atau constituent.

        sehingga baris program di vba:

        ’2. membaca atau set variable periode tiap Constituents

        Dim w(1 To 15) As Double ‘periode dari 15 Constituents pasut
        Dim pi As Double
        pi = 4 * Atn(1)
        w(1) = 2# * pi / 12.4206 ‘M2
        w(2) = 2# * pi / 12# ‘S2
        w(3) = 2# * pi / 12.6582 ‘N2
        w(4) = 2# * pi / 11.9673 ‘K2
        w(5) = 2# * pi / 23.9346 ‘K1
        w(6) = 2# * pi / 25.8194 ‘O1
        w(7) = 2# * pi / 24.0658 ‘P1
        w(8) = 2# * pi / 6.2103 ‘M4
        w(9) = 2# * pi / 6.1033 ‘MS4
        w(10) = 2# * pi / <>
        w(11) = 2# * pi / <>
        w(12) = 2# * pi / <>
        w(13) = 2# * pi / <>
        w(14) = 2# * pi / <>
        w(15) = 2# * pi / <>

        Perubahan ini akan berefek ke dimensi matrix A. Misal pengukuran selama 30 hari dengan interval 1 jam, maka banyaknya observasi level=30 x 24 x 1 = 720jam Matrix A (720 x 31), MatrixX(31 x 1) dan MatrixL (720 x 1)

    1. dengan data yang berbeda 9 kompenen tentu akan berbeda hasilnya. jadi data hasil perhitungan dengan data saya akan berbeda hasilnya jika diterapkan di tempat lain

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out /  Change )

Google+ photo

You are commenting using your Google+ account. Log Out /  Change )

Twitter picture

You are commenting using your Twitter account. Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. Log Out /  Change )

Connecting to %s