from rsf.proj import *
Fetch('marmvel.hh','marm')
Flow('vel','marmvel.hh',
'''
dd form=native |
put label1=Z label2=X label=Velocity unit=m/s |
window j1=2 j2=4 |
smooth rect1=20 rect2=20 repeat=1 | put d1=25 d2=25
''')
Result('mvel','vel',
'''
grey color=j allpos=y bias=1500 scalebar=y title="Marmousi Model"
barlabel=Velocity barreverse=y barunit=m/s
labelfat=6 titlefat=6 labelsz=9 titlesz=10
unit1=m unit2=m
''')
Flow('grad1','vel','igrad square=n adj=n | math output="input/25" ')
Flow('grad2','vel','transp| igrad square=n adj=n | transp | math output="input/25" ')
Result('grad1','grey scalebar=y color=j title="Velocity gradient in Z direction" ')
Result('grad2','grey scalebar=y color=j title="Velocity gradient in X direction" ')
time=1.8
nt=121
dt=0.015
dtt=0.0005
factor=dt/dtt
ntt=(nt-1)*factor+1
ktt=(ntt-1)/10+1
Flow('source1',None,
'''
spike n1=%d d1=%g k1=%d |
ricker1 frequency=15
'''%(ntt,dtt,ktt))
Flow('real','source1','math "output=0"')
Flow('imag','source1','envelope hilb=y | halfint | halfint | math output="input/2" ')
Flow('csource1','real imag','cmplx ${SOURCES[1]}')
Flow('csource','csource1','window j1=%d'% factor)
Flow('source','source1','window j1=%d'% factor)
Result('source','graph title="Source Wavelet" ')
Flow('fft','vel','rtoc | fft3 axis=1 pad=1 | fft3 axis=2 pad=1')
Flow('refl','vel','spike k1=150 k2=288 | smooth rect1=2 rect2=2 repeat=2')
Flow('right left','vel grad1 grad2 fft',
'''
cisolr2grad seed=2010 dt=%g grad1=${SOURCES[1]} grad2=${SOURCES[2]} fft=${SOURCES[3]} left=${TARGETS[1]} npk=50 eps=1e-4
''' % dt)
Flow('wave','csource refl left right',
'''
cfftwave2omp ref=${SOURCES[1]} left=${SOURCES[2]} right=${SOURCES[3]} verb=y cmplx=n
''')
Plot('wave',
'''
window j3=25 |
grey label2="Z" label1="X" title="Isotropic"
yreverse=y gainpanel=all
''')
Result('wavesnap','wave',
'''
window n3=1 min3=%g |
grey title="With Gradient"
yreverse=y gainpanel=all
label1=Z unit1=m label2=X unit2=m
labelfat=6 titlefat=6 labelsz=9 titlesz=10
''' % time)
Plot('wavecross','wave',
'''
window n3=1 min3=%g n1=1 min1=8000 |
graph label2="Amplitude" label1="X" title="Isotropic" min2=-0.0004 max2=0.0005
''' % time)
Flow('right1 left1','vel fft',
'''
cisolr2 seed=2010 dt=%g fft=${SOURCES[1]} left=${TARGETS[1]} npk=50 eps=1e-4
''' % dt)
Flow('wave1','csource refl left1 right1',
'''
cfftwave2 ref=${SOURCES[1]} left=${SOURCES[2]} right=${SOURCES[3]} verb=y cmplx=n
''')
Plot('wave1',
'''
window j3=250 |
grey label2="X" label1="Z" title="Isotropic"
yreverse=y gainpanel=all
''')
Result('wave1snap','wave1',
'''
window n3=1 min3=%g |
grey title="Without Gradient"
yreverse=y gainpanel=all
label1=Z unit1=m label2=X unit2=m
labelfat=6 titlefat=6 labelsz=9 titlesz=10
''' % time)
Plot('wavecross1','wave1',
'''
window n3=1 min3=%g n1=1 min1=8000 |
graph label2="Amplitude" label1="X" title="Isotropic" min2=-0.0004 max2=0.0005
''' % time)
nt=1201
dt=0.0015
dtt=0.0005
factor=dt/dtt
ntt=(nt-1)*factor+1
ktt=(ntt-1)/10+1
Flow('source1-b',None,
'''
spike n1=%d d1=%g k1=%d |
ricker1 frequency=15
'''%(ntt,dtt,ktt))
Flow('real-b','source1-b','math "output=0"')
Flow('imag-b','source1-b','envelope hilb=y | halfint | halfint | math output="input/2" ')
Flow('csource1-b','real-b imag-b','cmplx ${SOURCES[1]}')
Flow('csource-b','csource1-b','window j1=%d'% factor)
Flow('source-b','source1-b','window j1=%d'% factor)
Result('source-b','graph title="Source Wavelet" ')
Flow('right1-b left1-b','vel fft',
'''
cisolr2 seed=2010 dt=%g fft=${SOURCES[1]} left=${TARGETS[1]} npk=50 eps=1e-4
''' % dt)
Flow('wave1-b','csource-b refl left1-b right1-b',
'''
cfftwave2 ref=${SOURCES[1]} left=${SOURCES[2]} right=${SOURCES[3]} verb=y cmplx=n
''')
Plot('wave1-b',
'''
window j3=25 |
grey label2="X" label1="Z" title="Isotropic"
yreverse=y gainpanel=all
''')
Result('wave1snap-b','wave1-b',
'''
window n3=1 min3=%g |
grey title="Reference"
yreverse=y gainpanel=all
label1=Z unit1=m label2=X unit2=m
labelfat=6 titlefat=6 labelsz=9 titlesz=10
''' % time)
Plot('wavecross1-b','wave1-b',
'''
window n3=1 min3=%g n1=1 min1=8000 |
graph label2="Amplitude" label1="X" title="Isotropic" min2=-0.0004 max2=0.0005
''' % time)
Flow('last','wave','window n3=1 min3=%g | scale axis=3 ' % time)
Flow('last1','wave1','window n3=1 min3=%g | scale axis=3 '% time)
Flow('last2','wave1-b','window n3=1 min3=%g | scale axis=3 '% time)
Flow('diff','last last2','math truth=last2.rsf output="abs(input-truth)" ')
Flow('diff1','last1 last2','math truth=last2.rsf output="abs(input-truth)" ')
Result('mdiff','diff','grey maxval=0.15 minval=0 clip=0.15 scalebar=y label1=Z unit1=m label2=X unit2=m title="Error (with gradient)" labelfat=6 titlefat=6 labelsz=9 titlesz=10 barlabel=Normalized\ Amplitude')
Result('mdiff1','diff1','grey maxval=0.15 minval=0 clip=0.15 scalebar=y label1=Z unit1=m label2=X unit2=m title="Error (without gradient)" labelfat=6 titlefat=6 labelsz=9 titlesz=10 barlabel=Normalized\ Amplitude ')
pos=8800
Flow('trace','wave','window n1=1 min1=%d min3=%g n3=1 min2=4000 max2=10000 ' % (pos,time))
Flow('trace1','wave1','window n1=1 min1=%d min3=%g n3=1 min2=4000 max2=10000 ' % (pos,time))
Flow('trace2','wave1-b','window n1=1 min1=%d min3=%g n3=1 min2=4000 max2=10000 ' % (pos,time))
Flow('comp','trace trace1 trace2','cat axis=2 ${SOURCES[1:3]} | scale axis=1')
Result('mcompp','trace trace2','cat axis=2 ${SOURCES[1]} | scale axis=1 | graph dash=0,1 plotcol=6,7 plotfat=3 label1=X unit1=m label2="Normalized Amplitude" unit2= wanttitle=y title="With gradient" screenwd=14.5 screenht=5.5 labelfat=4 labelsz=5.5')
Result('mcompp1','trace1 trace2','cat axis=2 ${SOURCES[1]} | scale axis=1 | graph dash=4,1 plotcol=5,7 plotfat=3 label1=X unit1=m label2="Normalized Amplitude" unit2= wanttitle=y title="Without gradient" screenwd=14.5 screenht=5.5 labelfat=4 labelsz=5.5')
End() |