from rsfproj import *
import zomig,rwezo

par = {
    'nz':751,  'dz':0.004,'oz':0,     'lz':'z', 'uz':'km',                 # CC axis 1
    'nx':1500, 'dx':0.004,'ox':2,     'lx':'x', 'ux':'km',                 # CC axis 2
    'nt':550,  'dt':0.002,'ot':0,     'lt':'"\F10 t\F3 "','ut':'s',     # RC axis 1
    'ng':512,  'dg':0.350,'og':-89.6, 'lg':'"\F10 g\F3 "','ug':'"\^o\_"', # RC axis 2
    'nT':1250, 'dT':0.002,'oT':0, 'kT':1000,    # time
    'tmx':10,'eps':0.0000001, 'verb':'y','nrmax':1,
    'ow':1, 'nw': 275, 'incore':'y',
    'spx':750,
    'prefix':'M-'
    }

# compute parameters
rwezo.param(par)

Fetch('marmvel.hh','marm')

# true velocity
Flow('vel','marmvel.hh',
     '''
     dd form=native |
     window min2=2000 n2=1500 |
     smooth rect1=5 rect2=5 |
     scale rscale=0.001 |
     put label1=z label2=x
     o1=%(oz)g d1=%(dz)g
     o2=%(ox)g d2=%(dx)g
     ''' % par)
Plot('vel',rwezo.cgrey('allpos=y',par))

# smooth velocity for HWT
Flow('vsm','vel','smooth rect1=201 rect2=201')

# ------------------------------------------------------------
# wavefronts (by HWT)
Flow('hwt','vsm',
     '''
     hwt2d verb=n xsou=5 zsou=0 
     nt=%(nt)d ot=%(ot)g dt=%(dt)g
     ng=%(ng)d og=%(og)g dg=%(dg)g
     ''' % par)

# coordinate system: rays(g,t)
Flow('cos','hwt','window min2=%(ot)g n2=%(nt)d' % par)

# coordinate system plot
rwezo.cos('cos',40,40,'plotcol=0',par)
Result(par['prefix']+'cos' % par,'vel cos','Overlay')

# slowness
rwezo.slow('sloCC','sloRC','vel','cos',par)
Result(par['prefix']+'sloRC','sloRC','transp |' % par
       + rwezo.rgrey('pclip=100 allpos=y bias=0.18',par))

# ------------------------------------------------------------
# coordinate system coefficients (A,B,M=mask) and references (A,B)
rwezo.abm('abmRC','abrRC','sloRC','cos',par)
rwezo.abmplot('abmRC',par)

# ------------------------------------------------------------
# time data CC: (t,x)
Flow('datCC',None,
     '''
     spike nsp=1 mag=1
     n1=%(nT)d o1=%(oT)g d1=%(dT)g k1=%(kT)d
     n2=%(nx)d o2=%(ox)g d2=%(dx)g k2=%(spx)d l2=%(spx)d |
     ricker1 frequency=15 |
     smooth rect2=3 |
     put label1=t label2=x
     ''' % par )
Result(par['prefix']+'datCC' % par,'datCC','grey pclip=100 wanttitle=n')

# ------------------------------------------------------------
# frequency data RC: (g,w)
rwezo.frq('frqRC','frqCC','datCC','cos',par)

# ------------------------------------------------------------
# WEM image CC: (z,x)
zomig.image('imgCC','sloCC','frqCC',par)
Plot  ('imgCC','window | transp |' + rwezo.cgrey('pclip=99',par))
Result(par['prefix']+'imgCC',['imgCC','cos'],'Overlay')

# ------------------------------------------------------------
# RWE image
rwezo.mig('migCC','migRC','frqRC','abmRC','abrRC','cos',par)
rwezo.mod('modCC','modRC','migRC','abmRC','abrRC','cos',par)

# ------------------------------------------------------------
# PLOTS
rwezo.plots(par)

End()