Skip to content
Snippets Groups Projects
PlaneWave.scr 4.69 KiB
#!/bin/bash

export OMP_NUM_THREADS=4

# The commands used to generate the wavelet and model are:
#
# Model: (Warning: Overwrites tmp_cp.su & tmp_ro.su)
# makemod sizex=6000 sizez=2000 dx=$dx dz=$dx cp0=1900  ro0=1200 \
# 	orig=-3000,0 file_base=tmp.su verbose=2 \
# 	intt=def x=-3000,500,3000 z=195,195,195 poly=1 cp=1950 ro=3700 \
# 	intt=def x=-3000,3000 z=600,600 poly=0 cp=2050 ro=1750 \
# 	intt=def x=-3000,3000 z=680,680 poly=0 cp=2150 ro=2220 \
# 	intt=def x=-3000,3000 z=780,780 poly=0 cp=2230 ro=1700 \
# 	intt=def x=-3000,-2200,-1500,0,1300,2100,3000 z=520,580,680,840,680,600,500 poly=2 cp=2400 ro=2800 \
# 	intt=def x=-3000,0,3000 z=1110,1110,1110 poly=0 cp=2350 ro=1950 \
# 	intt=def x=-3000,3000 z=1180,1180 poly=0 cp=2480 ro=1820 \
# 	intt=def x=-3000,0,3000 z=1290,1290,1370 poly=0 cp=2490 ro=2000 \
# 	intt=def x=-3000,3000 z=1380,1380 poly=0 cp=2520 ro=2050 \
# 	intt=def x=-3000,3000 z=1480,1480 poly=0 cp=2550 ro=1850;
# rm -f tmp_ro.su; #We don't need the density model.
# sushw <tmp_cp.su key=tracl,tracf,trid,timbas,trwf a=1,1,0,0,0 b=1,1,0,0,0 >Synform_cp.su;
# rm -f tmp_cp.su;
#
#
# Source Wavelet: (Warning: Overwrites tmp.su)
# makewave nt=4001 dt=0.0005 fp=40 shift=1 file_out=tmp.su;
# sushw <tmp.su key=scalel,scalco,hcf,trwf,d2 a=1,1,0,0,0 >SrcWav_Ricker_40Hz_2s.su;
#
# NOTE: This requires the fdelmodc package from Jan-Willem
#       Thorbecke. You can find it at:
#       https://janth.home.xs4all.nl/Software/Software.html

# 0.3: Display Velocity Model
suximage <Synform_cp.su title=Synform_Velocity_Model &

# 0.4: Display Ricker Wavelet
suxwigb <SrcWav_Ricker_40Hz_2s.su title=40Hz_Ricker_Wavelet &

# 1.1: Generate Base Source Array
rm -f SrcArr.su; #Makes sure that we start off with a clean file
for i in {1..1601};do cat SrcWav_Ricker_40Hz_2s.su >> SrcArr.su; done

# 1.2: Generate Base Source Array
sushw <SrcArr.su key=tracl,tracr,fldr,tracf,scalco,sx a=1,1,1,1,-10,-20000 b=1,1,0,1,0,25 >tmp.su; mv -f tmp.su SrcArr.su;

# 1.3: Model Data
../../fdacrtmc        '#The RTM engine'\
	file_cp=Synform_cp.su  '#The input acoustic velocity model'\
	file_src=SrcArr.su     '#The input source array'\
	file_rcv=RcvArr.su     '#The output receiver data base filename'\
	top=2                  '#Top boundary is absorbing, others are by default'\
	npml=50                '#Number of absorbing layers'\
	rcv_top=1              '#Specifies to record along the top of the model'\
	rcv_p=1                '#Specifies to record pressure data'\
	rcv_write=1            '#Write out recorded data'\
	mig_mode=0             '#Specifies to only model data, not migrate data'\
	verbose=2;              #Verbosity level, >0 so that we see something

# 2.1: Create Direct-Wave Models
suwind <Synform_cp.su itmax=50 | sugain dt=1 scale=0.0 | sugain dt=1 bias=1900 >Direct_cp.su;

# 2.2: Model Direct-Wave Data
../../fdacrtmc file_cp=Direct_cp.su file_src=SrcArr.su\
	file_rcv=DirArr.su\
	top=2 npml=50 rcv_top=1 rcv_p=1 rcv_write=1 mig_mode=0 verbose=2;

# 2.3: Remove Direct Wave
sudiff RcvArr_rp.su DirArr_rp.su >RcvArr_ND.su; #Direct-Wave-Free Data

# 2.4: Taper Trace Ends
ntr=$(surange <RcvArr_ND.su key=tracl | sed -n '1p' | awk '{print $1}');      #We need to get the number of traces in the file
sutaper <RcvArr_ND.su ntr=${ntr} tend=50  >tmp.su; mv -f tmp.su RcvArr_ND.su; #This tapers the data at the end of the trace

# 2.5: Clean-up
rm -f Direct_cp.su Direct_ro.su RcvArr_rp.su DirArr_rp.su;

# 3.1: Migrate!
# Without compression this migration needs about 1.30 GiB of RAM
# With    compression this migration needs about 0.75 GiB of RAM at a tol.  of 1e-6
#                                                0.47 GiB of RAM at a tol.  of 1e-3
#                                                0.16 GiB of RAM at a prec. of 1e-3
../../fdacrtmc file_cp=Synform_cp.su file_src=SrcArr.su\
	file_rcv=RcvArr_ND.su    '#We now use the direct-wave free reflection data'\
	file_mig=Conventional.su '#Base filename for migrated image'\
	top=2 npml=50\
	mig_mode=1               '#This is the default (no decomposition) migration mode'\
	migdt=0.0025             '#This upsamples the cross-correlation, reducing memory & compute time'\
	migdx=5 migdz=5          '#This upsamples the final image,       reducing memory & compute time'\
	compress=1               '#This toggle compression on (1) or off (0).'\
	verbose=2;

# 3.2: Display Migrated Image
suximage <Conventional_mig.su title=Migrated_Image clip=2e5 &

# 3.3: Migrate Using The Hilbert Transform Imaging Condition
../../fdacrtmc file_cp=Synform_cp.su file_src=SrcArr.su file_rcv=RcvArr_ND.su\
	file_mig=Hilbert.su '#Base filename for migrated image'\
	top=2 npml=50 mig_mode=5 migdt=0.0025 migdx=5 migdz=5 compress=1 verbose=2;

# 3.4: Display Migrated Image
suximage <Hilbert_mig.su title=Hilbert_Image clip=2e5 &
exit;