Thursday, August 15, 2019

Tutorial VELEST


Tutorial VELEST 3.1 (Windows)

By Kris H.P David (krisdavid98@gmail.com)

1. Persiapan Data

A. File *.cnv (local earthquake data)
Isi dalam file ini adalah data katalog yang di buat dari script matlab GADtoVelest.m

Yang perlu disiapkan untuk di run oleh script matlab adalah:
·      File traveltime.txt berisikan data perkolom yaitu:

1. Tahun
2. Bulan
3. Tanggal
4. Jam
5. Menit
6. Stasiun (Maks 4 huruf, kemudian dijadikan kode berupa angka)
7. Traveltime Gel.P (Tp-T0)
8, Traveltime Gel.S (Ts-T0)
9. Nomor Urut Event

Gambar 1. Contoh file traveltime.txt

Catatan Penting : Yang dimasukkan ke dalam file traveltime.txt adalah angka saja, bukan dengan keterangan tahun, bulan, tanggal, dll. 
Event gempa yang digunakan pada Windows max 650 Event.

·      File lokasi.txt berisikan data perkolom yaitu:

1. Tahun
2. Bulan
3. Tanggal
4. Jam
5. Menit
6. Detik (Origin time (T0 dari GAD))
7. Latitude
8. Longitude
9. Depth (Km)
10. Nomor Urut Event 

Gambar 2. Contoh file lokasi.txt
Ikuti Format sesuai diatas agar tidak error, jangan lupa untuk menyamakan angka di depan dan di belakang koma ya. Gunakan Ms. Excel agar sesuai lalu masukkan ke notepad++ dalam format lokasi.txt, Kemudian di run matlabnya, cek hasil *.cnv. Buka hasil di dalam notepad++, tambahkan 3x enter di akhir. Pastikan spasi rapi dan sesuai karena software VELEST 3.1 sensitif.

Catatan Penting : Yang dimasukkan ke dalam file lokasi.txt adalah angka saja, bukan dengan keterangan tahun, bulan, tanggal, dll

Script Matlab: 

clc
%--- initiate data
%fid1=fopen(file_catalog);%file catalog menjadi file variabel
str=fopen('lokasi.txt','r');%baca file lokasi
%--- read data
lokasi=textscan(str,'%d %d %d %d %d %6.3f %7f %*5f %9f %*3f %f %d');
%ftell(str)
fclose(str);
tahun=lokasi{1};
bulan=lokasi{2};
hari=lokasi{3};
jam=lokasi{4};
menit=lokasi{5};
detik=lokasi{6};
lat=lokasi{7};
long=lokasi{8};
depth=lokasi{9};
ID=lokasi{10};
South='S';
East='E';
E=0.00;
Y=0;
k=1;%bilangan untuk membantu menentukan jumlah stasiun dalam file katalog yg akan kita buat (format velest)
h=0;%bilangan untuk menerangkan batas data maksimum travel time dalam 1 baris file .cnv adalah 6
fidP=fopen('datavelest2.cnv','wt');
%Scrip untuk mengetahui jumlah event
all_arr_time_100_kamojang=importdata('traveltime.txt');
event=0;
for p=1:length(all_arr_time_100_kamojang) %tergantung nama file yang mengandung catalog gempa
    if isnan(all_arr_time_100_kamojang(p,9))==0
        event=event+1;
    end;
end;

z = 1;
for i=1:event% i sampai jumlah event
    fprintf(fidP,'%2i%2i%2i %2i%2i%6.2f%8.4f%1c %8.4f%1c %6.2f   %4.2f      %1i\n',tahun(i),bulan(i),hari(i),jam(i),menit(i),detik(i),-lat(i),South,long(i),East,depth(i),E,Y);   
    if all_arr_time_100_kamojang(z,7)~=99.990 && all_arr_time_100_kamojang(z,7)> detik(i)
        fprintf(fidP,'%4i%1c%1i%6.2f',all_arr_time_100_kamojang(z,6),'P',1,all_arr_time_100_kamojang(z,7)-detik(i));
        h=h+1;
    else if all_arr_time_100_kamojang(z,7)~=99.990 && all_arr_time_100_kamojang(z,7)< detik(i)
            fprintf(fidP,'%4i%1c%1i%6.2f',all_arr_time_100_kamojang(z,6),'P',1,(all_arr_time_100_kamojang(z,7)+60)-detik(i));
            h=h+1;
        end;
    end;
    if all_arr_time_100_kamojang(z,8)~=99.990 && all_arr_time_100_kamojang(z,8)> detik(i)
        fprintf(fidP,'%4i%1c%1i%6.2f',all_arr_time_100_kamojang(z,6),'S',1,all_arr_time_100_kamojang(z,8)-detik(i));
        h=h+1;
    else if all_arr_time_100_kamojang(z,8)~=99.990 && all_arr_time_100_kamojang(z,8)< detik(i)
            fprintf(fidP,'%4i%1c%1i%6.2f',all_arr_time_100_kamojang(z,6),'S',1,(all_arr_time_100_kamojang(z,8)+60)-detik(i));
            h=h+1;
        end;
    end;
    z=z+1;
   
    while isnan(all_arr_time_100_kamojang(z,9))==1
        z;
        %k=k+1
        if all_arr_time_100_kamojang(z,7)~=99.990 && all_arr_time_100_kamojang(z,7)> detik(i)
        fprintf(fidP,'%4i%1c%1i%6.2f',all_arr_time_100_kamojang(z,6),'P',1,all_arr_time_100_kamojang(z,7)-detik(i));
        h=h+1
        else if all_arr_time_100_kamojang(z,7)~=99.990 && all_arr_time_100_kamojang(z,7)< detik(i)
                fprintf(fidP,'%4i%1c%1i%6.2f',all_arr_time_100_kamojang(z,6),'P',1,(all_arr_time_100_kamojang(z,7)+60)-detik(i));
                h=h+1;
            end;
        end;
        if mod(h,6)==0 && h~=0
            fprintf(fidP,'\n')
            h=0;
        end;
        if all_arr_time_100_kamojang(z,8)~=99.990 && all_arr_time_100_kamojang(z,8)> detik(i)
        fprintf(fidP,'%4i%1c%1i%6.2f',all_arr_time_100_kamojang(z,6),'S',1,all_arr_time_100_kamojang(z,8)-detik(i));
        h=h+1;
        else if all_arr_time_100_kamojang(z,8)~=99.990 && all_arr_time_100_kamojang(z,8)< detik(i)
                fprintf(fidP,'%4i%1c%1i%6.2f',all_arr_time_100_kamojang(z,6),'S',1,(all_arr_time_100_kamojang(z,8)+60)-detik(i));
                h=h+1
            end;
        end;
        if mod(h,6)==0 && h~=0
            fprintf(fidP,'\n')
            h=0;
        end;
       
        z=z+1;
    end;
    if h~=0
    fprintf(fidP,'\n')
    end;
    fprintf(fidP,'\n')
    h=0;
end;


Salin script di atas ke dalam notepad++, kemudian save dalam GADtoVelest.m
Contoh hasil *.cnv sebagai berikut:

Gambar 3.  Contoh hasil *.cnv

B. File stasiun.sta (stasiun list)
Berisikan file stasiun dengan pembobotan. Angka pembobotan ditentukan oleh jumlah event gempa yang di catat oleh stasiun. Tentukan jumlah event terbanyak yang direkam oleh stasiun, kemudian di bagi tiga, maka skala yang digunakan adalah 1 untuk bobot dengan perekaman sedikit, 2 untuk bobot perekaman sedang dan 3 untuk perekaman banyak. Perekaman tersebut dilihat dari file arrivaltime.txt

Contoh : stasiun KARA merekam 120 event, maka stasiun yang merekam event (x<=30) memiliki bobot 1, event (30<x<=60) memiliki bobot 2, event (60<x<=120) memiliki bobot 3.

Gambar 4. file Ms.Excel stasiun.sta

Yang dimasukkan adalah data seperti di bawah ini di Notepad++ dengan format *.sta 
Gambar 5. Contoh file stasiun.sta

Ditambahkan enter 1x setelah data. Pakai stasiun yang mencatat event gempa saja, jika tidak mencatat maka tidak perlu dimasukkan.

C. File velocity.mod (initial velocity model)
Berisi file kecepatan. Dapat di replace dari file yang saya sediakan. Dari kolom kiri ke kanan velocity, kedalaman (km), damping(ini coba2, kalo bingung kasil 01.0 aja.  Yang baris 2 dan 10 (di file contoh) itu ada tulisan 7 merupakan jumlah lapisan. Lalu pada kedalaman di velest itu minus jika di atas permukaan.
pada penelitian saya sebelumnya, model kecepatan yang saya gunakan pada gambar 6 di bawah adalah model kecepatan regional Jawa Tengah (Koulakov dkk., 2009).
Ini contohnya:

Gambar 6.  Contoh file velocity.mod

Jangan lupa di tambahkan enter 2x di akhir. CR LF itu merupakan tanda untuk Windows.
Kalau Mau lebih Jelas dapat di lihat di VELEST User’s Guide (e-book).

D. File velest.cmn (parameter control)

Berisikan parameter control dari velest. Formatnya sama seperti di VELEST User’s Guide (e-book).

Berikut contoh Penggunaan dari penelitian saya:

******* CONTROL-FILE FOR PROGRAM  V E L E S T  (28-SEPT1993) *******

***
*** ( all lines starting with  *  are ignored! )
*** ( where no filename is specified, 
***   leave the line BLANK. Do NOT delete!)
***
*** next line contains a title (printed on output):

CALAVERAS area7 1.10.93 EK startmodell vers. 1.1  

***      starting model 1.1 based on Castillo and Ellsworth 1993, JGR

***  olat       olon           icoordsystem      zshift   itrial   ztrial    ised
      -7.9948  -110.4586            0                      0       0        0.00       0
***
*** neqs   nshot   rotate
         137      0      0.0
***
*** isingle   iresolcalc
            0              0
***
*** dmax    itopo    zmin     veladj    zadj   lowveloclay
        46.00     0          0.05       0.01      0.02         1
***
*** nsp    swtfac   vpvs       nmod
          2      0.50     1.730          2
***
*** othet   xythet    zthet    vthet   stathet
          0.01    0.20      0.20      1.0        0.05
***
*** nsinv   nshcor   nshfix     iuseelev    iusestacorr
           1            0            0               1                   0
***
*** iturbo    icnvout   istaout   ismpout
             1             1              1             0
***
*** irayout   idrvout   ialeout   idspout   irflout   irfrout   iresout
              1             1              1             1              1           1            1
***
*** delmin   ittmax   invertratio
         0.010      2                  2
***
*** Modelfile:
velocity.mod
***
*** Stationfile:
station.sta
***
*** Seismofile:                                                                             
***
*** File with region names:
***
*** File with region coordinates:
***
*** File #1 with topo data:                                                                               
***
*** File #2 with topo data:                                                                            
***
*** DATA INPUT files:
***
*** File with Earthquake data:
datavelest2.cnv
***
*** File with Shot data:                                                                         
***
*** OUTPUT files:
***
*** Main print output file:
Output.OUT
***
*** File with single event locations:
***
*** File with final hypocenters in *.cnv format:
FinalHypo.CNV
***
*** File with new station corrections:
stationcorr.OUT
***
*** File with summary cards (e.g. for plotting):
***
*** File with raypoints:
***
*** File with derivatives:
***
*** File with ALEs:
***
*** File with Dirichlet spreads:
***
*** File with reflection points:
***
*** File with refraction points:
***
*** File with residuals:
VELOUT.RES
***
******* END OF THE CONTROL-FILE FOR PROGRAM  V E L E S T  *******


2.   Setelah semua proses ke-1 telah selesai, pastikan semua data telah tersedia, yaitu file cnv, sta, mod, dan cmn, kemudian running velest.exe. 


3.   Berdasarkan penelitian yang telah saya lakukan, koreksi/cek nilai-nilai berikut :
  • RMS, jika nilai RMS konvergen (menurun) dalam setiap iterasi maka hasilnya baik.
  • Final Hypocenter, lakukan plotting hasil relokasi tersebut (Latitude, Longitude) ke dalam Ms.Excel, apabila kurang baik (Hasil relokasi jauh berbeda dengan relokasi awal), maka coba-coba kembali pengolahan data kalian (terutama pada file *.cmn kalian.
  • Update Velocity, lihat hasil perubahan dari velocity apabila memiliki nilai yang jauh berbeda, maka cek kembali file mod dan velocity. 
  • GAP, berdasarkan penelitian saya nilai GAP yang baik adalah <180°
4. Berdasarkan pengalaman saya, ada beberapa masalah error yang sering terjadi:
  • File input format masih salah, Koreksi file *.cnv, kalo ga rapi atau ada yang kelebihan kekurangan spasi dia udah ga jalan), inget diakhir ditambah enter 3x
  • Harus ada beberapa enter diakhir file
  • Huruf besar atau kecil
  • Penamaan stasiun (pada penelitian saya mengubah menjadi angka)
  • Parameter kontrol. Kalo di output.out semua input udah ada disitu, menurut saya yang berpotensi error adalah parameter kontrol ini (*.cmn). biasanya beda penanganannya, jadi sering coba ubah-ubah. Intinya trial and error. 
5.  Selamat mencoba, Semoga tidak lama errornya

Refrensi 

Kissling, E. 1995. Program VELEST User’s GuideShort Introduction. Institute of Geophysics and Swiss Seismological Service, ETHZuerich, Switzerland.

Koulakov, I., Jakovlev, A., dan Luehr, B. (2009): Anisotropic structure beneath central Java from local earthquake tomography, Geochemistry Geophysics Geosystems, 10, 56–87.

CERITA AKUISISI DATA GEOLISTRIK

  Pengalaman untuk melakukan penelitian bersama dosen adalah hal yang sangat menyenangkan, banyak hal yang bisa di petik. Ini adalah cerit...