up [pdf]
from rsf.proj import*
from rsf.recipes.beg import server

#######################################################################
## parameters definition
#######################################################################

clip=0.5 #display percentage
nfw=9

#######################################################################
## module definition
#######################################################################
def Grey(data,other):
        Result(data,'grey %s label2=Trace unit2="" label1=Time unit1="s" title="" wherexlabel=b wanttitle=n  '%other)

def plot (other):
        return'''       
        grey label2=Trace unit2="" label1=Time unit1="s" title="" wherexlabel=b   %s'
        '''%other

#######################################################################
## data format conversion # segy->rsf
#######################################################################
for file in [
         '2-SOURCE-simulated_receiver_gather_SEGY','1-SOURCE-simulated_receiver_gather_SEGY']:
         if file == '2-SOURCE-simulated_receiver_gather_SEGY':
                file1='fairblended'
         elif file == '1-SOURCE-simulated_receiver_gather_SEGY':
                file1='fairunblended'

         Fetch(file+'.sgy','fairfield',server)

         Flow([file1,'t'+file1,file1+'.asc',file1+'.bin'],
                file+'.sgy',
                '''
                segyread
                tfile=${TARGETS[1]}
                hfile=${TARGETS[2]}
                bfile=${TARGETS[3]}
                ''')
         #Result(file1,file1,'grey clip=%g'%(clip))

Flow('fairblended1','fairblended','window n2=201 | put o2=0 ')
Flow('fairblended2','fairblended','window f2=201 | put o2=0 ') # This is what we want to use.
Flow('fairunblended1','fairunblended','window n2=201 | put o2=0 ')
Flow('fairunblended2','fairunblended','window f2=201 | put o2=0 ') # This is what we want to use.

#Result('fairblended1','grey title="Source1" clip=%g '%(clip))
#Result('fairblended2','grey title="Source2" clip=%g '%(clip))

#######################################################################
## read shottime 
#######################################################################
Flow('shottime1','tfairblended',
        '''
        dd type=float
        | headermath output="cdpx" | window n2=201 |transp
        ''')

Flow('shottime2','tfairblended',
        '''
        dd type=float
        | headermath output="cdpx" | window f2=201 |transp
        ''')

#######################################################################
## mute data 1 
#######################################################################
## 
Flow('fairblended1mutedown','fairblended1','mutter half=n t0=1 slope0=0.0224')
Flow('fairblended1muteup','fairblended1','mutter half=n t0=0.4 slope0=0.016 | add scale=-1,1 ${SOURCES[0]}')
Flow('fairblended1mutedir','fairblended1mutedown fairblended1muteup fairblended1','add scale=-1,-1,1 ${SOURCES[1:3]}')

Flow('fairblended1mutedisp','fairblended1mutedown fairblended1mutedir fairblended1muteup','cat axis=2 ${SOURCES[1:3]}')
Flow('fairblended1mute','fairblended1 fairblended1mutedir','add scale=1,-1 ${SOURCES[1]}')
Flow('fairblended1mute-dir','fairblended1mute fairblended1mutedir','cat axis=2 ${SOURCES[1]}')

#######################################################################
## filtering using median filter and space-varying median filter
#######################################################################
## using median filter
Flow('fairdeblended2mf','fairblended2','transp | mf nfw=%d | transp'%nfw)
#Flow('fairdeblended1mf','fairblended1mute fairblended1mutedir','transp | mf nfw=%d | transp | add scale=1,1 ${SOURCES[1]}'%nfw)
Flow('fairdeblended1mf','fairblended1','transp | mf nfw=%d | transp '%nfw)

## using svmf filter
Flow('fairdeblended2tsmf fairL2-tsmf','fairblended2','transp | tsmf L=${TARGETS[1]} nfw=%d| transp'%nfw)

Flow('fairdeblended2tsmf2','fairblended2',
                        '''
                        transp | tsmf nfw=%d ne=2499 |
                        tsmf nfw=%d ns=2500 ne=3000 N=150000 l1=2 l2=0 l3=6 l4=8 |
                        tsmf nfw=%d ns=3001 |
                         transp'''%(nfw,nfw,nfw+4)) # more detailed MF

Flow('fairsimi2','fairblended2 fairdeblended2mf','similarity other=${SOURCES[1]} rect1=5 rect2=5 | transp | scale axis=2')
Flow('fairdeblended2svmf fairL2','fairblended2 fairsimi2','transp |svmf L=${TARGETS[1]} similarity=${SOURCES[1]} nfw=%d lambda1=0.50 lambda2=0.75 lambda3=0.85 lambda4=0.95| transp'%nfw)

Flow('fairdeblended2svmf2','fairblended2 fairsimi2',
                        '''
                        transp | svmf nfw=%d ne=2499 similarity=${SOURCES[1]} lambda1=0.50 lambda2=0.75 lambda3=0.85 lambda4=0.95|
                        svmf nfw=%d ns=2500 ne=3000 l1=8 l2=6 l3=6 l4=8 similarity=${SOURCES[1]} lambda1=0.50 lambda2=0.75 lambda3=0.85 lambda4=0.95|
                        svmf nfw=%d ns=3001 similarity=${SOURCES[1]} lambda1=0.50 lambda2=0.75 lambda3=0.85 lambda4=0.95 |
                         transp'''%(nfw,nfw,nfw+4)) # more detailed MF
Result('fairsimi2','transp | grey color=j label2=Trace unit2="" label1=Time unit1="s" title="Signal reliability" wherexlabel=b wanttitle=y scalebar=y allpos=y wheretitle=t')



#Flow('fairdeblended1svmf fairL1','fairblended1mute fairblended1mutedir','transp | tsmf  L=${TARGETS[1]} nfw=%d l3=2 l4=4 N=10000000 | transp | add scale=1,1 ${SOURCES[1]}'%nfw)
Flow('fairdeblended1svmf fairL1','fairblended1','transp | tsmf nfw=%d l3=6 l4=8 L=${TARGETS[1]} | transp'%nfw)

Flow('fairdeblendedmf','fairdeblended1mf fairdeblended2mf','cat axis=2 ${SOURCES[1]}')
Flow('fairdeblendedsvmf','fairdeblended1svmf fairdeblended2svmf','cat axis=2 ${SOURCES[1]}')
Flow('noisemf','fairblended fairdeblendedmf','add scale=1,-1 ${SOURCES[1]}')
Flow('noisesvmf','fairblended fairdeblendedsvmf','add scale=1,-1 ${SOURCES[1]}')
Flow('fairdeblended1comp','fairdeblended1mf fairdeblended1svmf','cat axis=2 ${SOURCES[1]}')
Flow('fairdeblended2comp','fairdeblended2mf fairdeblended2svmf','cat axis=2 ${SOURCES[1]}')
Flow('noise1comp','diff1mf diff1svmf','cat axis=2 ${SOURCES[1]}')
Flow('noise2comp','diff2mf diff2svmf2','cat axis=2 ${SOURCES[1]}')
## noise section using mf
Flow('diff2mf','fairblended2 fairdeblended2mf','add scale=1,-1 ${SOURCES[1]}')
Flow('diff1mf','fairblended1 fairdeblended1mf','add scale=1,-1 ${SOURCES[1]}')

## noise section using svmf
Flow('diff2svmf','fairblended2 fairdeblended2svmf','add scale=1,-1 ${SOURCES[1]}')
Flow('diff2svmf2','fairblended2 fairdeblended2svmf2','add scale=1,-1 ${SOURCES[1]}')
Flow('diff1svmf','fairblended1 fairdeblended1svmf','add scale=1,-1 ${SOURCES[1]}')

## zoomed part
Flow('diff2mfzoom','diff2mf','window n1=550 f1=2000')      #n1=550
Flow('diff2svmfzoom','diff2svmf','window n1=550 f1=2000') #n1=550
Flow('diff1mfzoom','diff1mf','window n1=1125 f1=500')
Flow('diff1svmfzoom','diff1svmf','window n1=1125 f1=500')


#######################################################################
## Ploting
#######################################################################
Grey('fairblended','scalebar=n clip=%f'%clip)
Grey('fairunblended','scalebar=n clip=%f'%clip)
Grey('fairblended1','scalebar=n clip=%f'%clip)
Grey('fairblended2','scalebar=n clip=%f'%clip)
Grey('fairunblended1','scalebar=n clip=%f'%clip)
Grey('fairunblended2','scalebar=n clip=%f'%clip)

Grey('fairdeblended2mf','scalebar=n clip=%f'%clip)
Grey('fairdeblended2svmf','scalebar=n clip=%f '%clip)
Grey('fairdeblended2svmf2','scalebar=n clip=%f '%clip)
Grey('diff2mf','scalebar=n clip=%f '%clip)
Grey('diff2svmf','scalebar=n clip=%f '%clip)
Grey('diff2svmf2','scalebar=n clip=%f '%clip)

Grey('fairdeblended1mf','scalebar=n clip=%f'%clip)
Grey('fairdeblended1svmf','scalebar=n clip=%f '%clip)
Grey('diff1mf','scalebar=n clip=%f '%clip)
Grey('diff1svmf','scalebar=n clip=%f '%clip)
Grey('fairblended1mutedisp','clip=%g '%clip)
Grey('fairblended1mute-dir','clip=%g '%clip)
Grey('fairblended1mute','clip=%g '%clip)
Grey('fairdeblendedmf','clip=%g'%clip)
Grey('fairdeblendedsvmf','clip=%g'%clip)
Grey('noisemf','clip=%g'%clip)
Grey('noisesvmf','cliplip=%f'%clip)
Grey('fairdeblended1comp','clip=%g'%clip)
Grey('fairdeblended2comp','clip=%g'%clip)
Grey('noise1comp','clip=%f'%clip)
Grey('noise2comp','clip=%f'%clip)

## zoomed part ## zoom window 8s->10.2s, 8/0.004=2000 -> 10.2/0.004=2550
Grey('diff2mfzoom','scalebar=n clip=%f '%clip)
Grey('diff2svmfzoom','scalebar=n clip=%f '%clip)
Grey('diff1mfzoom','scalebar=n clip=%f '%clip)
Grey('diff1svmfzoom','scalebar=n clip=%f '%clip)

## Creating framebox
x2=0
y2=8.0
w=200
w2=2.2
x1=0
y1=2.0
w=200
w1=4.5

## frame2
Flow('frame2.asc',None,'echo %s n1=10 data_format=ascii_float in=$TARGET'% \
        string.join(map(str,(x2,y2,x2+w,y2,x2+w,y2+w2,x2,y2+w2,x2,y2))))
Plot('frame2','frame2.asc',
        '''
        dd type=complex form=native |
        graph min1=0 max1=200 min2=0 max2=16 pad=n plotfat=10 plotcol=3 
        wantaxis=n wanttitle=n yreverse=y
        ''')
## frame1
Flow('frame1.asc',None,'echo %s n1=10 data_format=ascii_float in=$TARGET'% \
        string.join(map(str,(x1,y1,x1+w,y1,x1+w,y1+w1,x1,y1+w1,x1,y1))))
Plot('frame1','frame1.asc',
        '''
        dd type=complex form=native |
        graph min1=0 max1=200 min2=0 max2=16 pad=n plotfat=10 plotcol=3 
        wantaxis=n wanttitle=n yreverse=y
        ''')

Result('diff1mfframe','./Fig/diff1mf.vpl frame1','Overlay')
Result('diff1svmfframe','./Fig/diff1svmf.vpl frame1','Overlay')
Result('diff2mfframe','./Fig/diff2mf.vpl frame2','Overlay')
Result('diff2svmfframe','./Fig/diff2svmf.vpl frame2','Overlay')

## writing to segy file
Flow('fairdeblendedmf.segy','fairdeblendedmf tfairblended','segywrite tfile=${SOURCES[1]}')
Flow('fairdeblendedsvmf.segy','fairdeblendedsvmf tfairblended','segywrite tfile=${SOURCES[1]}')
Flow('blendingnoisemf.segy','noisemf tfairblended','segywrite tfile=${SOURCES[1]}')
Flow('blendingnoisesvmf.segy','noisesvmf tfairblended','segywrite tfile=${SOURCES[1]}')

Result('compppt','fairblended fairdeblendedmf fairdeblendedsvmf','cat axis=2 ${SOURCES[1:3]} | grey label2=Trace unit2="" label1=Time unit1="s" title="" wherexlabel=b clip=0.5 wanttitle=n ')

Flow('fkfairunblended2','fairunblended2','cut f1=2501 | put d2=1 | fft1 | fft3 axis=2 pad=1')
Flow('fkfairblended2','fairblended2','cut f1=2501 |put d2=1 | fft1 | fft3 axis=2 pad1=1')
Flow('fkfairdeblended2mf','fairdeblended2mf','cut f1=2501 |put d2=1 | fft1 | fft3 axis=2 pad1=1')
Flow('fkfairdeblended2svmf','fairdeblended2svmf','cut f1=2501 | put d2=1 | fft1 | fft3 axis=2 pad1=1 ')
Result('fkfairunblended2','cabs | grey title="" label2="Normalized wavenumber" unit2="" color=j ')
Result('fkfairblended2','cabs | grey title="" label2="Normalized wavenumber" unit2="" color=j ')
Result('fkfairdeblended2mf','cabs | grey title="" label2="Normalized wavenumber" unit2="" color=j ')
Result('fkfairdeblended2svmf','cabs | grey title="" label2="Normalized wavenumber" unit2="" color=j ')

Result('fairL1','transp | grey label2=Trace unit2="" label1=Time unit1="s" title="" wherexlabel=b wanttitle=n scalebar=y color=j')
Result('fairL2','transp | grey label2=Trace unit2="" label1=Time unit1="s" title="Filter length" wherexlabel=b wanttitle=y wheretitle=t scalebar=y color=j allpos=y')

End()

sfsegyread
sfwindow
sfput
sfdd
sfheadermath
sftransp
sfmutter
sfadd
sfcat
sfmf
sftsmf
sfsimilarity
sfscale
sfsvmf
sfgrey
sfgraph
sfsegywrite
sfcut
sffft1
sffft3
sfcabs