Monday 3 June 2013

Membuat Peta Hyposenter Gempa Dengan GMT

Kali ini saya akan berbagi cara membuat peta hyposenter dengan GMT. Sebenarnya cara berikut tidak terbatas pada data gempa saja, tetapi bisa dilakukan pada semua data yang memiliki data elevasi . Sebelumnya buatlah data gempa seperti berikut :
 PDE    1973  01 02 005320.30  -9.85  117.43  66  5.5 mbGS    ... .......     
 PDE    1973  01 02 022709.20   1.03  126.21  61  5.4 mbGS    ... .......     
 PDE    1973  01 07 221702.60   5.68  127.30  79  5.1 mbGS    ... .......     
 PDE    1973  01 09 061425.30   6.98  126.15  51  5.3 mbGS    ... .......     
 PDE    1973  01 13 111020.40  -2.70  101.27 105  5.2 mbGS    ... .......      
.............................................................
Simpan dengan nama "gempa.xyz".
Kolom yang kita pakai adalah kolom ke 7 dan 6 untuk koordinat, kolom ke 8 untuk data elevasi, dan kolom 9 untuk klasifikasi warna berdasarkan magnitude. Anda dapat merubah format data diatas, tetapi anda juga harus merubah script "gawk" agar sesuai dengan kolom yang anda buat.


Selanjutnya buatlah file color palette berikut :
5       0 0 255        5.5    0 0 255   
5.5    255 255 0    6       255 255 0
6       255 0 0        9       255 0 0    
Simpan dengan nama "file.cpt".
Arti dari file diatas adalah untuk magnitude 5 - 5.5 kita beri warna biru, magnitude 5.5 - 6 warna kuning, dan magnitude 6 - 9 warna merah. Format warna diatas dalam RGB dan anda dapat merubahnya sesuai keinginan.

Terakhir, buatlah script berikut :

set F=height.ps
set R=92/100/0/6
set azimuth=165
set elevation=25

makecpt -Cglobe -Z > elev.cpt
grdview indo.nc -Jx2.5c -Jy2.5c -R%R%/-170000/6000 -Celev.cpt -JZ8c -Qi  -E%azimuth%/%elevation% -K > %F%
gawk "{print $7, $6, $8*(-1), $9}" gempa.xyz | psxyz -R%R%/-170/6 -Jx -Jy -JZ -Su0.1c -Cmag.cpt -W -O -E%azimuth%/%elevation% -B1::/1::/30:Km: -U"Membuat Peta Hyposenter Dengan GMT"  >> %F%
ps2raster %F% -Gc:\programs\gs\gs8.53\bin\gswin32c -Tj -P -Fhypo

Simpan dengan ekstensi ".bat".
Anda dapat merubah azimuth dan sudut elevasi untuk merubah arah sudut pandang anda. Pada perintah kedua, -R saya beri nilai -170000 untuk batas kedalaman dalam meter dan 6000 untuk batas ketinggian dalam meter. Sedangkan pada perintah ketiga, -R saya beri nilai -170 untuk batas kedalaman gempa dalam kilometer dan 6 untuk batas ketinggian dalam kilometer. Anda dapat menyesuaikannya sesuai data yang anda miliki.

Azimuth 165 dan sudut elevasi 25 :

Azimuth 180 dan sudut elevasi 90 :

Azimuth 147 dan sudut elevasi 0.1 (sudut elevasi tidak boleh sama dengan 0)



gawk "{print $7, $6, $8*(-1), $9}" gempa.xyz artinya adalah bujur terletak di kolom ke 7, lintang pada kolom ke 6, kedalaman pada kolom ke 8 (karena data saya bernilai positif sedangkan batas peta yang saya gunakan bernilai negatif, maka nilai kolom ke 8 saya kali dengan -1 agar menjadi negatif), dan magnitude pada kolom ke 9.

Sekian untuk kali ini. Selamat mencoba dan semoga bermanfaat. Salam Orang Indonesia...

5 comments:

  1. Bang, untuk file.nc itu selain Indonesia, misalnya wilayah Jepang nyari di mana ya?
    thanks.

    ReplyDelete
  2. scriptnya gak jalan bang, mungkin bisa dicheck ulang, dengan petunjuk yang sama, gak mampu jalanin..

    ReplyDelete
  3. untuk file tersebut bisa didownload di https://www.bodc.ac.uk/data/online_delivery/gebco/
    atau sumber lain juga bisa, asalkan format datanya netCDF

    untuk scriptnya mungkin ada keterangan errornya?

    ReplyDelete
  4. E:\>grdview indo.nc -Jx2.5c -Jy2.5c -R92/100/0/6/-170000/6000 -Celev.cpt
    -JZ8c -Qi -E165/25 -K 1>height.ps
    grdview: Warning: ╠╦‼ not a valid number and may not be decoded properly.
    grdview: GMT ERROR: Could not decode 2.5c, return NaN.
    grdview: GMT ERROR: Could not decode , return NaN.
    grdview: GMT SYNTAX ERROR -J option. Correct syntax:
    -Jy// OR -JY//
    or inch/degree, or use in inch
    E:\>gawk "{print $7, $6, $8*(-1), $9}" gempa.xyz | psxyz -R92/100/0/6/
    -170/6 -Jx -Jy -JZ -Su0.1c -Cmag.cpt -W -O -E165/25 -B1::/1::/30:Km: -U"Membuat P
    eta Hyposenter Dengan GMT" 1>>height.ps
    psxyz: Warning: ÿ█‼ not a valid number and may not be decoded properly.
    psxyz: GMT ERROR: Could not decode 2.5c, return NaN.
    psxyz: GMT ERROR: Could not decode , return NaN.
    psxyz: GMT SYNTAX ERROR -J option. Correct syntax:
    -Jy// OR -JY//
    or inch/degree, or use in inch
    gitu bang hasilnya..

    ReplyDelete
  5. bang height.ps nya drimana ya?

    ReplyDelete