Analisis Data Geofisika

Analisis Data Geofisika: Memahami Teori Inversi

Dr. Eng. Supriyanto, M.Sc Edisi I

Departemen Fisika-FMIPA Univeristas Indonesia 2007

Kata Pengantar

Buku Buku ini semula semula merupa merupakan kan diktat diktat perkuli perkuliaha ahan n mata mata kuliah kuliah Analisi Analisiss Data Data Geofisi Geofisika ka yang yang diberi diberikan kan kepada mahasiswa Geofisika pada Departemen Fisika, Fakultas MIPA, Universitas Indonesia. Acuan utama untuk edisi perdana ini adalah buku Geophysica Geophysicall Data Analysis: Analysis: Understandin Understanding g Inverse Problem Theory and Practice yang Practice yang ditulis oleh Max A. Meju dan diterbitkan oleh Society of Exploration Geophysicists pada tahun 1994. Semoga Semoga kebera keberadaa daan n buku buku ini dapat dapat memban membantu tu mahasi mahasiswa swa geofisi geofisika ka untuk untuk memilik memilikii kemampuan memformulasikan masalah, menyusun hipotesis, metode dan solusi sehingga mampu menyelesaikan masalah-masalah geofisika secara mandiri. Terima kasih yang tak terhingga ingin saya sampaikan kepada Dede Djuhana yang telah bersedia berbagi memberikan file format buku dalam L A TE X sehingga tampilan buku ini menjadi jauh lebih baik. Depok, 15 Maret 2007 Supriyanto S.

v

Kata Pengantar

Buku Buku ini semula semula merupa merupakan kan diktat diktat perkuli perkuliaha ahan n mata mata kuliah kuliah Analisi Analisiss Data Data Geofisi Geofisika ka yang yang diberi diberikan kan kepada mahasiswa Geofisika pada Departemen Fisika, Fakultas MIPA, Universitas Indonesia. Acuan utama untuk edisi perdana ini adalah buku Geophysica Geophysicall Data Analysis: Analysis: Understandin Understanding g Inverse Problem Theory and Practice yang Practice yang ditulis oleh Max A. Meju dan diterbitkan oleh Society of Exploration Geophysicists pada tahun 1994. Semoga Semoga kebera keberadaa daan n buku buku ini dapat dapat memban membantu tu mahasi mahasiswa swa geofisi geofisika ka untuk untuk memilik memilikii kemampuan memformulasikan masalah, menyusun hipotesis, metode dan solusi sehingga mampu menyelesaikan masalah-masalah geofisika secara mandiri. Terima kasih yang tak terhingga ingin saya sampaikan kepada Dede Djuhana yang telah bersedia berbagi memberikan file format buku dalam L A TE X sehingga tampilan buku ini menjadi jauh lebih baik. Depok, 15 Maret 2007 Supriyanto S.

v

Daftar Isi

Lembar Persembahan

i

Kata Pengantar

v

Daftar Isi

vii

Daftar Gambar

ix

Daftar Tabel

xi

1 Pendahuluan

1

1.1 Definisi dan Konsep Dasar . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

1

1.2 Proses geofisika . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

3

1.3 Eksplorasi geofisika dan inversi . . . . . . . . . . . . . . . . . . . . . . . . . . . .

3

1.4 Macam-macam data geofisika . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

3

1.5 Deskripsi proses geofisika: Model matematika ika . . . . . . . . . . . . . . . . . . . .

4

1.6 Diskritisasi dan linearisasi . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

5

2 Formulasi Masalah Inversi

9

2.1 Klasifikasi masalah inversi . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

9

2.2 Inversi Model Garis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

9

2.3 Inversi Model Parabola . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

15

2.4 Inversi Model Bidang . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

20

2.5 Contoh aplikasi . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

24

2.5.1

Menghitung gravitasi di planet X . . . . . . . . . . . . . . . . . . . . . . .

24

2.5. 2.5.2 2

Anal Analis isa a data data seis seismi mik k pada pada refle reflekt ktor or tung tungga gall hori horizo zont ntal al . . . . . . . . . . .

28

2.5. 2.5.3 3

Anal Analis isa a data data seis seismi mik k pada pada refle reflekt ktor or tung tungga gall miri miring ng . . . . . . . . . . . . .

29

2.6 Kesimpulan . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

30

3 Penyelesaian Masalah Overdetermined

31

3.1 Regresi linear sederhana . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

31

3.2 3.2 Meto Metode de least least square . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

31

3.3 3.3 Ap Apli lika kasi si regr regres esii line linear ar pada pada anal analis isis is data data seis seismi mik k refr refrak aksi si . . . . . . . . . . . . . .

34

3.4 Regresi linear dengan gan pendekatan matriks . . . . . . . . . . . . . . . . . . . . . .

36

3.5 Soal . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

39

vii

viii 4 Constrained Linear Least Squares Inversion 4.1 Inversi dengan informasi awal . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

41 41

4.1. 4.1.1 1

Memf emformu ormula lasi sika kan n pers persam amaa aan n ter terkons konstr trai ain n . . . . . . . . . . . . . . . . . .

43

4.1.2

Contoh aplikasi inversi terkonstr strain . . . . . . . . . . . . . . . . . . . . . .

43

4.2 Inversi dengan Smoothness . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

46

4.2.1

Formulasi masalah . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

46

4.2.2

Solusi masalah . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

47

Daftar Pustaka

49

Indeks

51

Daftar Gambar

1.1 Alur pemodelan inversi . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

2

1.2 Alur pemodelan forward . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

2

1.3 1.3 Alur Alur eksp eksper erim imen en lapa lapang ngan an dan dan eksp eksper erim imen en labo labora rato tori rium um . . . . . . . . . . . . . .

4

1.4 1.4 Konfig onfigur uras asii ele elektr ktroda oda pada ada meto etode Schl Schlum umbe berg rger er . . . . . . . . . . . . . . . . . .

6

2.1 Data Data observ observasi asi peruba perubahan han suhu suhu terhad terhadap ap kedalam kedalaman an dari dari permuk permukaan aan tanah tanah . . .

10

2.2 2.2 Hasi Hasill inver inversi si atas atas data data obse observ rvas asii peru peruba baha han n suhu suhu terh terhad adap ap keda kedala lama man n . . . . . .

14

2.3 Data Data observ observasi asi peruba perubahan han suhu suhu terhad terhadap ap kedalam kedalaman an dari dari permuk permukaan aan tanah tanah . . .

15

2.4 Hasil Hasil invers inversii atas atas data data observ observasi asi perubah perubahan an suhu terhad terhadap ap kedalama kedalaman. n. Tanda anda titik merah merah adalah adalah data observa observasi si sementara sementara kurva kurva biru adalah adalah kurva hasil hasil inversi inversi

18

2.5 2.5 Data Data obse observ rvas asii peru peruba baha han n nila nilaii terh terhad adap ap posi posisi si koor koordi dina natt x dan dan y . . . . . . . . .

20

2.6 2.6 Hasi Hasill inv inver ersi si berb berben entu tuk k bid bidan ang. g. Saya Saya sebu sebutt saj saja a seb sebag agai ai bida bidang ng inve invers rsii . . . . . .

23

2.7 Ini adalah hasil hasil inversi inversi yang sama, hanya saja saja bidang inversi inversi diputar diputar (di-rotasi) (di-rotasi) pada suatu arah sehingga kita bisa saksikan dan pastikan bahwa sebaran titik data observasi berada di-dekat berada di-dekat permukaan bidang inversi . . . . . . . . . . . . .

23

2.8 Grafik data pengukuran gerak batu . . . . . . . . . . . . . . . . . . . . . . . . . .

26

2.9 Grafik hasil inversi parabola . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

28

2.10 Reflektor mendatar pada kedalaman z. Kecepata Kecepatan n gelombang gelombang v dianggap konstan. S adalah sumber gelombang seismik dan R adalah penerima gelombang seismik. seismik. Jarak antara antara S dan R disebut offset offset (x). Sementara Sementara garis merah merah yang ada ada panah anahny nya a adal adalah ah lint lintas asan an gelo gelomb mban ang g seis seism mik. ik.

. . . . . . . . . . . . . . . . .

29

2.11 Reflektor miring dengan sudut kemiringan sebesar α. Kecepata Kecepatan n gelombang gelombang v dianggap konstan. S adalah adalah sumber gelombang seismik dan R adalah penerima gelomb gelombang ang seismik. seismik. Jarak Jarak antara antara S dan R disebu disebutt offset offset (x). Sement Sementara ara garis garis mera merah h yang yang ada ada pana panahn hnya ya adal adalah ah linta lintasa san n gelo gelomb mban ang g seis seismi mik. k. . . . . . . . . . .

30

3.1 3.1 Hasi Hasill plotting dat data obse obserrvasi vasi dala dalam m sumb sumbuu-x x dan dan sumb sumbuu-y y . . . . . . . . . . . . .

32

3.2 Contoh solusi regresi linear . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

34

ix

Bab 1

Pendahuluan

1.1 Definisi dan Konsep Dasar Dalam geofisika, kegiatan pengukuran lapangan selalu dilakukan berdasarkan prosedur yang sudah ditentukan. Kemudian, hasil pengukuran dicatat dan disajikan dalam bentuk tabel angkaangka pengukuran. Hasil pengukuran tersebut sudah barang tentu sangat tergantung pada kondisi dan sifat fisis batuan bawah permukaan. Tabel angka-angka itu selanjutnya disebut data observasi atau juga biasa disebut data lapangan . Kita berharap data eksperimen dapat memberi informasi sebanyak-banyaknya, tidak sekedar mengenai sifat fisis batuan saja, melainkan juga kondisi geometri batuan bawah permukaan dan posisi kedalaman batuan tersebut. Informasi itu hanya bisa kita dapat bila kita mengetahui hubungan antara sifat fisis batuan tersebut dan data observasinya. Penghubung dari keduanya hampir selalu berupa persamaan matematika atau kita menyebutnya sebagai model matematika. Maka dengan berdasarkan model matematika itulah, kita bisa mengekstrak parameter fisis batuan dari data observasi. Proses ini disebut proses inversi atau istilah asingnya disebut inverse modelling , lihat Gambar 2.9. Sementara proses kebalikannya dimana kita ingin mem-

peroleh data prediksi hasil pengukuran berdasarkan parameter fisis yang sudah diketahui, maka proses ini disebut proses forward atau forward modelling , lihat Gambar 1.2. Proses inversi adalah suatu proses pengolahan data lapangan yang melibatkan teknik penyelesaian matematika dan statistik untuk mendapatkan informasi yang berguna mengenai distribusi sifat fisis bawah permukaan. Di dalam proses inversi, kita melakukan analisis terhadap data lapangan dengan cara melakukan curve fitting (pencocokan kurva) antara model matematika dan data lapangan. Tujuan dari proses inversi adalah untuk mengestimasi parameter fisis batuan yang tidak diketahui sebelumnya (unknown parameter ). Proses inversi terbagi dalam level-level tertentu mulai dari yang paling sederhana seperti fitting garis untuk data seismik refraksi sampai kepada level yang rumit seperti tomografi akustik dan matching (pencocokan) kurva resistivity yang multidimensi. Contoh problem inversi dalam bidang geofisika adalah 1. Penentuan struktur bawah tanah 2. Estimasi parameter-parameter bahan tambang 1

2

BAB 1. PENDAHULUAN

Gambar 1.1: Alur pemodelan inversi

Gambar 1.2: Alur pemodelan forward

1.2. PROSES GEOFISIKA

3

3. Estimasi parameter-parameter akumulasi sumber energi 4. Penentuan lokasi gempa bumi berdasarkan waktu gelombang datang 5. Pemodelan respon lithospere untuk mengamati proses sedimentasi 6. Analisis sumur bor pada hidrogeologi

1.2 Proses geofisika Perambatan gelombang seismik, perambatan gelombang elektromagnetik di bawah tanah dan juga aliran muatan (arus listrik) ataupun arus fluida pada batuan berpori adalah contoh-contoh proses geofisika. Data lapangan tak lain merupakan refleksi dari kompleksitas sistem geofisika yang sedang diamati, yang dikontrol oleh distribusi parameter fisis batuan berikut struktur geologinya.

1.3 Eksplorasi geofisika dan inversi Tujuan utama dari kegiatan eksplorasi geofisika adalah untuk membuat model bawah permukaan bumi dengan mengandalkan data lapangan yang diukur bisa pada permukaan bumi atau di bawah permukaan bumi atau bisa juga di atas permukaan bumi dari ketinggian tertentu. Untuk mencapai tujuan ini, idealnya kegiatan survey atau pengukuran harus dilakukan secara terus menerus, berkelanjutan dan terintegrasi menggunakan sejumlah ragam metode geofisika. Seringkali –bahkan hampir pasti– terjadi beberapa kendala akan muncul dan tak bisa dihindari, seperti kehadiran noise pada data yang diukur. Ada juga kendala ketidaklengkapan data atau malah kurang alias tidak cukup. Namun demikian, dengan analisis data yang paling mungkin, kita berupaya memperoleh informasi yang relatif valid berdasarkan keterbatasan data yang kita miliki. Dalam melakukan analisis, sejumlah informasi mengenai kegiatan akuisisi data juga diperlukan, antara lain: berapakah nilai sampling rate yang optimal? Berapa jumlah data yang diperlukan? Berapa tingkat akurasi yang diinginkan? Selanjutnya –masih bagian dari proses analisis– model matematika yang cocok mesti ditentukan yang mana akan berperan ketika menghubungkan antara data lapangan dan distribusi parameter fisis yang hendak dicari. Setelah proses analisis dilalui, langkah berikutnya adalah membuat model bawah permukaan yang nantinya akan menjadi modal dasar interpretasi. Ujung dari rangkaian proses ini adalah penentuan lokasi pemboran untuk mengangkat sumber daya alam bahan tambang/mineral dan oil-gas ke permukaan. Kesalahan penentuan lokasi berdampak langsung pada kerugian meteril yang besar dan waktu yang terbuang percuma. Dari sini terlihat betapa pentingnya proses analisis apalagi bila segala keputusan diambil berdasarkan data eksperimen.

1.4 Macam-macam data geofisika Data geofisika bisa diperoleh dari pengukuran di lapangan atau bisa juga dari pengukuran di laboratorium. Gambar 1.3 memperlihatkan alur pengambilan data dari masing-masing pen-

4

BAB 1. PENDAHULUAN

Gambar 1.3: Alur eksperimen lapangan dan eksperimen laboratorium gukuran. Pada pengukuran lapangan, data geofisika yang terukur antara lain bisa berupa densitas, kecepatan gelombang seismik, modulus bulk, hambatan jenis batuan, permeabilitas batuan, suseptibilitas magnet dan lain sebagainya yang termasuk dalam besaran fisis sebagai karakteristik bawah permukaan bumi. Pada pengukuran di laboratorium, model lapisan bumi ataupun keberadaan anomali dalam skala kecil dapat dibuat dan diukur respon-nya sebagai data geofisika. Diharapkan hasil uji laboratorium tersebut bisa mewakili kondisi lapangan yang sesungguhnya yang dimensinya jauh lebih besar. Jika suatu pengukuran diulang berkali-kali, entah itu di lapangan maupun di laboratorium, seringkali kita temukan hasil pengukuran yang berubah-ubah, walaupun dengan variasi yang bisa ditolerir. Variasi ini umumnya disebabkan oleh kesalahan instrumen pengukuran (instrumental error ) atau bisa juga dikarenakan kesalahan manusia (human error ). Seluruh variasi ini bila di-plot kedalam histogram akan membentuk distribusi probabilistik.

1.5 Deskripsi proses geofisika: Model matematika Seluruh proses geofisika dapat dideskripsikan secara matematika. Sebagaimana yang telah disebutkan diawal, suatu formulasi yang bisa menjelaskan sistem geofisika disebut model. Namun perlu ditekankan juga bahwa istilah model memiliki ragam konotasi berbeda di kalangan geosaintis. Misalnya, orang geologi kerapkali menggunakan istilah model konseptual, atau istilah model fisik yang digunakan untuk menyebutkan hasil laboratorium, atau dalam catatan ini kita menggunakan istilah model matematika yang merupakan istilah umum dikalangan para ahli

1.6. DISKRITISASI DAN LINEARISASI

5

geofisika. Kebanyakan proses geofisika dapat dideskripsikan oleh persamaan integral berbentuk z

di =



K i (z) p(z)dz

(1.1)

0

dimana d i adalah respon atau data yang terukur, p(z) adalah suatu fungsi yang berkaitan dengan parameter fisis yang hendak dicari (misalnya: hambatan jenis, densitas, kecepatan, dan lain-lain) yang selanjutnya disebut parameter model, dan K i disebut data kernel. Data kernel menjelaskan hubungan antara data dan parameter model p(z). Parameter model (misalnya kecepatan, resistivitas dan densitas) bisa jadi merupakan fungsi yang kontinyu terhadap jarak atau posisi. Sebagai contoh, waktu tempuh t antara sumber gelombang seismik dengan penerimanya sepanjang lintasan L dalam medium, yang distribusi kecepatan gelombangnya kontinyu v(x, z), ditentukan oleh

t =

 L

1 dl v(x, z)

(1.2)

Deskripsi matematika terhadap sistem geofisika seperti contoh di atas disebut forward modelling . forward modelling digunakan untuk memprediksi data simulasi berdasarkan hipotesa

kondisi bawah permukaan. Data simulasi tersebut biasanya dinamakan data teoritik atau data sintetik atau data prediksi atau data kalkulasi. Cara seperti ini disebut pendekatan forward atau lebih dikenal sebagai pemodelan forward (Gambar 1.2).

1.6 Diskritisasi dan linearisasi Dalam banyak kasus, model bumi selalu fungsi kontinyu terhadap jarak dan kedalaman. Mari kita ambil kasus massa dan momen inersia bumi. Keduanya terkait dengan densitas bawah permukaan sesuai rumus-rumus berikut R

Massa = 4π

 

r 2 ρ(r)dr

(1.3)

(1.4)

0

Moment inersia =

8π 3

R

r 4 ρ(r)dr

0

dimana R adalah jejari bumi dan ρ(r) merupakan fungsi densitas terhadap jarak r . ρ(r) juga berhubungan dengan p(z) pada persamaan (1.1). Persamaan (1.4 dan 1.4) dapat dinyatakan dalam formulasi yang lebih umum yaitu R

di =



K i (r) p(r)dr

(1.5)

0

sama persis dengan persamaan (1.1). Integral ini relatif mudah dievaluasi secara komputasi dengan matematika diskrit. Pendekatan komputasi memungkinkan kita untuk menyederhanakan

ρ(r)dr menjadi m , sementara K i menjadi G i sehingga persamaan (1.5) dapat dinyatakan sebagai

di =



Gij m j

(1.6)

6

BAB 1. PENDAHULUAN

Gambar 1.4: Konfigurasi elektroda pada metode Schlumberger Ini adalah bentuk diskritisasi. Secara umum, memang pada kenyataannya ketika melakukan eksperimen di lapangan, data pengukuran maupun paremeter model selalu dibatasi pada inter val tertentu. Kita sering berasumsi bahwa bawah permukaan bumi terdiri dari lapisan-lapisan yang masing-masing memiliki sifat fisis atau parameter fisis p(z) yang seragam. Misalnya lapisan tertentu memiliki densitas sekian dan ketebalan sekian. Langkah praktis ini yang terkesan menyederhanakan obyek lapangan disebut langkah parameterisasi. Dalam kuliah ini, kita akan selalu memandang model yang diskrit dan juga parameter yang diskrit daripada model dan paremeter yang kontinyu. Sehingga proses inversi yang akan kita lakukan disebut sebagai teori inversi diskrit dan bukan teori inversi kontinyu. Dalam bentuk diskrit, persamaan (1.2) bisa dinyatakan sebagai p

ti =



j =1

Lij v j

(1.7)

Perlu dicatat disini bahwa waktu tempuh t tidak berbanding lurus dengan parameter model

v , melainkan berbanding terbalik. Hubungan ini dinamakan non-linear terhadap v . Namun demikian, jika kita mendefinisikan parameter model c = 1/v , dimana c adalah slowness gelombang seismik, maka problem ini bisa dinyatakan sebagai p

ti =



Lij c j

(1.8)

j =1

Hubungan ini disebut linear. Persamaan memenuhi bentuk d = Gm . Operasi transformasi seperti itu dinamakan linearisasi parameter. Dan proses menuju kesana dinamakan linearisasi. Sekarang mari kita lihat problem dari pengukuran resistivitas semu dengan metode Schlumberger untuk mengamati lapisan bawah permukaan yang diasumsikan terdiri dari dua lapisan. Formula model yang diturunkan oleh Parasnis, 1986 adalah



ρa (L) = ρ 1 1 + 2L2



∞

0

K (λ)J 1 (λL)λdλ



(1.9)

1.6. DISKRITISASI DAN LINEARISASI

7

dimana L = AB/2 adalah jarak masing-masing elektroda terhadap titik tengah, J 1 adalah fungsi Bessel orde 1 dan K (λ) adalah fungsi parameter (resistivitas masing-masing lapisan yaitu ρ1 dan ρ 2 serta ketebalan lapisan paling atas t ) dari sistem yang kita asumsikan. K (λ) dinyatakan sebagai (−2λt)

K (λ) =

−k1,2

(−2λt)

1 + −k1,2

dimana

k1,2 =

ρ1 − ρ2 ρ1 + ρ2

Kita bisa lihat bahwa persamaan (1.9) tidak bisa didekati dengan d = Gm sebagaimana yang dilakukan pada persamaan (1.2). Oleh karena itu persamaan resistivitas semu di atas disebut highly non-linear.

Bab 2

Formulasi Masalah Inversi

2.1 Klasifikasi masalah inversi Dalam masalah inversi, kita selalu berhubungan dengan parameter model ( M ) dan data (N ); yang mana jumlah dari masing-masing akan menentukan klasifikasi permasalahan inversi dan cara penyelesaiannya. Bila jumlah model parameter lebih sedikit dibandingkan data observasi (M < N ), maka permasalahan inversi ini disebut overdetermined. Umumnya masalah ini diselesaikan menggunakan pencocokan (best fit ) terhadap data observasi. Dalam kondisi yang lain dimana jumlah parameter yang ingin dicari (M ) lebih banyak dari pada jumlah datanya ( N ), maka masalah inversi ini disebut underdetermined. Dalam kasus ini terdapat sekian banyak model yang dapat sesuai kondisi datanya. Inilah yang disebut dengan masalah non-uniqness. Bagaimana cara untuk mendapatkan model yang paling mendekati kondisi bawah permukaan? Menurut Meju, 1994 persoalan ini bisa diselesaikan dengan model yang parameternya berbentuk fungsi kontinyu terhadap posisi. Kasus yang terakhir adalah ketika jumlah data sama atau hampir sama dengan jumlah parameter. Ini disebut evendetermined. Pada kasus ini model yang paling sederhana dapat diperoleh menggunakan metode inversi langsung. Pada bab ini, saya mencoba menyajikan dasar teknik inversi yang diaplikasikan pada model garis, model parabola dan model bidang. Uraian aplikasi tersebut diawali dari ketersediaan data observasi, lalu sejumlah parameter model ( unknown parameter) mesti dicari dengan teknik inversi. Mari kita mulai dari model garis.

2.2 Inversi Model Garis Secara teori, variasi temperatur bawah permukaan akan semakin meningkat ketika temperatur tersebut diukur semakin kedalam permukaan bumi. Misalnya telah dilakukan sebanyak sepuluh kali (N = 10) pengukuran temperatur (T i ) pada kedalaman yang berbeda beda ( zi ) sebagaimana ditunjukan datanya pada Tabel 2.6.

9

10

BAB 2. FORMULASI MASALAH INVERSI

Tabel 2.1: Data temperatur bawah permukaan tanah terhadap kedalaman Pengukuran ke-i Kedalaman (m) Temperatur ( O C )

z1 = 5 z2 = 16 z3 = 25 z4 = 40 z5 = 50 z6 = 60 z7 = 70 z8 = 80 z9 = 90 z10 = 100

1 2 3 4 5 6 7 8 9 10

T 1 = 35.4 T 2 = 50.1 T 3 = 77.3 T 4 = 92.3 T 5 = 137.6 T 6 = 147.0 T 7 = 180.8 T 8 = 182.7 T 9 = 188.5 T 10 = 223.2

Variasi Suhu vs Kedalaman 250

200

) s 150 u i c l e C ( u h u S100

50

0 0

20

40 60 Kedalaman (m)

80

100

Gambar 2.1: Data observasi perubahan suhu terhadap kedalaman dari permukaan tanah Source code untuk menggambar grafik tersebut dalam Matlab adalah 1

clc

2 clear 3 close 4 5

% Data observasi

6

z = [5 1 6 2 5 4 0 5 0 6 0 7 0 8 0 9 0 1 00];

7

T = [35.4 50.1 77.3 92.3 137.6 147.0 180.8 182.7 188.5 223.2];

8 9

% Plot data observasi

10 plot(z,T,’*r’); 11 grid; 12

xlabel(’Kedalaman (m)’);

13

ylabel(’Suhu (Celcius)’);

14

title(’\fontsize{14} Variasi Suhu vs Kedalaman’);

2.2. INVERSI MODEL GARIS

11

Lalu kita berasumsi bahwa variasi temperatur terhadap kedalaman ditentukan oleh rumus berikut ini:

m1 + m2 zi = T i

(2.1)

dimana m 1 dan m 2 adalah konstanta-konstanta yang akan dicari. Rumus di atas disebut model matematika. Sedangkan m1 dan m2 disebut parameter model atau biasa juga disebut unknown parameter. Pada model matematika di atas terdapat dua buah parameter model,

(M = 2). Sementara jumlah data observasi ada empat, (N = 10), yaitu nilai-nilai kedalaman, z i , dan temperatur, T i . Berdasarkan model tersebut, kita bisa menyatakan temperatur dan kedalaman masing-masing sebagai berikut:

m1 + m2 z1 = T 1 m1 + m2 z2 = T 2 m1 + m2 z3 = T 3 m1 + m2 z4 = T 4 m1 + m2 z5 = T 5 m1 + m2 z6 = T 6 m1 + m2 z7 = T 7 m1 + m2 z8 = T 8 m1 + m2 z9 = T 9 m1 + m2 z10 = T 10 Semua persamaan tersebut dapat dinyatakan dalam operasi matrik berikut ini:

       

1

z1

1

z2

1 1

z3 z4

1 1

z5 z6

1 1

z7 z8

1

z9

1 z10

                          T 1 T 2 T 3 T 4

m1 m2

=

T 5 T 6

(2.2)

T 7 T 8 T 9

T 10

Lalu ditulis secara singkat

Gm = d

(2.3)

dimana d adalah data yang dinyatakan dalam vektor kolom, m adalah model parameter, juga dinyatakan dalam vektor kolom, dan G disebut matrik kernel. Lantas bagaimana cara menda-

12


patkan nilai m 1 dan m 2 pada vektor kolom m? Manipulasi1 berikut ini bisa menjawabnya

GT Gm = G T d

(2.4)

dimana T disini maksudnya adalah tanda transpos matrik. Selanjutnya, untuk mendapatkan elemen-elemen m, diperlukan langkah-langkah perhitungan berikut ini:

GT Gm = GT d [GT G]

1

−

GT Gm = [GT G] m

= [GT G]

1

GT d

1

GT d

−

−

(2.5)

1. Tentukan transpos dari matrik kernel, yaitu GT

G =

      

1 1

z1

1

z2 z3

1 1

z4 z5

1 1

z6 z7

1 1

z8 z9

1 z10

      

⇒

GT =



1

1

1

1

1

1

1

1

1

1

z1 z2 z3 z4 z5 z6 z7 z8 z9 z10



2. Lakukan perkalian matriks G T G

GT G =



1 1 1 1 1 1 1 1 1 1 z1 z2 z3 z4 z5 z6 z7 z8 z9 z10

        

1

z1

1

z2

1 1

z3 z4

1 1

z5 z6

1 1

z7 z8

1

z9

1 z10

             =

N zi

zi zi2

dimana N = 10 atau sesuai dengan jumlah data observasi; sementara i = 1, 2, 3,..., 10.

1

Matrik G biasanya tidak berbentuk bujursangkar. Akibatnya tidak bisa dihitung nilai invers-nya. Dengan mengalikan matrik G dan transpose matrik G, maka akan diperoleh matrik bujursangkar

2.2. INVERSI MODEL GARIS

13

3. Kemudian tentukan pula G T d

GT d =



1 1 1 1 1 1 1 1 1 1 z1 z2 z3 z4 z5 z6 z7 z8 z9 z10

       

T 1 T 2 T 3 T 4 T 5 T 6 T 7 T 8 T 9 T 10

            =

T i zi T i

4. Dengan menggunakan hasil dari langkah 2 dan langkah 3, maka persamaan (2.4) dapat dinyatakan sebagai

GT Gm = GT d



         N zi

zi zi2

m1 m2

T i zi T i

=

(2.6)

Berdasarkan data observasi pada tabel di atas, diperoleh m

= [GT G]

1

−

GT d

m1 m2

m1 m2

=

=

N

zi

zi

1

−

           z2

T i

zi T i

i

10 536 536 38006

    

1

−

1314 88855

Operasi matriks di atas akan menghasilkan nilai m 1 = 24.9 dan m 2 = 2.0. Perintah di Matlab untuk menghitung elemen-elemen m, yaitu m=inv(G’*G)*G’*d

Secara lebih lengkap, source code Matlab untuk melakukan inversi data observasi adalah 1

clc

2 clear 3 close 4 5

% Data observasi

6

z = [5 1 6 2 5 4 0 5 0 6 0 7 0 8 0 9 0 1 00];

7

T = [35.4 50.1 77.3 92.3 137.6 147.0 180.8 182.7 188.5 223.2];

8

(2.7)

14 9



10 plot(z,T,’*r’); 11 grid; 12


13


14


15 16

% Membentuk matrik kernel G dan vektor d

17

n = length(z);

18

for k = 1 :n G(k,1) = 1;

19

G(k,2) = z(k);

20 21

end

22

d = T’;

23 24

% Perhitungan inversi dengan general least-squares

25

m = inv(G’*G)*G’*d;

26 27

% Plot hasil inversi (berupa garis least-squares)

28

hold on;

29

zz = 0:0.5:z(n);

30

TT = m(1) + m(2)*zz;

31 plot(zz,TT);


200

) s 150 u i c l e C ( u h u S100

50

0 0

20

40 60 Kedalaman (m)

80

100

Gambar 2.2: Hasil inversi atas data observasi perubahan suhu terhadap kedalaman

Demikianlah contoh aplikasi teknik inversi untuk menyelesaikan persoalan model garis. Anda bisa mengaplikasikan pada kasus lain, dengan syarat kasus yang anda tangani memiliki bentuk model yang sama dengan yang telah dikerjakan pada catatan ini, yaitu model garis:

y = m 1 + m2 x. Selanjutnya mari kita pelajari inversi model parabola.

2.3. INVERSI MODEL PARABOLA

15

2.3 Inversi Model Parabola Kembali kita ambil contoh variasi temperatur terhadap kedalaman dengan sedikit modifikasi data. Misalnya telah dilakukan sebanyak delapan kali (N = 8) pengukuran temperatur (T i ) pada kedalaman yang berbeda beda (zi ). Tabel pengukuran yang diperoleh adalah: Data observasi Tabel 2.2: Data temperatur bawah permukaan tanah terhadap kedalaman Pengukuran ke-i Kedalaman (m) Temperatur ( O C )

z1 = 5 z2 = 8 z3 = 14 z4 = 21 z5 = 30 z6 = 36 z7 = 45 z8 = 60

1 2 3 4 5 6 7 8

T 1 = T 2 = T 3 = T 4 = T 5 = T 6 = T 7 = T 8 =

20, 8 22, 6 25, 3 32, 7 41, 5 48, 2 63, 7 74, 6

tersebut selanjutnya di-plot ke dalam grafik variasi suhu terhadap kedalaman.


70

60

) s u i c l e C ( 50 u h u S

40

30

20 5

10

15

20

25 30 Kedalaman (m)

35

40

45

50

Gambar 2.3: Data observasi perubahan suhu terhadap kedalaman dari permukaan tanah Lalu kita berasumsi bahwa variasi temperatur terhadap kedalaman memenuhi model matematika berikut ini:

m1 + m2 zi + m3 zi2 = T i

(2.8)

dimana m1 , m2 dan m3 adalah unknown parameter. Jadi pada model di atas terdapat tiga buah model parameter, (M = 3) . Adapun yang berlaku sebagai data adalah nilai-nilai temperatur T 1 ,

T 2 ,..., dan T 8 . Berdasarkan model tersebut, kita bisa menyatakan temperatur dan kedalaman sebagai sistem persamaan simultan yang terdiri atas 8 persamaan (sesuai dengan jumlah data

16


observasi):

m1 + m2 z1 + m3 z12 = T 1 m1 + m2 z2 + m3 z22 = T 2 m1 + m2 z3 + m3 z32 = T 3 m1 + m2 z4 + m3 z42 = T 4 m1 + m2 z5 + m3 z52 = T 5 m1 + m2 z6 + m3 z62 = T 6 m1 + m2 z7 + m3 z72 = T 7 m1 + m2 z8 + m3 z82 = T 8 Semua persamaan tersebut dapat dinyatakan dalam operasi matrik berikut ini:

     

1 z1 z12 1 z2 z22 1 z3 z32 1 z4 z42 1 z5 z52 1 z6 z62 1 z7 z72 1 z8 z82

                      T 1 T 2 T 3

m1 m2 m3

=

T 4 T 5

(2.9)

T 6 T 7 T 8


Gm = d

(2.10)

dimana d adalah data yang dinyatakan dalam vektor kolom, m adalah model parameter, juga dinyatakan dalam vektor kolom, dan G disebut matrik kernel. Lantas bagaimana cara mendapatkan nilai m 1 , m 2 dan m 3 pada vektor kolom m? Manipulasi berikut ini bisa menjawabnya

GT Gm = GT d

(2.11)

dimana t disini maksudnya adalah tanda transpos matrik. Selanjutnya, untuk mendapatkan elemen-elemen m, diperlukan langkah-langkah perhitungan berikut ini:


17

1. Tentukan transpos dari matrik kernel, yaitu GT

G =

     

1 z1 1 z2 1 z3 1 z4 1 z5 1 z6 1 z7 1 z8

     

z12 z22 z32 z42 z52 z62 z72 z82

GT =

⇒

 

1

1

1

1

1

1

1

 

1

z1 z2 z3 z4 z5 z6 z7 z8 z12 z22 z32 z42 z52 z62 z72 z82

2. Tentukan G T G

GT G =

 

1

1

1

1

1

1

1

1


       

1 z1 1 z2 1 z3 1 z4 1 z5 1 z6 1 z7 1 z8

z12 z22 z32 z42 z52 z62 z72 z82

     

N

=

zi

zi2

zi2 zi3

zi3 zi4

         zi zi2

dimana N = 8 dan i = 1, 2, 3,..., 8.

3. Kemudian tentukan pula G T d

GT d =

 

1

1

1

1

1

1

1

1


       

T 1 T 2 T 3 T 4 T 5 T 6 T 7 T 8

     

=

     T i zi T i

zi2 T i

4. Sekarang persamaan (2.16) dapat dinyatakan sebagai

zi

zi2

m1

zi2 zi3

zi3 zi4

m2 m3

                  N

zi zi2

T i

=

zi T i zi2 T i

(2.12)

18

BAB 2. FORMULASI MASALAH INVERSI Berdasarkan data observasi pada tabel di atas, diperoleh

 

8

219

8547

219

8547

393423

8547 393423 19787859

       m1 m2

=

m3

349, 89 12894, 81 594915, 33

 

(2.13)

Matlab telah menyediakan sebuah baris perintah untuk menghitung elemen-elemen m, yaitu m=inv(G’*G)*G’*d

Sehingga operasi matriks di atas akan menghasilkan nilai m 1 = 21, m 2 = 0, 05 dan m 3 =

0, 02. Gabungan antara data observasi dan kurva hasil inversi diperlihatkan oleh Gambar 2.4 Script


70

60 ) s u i c 50 l e C ( u h 40 u S

30

20

10 0

10

20 30 Kedalaman (m)

40

50

Gambar 2.4: Hasil inversi atas data observasi perubahan suhu terhadap kedalaman. Tanda titik merah adalah data observasi sementara kurva biru adalah kurva hasil inversi Matlab yang lengkap untuk tujuan ini adalah 1

clc

2 clear 3 close 4 5

% Data observasi

6

z = [5 8 14 2 1 3 0 36 4 5 5 0];

7

T = [20.8 22.6 25.3 32.7 41.5 48.2 63.7 74.6];

8 9


10 plot(z,T,’*r’);


19

11 grid; 12


13


14


15 16


17

n = length(z);

18

for k = 1 :n

19

G(k,1) = 1;

20

G(k,2) = z(k);

21

G(k,3) = z(k).^2;

22

end

23

d = T’;

24 25


26


27 28


29

hold on;

30

zz = 0:0.5:z(n);

TT = m(1) + m(2)*zz + m(3)*zz.^2; 32 plot(zz,TT); 31

Demikianlah contoh aplikasi teknik inversi untuk menyelesaikan persoalan model parabola. Anda bisa mengaplikasikan pada kasus lain, dengan syarat kasus yang anda tangani memiliki bentuk model yang sama dengan yang telah dikerjakan pada catatan ini, yaitu model garis:

y = m1 + m 2 x + +m3 x2 . Selanjutnya mari kita pelajari inversi model bidang atau model 2dimensi (2-D).

20


2.4 Inversi Model Bidang Misalnya telah dilakukan sebanyak sepuluh kali (N = 10) pengukuran suatu nilai parameter pada suatu area sebagaimana ditunjukan datanya pada Tabel 2.3. Tabel 2.3: Distribusi suatu nilai dalam area tertentu Pengukuran ke-i X (m) Y (m) Nilai 1 2 3 4 5 6 7 8 9 10

2 5 7 4 1 3 6 9 8 4

3 6 2 7 8 9 4 1 5 5

10,6 23,5 27,3 20,8 11,1 18,9 25,4 33,5 33,2 24,1

Sebaran nilai terhadap X dan Y

35 30 25

i a l i N

20 15

10 10 8

5

Y (m)

0

2 0

4

10

6

X (m)

Gambar 2.5: Data observasi perubahan nilai terhadap posisi koordinat x dan y Source code untuk menggambar grafik tersebut dalam Matlab adalah 1

clc

2 clear 3 close 4 5

% Data observasi

6

x = [2 5 7 4 1 3 6 9 8 4];

7

y = [3 6 2 7 8 9 4 1 5 5];

8

nilai = [10.6 23.5 27.3 20.8 11.1 18.9 25.4 33.5 33.2 24.1];

2.4. INVERSI MODEL BIDANG

21

9 10


11 plot3(x,y,nilai,’*r’); 12 grid; 13

xlabel(’X (m)’);

14

ylabel(’Y (m)’);

15 zlabel(’Nilai’); 16

title(’\fontsize{14} Sebaran nilai terhadap X dan Y’);

Model matematika untuk 2-dimensi berikut ini digunakan untuk analisa data tersebut:

m1 + m2 xi + m3 yi = d i

(2.14)

dimana m1 , m 2 dan m3 merupakan unknown parameter yang akan dicari. Adapun yang berlaku sebagai data adalah d1 , d2 , d3 ,...,dN . Berdasarkan model matematika tersebut, kita bisa nyatakan

m1 + m2 x1 + m3 y1 = d1 m1 + m2 x2 + m3 y2 = d2 m1 + m2 x3 + m3 y3 = d3 .. .. .. .. .. . . . . . m1 + m2 xN + m3 yN = d N Semua persamaan tersebut dapat dinyatakan dalam operasi matrik berikut ini:

   

1 1

x1 x2

y1 y2

1 .. .

x3 .. .

y3 .. .

1 xN yN

                d1 d2

m1 m2 m3

=

d3 .. .

dN


Gm = d

(2.15)

dimana d adalah data yang dinyatakan dalam vektor kolom, m adalah unknown parameter, juga dinyatakan dalam vektor kolom, dan G disebut matrik kernel. Lantas bagaimana cara mendapatkan nilai m 1 , m 2 dan m 3 pada vektor kolom m? Sama seperti sebelumnya, kita harus membuat persamaan matriks berikut ini

GT Gm = GT d

(2.16)

dimana T disini maksudnya adalah tanda transpos matrik. Selanjutnya, untuk mendapatkan elemen-elemen m, diperlukan langkah-langkah perhitungan berikut ini:

22

BAB 2. FORMULASI MASALAH INVERSI 1. Tentukan transpos dari matrik kernel, yaitu GT

G =

   

1 1

x1 x2

y1 y2

1 .. .

x3 .. .

y3 .. .

1 xN yN

   

Gt =

⇒

 

1

1

1

···

1

x1 x2 x3 · · · xN y1 y2 y3 · · · yN

 

2. Tentukan G T G

T

G G =

 

1 1 1 ··· 1 x1 x2 x3 · · · xN y1

y2

y3 · · · yN

     

1

x1

y1

1 1 .. .

x2 x3 .. .

y2 y3 .. .

1 xN yN

   

=

         N xi yi

xi x2i

yi xi yi

xi yi

yi2

dimana N = jumlah data. dan i = 1, 2, 3,...,N . 3. Kemudian tentukan pula G T d

Gt d =

 

1

1

1

···

1

x1 x2 x3 · · · xN y1 y2 y3 · · · yN

     

d1 d2 d3 .. . dN

   

     di

=

xi di yi di

4. Sekarang, persamaan (2.16) dapat dinyatakan sebagai

                  N xi yi

xi x2i

yi xi yi

m1 m2

x i yi

yi2

m3

=

di xi di

(2.17)

yi di

5. Sampai disini, jika tersedia data observasi, maka anda tinggal memasukan data tersebut ke dalam persamaan di atas, sehingga nilai elemen-elemen m dapat dihitung dengan perintah matlab m=inv(G’*G)*G’*d

Langkah-langkah selanjutnya akan sama persis dengan catatan sebelumnya (model linear dan model parabola) Gabungan antara data observasi dan kurva hasil inversi diperlihatkan oleh Gambar 2.6

2.4. INVERSI MODEL BIDANG

23

Sebaran nilai terhadap X dan Y

40

30 i a l i 20 N

10

0 10 10 5

Y (m)

6 0

2 0

8

4 X (m)

Gambar 2.6: Hasil inversi berbentuk bidang. Saya sebut saja sebagai bidang inversi

Gambar 2.7: Ini adalah hasil inversi yang sama, hanya saja bidang inversi diputar (di-rotasi) pada suatu arah sehingga kita bisa saksikan dan pastikan bahwa sebaran titik data observasi berada di-dekat permukaan bidang inversi Script Matlab yang lengkap untuk tujuan ini adalah 1

clc

24


2 clear 3 close 4 5

% Data observasi

6

x = [2 5 7 4 1 3 6 9 8 4];

7

y = [3 6 2 7 8 9 4 1 5 5];

8

nilai = [10.6 23.5 27.3 20.8 11.1 18.9 25.4 33.5 33.2 24.1];

9 10


11 plot3(x,y,nilai,’*r’); 12 grid; 13

xlabel(’X (m)’);

14

ylabel(’Y (m)’);

15 zlabel(’Nilai’); 16

title(’\fontsize{14} Sebaran nilai terhadap X dan Y’);ndata = length(nilai);

17 18


19

ndata = length(x);

20

for k = 1:ndata

21

G(k,1) = 1;

22

G(k,2) = x(k); G(k,3) = y(k);

23 24

end

25

d = nilai’;

26 27


28


29 30


31

hold on;

32

[X,Y] = meshgrid(min(x):max(x),min(y):max(y));

Z = m (1) + X .*m(2) + Y.*m(3); 34 surf(X,Y,Z); 33

Anda bisa mengaplikasikan data pengukuran yang anda miliki, dengan syarat kasus yang anda tangani memiliki bentuk model yang sama dengan yang telah dikerjakan pada catatan ini, yaitu memiliki tiga buah model parameter yang tidak diketahui dalam bentuk persamaan bidang (atau 2-dimensi): d = m 1 + m2 x + m3 y .

2.5 Contoh aplikasi 2.5.1 Menghitung gravitasi di planet X Seorang astronot tiba di suatu planet yang tidak dikenal. Setibanya disana, ia segera mengeluarkan kamera otomatis, lalu melakukan ekperimen kinematika yaitu dengan melempar batu vertikal ke atas. Hasil foto-foto yang terekam dalam kamera otomatis adalah sebagai berikut Plot data pengukuran waktu vs ketinggian diperlihatkan pada Gambar 2.8. Anda diminta untuk membantu proses pengolahan data sehingga diperoleh nilai konstanta gravitasi di planet tersebut dan kecepatan awal batu. Jelas, ini adalah persoalan inversi, yaitu mencari unkown parameter (konstanta gravitasi dan kecepatan awal batu) dari data observasi (hasil foto gerak

sebuah batu). Langkah awal untuk memecahkan persoalan ini adalah dengan mengajukan asumsi model

2.5. CONTOH APLIKASI

25

Tabel 2.4: Data ketinggian terhadap waktu dari planet X Waktu (dt) Ketinggian (m) Waktu (dt) Ketinggian (m) 0,00 0,25 0,50 0,75 1,00 1,25 1,50 1,75 2,00 2,25 2,50

5,00 5,75 6,40 6,94 7,38 7,72 7,96 8,10 8,13 8,07 7,90

2,75 3,00 3,25 3,50 3,75 4,00 4,25 4,50 4,75 5,00

7,62 7,25 6,77 6,20 5,52 4,73 3,85 2,86 1,77 0,58

matematika, yang digali dari konsep-konsep fisika, yang kira-kira paling cocok dengan situasi pengambilan data observasi. Salah satu konsep dari fisika yang bisa diajukan adalah konsep tentang Gerak-Lurus-Berubah-Beraturan (GLBB), yang formulasinya seperti ini

1 ho + vo t − gt 2 = h 2 Berdasarkan tabel data observasi, ketinggian pada saat t = 0 adalah 5 m. Itu artinya h o = 5 m. Sehingga model matematika (formulasi GLBB) dapat dimodifikasi sedikit menjadi

1 vo t − gt 2 = h − ho 2

(2.18)

Selanjut, didefinisikan m 1 dan m 2 sebagai berikut

m1 = v o

1 m2 = − g 2

(2.19)

sehingga persamaan model GLBB menjadi

m1 ti + m2 t2i = h i − 5

(2.20)

dimana i menunjukkan data ke-i. Langkah berikutnya adalah menentukan nilai tiap-tiap elemen matrik kernel, yaitu dengan memasukan data observasi kedalam model matematika (persamaan (2.20))

m1 t1 + m2 t21 = h1 − 5 m1 t2 + m2 t22 = h2 − 5 m1 t3 + m2 t23 = h3 − 5 .. .. . = .. . . m1 t20 + m2 t220 = h20 − 5

26

BAB 2. FORMULASI MASALAH INVERSI 9

8

7

6 ) r e5 t e m ( i g g n4 i T

3

2

1

0 0

0.5

1

1.5

2

2.5 Waktu (detik)

3

3.5

4

4.5

5

Gambar 2.8: Grafik data pengukuran gerak batu

Semua persamaan tersebut dapat dinyatakan dalam operasi matrik berikut ini:

    

t1 t2

t21 t22

t3 t4 .. .

t23 t24 .. .

t19 t219 t20 t220

           m1 m2

=

h1 − 5 h2 − 5 h3 − 5 .. . h19 − 5 h20 − 5

    

Operasi matrik di atas memenuhi persamaan matrik

Gm = d Seperti yang sudah dipelajari pada bab ini, penyelesaian masalah inversi dimulai dari proses manipulasi persamaan matrik sehingga perkalian antara Gt dan G menghasilkan matriks bujursangkar

Gt Gm = G t d

Selanjutnya, untuk mendapatkan m 1 dan m 2 , prosedur inversi dilakukan satu-per-satu

(2.21)


27

1. Menentukan transpos matrik kernel, yaitu G t

G =

    

t1 t2 t3 t4 .. .

t21 t22 t23 t24 .. .

t19 t219 t20 t220

    

⇒

Gt =





t1 t2 t3 t4 . . . t19 t20 t21 t22 t23 t24 . . . t219 t220

2. Menentukan Gt G

t

G G =



     

t1 t2 t3 t4 . . . t19 t20 t21 t22 t23 t24 . . . t219 t220

t21 t22 t23 t24

t1 t2 t3 t4 .. .

.. .

t19 t219 t20 t220

           =

t2i

t3i

t3i

t4i

dimana N = 20 dan i = 1, 2,...,N . 3. Kemudian menentukan hasil perkalian Gt d

Gt d =



t1 t2 t3 t4 . . . t19 t20 t21 t22 t23 t24 . . . t219 t220

     

h1 h2 h3 h4 .. . h19 h20

        =

ti hi

t2i hi


           t2i

t3i

m1

t3i

t4i

m2

=

ti h i

t2i hi

(2.22)

Berdasarkan data observasi, diperoleh



179, 4 689, 1 689, 1 2822, 9

   m1 m2

=

273, 7 796, 3



Hasil operasi matriks ini dapat diselesaikan dengan satu baris statemen di matlab yaitu

28

BAB 2. FORMULASI MASALAH INVERSI 9

8

7

6 ) m (

5

n a i g g n i t 4 e K

3

2

1

0

0

0.5

1

1.5

2

2.5 Waktu (dt)

3

3.5

4

4.5

5

Gambar 2.9: Grafik hasil inversi parabola m=inv(G’*G)*G’*d

Hasil inversinya adalah nilai kecepatan awal yaitu saat batu dilempar ke atas adalah sebesar

m1 = v o = 3,2009 m/dt. Adapun percepatan gravitasi diperoleh dari m 2 dimana m 2 = − 12 g = -0,8169; maka disimpulkan nilai g adalah sebesar 1,6338 m/dt2 . Gambar 2.9 memperlihatkan kurva hasil inversi berserta sebaran titik data observasi. Garis berwarna biru merupakan garis kurva fitting hasil inversi parabola. Sedangkan bulatan berwarna merah adalah data pengukuran ketinggian (m) terhadap waktu (dt). Jelas terlihat bahwa garis kurva berwarna biru benar-benar cocok melewati semua titik data pengukuran. Ini menunjukkan tingkat akurasi yang sangat tinggi. Sehingga nilai kecepatan awal dan gravitasi hasil inversi cukup valid untuk menjelaskan gerak batu di planet X. 2.5.2 Analisa data seismik pada reflektor tunggal horizontal Suatu survei seismik dilakukan untuk mengetahui kedalaman sebuah reflektor mendatar sebagaimana tampak pada gambar Waktu tempuh gelombang (t), yang bergerak sesuai dengan lintasan warna merah, memenuhi model matematika berikut ini

4z 2 x 2 + 2 = t 2 2 v v

(2.23)

Data observasi yang berhasil dihimpun dari survei tersebut adalah Berdasarkan data tersebut, tuliskan langkah demi langkah proses inversi selengkap mungkin untuk menentukan:

• kecepatan gelombang seismik (v ) pada lapisan • kedalaman reflektor mendatar ( z ) terhadap permukaan ( surface)


29

S

R3

R2

R1

R4

R5

R6

R7

R8

Surface

z v Reflektor

Gambar 2.10: Reflektor mendatar pada kedalaman z. Kecepatan gelombang v dianggap konstan. S adalah sumber gelombang seismik dan R adalah penerima gelombang seismik. Jarak antara S dan R disebut offset (x). Sementara garis merah yang ada panahnya adalah lintasan gelombang seismik. Tabel 2.5: Data variasi offset ( x) dan travel time ( t) Receiver R ke- i Offset (x), meter Travel time (t), detik 1 2 3 4 5 6 7 8

60 80 100 120 140 160 180 200

0, 5147 0, 5151 0, 5155 0, 5161 0, 5167 0, 5175 0, 5183 0, 5192

Petunjuk: sederhanakan model matematika di atas menjadi

m1 + m2 x2 = t 2

(2.24)

Hasil yang benar adalah: kecepatan = 2797 m/dt ; kedalaman = 719 meter. Jika anda tidak mendapatkan kedua angka tersebut, maka proses inversi anda mungkin masih keliru!

2.5.3 Analisa data seismik pada reflektor tunggal miring Suatu survei seismik dilakukan untuk mengetahui keberadaan sebuah reflektor miring sebagaimana tampak pada gambar Waktu tempuh gelombang (t), yang bergerak sesuai dengan lintasan warna merah, memenuhi model matematika berikut ini

4z 2 4xz sin α x 2 t = 2 + + 2 v v2 v 2

(2.25)

Data observasi yang berhasil dihimpun dari survei tersebut adalah Berdasarkan data tersebut, tentukan:

• kecepatan gelombang seismik (v ) pada lapisan • kedalaman reflektor miring ( z ) terhadap permukaan ( surface) — jarak terdekat ke sumber gelombang seismik

• sudut kemiringan reflektor ( α)

30

BAB 2. FORMULASI MASALAH INVERSI x R2

R1

S

Surface

z

R3

R4

R5

v



R ef l ek to r

O

Gambar 2.11: Reflektor miring dengan sudut kemiringan sebesar α. Kecepatan gelombang v dianggap konstan. S adalah sumber gelombang seismik dan R adalah penerima gelombang seismik. Jarak antara S dan R disebut offset (x). Sementara garis merah yang ada panahnya adalah lintasan gelombang seismik. Tabel 2.6: Data variasi offset ( x) dan travel time ( t) Receiver R ke- i Offset (x), meter Travel time (t), detik 1 2 3 4 5

60 80 100 120 140

0, 4877 0, 4900 0, 4924 0, 4949 0, 4974

Petunjuk: sederhanakan model matematika di atas menjadi

m1 + m2 x + m3 x2 = t 2

(2.26)

2.6 Kesimpulan Dari sejumlah contoh pada bab ini, terlihat bahwa matrik kernel kerap kali berubah-ubah, sesuai dengan model matematika. Jadi, model matematika secara otomatis akan mempengaruhi bentuk rupa matrik kernelnya.

Bab 3

Penyelesaian Masalah Overdetermined

3.1 Regresi linear sederhana Jika suatu masalah inversi dapat direpresentasikan kedalam persamaan d = Gm , maka ia disebut linear. Kita dapat menjalankan prosedur yang sederhana untuk memperoleh nilai m dari data observasi. Dalam kenyataannya, tidak semua data observasi berhimpit dengan satu garis lurus. Jika kita mencoba melakukan fitting terhadap semua titik data observasi kepada satu garis, maka garis yang didapat disebut garis regresi. Misalnya, ada satu set data observasi yang ditulis sebagai (x1 , y1 ),(x2 , y2 ),...,(xn , yn ), garis regresi dinyatakan sebagai

y = a 0 + a1 x

(3.1)

dan setiap data memenuhi relasi berikut

yi = a 0 + a1 xi + ei

(3.2)

dimana e i disebut error, residual, atau sering juga disebut misfit atau kesalahan prediksi ( prediction error ). Garis regresi tidak akan berhimpit dengan setiap data observasi dan biasanya untuk kasus inversi seperti ini selalu overdetermined . Secara umum, tipe masalah inversi seperti ini diselesaikan dengan metode least squares . Dengan metode least squares , kita mencoba meminimalkan error , ei , dengan cara menentukan nilai a0 dan a1 sedemikian rupa sehingga diperoleh jumlah-kuadrat-error , (S ), yang minimal.

3.2

Metode least square

Diketahui data eksperimen tersaji pada Tabel 3.1 Lalu data tersebut di-plot dalam sumbu x dan y . Sekilas, kita bisa melihat bahwa data yang telah di-plot tersebut dapat didekati dengan sebuah persamaan garis, yaitu a1 xi + a 0 . Artinya, kita melakukan pendekatan secara linear, dimana fungsi pendekatan-nya adalah

P (xi ) = a 1 xi + a0 31

(3.3)

32

BAB 3. PENYELESAIAN MASALAH OVERDETERMINED

Tabel 3.1: Contoh data observasi yang dapat diolah oleh least squares xi y i xi yi 1 2 3 4 5

1,3 3,5 4,2 5,0 7,0

6 7 8 9 10

8,8 10,1 12,5 13,0 15,6

16 14 12 10

Y 8 6 4 2 0 1

2

3

4

5

6

7

8

9

10

X

Gambar 3.1: Hasil plotting data observasi dalam sumbu-x dan sumbu-y

Problemnya adalah berapakah nilai konstanta a1 dan a0 yang sedemikian rupa, sehingga posisi garis tersebut paling mendekati atau bahkan melalui titik-titik data yang telah di-plot di atas? Dengan kata lain, sebisa mungkin y i (pada persamaan (3.2)) sama dengan P (xi ) (pada persamaan (3.3)) atau dapat diformulasikan sebagai m



yi − P (xi ) = 0

(3.4)

yi − (a1 xi + a0 ) = 0

(3.5)

i=1

m

 i=1

dimana jumlah data, m = 10. Suku yang berada disebelah kiri dinamakan fungsi error (error function), yaitu m

E (a0 , a1 ) =



yi − (a1 xi + a0 )

(3.6)

i=1

Semua data yang diperoleh melalui eksperimen, fungsi error-nya tidak pernah bernilai nol. Jadi, tidak pernah didapatkan garis yang berhimpit dengan semua titik data ekperimen. Namun demikian, kita masih bisa berharap agar fungsi error menghasilkan suatu nilai, dimana nilai tersebut adalah nilai yang paling minimum atau paling mendekati nol. Harapan tersebut di wujudkan oleh metode least square dengan sedikit modifikasi pada fungsi error-nya sehingga

3.2. METODE LEAST SQUARE

33

menjadi m

E (a0 , a1 ) =



[yi − (a1 xi + a0 )]2

(3.7)

i=1

Agar fungsi error bisa mencapai nilai minimum, maka syarat yang harus dipenuhi adalah:

∂E (a0 , a1 ) =0 ∂a i

(3.8)

dimana i = 0 dan 1 , karena dalam kasus ini memang cuma ada a 0 dan a 1 . Maka mesti ada dua buah turunan yaitu: m

∂E (a0 , a1 ) ∂ = ∂a 0 ∂a 0

 

[yi − (a1 xi + a0 )]2 = 0

i=1

m

2

(yi − a1 xi − a0 )(−1) = 0

i=1

a0 .m + a1

m

m





xi =

i=1

yi

(3.9)

xi yi

(3.10)

i=1

dan

∂E (a0 , a1 ) ∂ = ∂a 1 ∂a 1

m

 

[yi − (a1 xi + a0 )]2 = 0

i=1

m

2

(yi − a1 xi − a0 )(−xi ) = 0

i=1

a0

m

m





xi + a1

i=1

m

2

xi

=

i=1

 i=1

Akhirnya persamaan (3.9) dan (3.10) dapat dicari solusinya berikut ini:

a0 = dan

m i=1

x2i

m m m i=1 yi − i=1 xi yi i=1 2 m m 2 i=1 xi − ( i=1 xi )

            m

a1 =

m i=1

m

m

xi yi − m 2 i=1 xi

m i=1

− (

xi

m i=1

m i=1 yi

xi )2

xi

(3.11)

(3.12)

Berdasarkan data ekperimen yang ditampilkan pada tabel diawal catatan ini, maka didapat:

a0 =

385(81) − 55(572, 4) = −0, 360 10(385) − (55)2

dan

a1 =

10(572, 4) − 55(81) = 1, 538 10(385) − (55)2

34


Jadi, fungsi pendekatan-nya, P (xi ), adalah

P (xi ) = 1, 538xi − 0, 360 Nilai a 0 dan a 1 disebut koefisien regresi. Lebih jauh lagi a 0 disebut intercept (titik perpotongan) terhadap sumbu y sedangkan a 1 adalah gradient atau slope (kemiringan garis). Gambar di bawah ini menampilkan solusi regresi linear tersebut berikut semua titik datanya 16

P(x) = 1.538*x − 0.36

14 12 10 8 6 4 2 0 −2 0

2

4

6

8

10

Gambar 3.2: Contoh solusi regresi linear Teknik diatas diterapkan secara rutin dalam analisis data geofisika, khususnya ketika kita mencoba meng-esktrak satu atau dua parameter model dari data observasi. Teknik ini disebut analisis regresi linear (linear regression analysis ) atau classical least squares fitting . Teknik ini pertama kali dipakai oleh Gauss pada tahun 1809. Teknik ini pada mulanya digunakan untuk mencari solusi dari masalah overdetermined namun pada perkembangannya teknik ini diterapkan juga pada underdetermined problem setelah dimodifikasi. Ketika kita ingin mendapatkan lebih dari 2 parameter maka teknik ini disebut analisis regresi multipel ( multiple regression analysis ).

3.3 Aplikasi regresi linear pada analisis data seismik refraksi Tabel 3.2 menampilkan data survei seismik refraksi yang dilakukan dengan jarak offset x i dan waktu tempuh gelombang dari source ke receiver dicatat dalam t i . Anggap saja persamaan travel time adalah

ti =

xi + T h v

(3.13)

Parameter xi dan ti berperan sebagai known parameter (parameter yang diketahui). Sementara kecepatan gelombang pada medium, v , dan waktu tempuh gelombang secara vertikal, T h bertindak sebagai unknowm parameter (parameter yang tidak diketahui). Metode regresi linear, atau yang juga akrab disebut least square, dilakukan untuk memecahkan unknown parameter tersebut, artinya metode regresi linear berupaya mencari nilai v dan T h . Pertama-tama proses

3.3. APLIKASI REGRESI LINEAR PADA ANALISIS DATA SEISMIK REFRAKSI

35

Tabel 3.2: Data seismik refraksi: waktu-datang gelombang, ti , dan jarak antara source dan receiver atau jarak offset, x i

xi (m) 2 4 6 8

Trace 1 2 3 4

ti (ms) 5,1 9,2 11,9 14,9

linearisasi dilakukan terhadap persamaan (3.13) sehingga menjadi

ti = a 0 + a1 xi

a0 = T h

dimana

dan

a1 =

1 v

(3.14)

Kesalahan (error ) diasumsikan hanya berasal dari cuplikan waktu gelombang datang. Penerapan metode regresi linear yang berusaha meminimalkan jumlah kuadrat dari error, ei =

ti − (a0 + a 1 xi ), dapat dinyatakan sesuai persamaan (3.11) dan (3.12). Dengan jumlah data, m = 4, maka 4 4 4 4 2 i=1 xi i=1 ti − i=1 xi ti i=1 xi a0 = = 2, 25 2 4 4 2 4 i=1 xi − i=1 xi dan

              

a1 =

4 i=1 xi ti −

4

4

4 2 i=1 xi

4 i=1 xi

4 i=1 ti 2

4 i=1 xi

−

= 1, 605

Dengan demikian kecepatan gelombang seismik pada medium tersebut 1 adalah v = 1/a1 =

623, 053 m/dt. Adapun standard error χ2a dan χ 2a ditentukan oleh rumus berikut 1

χ2a

1

2

χa

0

0

= m = χ

χ2 D

2

x2 D



(3.15)

(3.16)

dimana m

D = m

2

xi

−

i=1

χ2 =

2

m

    

1 m − 2

xi

i=1

m

E i2

i=1

Sebagai catatan tambahan, χ2 adalah nilai deviasi rms (root mean square ) dari data ti terhadap garis regresi hasil analisis ( a0 + a1 xi ) dengan faktor (m − 2) karena dalam masalah ini hanya dicari 2 parameter model (a0 dan a 1 ). 1

Nilai kecepatan yang diperoleh adalah nilai hasil analisis regresi linear. Namun itu bukan nilai yang mutlak, karena bisa jadi, pengukuran kecepatan gelombang menggunakan (misalnya) sonic-log memperoleh hasil yang berbeda, walaupun perbedaan tersebut tidak boleh terlalu jauh atau masih bisa diterima

36


3.4 Regresi linear dengan pendekatan matriks Metode regresi linear dapat didekati dengan operasi matriks. Caranya sama persis dengan yang telah dibahas pada Bab 2. Pendekatan matrik ini dimulai dari upaya membuat persamaan matrik d = Gm dari persamaan model (persamaan (3.13)). Unknown parameter , yang hendak dipecahkan, terkandung pada elemen-elemen vektor m . Jika data yang kita miliki sangat ideal dalam arti tidak ada error sama sekali, maka m bisa diperoleh sebagai berikut

m = G

1

−

d

Akan tetapi, pada kenyataannya semua data pengukuran pasti memiliki error yang besarnya relatif bervariasi. Karenanya, data observasi tak akan pernah fit secara sempurna dengan model. Itu artinya keberadaan misfit tidak pernah bisa dihindari. Konsekuensinya, misfit (baca: error) tersebut mesti disertakan pada d = Gm

d = Gm + ei

(3.17)

dan selanjutnya, solusi regresi linear diupayakan dengan cara meminimalkan jumlah kuadrat dari error, e i . Cara ini tak lain berupaya untuk memperoleh misfit terkecil yaitu jarak perbedaan terkecil antara data survei dan model. Dalam formulasi matematika, kuadrat error tersebut dinyatakan dengan

q = e T e = (d − Gm)T (d − Gm)

(3.18)

Dimana simbol T maksudnya adalah operasi transpose. Agar kuadrat error di atas menghasilkan nilai minimal, maka persamaan (3.18), diturunkan terhadap m dan hasilnya harus sama dengan nol

∂ dT d − dT Gm − mT GT d + mT GT Gm ∂q = =0 ∂m j ∂m j





sehingga

−dT G − GT d + GT Gm + mT GT G = 0 akhirnya diperoleh

2GT Gm = 2GT d GT Gm = GT d

(3.19)

Persamaan (3.19) disebut persamaan normal2 . Dengan persamaan normal, estimasi unknown parameter yang terkandung pada vektor m ditentukan oleh

m = GT G

 

1

−

GT d

(3.20)

Persamaan (3.20) disebut unconstrained least squares terhadap masalah inversi d = Gm. Bagian GT G 2

 

1

−

GT dinamakan Generalized Inverse yang mengolah data d untuk memperoleh

Sebetulnya kita sudah menuliskannya di Bab 2 sebagai persamaan (2.4), namun tanpa penurunan rumus

3.4. REGRESI LINEAR DENGAN PENDEKATAN MATRIKS

37

parameter model m. Untuk menyelesaikan persamaan (3.20) dengan operasi matriks secara numerik atau komputasi bisa menggunakan beberapa metode, diantaranya metode Eliminasi Gauss, LU-Decomposition, Iterasi Gauss-Seidel, dan Singular Value Decomposition. 3.4.0.1

Aplikasi pada data survei seismik refraksi

Kembali kita tampilkan Tabel 3.3 merupakan data survei seismik refraksi yang dilakukan dengan jarak offset x i , dan persamaan travel time adalah

ti =

xi + T h v

Tabel 3.3: Data seismik refraksi: waktu-datang gelombang (ti ) pada empat posisi geophone ( xi )

xi (m) 2 4 6 8

Trace 1 2 3 4

ti (ms) 5,1 9,2 11,9 14,9

Pertama-tama proses linearisasi dilakukan terhadap persamaan travel time sehingga menjadi

ti = a 0 + a1 xi

dimana

a0 = T h

dan

a1 =

1 v

Berdasarkan linearisasi tersebut, kita bisa menyusun sistem persamaan linear sebagai berikut:

t1 = a 0 + a1 x1 t2 = a 0 + a1 x2 t3 = a 0 + a1 x3 t4 = a 0 + a1 x4 Semua persamaan tersebut dapat dinyatakan dalam operasi matrik ( d = Gm ) berikut ini:

  

t1 t2 t3 t4

  

=

  

1 x1 1 x2 1 x3 1 x4

    a0 a1

dimana

d =

  

t1 t2 t3 t4

  

G =

  

1 x1 1 x2 1 x3 1 x4

  

dan

m =

  a0 a1

38


Selanjutnya, unknown parameter, a0 dan a1 , dipecahkan dengan menggunakan persamaan (3.20). Untuk menuju kesana, langkah-langkah penyelesaian berikutnya adalah 1. Menentukan transpos dari matrik G, yaitu G T

G =

  

1 x1 1 x2 1 x3 1 x4

  

⇒

GT =



1

1

1

1

x1 x2 x3 x4



2. Melakukan perkalian matriks G T G

GT G =



1 1 1 1 x1 x2 x3 x4

  

       

1 x1 1 x2 1 x3

N xi

=

1 x4

xi x2i

dimana N = 4 dan i = 1, 2, 3, 4. 3. Kemudian menentukan perkalian GT d

GT d =



1 1 1 1 x1 x2 x3 x4

  

      

t1 t2 t3

ti xi ti

=

t4


   a0 a1

=

1

−

       N xi

xi x2i

ti xi ti

(3.21)

Berdasarkan data survei pada Tabel 3.3 di atas, diperoleh

   a0 a1

=

4

20

20 120

1

−

 

41, 1 237, 6



(3.22)

Operasi matriks di atas akan menghasilkan nilai a0 = 2, 25 dan a1 = 1, 605. Dengan demikian v = 1/a1 = 623, 053 m/dt

3.5. SOAL

39

3.5 Soal Diketahui data temperatur borehole sebagai berikut. Tentukan slope dan intercept pada z-axis dan kemudian perkirakan temperatur pada kedalaman 390m. Depth, z(m)

Temp, t(o C)

30 70 180 250 300

25,0 26,2 29,7 34,3 35,5

Bab 4

Constrained Linear Least Squares Inversion

Geofisika merupakan ilmu yang bertujuan untuk memahami kondisi bawah permukaan bumi dengan menggunakan prinsip-prinsip ilmu fisika. Karena letaknya yang berada dibawah permukaan bumi, maka obyek geofisika tidak bisa diamati dengan mata telanjang. Sebagai gantinya, mau tidak mau obyek tersebut dipelajari dengan mengamati respon bawah permukaan bumi ketika suatu gangguan fisis diberikan padanya. Selanjutnya, respon tersebut dianalisis hingga diperoleh solusi yang nantinya siap menjadi acuan interpretasi. Masalahnya, pada kebanyakan pengamatan geofisika, sangat mungkin mendapatkan solusi yang berbeda-beda, yang semuanya itu bisa saja dipakai untuk menjelaskan data observasi. Inilah yang disebut dengan istilah non-uniqueness problem atau non-unik. Namun demikian pada akhirnya, kita harus memilih satu buah solusi yang paling baik menurut kita. Karena kita sadar bahwa obyek yang sedang kita incar dibawah sana, pastilah hanya berada dalam satu kondisi tertentu yang unik. Untuk mendapatkan solusi yang unik, kita harus menambahkan sejumlah informasi yang sebelumnya tidak ada pada persamaan least square d = Gm. Informasi tambahan ini disebut a priori information, yang selanjutnya akan digunakan untuk meng-constrain solusi sehingga diperoleh solusi yang dianggap paling tepat untuk menggambarkan kondisi bawah permukaan. A priori information , atau yang saya indonesiakan menjadi informasi awal, ini didapat dari data

geofisika yang lainnya, atau bisa juga dari data borehole, atau juga bersumber dari data geologi.

4.1 Inversi dengan informasi awal Kita dapat menambahkan informasi awal ( h) kepada parameter model dalam suatu proses in versi. Secara umum, informasi awal (h) tersebut diharapkan membantu pemodelan inversi sehingga diperoleh hasil yang unik dari sejumlah kemungkinan solusi. Sekali lagi, proses ini disebut meng-constrain, sementara solusi yang diperoleh disebut solusi ter-constrain. Constrain 41

42

BAB 4. CONSTRAINED LINEAR LEAST SQUARES INVERSION

terhadap suatu data dirumuskan sebagai berikut

Dm = h

(4.1)

Dimana D adalah matrik -dengan seluruh elemen selain diagonal bernilai nol- yang beroperasi pada parameter model m sedemikian rupa sehingga parameter m "dipaksa" untuk sama persis dengan informasi awal (h). Kalau kita menyelesaikan persamaan (4.1), berarti kita telah melakukan apa yang disebut dengan linear equality constraints. Formulasi matematika-nya adalah sebagai berikut

φ = (d − Gm)T (d − Gm) + β 2 (Dm − h)T (Dm − h)

(4.2)

Untuk mendapatkan error minimum maka turunan φ terhadap parameter model m dibuat sama dengan nol

2GT Gm − 2GT d + 2β 2 DT Dm − 2β 2 DT h = 0 sehingga diperoleh persamaan normal

(GT G + β 2 DT D)m = G T d + β 2 DT h

(4.3)

Ketika D adalah matrik identitas, maka persamaan normal berubah menjadi

(GT G + β 2 I )m = (GT d + β 2 h)

(4.4)

Dari sini solusi ter-constrain didapat sebagai berikut

m ˆc = (GT G + β 2 I )

1

−

(GT d + β 2 h)

(4.5)

Formula ini dinamakan inversi linear terkonstrain atau disebut juga the biased linear estimation technique. Keuntungannya adalah formula ini dapat membantu mendapatkan satu solusi

yang unik dari sejumlah solusi yang mungkin pada masalah overdetermined, dimana didalamnya terdapat ketidakpastian sebagai akibat dari kesalahan pengukuran ( observational errors ) dan ketidakpastian (uncertainties). Parameter β ditentukan secara coba-coba (trial and error), namun biasanya ia memiliki nilai antara nol dan satu. Konstanta β disebut faktor pengali undetermined atau disebut juga faktor pengali Lagrange. Sehingga metode ini disebut Lagrange multiplier method (metode pengali Lagrange).

4.1. INVERSI DENGAN INFORMASI AWAL

43

4.1.1 Memformulasikan persamaan terkonstrain Persamaan Dm = h secara umum memiliki bentuk

  

1 .. .

1 .. .

.. .

1

   

  

m1 m2 .. . m p

=

  

h1 h2 .. . h p

  

(4.6)

Namun demikian, persamaan matrik di atas dapat dimodifikasi sesuai kebutuhan. Misalnya jika informasi awal yang diketahui hanya ada satu, modifikasinya menjadi



1 0 ... 0

   

m1 m2 .. . m p

    =

hknown



(4.7)

Jika pada kasus lain, kita punya informasi awal yaitu parameter pertama dan parameter ke empat, maka persamaan matrik terkonstrain menjadi

  

1 0 0 1

   

m1 m2 m3 m4

  

=

  

h1 0 0 h4

  

(4.8)

4.1.2 Contoh aplikasi inversi terkonstrain Contoh 1: Least squares garis terkonstrain Sekarang kita ingin menerapkan inversi terkonstrain pada pengolahan data first arrivals dari seismik refraksi. Persamaan least square garis adalah

ti = m 1 + m2 xi

(4.9)

dengan parameter model yang terdiri atas m1 dan m2 . Sementara data lapangan merupakan pasangan dari jarak offset xi dan waktu datang gelombang ( first arrival time) ti . Jika ada 4 data lapangan, sistem persamaan linear yang bisa disusun adalah

t1 = m 1 + m2 x1 t2 = m 1 + m2 x2 t3 = m 1 + m2 x3 t4 = m 1 + m2 x4 Sekarang kita berasumsi memiliki informasi konstrain dari kegiatan explorasi sebelumnya bah wa solusi least square harus melewati titik koordinat (xc , tc ). Jadi kita harus meng-konstrain

44


solusi least square untuk mengakomodasi informasi awal tersebut.



1 xc

    m1 m2

=

tc

(4.10)

Persamaan matrik di atas harus diintegrasikan dengan d = Gm sehingga solusi akhir merupakan solusi terkonstrain yang kita harapkan lebih akurat dibandingkan jika tidak terkonstrain.

Tentu anda masih ingat pada least square tidak terkonstrain, dimana dikenal dua komponen berikut G T G dan G T d.

GT G =



    n xi

dan

GT d =

xi x2i

   ti xi ti

Agar menjadi terkonstrain kedua komponen itu mesti dimodifikasi sedemikian rupa sehingga menjadi

(GT G + β 2 I ) =

     n xi 1

xi 1 x2i xc

xc

0

 

(4.11)

dan

(GT d + β 2 h) =

     ti xi ti

(4.12)

tc

Sebelum dilanjut, mari kita tengok lagi persamaan normal (4.4) secara utuh, yaitu

(GT G + β 2 I )m = (GT d + β 2 h) Sekarang mari kita masukkan kedua komponen kedalam persamaan (4.4)

     n

xi 1

xi

1

x2i xc xc 0

           m1 m2 β

ti

=

xi ti tc

(4.13)

dimana nilai β mesti kita tentukan secara coba-coba namun dalam batas antara nol dan satu.

4.1. INVERSI DENGAN INFORMASI AWAL

45

Tabel 4.1: Data seismik refraksi: waktu-datang gelombang (ti ) pada empat posisi geophone ( xi ) Trace 1 2 3 4

xi (m) 2 4 6 8

ti (ms) 5,1 9,2 11,9 14,9

Akhirnya, solusi terkonstrain terhadap inversi garis yang melewati (xc , tc ) adalah

         m1

mˆc =

m2 β

n

=

xi

1

x2i xc xc 0

xi 1

1

−

 

     ti

(4.14)

xi ti tc

Sekarang, kita melangkah pada kasus nyata. Tabel 4.1 menunjukkan data pengukuran seismik refraksi. Tentukan parameter model untuk persamaan garis yang melewati titik ( xc = 8, yc =

14, 9)! Langkah pertama kita hitung komponen matrik

     n

(GT G + β 2 I ) =

xi 1

xi

1

x2i xc xc 0

    =

4

20

1

20 120 8 1 8 0

 

(4.15)

Kemudian kita hitung komponen matrik yang lain

       ti

(GT d + β 2 h) =

xi ti tc

=

41, 1 237, 6 14, 9

 

(4.16)

Problem ini dapat diselesaikan secara numerik dengan metode Eliminasi Gauss. Solusi yang diperoleh adalah m 1 = 2, 3857, m 2 = 1, 5643 dan β = 0, 2714. Berikut ini adalah script lengkapnya dalam Matlab clear all clc; a(1,1)=4; a(1,2)=20; a(1,3)=1; a(1,4)=41.1; a(2,1)=20; a(2,2)=120; a(2,3)=8; a(2,4)=237.6; a(3,1)=1; a(3,2)=8;

46


a(3,3)=0; a(3,4)=14.9; a n=3; % n = jumlah baris % berikut ini proses triangularisasi ----------for j=1:(n-1) kk=j+1; for k=kk:n m(k,j)=a(k,j)/a(j,j) for i=1:(n+1) a(k,i)=a(k,i)-m(k,j)*a(j,i); end end end % akhir dari proses triangularisasi ------------a x(n,1)=a(n,n+1)/a(n,n); % berikut ini proses substitusi mundur ---------for k=1:n-1 i=n-1-k+1; j=i+1; sum=0; for k=j:n sum = sum + a(i,k) * x(k,1); end x(i,1)=(a(i,n+1)-sum)/a(i,i); end % akhir dari proses substitusi mundur ------------x

4.2 Inversi dengan Smoothness Cara yang paling efektif untuk menginversi data yang terbatas adalah dengan menentukan konstrain sehingan solusi yang diinginkan menjadi smooth (halus). Tingkatan smooth dapat diukur berdasarkan kondisi fisis atau kondisi geologi.

4.2.1 Formulasi masalah Mari kita pelajari bagaimana suatu masalah dapat diformulasikan menuju solusi yang smooth. Jika diinginkan parameter model bervariasi dengan jarak selisih yang kecil, maka lakukanlah proses minimalisasi perbedaan paramter yang berdekatan ( m1 −m2 ), ( m2 −m3 ),..., (m p

1 −m p ).

−

4.2. INVERSI DENGAN SMOOTHNESS

47

Perbedaan ini dinyatakan sebagai persamaan konstrain Dm = h

  

1 −1 1

−1 1 −1

   

m1 m2 .. . m p

  

=

  

0 0 0 0

  

(4.17)

dimana D adalah operator selisih yang bertindak sebagai matrik smoothness dan Dm disebut penghalus (flatness) solusi vektor parameter model m . Jika parameter model tidak bervariasi secara smooth terhadap posisi, maka gunakanlah persamaan konstrain berbentuk

  

1 .. .

1 .. .

.. .

1

    

m1 m2 .. . m p

  

=

  

0 0 0 0

  

(4.18)

Dalam kasus ini, D adalah matrik identitas dengan dimensi p × p dan dimensi h adalah p × 1. Operasi ini akan mendorong proses inversi menuju kondisi stabil. Untuk mendapatkan solusi yang smooth, kita gunakan ukuran selisih seperti persamaan (4.2), dinyatakan sebagai

q 2 (m) = (Dm − h)T (Dm − h) = m T DT Dm = m T Hm

(4.19)

dimana H = D T D. Kita nyatakan mengenai masalah terkonstrain adalah: Dimulai dari data lapangan yang tidak komplit, tidak lengkap, tidak cukup, kita mencari seluruh kemungkinan solusi dengan residual

q 1 = |d − Gm|2 dan solusi yang paling smoothness dengan judgement dari ukuran q 2 (m). Secara matematik, pernyataan di atas memiliki maksud: meminimalkan q 2 = m T Hm dibawah kondisi |d − Gm|2 = q 1 atau secara umum |d − Gm|2 ≤ q T , dimana q T adalah nilai toleransi maksimum dari residual atau misfit. Masalah konstrain membutuhkan minimalisasi d − Gm2 dan q 2 (m) secara bersamaan,

φ = (d − Gm)T (d − Gm) + β 2 (mT DT Dm)

(4.20)

4.2.2 Solusi masalah Untuk mendapatkan solusi parameter model, perlu dilakukan minimalisasi terhadap persamaan (4.20),

∂ dT d − mT GT d − dT Gm + mT GT Gm + β 2 mT Hm =0 ∂m j





sehingga

GT G + β 2 H m = G T d





48


Ini adalah persamaan normal yang baru. Dan akhirnya solusi smoothness diturunkan sebagai berikut

ms = GT G + β 2 H Dan bila D = I ,

 

 

ms = GT G + β 2 I

1

GT d

(4.21)

1

GT d

(4.22)

−

−

Persamaan (4.22) lebih populer disebut Damped Least Squares solution atau solusi Least Square Teredam. Nama lainya yang juga cukup terkenal adalah inversi Marquardt.

Daftar Pustaka

[1] Meju, A Max., Geophysical Data Analysis: Understanding Inverse Problem Theory and Practice, (1994), Society of Exploration Geophysicists (SEG)

49

Analisis Data Geofisika

Recommend Documents