# Author: G. Hennenfent
#         Seismic Laboratory for Imaging and Modeling
#         Department of Earch & Ocean Sciences
#         The University of British Columbia
#         
# Date  : August, 07

from rsfproj import *
import os

# define customized plotting function
def mygrey(title,label1="Time",label2="Offset",unit1="s",unit2="m"):
    return '''
    sfgrey title="%s" clip=15
    label1="%s" unit1="%s" label2="%s" unit2="%s" crowd1=.5
    '''%(title,label1,unit1,label2,unit2)
#crowd1=.5 crowd2=.75
def myfkgrey(title,label1="Frequency",label2="Wavenumber",unit1="Hz",unit2="1/m"):
    return '''
    sfgrey title="%s" clip=15
    label1="%s" unit1="%s" label2="%s" unit2="%s"
    '''%(title,label1,unit1,label2,unit2)
#crowd1=.5 crowd2=.25

# define SLIM FTP server information
FTPserver = {
    'server': os.environ.get('SLIM_SERVER'),
    'login': os.environ.get('SLIM_LOGIN'),
    'password':os.environ.get('SLIM_PASSWORD')}

dtype = 'real'
src = 'Chevron'

loc = os.path.join(dtype,src)
dataset = 'L1318P1_SP868C1_LACONA_KAA625_Gilles.rsf'

# fetch data from FTP server
Fetch('L1318P1_SP868C1_LACONA_KAA625_Gilles.rsf',loc,FTPserver)

Flow('model',dataset,
     '''sfdd form=native | sfwindow f2=539 n2=255 min1=2 | sfcostaper
     nw1=30 nw2=30 | sfbandpass fhi=125 | sfwindow j1=2''')
Result('model',mygrey("") )

Flow('fkmodel','model','sfspectra2 | sfwindow')
Result('fkmodel',myfkgrey("",label1="Frequency",
                          label2="Wavenumber",unit1="Hz",unit2="1/m") )

jump = 3
name = 'rdata_subREG'
fkname = 'fkrdata_subREG'
Flow([name,'mask_subREG2D'],'model',
     'sfwindow j2=%d | sflpad mask=${TARGETS[1]} jump=%d '%(jump,jump) ) #| sfwindow n2=481
mask = 'mask_subREG'
Flow(mask,'mask_subREG2D',
     'sftransp | sfwindow n2=1 | sfput n3=1 | sfdd type=int')

Flow(fkname,name,'sfspectra2 | sfwindow')
Result(name,mygrey("") )
Result(fkname,myfkgrey("",
                     label1="Frequency",label2="Wavenumber",
                     unit1="Hz",unit2="1/m") )
name2 = 'rdata_subIRREG'
fkname2 = 'fkrdata_subIRREG'
mask2 = 'mask_subIRREG'
Flow(mask2,'model',
     'sfjitter perc=%f seed=1982 | sfdd type=int'%(1-1./jump) )
Flow(name2,['model',mask2],
     'sfheadercut mask=${SOURCES[1]}')
Flow(fkname2,name2,'sfspectra2 | sfwindow')
Result(name2,mygrey("") )
Result(fkname2,myfkgrey("",
                        label1="Frequency",label2="Wavenumber",
                        unit1="Hz",unit2="1/m") )

# INTERPOLATION BY CRSI
# curvelet transform parameters
nbs = 6         # nb of scales
nba = 64        # nb of angles at 2nd coarsest scale
ac = 1          # all curvelets
# solver parameters
nbout = 50      # nb of optimization subproblems (SPs)
nbin = 5        # nb of itr used to approx. solution of SP
# sparsity constraint parameters
thrmax = .01    # perc of curvelet coeff. kept for 1st SP
thrmin = .99    # perc of curvelet coeff. kept for last SP
# angular constraint parameters
angmaxL = 90    # max dip (degres) for negative wave numbers
angmaxR = 90    # max dip (degres) for positive wave numbers

# look for interpolation.py in your path (see SLIMpy installation
# instructions)
interp = WhereIs('interpolation.py')

for mode in ['REG','IRREG']:
    # interpolate
    res = 'rcrsi_sub%(mode)s'%vars()
    Flow(res,['rdata_sub%(mode)s'%vars(),'mask_sub%(mode)s'%vars(),interp],
         interp + '''
         -q mask=${SOURCES[1]} data=$SOURCE output=$TARGET
         solverparams=%(nbout)d,%(nbin)d transparams=%(nbs)d,%(nba)d,%(ac)d
         thrparams=%(thrmax)f,%(thrmin)f angconst=%(angmaxL)f,%(angmaxR)f
         '''%vars(),
         stdin=0,stdout=-1)
    Result(res,mygrey(""))

    diff = 'df'+res
    Flow(diff,['model',res],'sfmath x=${SOURCES[1]} output="4*(x-input)" ')
    Result(diff,mygrey(""))

    fkres = 'fk'+res
    Flow(fkres,res,'sfspectra2 | sfwindow')
    Result(fkres,myfkgrey("",
                          label1="Frequency",label2="Wavenumber",
                          unit1="Hz",unit2="1/m") )

End()