from rsf.proj import *
import math

## Plot font and screen ratio, width, height for uniform ppt figs.
p1=0.7
sr=1.0
sw=7.0
sh=10.0
xll=2.0
fat=2

#####
##Synthetic CMP Parameters
nx=101
delx=0.0125*2
ox=-1.25
ny=nx
dely=delx
oy=ox
nt=1001
dt=0.004
fw=10.0
####

################################################################################################################
########  4 Reflection Events ################################################

#111111111 Make a layer: z = z0 + x*dx + y*dy

dx=0.4
dy=0.1
dt=.004
z0=0.8
##z0=0.75

cg = 1/math.sqrt(1+dx*dx+dy*dy)
ca = -dx*cg
cb = -dy*cg
d = z0*cg

mx = 0
my = 0

D = d - mx*ca - my*cb

v0 = 2.5
t0 = 2*D/v0
it0=int(t0/dt)
it1=it0
#print it1
Flow('spike1',None,
     '''
     spike n1=%d k1=%d nsp=1 | ricker1 frequency=%g |
     spray axis=2 n=%d o=%g d=%g |
     spray axis=3 n=%d o=%g d=%g
     ''' % (nt,it1,fw,ny,oy,dely,nx,ox,delx))

wx = 1-ca*ca
wy = 1-cb*cb
wxy = - ca*cb

Flow('vel1.asc',None,'echo %g %g %g n1=3 data_format=ascii_float in=$TARGET' % (wx,wy,wxy))
Flow('vel1','vel1.asc','dd form=native | scale dscale=%g | spray axis=1 n=%d' % (1/(v0*v0),nt),local=1)


Flow('cmp1','spike1 vel1','inmo3 velocity=${SOURCES[1]}',local=1)

#222222222222 Make another layer: z = z0 + x*dx + y*dy

dx=0.7
dy=0.41
dt=.004
z0=0.85

cg = 1/math.sqrt(1+dx*dx+dy*dy)
ca = -dx*cg
cb = -dy*cg
d = z0*cg

mx = 0
my = 0

D = d - mx*ca - my*cb

v0 = 1.7
t0 = 2*D/v0
it0=int(t0/dt)
it2=it0
#print it2
Flow('spike2',None,
     '''
     spike n1=%d k1=%d nsp=1 | ricker1 frequency=%g |
     spray axis=2 n=%d o=%g d=%g |
     spray axis=3 n=%d o=%g d=%g
     ''' % (nt,it2,fw,ny,oy,dely,nx,ox,delx))

wx = 1-ca*ca
wy = 1-cb*cb
wxy = - ca*cb

Flow('vel2.asc',None,'echo %g %g %g n1=3 data_format=ascii_float in=$TARGET' % (wx,wy,wxy))
Flow('vel2','vel2.asc','dd form=native | scale dscale=%g | spray axis=1 n=%d' % (1/(v0*v0),nt),local=1)


Flow('cmp2','spike2 vel2','inmo3 velocity=${SOURCES[1]}',local=1)


#3333333333333 Make yet another layer: z = z0 + x*dx + y*dy

dx=0.1
dy=0.9
dt=.004
z0=1.1

cg = 1/math.sqrt(1+dx*dx+dy*dy)
ca = -dx*cg
cb = -dy*cg
d = z0*cg

mx = 0
my = 0

D = d - mx*ca - my*cb

v0 = 1.75
t0 = 2*D/v0
it0=int(t0/dt)
it3=it0
#print it3
Flow('spike3',None,
     '''
     spike n1=%d k1=%d nsp=1 | ricker1 frequency=%g |
     spray axis=2 n=%d o=%g d=%g |
     spray axis=3 n=%d o=%g d=%g
     ''' % (nt,it3,fw,ny,oy,dely,nx,ox,delx))

wx = 1-ca*ca
wy = 1-cb*cb
wxy = - ca*cb

Flow('vel3.asc',None,'echo %g %g %g n1=3 data_format=ascii_float in=$TARGET' % (wx,wy,wxy))
Flow('vel3','vel3.asc','dd form=native | scale dscale=%g | spray axis=1 n=%d' % (1/(v0*v0),nt),local=1)

Flow('cmp3','spike3 vel3','inmo3 velocity=${SOURCES[1]}',local=1)

#444444444444 Make a 4th layer: z = z0 + x*dx + y*dy

dx=0.2
dy=0.1
dt=.004
z0=1.3

cg = 1/math.sqrt(1+dx*dx+dy*dy)
ca = -dx*cg
cb = -dy*cg
d = z0*cg

mx = 0
my = 0

D = d - mx*ca - my*cb

v0 = 2.0
t0 = 2*D/v0
it0=int(t0/dt)
it4=it0
#print it4
Flow('spike4',None,
     '''
     spike n1=%d k1=%d nsp=1 | ricker1 frequency=%g |
     spray axis=2 n=%d o=%g d=%g |
     spray axis=3 n=%d o=%g d=%g
     ''' % (nt,it4,fw,ny,oy,dely,nx,ox,delx))

wx = 1-ca*ca
wy = 1-cb*cb
wxy = - ca*cb

Flow('vel4.asc',None,'echo %g %g %g n1=3 data_format=ascii_float in=$TARGET' % (wx,wy,wxy))
Flow('vel4','vel4.asc','dd form=native | scale dscale=%g | spray axis=1 n=%d' % (1/(v0*v0),nt),local=1)

Flow('cmp4','spike4 vel4','inmo3 velocity=${SOURCES[1]}',local=1)

## Add events to create CMP.
Flow('cmp0','cmp1 cmp2 cmp3 cmp4',
     '''
     math cmp2=${SOURCES[1]} cmp3=${SOURCES[2]} cmp4=${SOURCES[3]} output="input+cmp2+cmp3+cmp4" |
     window max1=1.5 min1=0.5 j1=2 | put o1=0
     ''',local=1)

Flow('cmp','cmp0','noise seed=1 range=0.3 ')
Result('cmp0','cmp0','byte clip=0.3| grey3   frame1=100 frame2=50 frame3=50 flat=y wanttitle=n label2=X unit2=km label3=Y unit3=km ')
Result('cmp','cmp','byte clip=0.3| grey3   frame1=100 frame2=50 frame3=50 flat=y wanttitle=n label2=X unit2=km label3=Y unit3=km ')

Flow('fcmp','cmp','fft1 | transp plane=13 memsize=1000  | transp plane=12 memsize=1000')
"""
     shiftd2 ns=2 | transp plane=13 memsize=1000 | shiftd2 ns=2 |
     transp plane=24 memsize=1000 | transp plane=12 | put n3=25 n4=1 o3=0 d3=1 | 
     transp plane=45 memsize=1000 | transp plane=34 memsize=1000
"""
Flow('shifts','fcmp',
     '''
     cshifts2 ns1=2 ns2=2 | transp plane=34 memsize=300
     ''')
#Flow('sh1','shiftsa','window f4=13')
#Flow('shifts','shiftsa sh1','window n4=12 | cat ${SOURCES[1]} axis=4')
Flow('flt pre','shifts fcmp',
     'clpf match=${SOURCES[1]} pred=${TARGETS[1]} rect1=7 rect2=7 niter=30')

Flow('tpre','pre','transp plane=13 memsize=1000| transp plane=23 memsize=1000 | fft1 inv=y')

Result('tpre','tpre','put o1=0 |byte clip=0.3| grey3   frame1=100 frame2=50 frame3=50 flat=y wanttitle=n label2=X unit2=km label3=Y unit3=km')




#################
## 2D RNA
################

Flow('shifts2d','fcmp','shiftd1 ns=3 |transp plane=25 | window squeeze=y')

Flow('flt2d pre2d','shifts2d fcmp',
     'clpf match=${SOURCES[1]} pred=${TARGETS[1]} rect1=7 niter=400')

Flow('tpre2d','pre2d','transp plane=13 memsize=1000| transp plane=23 memsize=1000 | fft1 inv=y')
Result('tpre2d','tpre2d','put o1=0 |byte clip=0.3| grey3   frame1=100 frame2=50 frame3=50 flat=y wanttitle=n label2=X unit2=km label3=Y unit3=km')



End()