from rsf.proj import *
import rsf.recipes.pmig as pmig

Fetch('beinew.HH','midpts')
Flow('bei','beinew.HH',
     '''
     dd form=native |
     put
     label1=Time unit1=s
     label2=Offset unit2=km
     label3=Midpoint unit3=km
     ''')

Flow('cdp','bei','window n3=1')
Plot('cdp','grey title="(a) Input" ')
Plot('cdpw','cdp',
     '''
     wiggle transp=y yreverse=y poly=y 
     wherexlabel=t wheretitle=b title="(a) Input" 
     ''')

d2=0.0335
d1=0.004

def dip(slow2):
    return '''
    math output="(%g*x2/x1)" | mutter half=n v0=1.4 tp=0.4
    ''' % (slow2*d2/d1)

Flow('dip0','cdp',dip(1./3.-1./9.))
Flow('dip1','cdp',dip(1./9.))

#Flow('mask','dip0',
#     'math output=1 | nmostretch half=n v0=3 | mutter half=n v0=1.4 tp=0.4')

Flow('dip2','cdp dip0 dip1',
     '''
     nmostretch half=n v0=3 |
     twodip2 dip1=${SOURCES[1]} q0=0
     order=3 niter=10 eps=200 lam=40 |
     window n3=1 |
     math ref=${SOURCES[2]} output="input+ref" 
     ''')

Flow('dip3',['cdp','dip2'],
     '''
     dip idip=${SOURCES[1]} order=3
     niter=10 liter=100 rect1=50 rect2=20 pmin=0
     ''')
Plot('dip','dip3','grey title="(b) Slope" allpos=y color=j scalebar=y')

Result('pint','cdp dip','SideBySideAniso')

Flow('pwd','cdp dip3','pwd order=3 dip=${SOURCES[1]}')
Plot('pwd','grey title=Residual clip=3.9796e+06')

Result('mpwd','cdp dip pwd','SideBySideAniso',vppen='txscale=1.5')

Flow('txt','cdp dip3',
     '''
     noise rep=y seed=2004 |
     trismooth2 rect2=20 dip=${SOURCES[1]} |
     mutter v0=1.4
     ''')
Plot('txt','grey title=Texture')

Flow('txt2','cdp dip3',
     '''
     noise rep=y seed=2004 |
     pwdsmooth2 rect2=20 dip=${SOURCES[1]} |
     mutter v0=1.4
     ''')
Plot('txt2','grey title=Texture')

Result('mtxt','cdp dip txt','SideBySideAniso',vppen='txscale=1.5')

Flow('pnmo vel','cdp dip3','pnmo dip=${SOURCES[1]} vel=${TARGETS[1]}')
Plot('pnmo','grey title="(c) Oriented NMO" ')
Plot('pnmow','pnmo',
     '''
     wiggle transp=y yreverse=y poly=y 
     wherexlabel=t wheretitle=b title="(c) Oriented NMO" 
     ''')

Result('pnmo','cdp dip pnmo','SideBySideAniso',vppen='txscale=1.5')
Result('pnmow','cdpw dip pnmow','SideBySideAniso',vppen='txscale=1.5')

Flow('dip05','dip3','window min2=0.5')

Flow('vtr','cdp dip05',
     '''
     window min2=0.5 |
     pveltran v0=1.5 dv=0.01 nv=151 dip=${SOURCES[1]}
     ''')

Flow('vpk','vtr',
     '''
     envelope | scale axis=2 | pick rect1=20 | 
     window min1=1.2 | remap1 pattern=$SOURCE order=2
     ''')
Plot('vpk',
     '''
     graph transp=y yreverse=y pad=n min2=1.495 max2=3.005 
     plotcol=5 wanttitle=n wantaxis=n
     ''')

Plot('vtr',
     '''
     envelope | 
     grey allpos=y label2=Velocity unit2=km/s
     title="(b) Oriented Velocity Map" 
     ''')
Plot('vtr2','vtr vpk','Overlay')

Flow('vsc','cdp','vscan v0=1.5 dv=0.01 nv=151')
Plot('vsc',
     '''
     envelope | 
     grey allpos=y label2=Velocity unit2=km/s
     title="(a) Velocity Scan" 
     ''')
Plot('vsc2','vsc vpk','Overlay')

Result('bvsc','vsc2 vtr2','SideBySideAniso',vppen='txscale=1.5')

Flow('dips05','dip05','smooth rect1=100')

Flow('dipt05','dips05','deriv')
Flow('vint','cdp dips05 dipt05',
     '''
     window min2=0.5 |
     pveltran v0=1.5 dv=0.01 nv=151 dip=${SOURCES[1]}
     interval=y dipt=${SOURCES[2]}
     ''')

Result('bvint','vint',
     '''
     envelope | 
     grey allpos=y label2="Interval Velocity" unit2=km/s
     wanttitle=n 
     ''')

Flow('dips0','dip3','spray axis=3 n=250 d=0.0335 o=7.705')
Flow('dips','bei dips0',
     'dip idip=${SOURCES[1]} rect1=25 rect2=10 rect3=25')

Flow('hdip','dips','window n4=1')
Flow('ydip','dips','window f4=1')

Flow('pwds','bei dips','pwd dip=${SOURCES[1]}')

pmig.pmig('bei','hdip','ydip',500,125,20)

End()