from rsf.proj import *

dtt = 0.0005
nt = 1001
dt = 0.001
factor = dt/dtt

Flow('source1',None,
     '''
     spike n1=%d d1=%g k1=200 | 
     ricker1 frequency=50
     '''%(nt*factor,dtt))
Flow('real','source1','math "output=0"')
Flow('imag','source1','envelope hilb=y order=500 | 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" ')

par = {
    'dt' : dt,
    'w0' : 1500,
    'nz' : 201,
    'nx' : 201,
    'dz' : 8,
    'dx' : 8,
    'kz' : 100,
    'kx' : 101,
    'top': 131,
    'bot': 70,
    'topv' : 1800,
    'botv' : 3600,
    'topq' : 30,
    'botq' : 100
}
Flow('topv1',None,'spike n1=%(top)d n2=%(nx)d d1=%(dz)g d2=%(dx)g mag=%(topv)g' %par)
Flow('botv1',None,'spike n1=%(bot)d n2=%(nx)d d1=%(dz)g d2=%(dx)g mag=%(botv)g' %par)
Flow('vel1','topv1 botv1','cat axis=1 ${SOURCES[1]} | smooth rect1=3')
Flow('topq1',None,'spike n1=%(top)d n2=%(nx)d d1=%(dz)g d2=%(dx)g mag=%(topq)g' %par)
Flow('botq1',None,'spike n1=%(bot)d n2=%(nx)d d1=%(dz)g d2=%(dx)g mag=%(topq)g' %par)
Flow('dbq1','topq1 botq1','cat axis=1 ${SOURCES[1]} | smooth rect1=3')

Flow('botv2',None,'spike n1=%(bot)d n2=%(nx)d d1=%(dz)g d2=%(dx)g mag=%(topv)g' %par)
Flow('vel2','topv1 botv2','cat axis=1 ${SOURCES[1]} | smooth rect1=3')
Flow('botq2',None,'spike n1=%(bot)d n2=%(nx)d d1=%(dz)g d2=%(dx)g mag=%(botq)g' %par)
Flow('dbq2','topq1 botq2','cat axis=1 ${SOURCES[1]} | smooth rect1=3')

Flow('vel3','topv1 botv1','cat axis=1 ${SOURCES[1]} | smooth rect1=3')
Flow('dbq3','topq1 botq2','cat axis=1 ${SOURCES[1]} | smooth rect1=3')

Flow('refl',None,
     '''
     spike n1=%(nz)d n2=%(nx)d d1=%(dz)g d2=%(dx)g o1=0 o2=0 unit1=m unit2=m k1=%(kz)d k2=%(kx)d |
     smooth rect1=2 rect2=2 repeat=1
     '''% par)

Flow('fft','vel1','rtoc | fft3 axis=1 pad=1 | fft3 axis=2 pad=1')

Result('vel1','grey label1=Z unit1=m label2=X unit2=m allpos=y polarity=y bias=1800 title=Velocity screenratio=1')

for m in range (1,4):
    vel = 'vel%d' %m
    dbq = 'dbq%d' %m
    right = 'right-%d' %m
    left  = 'left-%d' %m
    wave = 'wave-%d' %m
    snap = 'snap-%d' %m
    Flow([right,left],[vel,'fft',dbq],
         '''
         zfraclr2 seed=2010 npk=10 eps=0.00001 dt=%(dt)g mode=0 rev=n sign=0
         fft=${SOURCES[1]} q=${SOURCES[2]}
         left=${TARGETS[1]} w0=%(w0)g
         ''' % par)

    Flow(wave,['csource','refl',left,right],
         '''
         cfftwave2 ref=${SOURCES[1]} left=${SOURCES[2]} right=${SOURCES[3]} verb=y
         cmplx=n
         ''')
    Flow(snap,wave,
           '''
           window n3=1 f3=330
           ''')
    Result(snap,
           '''
           grey wanttitle=n label1=Z label2=X gainpanel=all screenratio=1 clip=0.01 labelfat=3 titlefat=3
           ''')
    Result(wave,
           '''
           window j3=20 | grey wanttitle=n label1=Z label2=X gainpanel=all
           ''')

# acoustic

Flow(['right-0','left-0'],['vel3','fft','dbq3'],
     '''
     zfraclr2 seed=2010 npk=10 eps=0.00001 dt=%(dt)g mode=3 rev=n sign=0
     fft=${SOURCES[1]} q=${SOURCES[2]}
     left=${TARGETS[1]} w0=%(w0)g
     ''' % par)

Flow('wave-0',['csource','refl','left-0','right-0'],
     '''
     cfftwave2 ref=${SOURCES[1]} left=${SOURCES[2]} right=${SOURCES[3]} verb=y
     cmplx=n
     ''')
Flow('snap-0','wave-0',
       '''
       window n3=1 f3=330
       ''')
Result('snap-0',
       '''
       grey wanttitle=n label1=Z label2=X gainpanel=all screenratio=1 clip=0.01 labelfat=3 titlefat=3
       ''')
Result('wave-0',
       '''
       window j3=20 | grey wanttitle=n label1=Z label2=X gainpanel=all
       ''')

# compare with average
Flow('gamma','dbq3','math output="atan(1/input)/3.14159" ')

Flow(['right-4','left-4'],['vel3','fft','dbq3'],
     '''
     zfraclr2 seed=2010 npk=10 eps=0.00001 dt=%(dt)g mode=0 rev=n sign=0
     fft=${SOURCES[1]} q=${SOURCES[2]}
     left=${TARGETS[1]} w0=%(w0)g
     avg=y gamma=0.00799336
     ''' % par)
# gamma=0.00614674
#gamma=0.00318299

Flow('wave-4',['csource','refl','left-4','right-4'],
     '''
     cfftwave2 ref=${SOURCES[1]} left=${SOURCES[2]} right=${SOURCES[3]} verb=y
     cmplx=n
     ''')
Flow('snap-4','wave-4',
       '''
       window n3=1 f3=330
       ''')
Result('snap-4',
       '''
       grey wanttitle=n label1=Z label2=X gainpanel=all screenratio=1 clip=0.01 labelfat=3 titlefat=3
       ''')
Result('wave-4',
       '''
       window j3=20 | grey wanttitle=n label1=Z label2=X gainpanel=all
       ''')

Flow('diff','snap-3 snap-4','add ${SOURCES[1]} scale=1,-1')
Result('diff','grey wanttitle=n label1=Z label2=X gainpanel=all screenratio=1 pclip=99 labelfat=3 titlefat=3')

Flow('cross-3','snap-3','window f2=101 n2=1')
Flow('cross-4','snap-4','window f2=101 n2=1')
Flow('cross-d','diff','window f2=101 n2=1')
Flow('cross','cross-3 cross-4 cross-d','cat axis=2 ${SOURCES[1:3]}')
Result('cross','graph dash=0,4,1 plotcol=6,5,7 plotfat=2 label1=Distance unit1=m label2="Amplitude" unit2= wanttitle=n screenwd=14.5 screenht=5.5 labelfat=4 labelsz=5.5 labelfat=3 titlefat=3')


End()