include $(SEPINC)/SEP.top


F90LIBS =-lsupersetf90 -lsepgeef90 -lsuperset -lsepauxf90 -lsepmathf90 -lsep2df90 -lsep3df90 -lsep3d -lsepf90 -lsep


RESULTSER =
RESULTCR = 
RESULTSNR = 

#UF90FLAGS = -O0 -g -traceback -fpe:0 -check all
UF90FLAGS = -O3

BINDIR = ./bin
SRCDIR = ./src
RESDIR = ./Fig
OBJDIR = ./obj

R = ${RESDIR}
B = ${BINDIR}
S = ${SRCDIR}
O = ${OBJDIR}

HTMLDIR = .

default1: $R/latSpike.v $R/data1.v $R/mod1.v $R/filt2anew.v $R/filt2bnew.v $R/filt2aold.v $R/filt2bold.v $R/data2.v $R/mod2new.v $R/mod2old.v $R/data2.v $R/mod2new.v $R/mod2old.v $R/syn3.v $R/data3.v $R/mod3old.v $R/mod3new.v $R/filt3anew.v $R/filt3bnew.v $R/filt3aold.v $R/filt3bold.v $R/data4.v $R/mod4new.v $R/mod4old.v $R/filt4anew.v $R/filt4bnew.v $R/filt4aold.v $R/filt4bold.v $R/wavelet2old.v $R/wavelet2new.v $R/wavelet4new.v $R/wavelet4old.v $R/mod4newwavelet.v $R/mod4oldwavelet.v pdf_file



default:
	make -k -i default1; make -k -i default1

pdf_file:
	make -f pdf_makefile

${R}/%.pdf: ${R}/%.ps
	ps2pdf \
	-dDownsampleColorImages=false \
	-dDownsampleGrayImages=false \
	-dDownsampleMonoImages=false \
	-dColorImageFilter=/FlateEncode \
	-dGrayImageFilter=/FlateEncode \
	-dMonoImageFilter=/CCITTFaxEncode \
	-dAutoFilterColorImages=false \
	-dAutoFilterGrayImages=false \
	-dAutoFilterMonoImages=false \
	-dEPSCrop ${R}/$*.ps $@

dn = >/dev/null
grey = Grey $(dn)
FAT = 3
DATDIR = ./data

############################## alternating mixed-phase Blind Deconvolution ####################

clipdec = newclip=1 crowd=0.80 #wantscalebar=1
DECON_LABEL = label1='time(s)' label2='cmp_x(km)' label3=' ' title='' labelfat=1 labelsz=8

GRAPH_LABEL = label1='time(s)' label2=' ' label3=' ' title=' ' labelfat=1 labelstube z=4 screenratio=0.3 screenwd=12.0 screenht=3.6

data-decon = data/decon#
buf-trc = buf-trc#


#Very simple case, invert single min-phase wavelet
latSpike.H:
	Spike n1=256 k1=3 nsp=1 > $@

$R/latSpike.v:latSpike.H
	Graph< $<  out=$@ >/dev/null

############################## Decon on commmon offset data
COF_DECON_LABEL := label1='time(s)' label2='cmp_x(kft)' title='  ' labelfat=1 labelsz=9


##############Sythetic 1########
data-sythetic_1.H:
	Spike n1=50 k1=24,25,26 mag=3,7,2 nsp=3 >$@



mod-sythetic_1-decon_linear.H: $B/Sparse_decon_mixed_linear.x data-sythetic_1.H
	<data-sythetic_1.H $< max_overshot_shrink_cnt=10 \
		norm_dat=hybrid qt_period=10 filt_length=20\
		r_th_d=0.1  filter_a=linear_a.H filter_b=linear_b.H \
		niter_nl=1000 niter_lin=2 model_b=modtrc_b.H \
		wavelet=latSpike.H conv_res=wztrc_a_b_rev_conv.H >result.H
	< result.H Window f3=999 n3=1 >$@


$R/data1.v: data-sythetic_1.H 
	< $< Dots title=' ' label1=' ' label2=' ' out=$@ >/dev/null
$R/mod1.v: mod-sythetic_1-decon_linear.H
	< $< Dots title=' ' label1=' ' label2=' ' out=$@ >/dev/null



###==============Synthetic 2D data
###============sym wavelet
wavelet-sym.H: $B/Sparse_decon.x
	Spike nsp=1 k1=1 n1=128 > l.H 	
	$< <l.H filt=${data-decon}/wavelet-min.HH inv=0 > tmp.H
	< ${data-decon}/wavelet-min.HH Reverse > tmp1.H
	$< < tmp.H filt=tmp1.H inv=0 > $@

mod-symwavlet.H fita-symwavlet.H fitb-symwavlet.H:latSpike.H $B/Sparse_decon_mixed_nl.x $B/Sparse_decon.x wavelet-sym.H
	< $< $B/Sparse_decon.x inv=0 filt=wavelet-sym.H > data.H
	< data.H bin/Sparse_decon_mixed_nl.x  max_overshot_shrink_cnt=5 \
		norm_dat=hybrid qt_period=10 filt_length=32 \
		qnt_d=0.95  filter_a=fita_symwavlet.H filter_b=fitb_symwavlet.H \
		niter_nl=10 niter_lin=500 >dump.H model_b=mod_b_12.H \
		wavelet=latSpike.H conv_res=convRes12.H
	< dump.H Window Window f1=31 f3=9 n3=1 | Pad end1=31>mod-symwavlet.H
	< fita_symwavlet.H Window n3=1 f3=9 >fita-symwavlet.H
	< fitb_symwavlet.H Window n3=1 f3=9 >fitb-symwavlet.H

mod-symwavlet_linear.H fita-symwavlet_linear.H fitb-symwavlet_linear.H:latSpike.H $B/Sparse_decon_mixed_linear.x $B/Sparse_decon.x wavelet-sym.H
	< $< $B/Sparse_decon.x inv=0 filt=wavelet-sym.H > data.H
	< data.H bin/Sparse_decon_mixed_linear.x  max_overshot_shrink_cnt=5 \
		norm_dat=hybrid qt_period=10 filt_length=32 \
		qnt_d=0.95  filter_a=fita_symwavlet_linear.H filter_b=fitb_symwavlet_linear.H \
		niter_nl=1000 niter_lin=2 >dump.H model_b=mod_b_12.H \
		wavelet=latSpike.H conv_res=convRes12.H
	< dump.H Window n3=1 f3=999 >mod-symwavlet_linear.H
	< fita_symwavlet_linear.H Window n3=1 f3=999 >fita-symwavlet_linear.H
	< fitb_symwavlet_linear.H Window n3=1 f3=999 >fitb-symwavlet_linear.H


$R/filt2anew.v:  fita-symwavlet_linear.H
	< $< Graph  title=" " label2="Amp " label1="time(sec)" wherexlabel=t screenwd=13.5 screenht=3.8 out=$@ ${dn}
$R/filt2bnew.v:fitb-symwavlet_linear.H
	cp fitb-symwavlet_linear.H fitb-symwavlet_linear_cp.H
	echo o1=-0.124>> fitb-symwavlet_linear_cp.H
	Reverse <fitb-symwavlet_linear_cp.H which=01 opt=n |Graph  title=" " label2="Amp " label1="time(sec)" wherexlabel=t screenwd=13.5 screenht=3.8  out=$@ ${dn}

$R/filt2aold.v:  fita-symwavlet.H
	< $< Graph  title=" " label2="Amp " label1="time(sec)" wherexlabel=t screenwd=13.5 screenht=3.8 out=$@ ${dn}
$R/filt2bold.v:  fitb-symwavlet.H
	< $< Graph  title=" " label2="Amp " label1="time(sec)" wherexlabel=t screenwd=13.5 screenht=3.8  out=$@ ${dn}
$R/filt2.v:$R/filt2anew.v $R/filt2bnew.v $R/filt2aold.v $R/filt2bold.v
	vp_OverUnderAniso $R/filt2anew.v $R/filt2bnew.v $R/filt2aold.v $R/filt2bold.v>$@

$R/data2.v:  wavelet-sym.H
	< $< Graph label2="Amp " title=" " label2="Amp " label1="time(sec)" wherexlabel=t screenwd=13.5 screenht=3.8 min1=0 max1=0.5 out=$@ ${dn}

$R/mod2new.v:mod-symwavlet_linear.H
	Pad < $< end1=2|Window f1=2 >mod-symwavlet_linear1.H
	echo o1=0 >>mod-symwavlet_linear1.H
	< mod-symwavlet_linear1.H Graph title=" " label2="Amp " label1="time(sec)" wherexlabel=t screenwd=13.5 screenht=3.8 min1=0 max1=0.5 out=$@ ${dn}
$R/mod2old.v:mod-symwavlet.H
	cp mod-symwavlet.H mod-symwavlet_cp.H
	echo o1=0 >> mod-symwavlet_cp.H
	< mod-symwavlet_cp.H Graph title=" " label2="Amp " label1="time(sec)" wherexlabel=t screenwd=13.5 screenht=3.8 min1=0  max1=0.5 out=$@ ${dn}

wavelet_symwavelet.H: $B/bubble.x fita-symwavlet.H fitb-symwavlet.H
	Reverse< fitb-symwavlet.H which=1 opt=ignore >fitb-symwavlet_reverse.H
	src/Conv_c/bin/Conv_filter_c.x < fita-symwavlet.H >filter_ab_symwavelet.H\
		filter=fitb-symwavlet_reverse.H
	<filter_ab_symwavelet.H Pad beg1=0 n1out=128 | $< th=0 >$@

wavelet_symwavelet_linear.H: $B/bubble.x fita-symwavlet_linear.H fitb-symwavlet_linear.H
	 src/Conv_c/bin/Conv_filter_c.x < fita-symwavlet_linear.H >filter_ab_symwavelet_linear.H\
		filter=fitb-symwavlet_linear.H
	<filter_ab_symwavelet_linear.H Pad beg1=0 n1out=128 | $< th=0.0  >$@
$R/wavelet2old.v:wavelet_symwavelet.H src/shift/bin/shift.x
	src/shift/bin/shift.x < wavelet_symwavelet.H > wavelet_symwavelet_cp.H
	echo d1=0.004 >> wavelet_symwavelet_cp.H
	Graph <wavelet_symwavelet_cp.H title=" " label2="Amp " label1="time(sec)" wherexlabel=t screenwd=13.5 screenht=3.8 out=$@ ${dn}

$R/wavelet2new.v:wavelet_symwavelet_linear.H src/shift/bin/shift.x
	src/shift/bin/shift.x < wavelet_symwavelet_linear.H > wavelet_symwavelet_linear_cp.H
	echo d1=0.004 >> wavelet_symwavelet_linear_cp.H
	Graph <wavelet_symwavelet_linear_cp.H title=" " label2="Amp " label1="time(sec)" wherexlabel=t screenwd=13.5 screenht=3.8 out=$@ ${dn}

###Reflection======================================
refl-jon.H: ${data-decon}/refl-jon.HH
	Interp <$< d2out=0.064 type=0 | Add scale=100 > $@


data-syn-sym-2d.H: wavelet-sym.H $B/Sparse_decon.x refl-jon.H
	< refl-jon.H Pad end1=100 | $B/Sparse_decon.x inv=0 filt=$< >$@

     #Convolve Sym wavele
mod-sym-2d.H:$B/Sparse_decon_mixed_nl.x data-syn-sym-2d.H latSpike.H
	<data-syn-sym-2d.H $< max_overshot_shrink_cnt=5 \
		norm_dat=hybrid qt_period=10 filt_length=48 \
		r_th_d=0.1  filter_a=fit_a_sym.H filter_b=fit_b_sym.H \
		niter_nl=20 niter_lin=100 model_b=mod_b_sym.H \
		wavelet=latSpike.H conv_res=convRes_sym.H >mod-sym-2d-temp.H
	<mod-sym-2d-temp.H Window f1=47 f3=19 n3=1 | Pad end1=47 >mod-sym-2d.H



mod-sym-2d-decon_linear.H: $B/Sparse_decon_mixed_linear.x data-syn-sym-2d.H latSpike.H
	<data-syn-sym-2d.H $< max_overshot_shrink_cnt=5 \
		norm_dat=hybrid qt_period=10 filt_length=48\
		r_th_d=0.1 filter_a=fit_a_sym_linear.H filter_b=fit_b_sym_linear.H \
		niter_nl=1000 niter_lin=2 model_b=modtrc_b.H \
		wavelet=latSpike.H conv_res=wztrc_a_b_rev_conv.H >result33_linear.H
	< result33_linear.H  Window f3=999 n3=1 >$@ 

$R/syn3.v: refl-jon.H
	< $< Grey ${clipdec} ${DECON_LABEL} title=" " out=$@ ${dn}
$R/data3.v: data-syn-sym-2d.H
	< $< Grey ${clipdec} ${DECON_LABEL} title=" " out=$@ ${dn}

$R/mod3old.v: data/mod-sym-2d.H
	cp data/mod-sym-2d.H mod-sym-2d_cp.H
	echo o1=0 >> mod-sym-2d_cp.H
	< mod-sym-2d_cp.H Grey ${clipdec} ${DECON_LABEL} title="  " out=$@ ${dn}
$R/mod3new.v:mod-sym-2d-decon_linear.H
	< $< Grey ${clipdec} ${DECON_LABEL} title=" " out=$@ ${dn}


#=============================filters
fit_a_sym_linear_reshape.H:fit_a_sym_linear.H
	Reshape < $< reshape=1,3 |Window f2=999 n2=1> $@

fit_b_sym_linear_reshape.H:fit_b_sym_linear.H
	Reshape < $< reshape=1,3 |Window f2=999 n2=1> $@

$R/filt3anew.v: fit_a_sym_linear.H
	< $< Window f3=999 n3=1 |Graph label2="" title="  " out=$@ ${dn}
$R/filt3bnew.v:fit_b_sym_linear.H
	< $< Window f3=999 n3=1>fit_b_sym_linear_window.H
	Reverse <fit_b_sym_linear_window.H which=01 opt=ignore |Graph label2="" title=" " out=$@ ${dn}

fit_a_sym_reshape.H:fit_a_sym.H
	Reshape < $< reshape=1,3 |Window f2=19 n2=1> $@
fit_b_sym_reshape.H:fit_b_sym.H
	Reverse  < $<  which=01 opt=ignore |Reshape reshape=1,3 |Window f2=19 n2=1> $@

$R/filt3aold.v: fit_a_sym.H
	< $< Window f3=19 n3=1 |Graph label2="" title=" " out=$@ ${dn}
$R/filt3bold.v:fit_b_sym.H
	< $< Window f3=19 n3=1|Graph label2="" title=" " out=$@ ${dn}


#=============================



#########Real Data

data-COFnew.H: data/decon/data-COFnew.HH
	< $< Window min2=7200 max2=20000 > junk.H
	< junk.H Add scale=0.00125 > $@
## wavelet_init item Only provide the length of the desired wavelet


FILTLEN=32
THRESHOLD=0.1


mod-COFnew-decon.H: $B/Sparse_decon_mixed_nl.x data-COFnew.H latSpike.H
	<data-COFnew.H $< max_overshot_shrink_cnt=10 \
		norm_dat=hybrid qt_period=10 filt_length=32 \
		r_th_d=0.1  filter_a=fit_a_COF.H filter_b=fit_b_COF.H \
		niter_nl=8 niter_lin=100 model_b=modtrc_b.H \
		wavelet=latSpike.H conv_res=wztrc_a_b_rev_conv.H >result_1.H
	< result_1.H Window f1=31 f3=7 n3=1 | Pad end1=31 > mod-COFnew-decon.H


mod-COFnew-decon_linear.H: $B/Sparse_decon_mixed_linear.x data-COFnew.H latSpike.H
	<data-COFnew.H $< max_overshot_shrink_cnt=10 \
		norm_dat=hybrid qt_period=10 filt_length=32 \
		r_th_d=0.1  filter_a=fit_a_COF_linear.H filter_b=fit_b_COF_linear.H \
		niter_nl=300 niter_lin=2 model_b=modtrc_b.H \
		wavelet=latSpike.H conv_res=wztrc_a_b_rev_conv.H >result_2.H
	< result_2.H Window f1=0 f3=299 n3=1 >$@

mod-COFnew-decon_wavelet.H: $B/Sparse_decon_mixed_nl.x data-COFnew.H latSpike.H
	<data-COFnew.H $< max_overshot_shrink_cnt=10 \
		norm_dat=hybrid qt_period=10 filt_length=100 \
		r_th_d=0.1  filter_a=fit_a_COF_wavelet.H filter_b=fit_b_COF_wavelet.H \
		niter_nl=8 niter_lin=100 model_b=modtrc_b.H \
		wavelet=latSpike.H conv_res=wztrc_a_b_rev_conv.H >result_3.H
	< result_3.H Window f1=31 f3=7 n3=1 | Pad end1=31 > mod-COFnew-decon_wavelet.H


mod-COFnew-decon_linear_wavelet.H: $B/Sparse_decon_mixed_linear.x data-COFnew.H latSpike.H
	<data-COFnew.H $< max_overshot_shrink_cnt=10 \
		norm_dat=hybrid qt_period=10 filt_length=100 \
		r_th_d=0.1  filter_a=fit_a_COF_linear_wavelet.H filter_b=fit_b_COF_linear_wavelet.H \
		niter_nl=300 niter_lin=2 model_b=modtrc_b.H \
		wavelet=latSpike.H conv_res=wztrc_a_b_rev_conv.H >result_4.H
	< result_4.H Window f1=0 f3=299 n3=1 >$@


wavelet_COF.H: $B/bubble.x fit_a_COF_wavelet.H fit_b_COF_wavelet.H
	Reverse< fit_b_COF_wavelet.H which=1 opt=ignore|Window f3=7 n3=1 >fit_b_COF_wavelet_reverse.H
	<fit_a_COF_wavelet.H Window f3=7 n3=1 >fit_a_COF_wavelet_window.H
	src/Conv_c/bin/Conv_filter_c.x < fit_a_COF_wavelet_window.H >filter_ab_COF_wavelet.H\
		filter=fit_b_COF_wavelet_reverse.H
	<filter_ab_COF_wavelet.H Pad beg1=0 n1out=256 | $< th=0 >$@

wavelet_COF_linear.H: $B/bubble.x fit_a_COF_linear_wavelet.H fit_b_COF_linear_wavelet.H
	<fit_b_COF_linear_wavelet.H Window f3=299 n3=1 >fit_b_COF_linear_wavelet_window.H
	<fit_a_COF_linear_wavelet.H Window f3=299 n3=1 >fit_a_COF_linear_wavelet_window.H
	src/Conv_c/bin/Conv_filter_c.x < fit_a_COF_linear_wavelet_window.H >filter_ab_COF_linear_wavelet.H\
		filter=fit_b_COF_linear_wavelet_window.H
	<filter_ab_COF_linear_wavelet.H Pad beg1=0 n1out=256 | $< th=0 >$@

$R/wavelet4old.v:wavelet_COF.H src/shift/bin/shift.x
	src/shift/bin/shift.x < $< >wavelet_COF_cp.H
	echo d1=0.004 >> wavelet_COF_cp.H
	<wavelet_COF_cp.H Graph   title=" " label2="Amp " label1="time(sec)" wherexlabel=t screenwd=13.5 screenht=3 max2=2 min2=-2 out=$@ ${dn}

$R/wavelet4new.v:wavelet_COF_linear.H src/shift/bin/shift.x
	src/shift/bin/shift.x < $< >wavelet_COF_linear_cp.H
	echo d1=0.004 >> wavelet_COF_linear_cp.H
	<wavelet_COF_linear_cp.H Graph   title=" " label2="Amp " label1="time(sec)" wherexlabel=t screenwd=13.5 screenht=3 max2=2 min2=-2 out=$@ ${dn}


bc = -1
ec = 1
#ebclip = newclip=1 bclip=${bc} eclip=${ec}
ebclip = newclip=1 bpclip=1 epclip=99




$R/data4.v : data-COFnew.H
	echo o1=1.60 o2=7.2 d2=0.02667 >> $<
	< $< Add scale=-1|  Grey  ${ebclip} ${COF_DECON_LABEL}  title=" "  min1=1.6 max1=2.7 out=$@ ${dn}
$R/mod4new.v: mod-COFnew-decon_linear.H
#	echo o1=1.60 o2=7.2 d2=0.02667 >> $<
	< $< Add scale=-1| Grey title=" " ${ebclip} ${COF_DECON_LABEL}  title=" " min1=1.6 out=$@ ${dn}
$R/mod4old.v: mod-COFnew-decon.H
	echo o1=1.60 o2=7.2 d2=0.02667 >> $<
	< $< Add scale=-1| Grey  ${ebclip} ${COF_DECON_LABEL} title="  "  min1=1.6  out=$@ ${dn}

$R/filt4aold.v: fit_a_COF.H
	< $< Window f3=7 n3=1 |Graph  title=" " label2="Amp " label1="time(sec)" wherexlabel=t screenwd=13.5 screenht=3.8 out=$@ ${dn}
$R/filt4bold.v:fit_b_COF.H
	< $< Window f3=7 n3=1|Graph title=" " label2="Amp " label1="time(sec)" wherexlabel=t screenwd=13.5 screenht=3.8 out=$@ ${dn}
$R/filt4anew.v: fit_a_COF_linear.H
	< $< Window f3=299 n3=1 |Graph title=" " label2="Amp " label1="time(sec)" wherexlabel=t screenwd=13.5 screenht=4 out=$@ ${dn}
$R/filt4bnew.v:fit_b_COF_linear.H
	< $< Window f3=299 n3=1|Reverse which=01 opt=ignore |Graph title=" " label2="Amp " label1="time(sec)" wherexlabel=t screenwd=13.5 screenht=3.9 out=$@ ${dn}

$R/mod4newwavelet.v:mod-COFnew-decon_linear_wavelet.H
#	echo o1=1.60 o2=7.2 d2=0.02667 >> $<
	< $< Add scale=-1| Grey title=" " ${ebclip} ${COF_DECON_LABEL}  title=" " min1=1.6 out=$@ ${dn}
$R/mod4oldwavelet.v: mod-COFnew-decon_wavelet.H
	echo o1=1.60 o2=7.2 d2=0.02667 >> $<
	< $< Add scale=-1| Grey  ${ebclip} ${COF_DECON_LABEL} title="  "  min1=1.6  out=$@ ${dn}


.PHONY: view

view:


cleanfig:
	-rm -f fig/*.v
	-rm -f fig/*.ps
	
clean: jclean
	-rm ${BINDIR}/*.* -f
	-rm ${OBJDIR}/*.* -f
	-rm *.mod *.o *.H -f
	-rm $R/*.v  -f
	-rm $R/*.ps -f
	-rm $R/*.pdf -f
#	@(cd $(srcdir);${MAKE} clean)

buildbin:  #Re build the binaries
	-rm ${BINDIR}/*.* -f
	-rm ${OBJDIR}/*.* -f
	-rm *.mod

#Type 'make' to build all the figures
#And then type 'scons' to build the report. 

include $(SEPINC)/SEP.bottom


