#define FIGDIR ../Fig
#define BINDIR ../../Bin/MTYPE
#include <SEP.defs>
#define NULL head=/dev/null
#define FIGLIST derivatives waves2D
#define GRAPHMAXMIN min1=0 min2=0 max1=3.1415 max2=3.1415 labelsz=13 titlesz=16 fat=4
#define LABEL     label1=" "                    label2="Normalized Amplitude" 
#define LABELBOT  label1="Normalized Frequency" label2="Normalized Amplitude" 
#define TITLE     title=" "
#define TITLETOP  title="Spectrum"
#define PSTEXOPTS txscale=2

FIGDIR/derivatives.v: centered.v stagg.v centlim.v holberg.v
	vp_OverUnderAniso centered.v stagg.v centlim.v holberg.v >FIGDIR/derivatives.v

%1.dots.%2.v: %1.H
	 Dots < %1.H  title="%2" strings=1 dots=2 out=%1.dots.%2.v NULL labelsz=15

centered.spec.Centered.v: centered.H
	Pad <centered.H n1=128 >jpad.H
	Spectrum  <jpad.H >jspec.H
	echo "0 0 3.1425 3.1415" >jj.ascii
	atoF <jj.ascii >janal.H
	echo "esize=8 n1=2" >> janal.H
	Graph <janal.H out=janal.v	 NULL GRAPHMAXMIN title=" "
	Graph <jspec.H out=jcentered.spec.Centered.v    NULL GRAPHMAXMIN TITLETOP  LABEL
	vp_Overlay janal.v jcentered.spec.Centered.v >centered.spec.Centered.v
	
holberg.spec.Holberg.v: holberg.H
	Pad <holberg.H n1=128 >jpad.H
	Spectrum  <jpad.H >jspec.H
	echo "0 0 3.1425 3.1415" >jj.ascii
	atoF <jj.ascii >janal.H
	echo "esize=8 n1=2" >> janal.H
	Graph <janal.H out=janal.v	 NULL GRAPHMAXMIN title=" "
	Graph <jspec.H out=jholberg.spec.Holberg.v    NULL GRAPHMAXMIN TITLE LABELBOT
	vp_Overlay janal.v jholberg.spec.Holberg.v >holberg.spec.Holberg.v
	
%1.spec.%2.v: %1.H
	Pad <%1.H n1=128 >jpad.H
	Spectrum  <jpad.H >jspec.H
	echo "0 0 3.1425 3.1415" >jj.ascii
	atoF <jj.ascii >janal.H
	echo "esize=8 n1=2" >> janal.H
	Graph <janal.H out=janal.v	 NULL GRAPHMAXMIN title=" "
	Graph <jspec.H out=j%1.spec.%2.v    NULL GRAPHMAXMIN TITLE LABEL
	vp_Overlay janal.v j%1.spec.%2.v >%1.spec.%2.v
	
centered.v: centered.dots.Centered.v centered.spec.Centered.v
	vp_SideBySideAniso centered.dots.Centered.v centered.spec.Centered.v > centered.v
	
stagg.v: stagg.dots.Staggered.v stagg.spec.Staggered.v
	vp_SideBySideAniso stagg.dots.Staggered.v stagg.spec.Staggered.v > stagg.v

holberg.v: holberg.dots.Holberg.v holberg.spec.Holberg.v
	vp_SideBySideAniso holberg.dots.Holberg.v holberg.spec.Holberg.v > holberg.v
	
centlim.v: centlim.dots.Central_Limit.v centlim.spec.Central_Limit.v
	vp_SideBySideAniso centlim.dots.Central_Limit.v centlim.spec.Central_Limit.v > centlim.v

/* create the finite difference stencils */

centered.H: if not exist centered.H /*    very poor  */
        echo "0 0 0 0 0 0 -1. 0  1. 0 0 0 0 0 0 0 " >centered.ascii.H 
	atoF <centered.ascii.H >centered.H
	echo "d1=0.15915494  o1=-1.352817" >> centered.H  /* pi/128  */


stagg.H:    if not exist stagg.H /*    better  */
        echo "0 0 0 0 0 0 0 -1. 1. 0 0 0 0 0 0 0" >stagg.ascii.H
	atoF < stagg.ascii.H >stagg.H
	echo "d1=0.15915494  o1=-1.352817" >> stagg.H

centlim.H: if not exist centlim.H /* compromise Francis' central limit   */ 
	echo "0 0 0 0 0.0006975446 -0.0095703125 0.079752604 -1.1962891 1.1962891 -0.079752604 0.0095703125 0.0006975446 0 0 0 0" >centlim.ascii.H 
	atoF <centlim.ascii.H  >centlim.H
	echo "d1=0.15915494  o1=-1.352817" >> centlim.H  /* pi/128  */

holberg.H: if not exist holberg.H /*      Holberg's coefficients        */
	echo "0.1326224603581672E-04 -0.1260422755482438E-03 0.6766549997685445E-03 -0.2713077307504363E-02 0.9277622829117149E-02 -0.3045387209576234E-01 0.1178572129899180E+00 -0.1247830531633298E+01 0.1247830531633298E+01 -0.1178572129899180E+00  0.3045387209576234E-01 -0.9277622829117149E-02  0.2713077307504363E-02 -0.6766549997685445E-03  0.1260422755482438E-03 -0.1326224603581672E-04" >holberg.ascii.H
	 atoF <holberg.ascii.H >holberg.H
	echo "d1=0.15915494  o1=-1.352817" >> holberg.H  /* pi/128  */



FIGDIR/waves2D.ps:
	pstexpen FIGDIR/waves2D.v FIGDIR/waves2D.ps invras=n color=COLOR fat=1 


FIGDIR/waves2D.v: zsnap.ST.3.v zsnap.HB.3.v zsnap.CL.3.v circle.v
	vppen erase=once zsnap.ST.3.v circle.v |vppen align=lb >j1.v
	vppen erase=once zsnap.CL.3.v  |vppen align=lb >j2.v
	vppen erase=once zsnap.HB.3.v  |vppen align=lb >j3.v
	vp_SideBySideIso j1.v j2.v j3.v >FIGDIR/waves2D.v
	vppen j1.v j2.v j3.v >FIGDIR/waves2D.v3


circle.v: BINDIR/Circle
        Window f3=3 n3=1 <zsnap.HB.H >jj
        BINDIR/Circle <jj radius=3.2 x0=4.10 z0=4.10 >junk half=0
        Graph <junk yreverse=1 out=circle.v head=/dev/null xll=0 yll=0 min1=0 min2=0 max2=8.28 max1=8.28 xinch=3 yinch=3 title=" " plotfat=10 plotcol=7 dash=0  wantaxis=n label1=" " label2=" "
        Graph <junk yreverse=1 out=circle.v head=/dev/null xll=0 yll=0 min1=0 min2=0 max2=8.28 max1=8.28 xinch=3 yinch=3 title=" " dash=0 plotfat=5  wantaxis=n label1=" " label2=" "


#define SNAPSIZE xll=0 yll=0 min1=0 min2=0 max2=8.28 max1=8.28 xinch=3 yinch=3 label1="Depth [km]" label2="Distance [km]"

zsnap.CE.3.v: zsnap.CE.3.T
	Ta2vplot <zsnap.CE.3.T  title="Centered" SNAPSIZE out=zsnap.CE.3.v NULL

zsnap.ST.3.v: zsnap.ST.3.T
	Ta2vplot <zsnap.ST.3.T  title="Staggered" SNAPSIZE out=zsnap.ST.3.v NULL

zsnap.HB.3.v: zsnap.HB.3.T
	Ta2vplot <zsnap.HB.3.T  title="Holberg" SNAPSIZE out=zsnap.HB.3.v NULL

zsnap.CL.3.v: zsnap.CL.3.T
	Ta2vplot <zsnap.CL.3.T  title="Central Limit" SNAPSIZE out=zsnap.CL.3.v NULL

%.3.T: %.3.H
        Byte <%.3.H >%.3.T

%.3.H: %.H
	Window f3=3 n3=1 <%.H >%.3.H


/* FD model in 2D */


zsnap.CL.H acoustCL& : centlim.H homoac.H wavelet.vert.H  BINDIR/Ultiall 
        time BINDIR/Ultiall  type=1 modeling=1 moduli=homoac.H par=parfile  <wavelet.vert.H >junk ncoef=16 coeff=centlim.H zsnap=zsnap.CL.H zseis=zseis.CL.H

zsnap.HB.H acoustHB& : holberg.H homoac.H wavelet.vert.H  BINDIR/Ultiall 
        time BINDIR/Ultiall  type=1 modeling=1 moduli=homoac.H par=parfile  <wavelet.vert.H >junk ncoef=16 coeff=holberg.H zsnap=zsnap.HB.H zseis=zseis.HB.H

zsnap.ST.H acoustST& : stagg.H homoac.H wavelet.vert.H  BINDIR/Ultiall 
        time BINDIR/Ultiall  type=1 modeling=1 moduli=homoac.H par=parfile  <wavelet.vert.H >junk ncoef=16 coeff=stagg.H zsnap=zsnap.ST.H zseis=zseis.ST.H

zsnap.CE.H acoustCE& : centered.H homoac.H wavelet.vert.H  BINDIR/Ultiall 
	Wavelet n1=128 d1=0.001 domain=time wavelet=ricker1 tdelay=0.080 fund=15 fhigh=30. | Scale dscale=1000000. >wavelet.vert.ST.H
        time BINDIR/Ultiall  type=1 modeling=1 moduli=homoac.H par=parfile  <wavelet.vert.ST.H >junk ncoef=16 coeff=centered.H zsnap=zsnap.CE.H zseis=zseis.CE.H stagger=0



wavelet.vert.H: 
        Wavelet n1=128 d1=0.004 domain=time wavelet=ricker1 tdelay=0.080 fund=20 fhigh=30. | Scale dscale=1000000. >junk
        Spike n1=128 n2=1 n3=1 d1=0.004 k1=-1 k2=-1 k3=-1 >junk2
        Merge junk junk2 space=n axis=2 >wavelet.vert.H


homoiso.H: BINDIR/Geninput
        BINDIR/Geninput nx=208 nz=208 dx=40. dz=40.  symmetry=2 vp=4.0 vs=2.3  rho=2.5 >homoiso.H

homoac.H: homoiso.H
        Window n3=1 <homoiso.H >junk.c11
        Window f3=4 <homoiso.H >junk.rho
        Merge junk.c11 junk.rho axis=3 space=n >homoac.H



BINDIR/%:  if not {{ cd ../../Bin ; cake -q MTYPE/% }}
	( cd ../../Bin ; cake MTYPE/% )

clean&: jclean 
	RM_CMD *.v jj* *.H *.T

#include <SEP.obj.rules>
#include <SEP.prog.rules>
#include <SEP.idoc.rules>
