from rsf.proj import *
from rsf.recipes.beg import server as private


models='''uniform barrier wedge twolayers waveguid'''

models=Split(models)

size=160
freq=16

half=size/2
n=size-1
d=1.0/size
slice1=0.7*size
slice2=0.55*size


## generate source files 
Flow('OneShot.asc',None,'''echo '0.5 0.5 0.1' in=$TARGET data_format=ascii_float n1=3 n2=1''')
Flow('OneShot.rsf','OneShot.asc','dd form=native')

Flow('ThreeShot.asc',None,'''echo '0.5 0.5 0.1 0.25 0.25 0.1 0.75 0.75 0.5' in=$TARGET data_format=ascii_float n1=3 n2=3''')
Flow('ThreeShot.rsf','ThreeShot.asc','dd form=native')


## generate models 
Flow(models[0],None,
     '''
     math 
     n1=%d n2=%d n3=%d 
     d1=%g d2=%g d3=%g 
     label1=X label2=Y label3=Z 
     unit1=km unit2=km unit3=km 
     label=Velocity unit=km/s output=1
     ''' % (n,n,n,d,d,d))


## 1. High-velocity separator
Flow('separ',models[0],'cut max3=0.75 min2=0.25 max2=0.3')
Flow(models[1],'separ','math output="1e10*(1-input) + input"')


## 2. Wedge
Flow('mask',models[0],'math output="0.8-0.2*x2-x3" | mask min=0')
Flow('wedge1','mask','dd type=float | math output="1.5*input+3*(1-input)" ')
Flow('mask2','wedge1','math output="0.4+0.1*x2-x3" | mask min=0')
Flow(models[2],'mask2 wedge1','dd type=float | math other=${SOURCES[1]} output="2*input+other*(1-input)"')


## 3. Two decreasing layers
Flow('top',models[0],'window n2=%d | math output=4' % half)
Flow('bot',models[0],'window f2=%d | put o2=0 | math output=1' % half)
Flow(models[3],'top bot','cat axis=2 ${SOURCES[1]}')


## 4. Wave guide
Flow(models[4],models[0],'math output="1.25 * (1.0 - 0.4*exp(-32.*((x1-0.5)^2+(x2-0.5)^2)))"')

for model in models[1:5]:

    Result(model,'''byte gainpanel=all mean=y | 
    grey3 frame1=%d frame2=%d frame3=%d color=j title=%s screenratio=1
    ''' % (half,half,half,model.capitalize()))

    if model == 'barrier' or model== 'twolayers':
        source='ThreeShot.rsf'
    else:
        source='OneShot.rsf'

    if model == 'waveguid':
        slice=slice2
    else:
        slice=slice1

    solve='solution-'+model

    Flow(solve,[model, source],
         '''
         mpipsp velocity=${SOURCES[0]} source=${SOURCES[1]}
         freq=%g solution=$TARGET
         ''' % freq)

    Plot(solve,
           '''
           real | window n1=1 f1=%d | transp |
           grey screenratio=1 mean=y
               scalebar=y wanttitle=n 
           ''' % slice)

Result('solution-models','solution-barrier solution-wedge solution-twolayers solution-waveguid','TwoColumns')

End()