Advanced Seismic Unix Processing Part-4

Waveform Modelling

Seismic Unix memiliki fungsi sufdmod2 yang berguna untuk melakukan Waveform Modelling.

sufdmod2 memerlukan model kecepatan dalam domain kedalaman dengan format "C-style binary floating point numbers".

Berikut ini adalah script untuk melakukan Waveform Modelling untuk model kecepatan Marmousi.

Oleh karena sufdmod2 memerlukan memory yang cukup excessive, maka langkah pertama yang saya lakukan adalah mereduksi jumlah sampel (suresamp) velocity.su menjadi 25% (rf=0.25). Anda dapat menggunakan surange untuk mengevaluasi jumlah sampel sebelum dan setelah suresamp:

suresamp < velocity.su > velocity1.su rf=0.25

Selanjutnya melakukan suwind dengan j=4 untuk mereduksi jumlah trace menjadi 25%:
suwind < velocity1.su key=tracl j=4 > velocity2.su
mv velocity2.su velocity.su

Mengubah format dari su ke binary:
sustrip < velocity.su > velocity.h@ ftn=0

Konversi dari time ke depth domain:
suttoz < velocity.su > velocitydepth.su vfile=velocity.h@

Membatasi depth dengan maksimum depth=2400m:
suwind < velocitydepth.su > velocitydepth1.su tmax=2400
mv velocitydepth1.su velocitydepth.su

Mengubah format dari su ke binary untuk model kecepatan dalam domain kedalaman:
sustrip < velocitydepth.su > velocitydepth.h@ ftn=0

Ketika anda melakukan sustrip anda akan memperoleh informasi parameter berikut:
n1=348 n2=575 nt=348 ntr=575 ns=348


Gunakan text editor untuk membuat kode Waveform Modelling:

#! /bin/sh
##gunakan surange velocitydepth.su untuk memperoleh informasi d1 dan d2
n1=348
n2=575
d1=6.912000
nt=348
ntr=575
d2=9.216
ns=348
velfile=velocitydepth.h@

ximage < $velfile wbox=400 hbox=200 n1=$n1 n2=$n2 d1=$d1 d2=$d2 \
wbox=800 hbox=400 legend=1 title="Marmousi model" label1="Depth (m)" cmap=hsv2 \
label2="distance m " units="m/s" &

# input sufdmod2
xs=2000 zs=10 hsz=10 vsx=1000 verbose=2
vsfile="vsfile.su" ssfile="ssfile.su" hsfile="hsfile.su"
tmax=2.5 mt=18
label1="Depth m"
label2="Distance m"

sufdmod2 < $velfile nz=$n1 dz=$d1 nx=$n2 dx=$d2 verbose=1 \
xs=$xs zs=$zs hsz=$hsz vsx=$vsx dt=0.002 hsfile=$hsfile fmax=70 mono=0 \
vsfile=$vsfile ssfile=$ssfile verbose=$verbose \
tmax=$tmax abs=1,1,1,1 mt=$mt |
suxmovie clip=1.0 fframe=70 dframe=10 title=" %g ms" \
label1=$label1 label2=$label2 \
n1=$n1 d1=$d1 f1=$f1 n2=$n2 d2=$d2 f2=$f2 \
loop=1 sleep=20000 -geometry 800x400 &
exit 0

Gambar di bawah ini menunjukkan model Marmousi yang telah di-downsampling dengan parameter di atas:

Gambar di bawah ini menunjukkan snapshot untuk wavefield pada Model Marmousi dengan posisi source pada X=1000m dan depth=10m.


Gambar di bawah ini menunjukkan animasi wavefield pada Model Marmousi:

Advanced Seismic Unix Processing Part-3

PreStack Time Migration

Pada bagian ini saya akan menunjukkan bagaimana caranya melakukan PreStack Time Migration (PreSTM) untuk data Marmousi dengan Seismic Unix. Seismic Unix menyediakan fungsi suktmig2d untuk melakukan PreSTM dengan metoda Kirchhoff. Data input untuk PreSTM adalah RMS Velocity atau Average velocity atau Stacking Velocity serta data Seismic pada domain Common Offset.

Marilah kita gunakan data interval velocity dan data seismic yang sebelumnya saya gunakan di sini.

Oleh karena data kecepatan yang kita miliki adalah interval velocity, maka langkah pertama adalah melakukan konversi dari interval velocity menjadi Average Velocity.

Gunakan perintah berikut untuk mengkonversi format su menjadi ascii:
suascii < velocity.su bare=1 > velocity.ascii

Selanjutnya buatlah kode octave dengan gedit dan berilah nama rmsvelgen.m:

clear; clc
nocmp=574; %%no cmp veocity harus sama dengan no cmp seismic
nosamp=751; %%no samples
nocmporig=2301; %%no cmp origial pada model kecepatan (kita harus melakuka desimasi)
data1=load('velocity.ascii');
data2=reshape(data1,nosamp,[]);
ind=[4:4:nocmporig-4]; %% desimasi cmp
data3=data2(:,ind);

for i=1:nocmp
trace{i}=data3(:,i);
end

for j=1:nocmp
trace{j}=data(:,j);
trace{j}(1)=trace{j}(1);
for i=1:nosamp-1
trace{j}(i+1)=trace{j}(i)+trace{j}(i+1);
end
end

for j=1:nocmp
for i=1:nosamp-1
traceb{j}(i)=trace{j}(i)/i;
end
end

vel2=cell2mat(traceb');
vel3=vel2';
save -ascii velrms.txt vel3

Run kode octave di atas dengan mengetik rmsvelgen pada terminal octave.

Pada terminal linux, ketiklah perintah berikut untuk mengkonversi dari ascii ke binary:
a2b < velrms.txt > velrms.bin

Convert dari binary ke su, flip dan tampilkan lalu QC:
suaddhead ns=574 < velrms.bin | sushw key=dt a=4000 > velrms.su
suflip flip=0 < velrms.su > velrms1.su
mv velrms1.su velrms.su
suximage < velrms.su legend=1 &

Ubah menjadi binary kembali, lalu smoothing, batasi hingga t=2.75s (menghindari noise) dan QC dengan menampilkannya:

sustrip < velrms.su > velrms.h@ ftn=0
smooth2 < velrms.h@ n1=750 n2=574 r1=20 r2=20 > smoothvelrms
suaddhead ns=750 < smoothvelrms | sushw key=dt a=4000 > smoothvelrms.su
suwind < smoothvelrms.su > smoothvelrms1.su tmax=2.75
mv smoothvelrms1.su smoothvelrms.su
suximage < smoothvelrms.su legend=1 &

Setelah window ximage aktif, tekan huruf 'h' pada keyboard untuk mengubah skala warna.

Ubah kembali ke binary:
sustrip < smoothvelrms.su > smoothvelrms.h@ ftn=0

Untuk data seismik, lakukan sorting untuk membuat domain Common Offset dan membatasinya hingga t=2.75s

susort offset gx < data.su > data_co.su
suwind < data_co.su > data_co1.su tmax=2.75
mv data_co1.su data_co.su

QC Offset -800
suwind < data_co.su key=offset min=-800 max=-800 | suximage perc=90 title="offset -800m"


PreSTM dilakukan offset by offset, sehingga buatlah kode Bourne Shell untuk menjalankan PreSTM. Lakukan surange untuk mempelajari parameter-parameternya.

#!/bin/sh
##set -x
offsetmin=-2575
offsetmax=-200
offdist=25
filein=data_co.su
vfile=smoothvelrms.h@
#################################################
n=$offsetmin
while [ $n -le $offsetmax ]
do
suwind < $filein key=offset min=$n max=$n | suktmig2d vfile=$vfile dx=25 > data_pstm_co$n.su
n=$((n+$offdist))
done

Setelah kode Bourne Shell tersebut dijalankan, maka anda akan memperoleh hasil PreSTM untuk masing-masing offset.

Buatlah Kode Bourne Shell Berikut untuk melakukan penggabungan file Common Offset:

#!/bin/sh
##set -x
offsetmin=-2575
offsetmax=-200
offdist=25
fileout=data_pstm_co.su
echo "#!/bin/sh" >> gabung
echo "##set -x" >> gabung
echo "cat \\" >> gabung
n=$offsetmin
while [ $n -le $((offsetmax)) ]
do
echo "data_pstm_co$n.su \\" >> gabung
n=$((n+$offdist))
done
echo "| sushw key=tracl a=1 > $fileout" >> gabung

Setelah selesai jalankan sh gabung, sehingga diperoleh file data_pstm_co.su

Lakukan sorting data sehingga diperoleh cdp gather:
susort cdp offset < data_pstm_co.su > data_pstm_cdp.su

Lakukan QC sehingga diperoleh gambar berikut:
suwind < data_pstm_cdp.su key=cdp min=100 max=103 | suxwigb perc=90 &

Lakukan Mute:
sumute < data_pstm_cdp.su > data_pstm_cdp_mute.su xmute=-2575,-200 tmute=1.1,0.1 key=offset

Lakukan QC:
suwind < data_pstm_cdp_mute.su key=cdp min=100 max=103 | suxwigb perc=90 &

Lakukan Stacking untuk memperoleh hasil Stack PreSTM:
sustack < data_pstm_cdp_mute.su | suximage perc=90 &

Analisa kecepatan biasanya dilakukan pada interval CMP tertentu serta sampling waktu yang tidak beraturan. Untuk mengatasi hal ini, saya membuat kode Octave untuk melakukan interpolasi CMP serta regularisasi time sampling, sehingga hasilnya dapat digunakan untuk melakukan PreSTM.

clear; clc
nocmp=22; %%% katakanlah pada model Marmousi kita hanya melakukan analisa kecepatan untuk 22 cmp
cmpint=20; %%% dengan interval cmp=20
ts=0.004; %%%sampling waktu
tmin=0; %%%waktu t0
tmax=3-ts; %%%waktu maksimum pada data Marmousi awal
cmp1=100; %%%CMP pertama dari analisa kecepatan
cmp2=520; %%%CMP terakhir dari analisa kecepatan
cmp1x=1; %%%CMP awal pada seismik
cmp2x=574; %%%CMP akhir pada seismik
ns=750; %%%jumlah sample pada data seismik

%File dari hasil analisa kecepatan atau Velocity Picker,
%hapus 2 baris pertama, tnmo=, vnmo=, dan \ (spasi \).
%dengan demikian filenya bisa dibaca oleh octave

fid=fopen('vpick.txt','rt');
for i=1:(nocmp*2)
line=fgets(fid);
tv{i}=str2num(line);
end

for i=1:nocmp
time{i}=[tmin,tv{i*2-1},tmax];
end

for i=1:nocmp
a=size(tv{i*2})(1,1); b=size(tv{i*2})(1,2);
vel{i}=[tv{i*2}(a),tv{i*2},tv{i*2}(b)];
end

%interpolasi dan regularisasi
tnew=[tmin:ts:tmax];

for i=1:nocmp
vnew{i}=interp1(time{i},vel{i},tnew);
end

vall=cell2mat(vnew');
vel1=vall';

velearly=repmat(vel1(:,cmp1x),1,(cmp1-cmp1x));
vellast=repmat(vel1(:,nocmp),1,(cmp2x-cmp2));

for i=1:ns
valhornew{i}=interp1([cmp1:cmpint:cmp2]', vall(:,i),[cmp1:cmp2]');
end

velmid=cell2mat(valhornew);
vel2=velmid';

velocity=[velearly,vel2,vellast];
save -ascii velocity.txt velocity

Lakukan proses di atas untuk melakukan smoothing dan konversi ke binary.
Gambar berikut adalah RMS Velocity dengan interpolasi dari hasil analisa kecepatan.

Interactive Velocity Analysis SU

Analisa kecepatan lebih nyaman dilakukan secara interaktif. Jika anda bekerja dengan Seismic Unix, anda dapat menggunakan Velocity Picker 1.1 yang bisa didownload secara gratis.

Berikut adalah cara untuk mendownload, menginstall dan mempergunakan Velocity Picker 1.1:

• Download Velocity Picker 1.1 di sini , sehingga anda memperoleh file yang bernama: VelocityPicker V1.1.zip
• Oleh karena file tersebut mengandung spasi, rename file tersebut menjadi VelocityPicker (contoh)
• Pada terminal linux, ketik unzip VelocityPicker . Jika belum memiliki fungsi unzip, ketik sudo apt-get install unzip.
• Setelah diunzip, semua file akan tersimpan dalam directory yang bernama dist
• Masuklah ke dalam directory dist, lalu ketik java -jar "VelocityPicker_V2.jar"
  
• Sebuah windows Velocity Picker akan muncul
• File > Open Trace File (CMP Gathers dengan format SU)
• Nomor CMP yang akan dianalisa, dapat dipilih pada scroll box
• Pick pada puncak semblance dengan MB1
  
• Setelah selesai picking, Save Velocity Field (SU Format)
• File yang telah disave tersebut, selanjutnya anda copy dan paste di Bourne Shell nmo.sh
• Untuk mempelajari tentang fitur aplikasi ini, lihatlah Help > Instructions

Berikut adalah snapshot dari aplikasi Velocity Picker 1.1. yang saya gunakan untuk analisa CMP gathers Marmousi:

Advanced Seismic Unix Processing Part-2

PreStack Depth Migration

Seismic Unix (SU) memiliki kapabilitas untuk melakukan PreStack Depth Migration (PreSDM). Fungsi dalam SU untuk melakukan PreSDM adalah sukdmig2d, dimana metoda yang digunakannya adalah metoda Kirchhoff.

Untuk mendemonstrasikan bagaimana caranya melakukan PreSDM, saya akan menggunakan shot gathers serta model kecepatan Marmousi.

• Data Marmousi bisa didownload di Gzipped tar file (46 MB)
• Ekstrak data yang telah didownload dengan perintah tar -xvvf ieee.tar.gz
• Semua file yang telah diekstrak akan berada di folder ieee, file-file tersebut adalah: data.segy (seismic shot gathers), velocity.segy dan density.segy
• Lakukan konversi dari format segy ke su dengan perintah:
• segyread tape=data.segy conv=0 endian=0 > data.su
• segyread tape=velocity.segy conv=0 endian=0 > velocity.su

• Coba tampilkan model kecepatan dengan perintah suximage < velocity.su maka anda akan memperoleh model kecepatan yang terotasi 90 derajat.
• Lakukan rotasi -90 derajat untuk model kecepatan tersebut dengan perintah:
• suflip flip=0 < velocity.su > velocity1.su
• mv velocity1.su velocity.su
  
• Sekarang coba tampilkan kembali, lalu tekan huruf 'r' untuk mengubah skala warna. sehingga anda akan memperoleh gambar sbb (Model Marmousi dibuat berdasarkan tipikal model geologi di West Africa-Angola):

Sekarang lihatlah data seismik, lakukan sorting dari shot gather ke CMP gather dengan teknik yang saya jelaskan pada bagian ini, lakukan analisa kecepatan, koreksi NMO dan stack. Sehingga anda memperoleh gambar sbb:

• Gambar di atas merupakan stack dari Model Marmousi dengan menggunakan stacking velocity. Terlihat dominasi difraksi akibat zona yang kompleks pada Model Marmousi, demikian juga dengan pencitraan yang kurang baik pada bagian yang lebih dalam. Pada tahapan berikutnya kita akan berupaya untuk memperbaiki citra tersebut dengan menggunakan PreSDM.
• PreSDM memiliki tiga tahapan utama: smoothing (smooth2) model kecepatan, raytracing (rayt2d) dan depth migration (sukdmig2d). Dengan data input: shot gathers (domain time) dengan format su, velocity dengan format "C-style binary floating point numbers". Data output adalah: Stack dalam domain depth.
• Ubahlah format velocity dari su ke "C-style binary floating point numbers" dengan perintah: sustrip < velocity.su > velocity.h@ ftn=0
• Mengingat data seismik memiliki noise di atas 2.9s, maka kita harus menghilangkannya dengan perintah:
  

• suwind < data.su > data1.su tmax=2.9
• mv data1.su data.su

• Dengan gedit, copy-lah kode Bourne Shell berikut lalu save dengan nama presdm.sh:

#!/bin/sh
# Parameter Model Kecepatan
nz=751
nx=2301
dz=4
dx=4
# Parameter Seismik
nt=750
ntr=96
dt=4000
ifile=data.su #data input seismik shot gathers
ofile=datamig.su # data output stack depth migrated
tfile=tfile # traveltimes
vfile=velocity.h@ # input kecepatan dengan format C floating point

# smoothing model kecepatan
smooth2 < $vfile n1=$nz n2=$nx r1=20 r2=20 >smoothvel

#raytracing untuk meghasilkan tfile
rayt2d <smoothvel dt=0.004 nt=751 dz=$dz nz=$nz dx=$dx nx=$nx \
fxo=0 dxo=25 nxo=369 fxs=0 dxs=100 nxs=93 >$tfile

#migrasi
sukdmig2d infile=$ifile datain=$ifile outfile=$ofile dataout=$ofile \
ttfile=$tfile fzt=0 dzt=4 nzt=751 fxt=0 dxt=25 nxt=369 \
fs=0 ns=93 ds=100 nzo=751 dzo=4 dxm=25 mtr=1

• Pada terminal linux ketiklah sh presdm.sh
• Setelah selesai (Intel Atom, 1GB MEM: kurang lebih 60 menit), tampilkan hasilnya dengan suximage < datamig.su, sehingga diperoleh hasil sbb:
  

• Bandingkan hasilnya dengan stack di atas, menakjubkan?

Advanced Seismic Unix Processing Part-1

Radon Demultiple

Pada bagian ini, saya akan menunjukkan bagaimana caranya melakukan eliminasi multiple melalui Transformasi Radon.

Eliminasi multiple dengan teknik Radon atau dikenal juga dengan Radon Demultiple memiliki 4 tahapan utama:

1. Koreksi NMO (CMP gather)
2. Transformasi ke domain τ-slope index
   
3. Desain Mute
4. Inverse τ-slope index
   

Setelah koreksi NMO, karakteristik multiple akan berbeda dengan reflektor. Jika nilai kecepatan NMO yang dipilih tepat maka reflektor akan menjadi datar (flat) sedangkan multiple masih melengkung ke arah bawah (karena nilai kecepatannya lebih rendah). Jika ditransformasi ke dalam domain τ-slope index (baca: sumbu y adalah waktu dan x adalah kemiringan event), maka reflektor yang datar akan memiliki slope index bernilai sekitar nol dan multiple akan bernilai positif. Sedangkan untuk event yang melengkung ke atas (over NMO-ed, difraksi, side swipe) akan memiliki nilai slope index negatif.

Gambar di atas (kiri) menunjukkan CMP sintetik dengan reflektor pada waktu 0.2s dan reflektor yang bersamaan dengan multiple pada 0.4s, sedangkan gambar sebelah kanan adalah CMP dalam domain τ-slope index. Perhatikan reflektor yang flat pada waktu 0.2s dan 0.4s memiliki nilai slope index di sekitar nol, sedangkan multiple pada 0.4s memiliki nilai slope index sekitar 100.

Berikut adalah tahapan untuk menghasilkan gambar di atas:

1. Run kode Bourne Shell berikut:
   

#!/bin/sh
npl=3 # number of planes
nt=200 # number of time samples
ntr=100 # number of traces
taper=0 # no end-of-plane taper = 1 taper planes to zero at the end
offset=400 # offset
dt=0.004 # time sample interval in seconds
suplane npl=$npl nt=$nt ntr=$ntr taper=$taper offset=$offset \
dip1=0 len1=128 ct1=100 cx1=50 \
dip2=1 len2=128 ct2=115 cx2=50 \
dip3=0 len3=128 ct3=50 cx3=50 > planemul.su

2. Tampilkan:

sushw < planemul.su key=offset a=0 b=10 | sufilter f=0,5,50,60 > planemul_fil.su
suximage < planemul_fil.su title="Sebelum Radon Demultiple" key=offset label1="time (s)" label2="offset (meters)" wbox=500 hbox=500 &

3. Lakukan Tranformasi Radon dan tampilkan:

suradon < planemul_fil.su igopt=3 interp=4 choose=0 depthref=1000 interoff=0 offref=1190 pmin=-300 pmax=300 > radon.su
suximage < radon.su perc=99 label1="tau (s)" label2="slope index" title="suplane data Radon transformed" wbox=500 hbox=500 &

Gambar di atas adalah CMP dalam domain τ-slope index, sebelum (kiri) dan setelah mute.

Untuk melakukan Mute:

sumute < radon.su key=f2 \
xmute=-300,30,31,300 \
tmute=0.8,0.79,0.01,0 \
mode=1 ntaper=10 > radon_muted.su

suximage < radon_muted.su label1="tau (s)" label2="slope index" title="suplane data Radon transformed-muted" wbox=500 hbox=500 perc=99 &

Berikut adalah perintah untuk melakukan inversi τ-slope index setelah dimute serta menampilkan hasilnya:

suradon < radon_muted.su igopt=3 interp=4 choose=4 depthref=1000 interoff=0 offref=1190 pmin=-300 pmax=300 | suximage title="Setelah Radon Demultiple" key=offset label1="time (s)" label2="offset (meters)" wbox=500 hbox=500 &

Dari gambar di atas terlihat bahwa, multiple pada 0.4s telah dieliminasi dengan sempurna.

Gambar di bawah ini menunjukkan muting untuk data real pada domain τ-slope index . Kiri (sebelum mute) kanan (setelah mute).

Pemilihan batas pengambilan mute, harus dilakukan dengan hati-hati, jangan sampai melukai sinyal reflektor, mengingat berbagai hal seperti error pemilihan nilai kecepatan NMO, anisotropic, dll. Jadi, sebaiknya kita melakukan test dengan mengembalikannya ke domain waktu-offset, lalu lihatlah perilaku reflektornya. Gambar di bawah ini menunjukkan perbandingan CMP gather sebelum (kiri) dan setelah Radon Demultiple(kanan).

Seismic Processing with Seismic Unix - Part9

Tahap sebelumnya

Tahap 18:

Pada tahapan ini, saya akan menunjukkan bagaimana melakukan PoSTM (Post Stack Time Migration) dengan menggunakan Seismic Unix.

Seismic Unix menyediakan fungsi migrasi dengan menggunakan beberapa metodologi diantaranya Stolt Migration, Gazdag atau Phase-Shift Migration, Claerbout's Migration, dll. Teknik-teknik migrasi tersebut memiliki kelebihan dan kekurangannya masing-masing.

Saya akan menerapkan PoSTM pada data telah distack dengan menggunakan velocity analysis dan juga telah mengalami koreksi statik. Saya memilih Metoda Stolt, karena metoda ini sangat cepat dan cukup robust.

Sebelumnya, saya akan menerapkan mute dan tapering (walaupun pada sustolt sendiri diterapkan tapering) sehingga diperoleh ujung lintasan kiri, kanan dan atas yang lebih gradual. Hal ini penting dilakukan karena perbedaan amplitudo yang tiba-tiba akan menghasilkan migration artifacts.

Setelah itu migrasi diterapkan dengan menggunakan kecepatan RMS (yang saya ambil dari CMP ke 1000) dari hasil velocity analysis sebelumnya.

Berikut ini adalah kode Bourne Shell untuk melakukan mute, menampilkan stack sebelum migrasi, melakukan Stolt Migration dan sekaligus menampilkan hasilnya:

#!/bin/sh

sumute < Line001_stack_vel2_elev_rstat.su key=tracl xmute=1,150,1132,1132 tmute=3.0,0.2,0.2,3.0 ntaper=50 > Line001_stack_vel2_elev_rstat_mute.su

suximage < Line001_stack_vel2_elev_rstat_mute.su key=cdp perc=90 title='Before Stolt Migration' &

time=0.0187891,0.494781,0.914405,1.37787,1.94781,2.90605 vels=1992.35,2211.92,2488.77,2765.61,2975.64,3319.31

sustolt < Line001_stack_vel2_elev_rstat_mute.su \
cdpmin=1 cdpmax=1282 dxcdp=50 \
tmig=$time vmig=$vels \
smig=0.6 vscale=1 lstaper=50 lbtaper=50 | suximage title="After Stolt Migration" key=cdp min=100 max=1100 perc=90 verbose=0 &

exit

Copy-lah code di atas dengan menggunakan text editor (vi, pico, gedit, atau nedit) lalu save dengan nama tertentu katakanlah stoltmig. Lalu ketik sh stoltmig

Berikut adalah stack sebelum migrasi (untuk mengubah skala warna, tekan hurup r pada ximage yang aktif, untuk mengeksplorasi lebih jauh ketik ximage pada terminal linux).

Setelah Migrasi:

Animasi (klik untuk memperbesar)