from rsf.proj import *

Fetch('velan.HH','dix')
Flow('velan','velan.HH','dd form=native')

Flow('pick','velan',
     '''
     mutter half=n x0=1500 v0=800 inner=y |
     cut max1=0.2 | scale axis=2 |
     pick rect1=20 rect2=20 vel0=1500
     ''')

Result('pick',
       '''
       grey color=j mean=y scalebar=y title="Picked Velocity"
       unit1=s unit2=m barunit=m/s barreverse=y
       ''')

Flow('slice','velan pick','slice pick=${SOURCES[1]}')

layers = []
spline = []
for lay in range(1,7):
    layer = 'layer%d' % lay
    Fetch(layer+'.txt','dix')
    Flow(layer,layer+'.txt',
         '''
         echo in=$SOURCE n1=6 n2=%d data_format=ascii_float |
         dd form=native |
         window n1=2 | reverse which=1          
         ''' % (19,19,15,16,18,7)[lay-1],stdin=0)
    Plot(layer,
         '''
         dd type=complex | window | 
         graph min2=0 max2=2.556 min1=13593 max1=20553 yreverse=y title=Horizons
         ''')
    layers.append(layer)
    Flow(layer+'s',layer,'transp | linear n1=175 d1=40 o1=13593 rect=3')
    spline.append(layer+'s')

Result('layers',layers,'Overlay')

Flow('spline',spline,'cat axis=2 ${SOURCES[1:%d]}' % len(spline))
Plot('spline','graph min2=0 max2=2.1 scalebar=y yreverse=y pad=n wanttitle=n wantaxis=n plotcol=7')
Result('spline','graph min2=0 max2=2.556 yreverse=y pad=n title=Horizons')

Flow('dips','spline',
     '''
     igrad | scale dscale=%g | clip2 lower=-2 upper=2 |
     cat axis=3 $SOURCE order=2,1 | transp | transp plane=23 |
     linear n1=640 d1=0.004 o1=0 rect=20 niter=100 | window 
     ''' % (1/0.004))

Result('dips','grey color=j scalebar=y title="Reflector Dip" ')

Flow('refl','spline','unif2 v00=1,2,3,4,5,6,7 n1=640 d1=0.004 o1=0 | ai2refl | ricker1 frequency=10')

Flow('dip','refl','dip rect1=50 rect2=5 order=3')

Result('dip','grey color=j scalebar=y title="Reflector Dip" ')

Flow('vint','pick slice dip',
     '''
     dixshape dip=${SOURCES[2]} order=3
     weight=${SOURCES[1]} 
     niter=100 verb=y rect1=20 rect2=5 
     ''')

Plot('vint',
     '''
     window max1=1.75 |
     grey color=j mean=y scalebar=y title="Interval Velocity"
     unit1=s unit2=m barunit=m/s barreverse=y
     ''')

Result('vint','vint spline','Overlay')

vints = []

for var in (0,2,4,8):
    prior = 'prior%d' % var

    Flow(prior,'dip vint',
         '''
         noise seed=2016 rep=y var=%g | 
         pwsmooth dip=${SOURCES[0]} ns=5 order=3 | smooth rect1=5 | 
         add ${SOURCES[1]}
         ''' % (var*1e4))

    vint = 'vint%d' % var

    Flow(vint,['pick','slice','dip',prior],
         '''
         dixshape dip=${SOURCES[2]} order=3
         weight=${SOURCES[1]} prior=${SOURCES[3]}
         niter=100 verb=y rect1=20 rect2=5 
         ''')

    Plot(vint,
         '''
         window max1=1.75 |
         grey color=j mean=y wanttitle=n
         unit1=s unit2=m 
         ''')

    vints.append(vint)

Result('scale',vints,'TwoRows')

vints = []

for seed in range(2016,2020):
    prior = 'dprior%d' % var

    Flow(prior,'dip vint',
         '''
         noise seed=%d rep=y var=2e4 | 
         pwsmooth dip=${SOURCES[0]} ns=5 order=3 | smooth rect1=5 | 
         add ${SOURCES[1]}
         ''' % (seed))

    vint = 'dvint%d' % var

    Flow(vint,['pick','slice','dip',prior],
         '''
         dixshape dip=${SOURCES[2]} order=3
         weight=${SOURCES[1]} prior=${SOURCES[3]}
         niter=100 verb=y rect1=20 rect2=5 
         ''')

    Plot(vint,
         '''
         window max1=1.75 |
         grey color=j mean=y wanttitle=n
         unit1=s unit2=m 
         ''')

    vints.append(vint)

Result('dix-real',vints,'TwoRows')

End()