## 
 # elastic modeling; wavefield separation operators
 # 2nd, 4th, 6th, and 8th orders accuracy
 ##
from rsf.proj import *
import sys
import fdmod,fdd,spk,stiff


# ------------------------------------------------------------
par = {
    'nx':600, 'ox':0, 'dx':0.002,  'lx':'x', 'ux':'km',
    'nz':600, 'oz':0, 'dz':0.002,  'lz':'z', 'uz':'km',
    'nt':1601,'ot':0, 'dt':0.0002, 'lt':'t', 'ut':'s',
    'kt':150,
    'jsnap':200,
    'height':10,
    'nb':0,
    'frq':75,
    'ratio':1
    }
fdmod.param(par)

nframe=8

# ------------------------------------------------------------
# model parameters
# ------------------------------------------------------------
Flow('zero',None,
     '''
     spike nsp=1 mag=0.0
     n1=%(nz)d o1=%(oz)g d1=%(dz)g 
     n2=%(nx)d o2=%(ox)g d2=%(dx)g |
     put label1=%(lz)s label2=%(lx)s unit1=%(uz)s unit2=%(ux)s
     ''' % par)

Flow('vp',     'zero','add add=3.00')
Flow('vs',     'zero','add add=1.50')
Flow('ro',     'zero','add add=1.00')
Flow('epsilon','zero','add add=+0.25')
Flow('delta',  'zero','add add=-0.29')
Flow('nu',     'zero','add add=0.00')

stiff.iso2d('cI','vp','vs','ro',par)
stiff.tti2d('cA','vp','vs','ro','epsilon','delta','nu',par)

# ------------------------------------------------------------
# derivative operators

fdd.derivatives(par)
spk.delt('spk',64,64,par)

for order in range(2,10,2):
#order=[2,4,6,8]
    tag=str(order)

    fdd.separatorD('dzI'+tag,'dxI'+tag,'spk','cI',
                   'y','x','sine',1.0, order,5,5,par)
    fdd.separatorD('dzM'+tag,'dxM'+tag,'spk','cA',
                   'y','x','sine',1.0, order,25,25,par)

    fdd.oneirST('iop'+tag,'dzI'+tag,'dxI'+tag,10,10,'pclip=99.5',par)
    fdd.oneirST('mop'+tag,'dzM'+tag,'dxM'+tag,10,10,'pclip=100 titlesz=15',par)
    fdd.oneirST('mop5-'+tag,'dzM'+tag,'dxM'+tag,5,5,'pclip=99',par)



# the following studies the difference between 8th and 2nd order operators
Flow('ax8','mop8-dxM8','scale axis=123')
Flow('ax2','mop2-dxM2','scale axis=123')

Flow('az8','mop8-dzM8','scale axis=123')
Flow('az2','mop2-dzM2','scale axis=123')

Flow('diffaz',['az8','az2'],' add scale=1,-1 ${SOURCES[1]}')
Flow('diffax',['ax8','ax2'],' add scale=1,-1 ${SOURCES[1]}')

fdd.oneirST('diff','diffaz','diffax',5,5,'clip=5 color=Ic',par)


End()