from rsf.proj import *
import string

csv = 'qsiwell2.csv'

Fetch(csv,'1504_Seismic_petrophysics_1',
      server='https://raw.githubusercontent.com',
      top='seg/tutorials/master')

Flow('well',csv,'csv2rsf header=y')

# Examine with < well.rsf sfheaderattr segy=n

# Limit depth

Flow('depthmask','well','window n1=1 | mask min=2100 max=2400')
Flow('logs','well depthmask','headerwindow mask=${SOURCES[1]}')

# Calculate litho-facies log

for case in ('VSH','SW','PHI','IP','VPVS','DEPTH'):
    Flow(case,'logs','headermath output=%s segy=n' % case)

Flow('bsmask1','VSH','mask max=0.2')
Flow('bsmask2','SW', 'mask min=0.9')
Flow('brine-sand-mask','bsmask1 bsmask2','mul ${SOURCES[1]}')

Flow('oil-sand-mask','bsmask2 bsmask1','add add=-1 scale=-1 | mul ${SOURCES[1]}')

Flow('shale-mask','brine-sand-mask oil-sand-mask','add ${SOURCES[1]} | add add=-1 scale=-1')

hists = []
for case in ('brine-sand','oil-sand','shale'):
    Flow(case,['logs',case+'-mask'],'headerwindow mask=${SOURCES[1]}')
    hist = case+'-hist'
    Plot(hist,case,
         '''
         window n1=1 f1=10 | histogram n1=51 o1=4200 d1=80 | dd type=float | 
         graph plotcol=1 title="%s" label1=IP unit1="m/s*g/cc" label2=Occurence unit2=
         ''' % ' '.join(map(string.capitalize,case.split('-'))))
    hists.append(hist)

Result('histogram',hists,'SideBySideIso')

for case in ('VSH','SW','PHI','IP','VPVS'):
    graph = ''
    if case in ('VSH','SW','PHI'):
        graph += ' min2=-0.1 max2=1.1 label2=Vcl/phi/Sw unit2='
    if case == 'VSH':
        graph += ' plotcol=4'
    if case == 'SW':
        graph += ' plotcol=1 wantaxis=n wanttitle=n'
    if case == 'PHI':
        graph += ' plotcol=5 wantaxis=n wanttitle=n'
    if case == 'IP':
        graph += ' label2=Ip unit2="m/s*g/cc" '
    if case == 'VPVS':
        graph += ' label2=Vp/Vs unit2='
    Plot(case,['DEPTH',case],
         '''
         cmplx ${SOURCES[1]} | window | 
         graph transp=y yreverse=y label1=Depth unit1= wanttitle=n plotcol=7 labelsz=12
         ''' + graph)

Plot('three','VSH SW PHI','Overlay')

# Litho-Facies Log

Flow('LFC','brine-sand-mask oil-sand-mask shale-mask',
     'add scale=1,2,4 ${SOURCES[1:3]} | dd type=float | window')

# Make a special colormap

blue =   [0,0,1]
green =  [0,1,0]
red =    [1,0,0]
brown =  [0.6,0.4,0.2]

Flow('lfc.asc',None,
     '''
     echo %s n1=3 n2=256 in=$TARGET data_format=ascii_float
     ''' % ' '.join(map(str,blue*64+green*64+red*64+brown*64)))
Flow('lfc.csv','lfc.asc','disfil col=3 format="%g," number=n | sed s/,$//')

Plot('LFC','LFC lfc.csv',
     '''
     grey color=${SOURCES[1]} pclip=100 allpos=y bias=1 
     wantaxis=n title=LFC wherexlabel=b titlesz=15
     ''')

Result('four','three IP VPVS LFC','SideBySideAniso')

Flow('lfc-depth','LFC','byte pclip=100 allpos=y bias=1')

Result('crossplot','IP VPVS lfc-depth lfc.csv',
       '''
       cmplx ${SOURCES[1]} | window |
       graph symbol=* depth=${SOURCES[2]} color=${SOURCES[3]} grid=y 
       label1=Ip unit1="m/s*g/cc" label2=Vp/Vs unit2= wanttitle=n screenratio=1
       ''')


End()