up [pdf]
from rsf.proj import *
from math import *
from rsf.recipes.beg import server as private
from rsf.prog import RSFROOT
ns=80

def Grey(data,other):
        Result(data,
           '''
                        grey transp=y color=j 
                        unit1=s unit2= allpos=y scalebar=y 
                        parallel2=n labelsz=6  %s'''%(other))#screenratio=1

def Grey3(data,other):
        Result(data,
           '''
                        put d3=1 |byte pclip=99 |  transp plane=23 memsize=1000 | grey3 flat=n transp=y color=j 
                        unit1=s unit2= allpos=y  frame1=100 frame2=1 frame3=20 label2=Trace label1=Time
                        label3=Frequency unit1=s unit2= unit3=Hz framelabelcol=VP_BLACK  
                        parallel2=n labelsz=6 screenratio=1 point1=0.8 point2=0.8 %s'''%(other))

def Grey3tfl(data,other):
        Result(data,
           '''
                        put d3=1 |  transp plane=23 memsize=1000|byte pclip=99 | grey3 flat=n transp=y color=j 
                        unit1=s unit2= allpos=y  frame1=100 frame2=1 frame3=20 label2=Trace label1=Time
                        label3=Frequency unit1=s unit2= unit3=Hz framelabelcol=VP_BLACK  
                        parallel2=n labelsz=6 screenratio=1 point1=0.8 point2=0.8 %s'''%(other))

def Greytf(data1,data2,other):
        Result(data1,data2,
          '''
       window n3=1 f3=%d max2=180 | scale axis=2 | 
       grey color=j allpos=y wanttitle=y screenratio=1.5 label2=Frequency unit2=Hz label1=Time unit1=s %s
       '''%(ns,other))

def Greyf(data1,data2,fmin,other):
        Result(data1,data2,
       '''
       window n2=1 min2=%g | put d2=1 | 
       sfgrey color=j allpos=y scalebar=n title="%gHZ" label2="Trace" unit2= label1=Time unit1=s pclip=99 %s
       '''%(fmin,fmin,other))

#-----read data from segy file
Fetch('ln472_old.sgy','guochang',private)

Flow('old','ln472_old.sgy',
     '''
     segyread bfile=/dev/null hfile=/dev/null read=data |  
     put n1=751 n2=201 o1=0 d1=0.002 02=1 d2=0.01 label1="Time" unit1="s" 
     label2="X" unit2="km" | window min1=0
     ''')




def ref(trace):
    out = 'ref%d' % trace
    Flow(out+'.asc',None,
         '''
         echo %d 0 %d 1 n1=4 in=$TARGET data_format=ascii_float
         ''' % (trace,trace))
    Plot(out,out+'.asc',
         '''
         dd form=native type=complex | 
         graph min1=0 max1=201 min2=0 title="" wantaxis=n scalebar=n pad=n plotfat=4 dash=2
         ''')
    return out


Plot('old','old','grey title=" " scalebar=n label2="Distance" label1=Time unit1=s ')

ref(80)

Result('old-1',['old', 'ref80'],'Overlay')


Result('spectra','old','spectra all=y | scale axis=1 | graph pad=n label2="Scaled amplitude" wanttitle=n unit2="" plotfat=4')

############################### LTFT

Flow('ltft','old',
     '''
     timefreq rect=8 verb=n nw=377 dw=0.665779 niter=100
     ''')

#Result('tf-l','ltft',
#       '''
#       window n3=1 f3=%d max2=180 | scale axis=2 | 
#       grey color=j allpos=y scalebar=y wanttitle=n screenratio=1.5 
#       ''' % ns)

# Average frequency
Flow('num-l','ltft','window max2=110 | math output="input*x2" | stack axis=2 ')

Flow('num-l-1','num-l','window f1=200 | smooth rect1=1 rect2=1')
Flow('num-l-2','num-l','window n1=200 ')
Flow('num-l-a','num-l-2 num-l-1','cat ${SOURCES[1]} axis=1 o=0 d=0.002')

Flow('den-l','ltft','window max2=110 | stack axis=2')
Flow('lcf-l','num-l den-l','divn den=${SOURCES[1]} rect1=1  rect2=1')

Result('lcf-l','lcf-l',
       '''
        grey allpos=y scalebar=y  clip=86  barunit=Hz  color=j title="" label2="Distance" label1=Time unit1=s
       ''')




############## ST
Flow('st','old','st  | math output="abs(input)" | real ')
#Result('tf-s','st',
#       '''
#       window n3=1 f3=%d max2=180| scale axis=2 | 
#       grey color=j allpos=y scalebar=y wanttitle=n screenratio=1.5
#       ''' % ns)


########################################################################
# INITIALIZATION for Matlab
########################################################################
matlab         = WhereIs('matlab')
matROOT = '../Matfun'
matfun = 'lowf'
matlabpath = os.environ.get('MATLABPATH',os.path.join(RSFROOT,'lib'))

if not matlab:
    sys.stderr.write('\nCannot find Matlab.\n')
    sys.exit(1)

############################################################
## process field data
############################################################
Flow('tfsswt1 tfsswt2 tfsswt3 tfsswt4 tfsswt5 tfsswt6 tfsswt',[os.path.join(matROOT,matfun+'.m'),'old'],
     '''MATLABPATH=%(matlabpath)s %(matlab)s 
     -nosplash -nojvm -r "addpath %(matROOT)s;%(matfun)s('${SOURCES[1]}','${TARGET}','${TARGETS[1]}','${TARGETS[2]}','${TARGETS[3]}','${TARGETS[4]}','${TARGETS[5]}','${TARGETS[6]}');quit"
     '''%vars(),stdin=0,stdout=-1)

nt=751
nx=201
dt=0.002
dx=1
nf=126
df=(1/2/dt)/(nf-1);

Flow('lowf-tfsswt1','tfsswt1','put o1=%g d1=%g n1=%d o2=%g d2=%g n2=%d'%(0,dt,nt,0,dx,nx))
Flow('lowf-tfsswt2','tfsswt2','put o1=%g d1=%g n1=%d o2=%g d2=%g n2=%d'%(0,dt,nt,0,dx,nx))
Flow('lowf-tfsswt3','tfsswt3','put o1=%g d1=%g n1=%d o2=%g d2=%g n2=%d'%(0,dt,nt,0,dx,nx))
Flow('lowf-tfsswt4','tfsswt4','put o1=%g d1=%g n1=%d o2=%g d2=%g n2=%d'%(0,dt,nt,0,dx,nx))
Flow('lowf-tfsswt5','tfsswt5','put o1=%g d1=%g n1=%d o2=%g d2=%g n2=%d'%(0,dt,nt,0,dx,nx))
Flow('lowf-tfsswt6','tfsswt6','put o1=%g d1=%g n1=%d o2=%g d2=%g n2=%d'%(0,dt,nt,0,dx,nx))
Flow('lowf-tfsswt','tfsswt','put o1=%g d1=%g n1=%d o2=%g d2=%g n2=%d o3=%g d3=%g n3=%d'%(0,dt,nt,0,df,nf,0,dx,nx))

Grey('lowf-tfsswt1','title="20 Hz"')
Grey('lowf-tfsswt2','title="40 Hz"')
Grey('lowf-tfsswt3','title="60 Hz"')
Grey('lowf-tfsswt4','title="80 Hz"')
Grey('lowf-tfsswt5','title="100 Hz"')
Grey('lowf-tfsswt6','title="120 Hz"')

Grey3('lowf-tfsswt','title="SSWT"')
Grey3tfl('ltft','title="TFL" frame3=60')
Grey3('st','title="ST" frame3=60')


Greytf('tf-l','ltft','title="Trace 80"')
Greytf('tf-s','st','title="Trace 80"')
Greytf('tf-ss','lowf-tfsswt','title="Trace 80"')

Greyf('l1','ltft',20,'')
Greyf('l2','ltft',40,'')
Greyf('l3','ltft',60,'')

Greyf('s1','st',20,'')
Greyf('s2','st',40,'')
Greyf('s3','st',60,'')

Greyf('ss1','lowf-tfsswt',20,'')
Greyf('ss2','lowf-tfsswt',40,'')
Greyf('ss3','lowf-tfsswt',60,'')

## Creating framebox1
x=415
y=1.05
w=50
w1=0.2

Flow('frame1.asc',None,'echo %s n1=10 data_format=ascii_float in=$TARGET'% \
        string.join(map(str,(x,y,x+w,y,x+w,y+w1,x,y+w1,x,y))))
Plot('frame1','frame1.asc',
        '''
        dd type=complex form=native |
        graph min1=0 max1=471 min2=0 max2=2.048 pad=n plotfat=15 plotcol=2 
        wantaxis=n wanttitle=n yreverse=y 
        ''')

## Creating framebox2
x=402
y=1.35
w=50
w1=0.2

Flow('frame2.asc',None,'echo %s n1=10 data_format=ascii_float in=$TARGET'% \
        string.join(map(str,(x,y,x+w,y,x+w,y+w1,x,y+w1,x,y))))
Plot('frame2','frame2.asc',
        '''
        dd type=complex form=native |
        graph min1=0 max1=471 min2=0 max2=2.048 pad=n plotfat=15 plotcol=3 
        wantaxis=n wanttitle=n yreverse=y 
        ''')

## Creating framebox3
x=22
y=0.95
w=30
w1=0.2

Flow('frame3.asc',None,'echo %s n1=10 data_format=ascii_float in=$TARGET'% \
        string.join(map(str,(x,y,x+w,y,x+w,y+w1,x,y+w1,x,y))))
Plot('frame3','frame3.asc',
        '''
        dd type=complex form=native |
        graph min1=0 max1=471 min2=0 max2=2.048 pad=n plotfat=15 plotcol=4 
        wantaxis=n wanttitle=n yreverse=y 
        ''')

## Creating framebox4
x=58
y=1.45
w=40
w1=0.2

Flow('frame4.asc',None,'echo %s n1=10 data_format=ascii_float in=$TARGET'% \
        string.join(map(str,(x,y,x+w,y,x+w,y+w1,x,y+w1,x,y))))
Plot('frame4','frame4.asc',
        '''
        dd type=complex form=native |
        graph min1=0 max1=471 min2=0 max2=2.048 pad=n plotfat=15 plotcol=7 
        wantaxis=n wanttitle=n yreverse=y 
        ''')

Result('ss1-0','Fig/ss1.vpl frame1 frame2 frame3 frame4','Overlay')
Result('ss2-0','Fig/ss2.vpl frame1 frame2 frame3 frame4','Overlay')
Result('ss3-0','Fig/ss3.vpl frame1 frame2 frame3 frame4','Overlay')
Result('l1-0','Fig/l1.vpl frame1 frame2 frame3 frame4','Overlay')
Result('l2-0','Fig/l2.vpl frame1 frame2 frame3 frame4','Overlay')
Result('l3-0','Fig/l3.vpl frame1 frame2 frame3 frame4','Overlay')
Result('s1-0','Fig/s1.vpl frame1 frame2 frame3 frame4','Overlay')
Result('s2-0','Fig/s2.vpl frame1 frame2 frame3 frame4','Overlay')
Result('s3-0','Fig/s3.vpl frame1 frame2 frame3 frame4','Overlay')



#Flow('ss1','lowf-tfsswt','window min2=15 n2=1')
#Flow('ss2','lowf-tfsswt','window min2=31 n2=1')
#Flow('ss3','lowf-tfsswt','window min2=70 n2=1')


Plot('label1',None,
        '''
        box x0=5.8 y0=7.1 label="Gas?" xt=-0.5 yt=0.5 length=0.8 
        ''')

Plot('label2',None,
        '''
        box x0=5.9 y0=4.5 label="Gas?" xt=-0.5 yt=0.5 length=0.8 
        ''')




#Result('lowf','Fig/lowf0.vpl label1 label2','Overlay')
#Result('higf','Fig/higf0.vpl label1 label2','Overlay')




End()

sfsegyread
sfput
sfwindow
sfgrey
sfdd
sfgraph
sfspectra
sfscale
sftimefreq
sfmath
sfstack
sfsmooth
sfcat
sfdivn
sfst
sfreal
sfbyte
sftransp
sfgrey3
sfbox