#! /bin/sh

if test "$PSPREVIEWER" = ""
then
echo "Environment variable PSPREVIEWER not set" ; exit 1
fi

# model number
model=3

####  find the right directory
datadir=../datadir/

####  aux files
modelfile=${datadir}model${model}.data
modelpsfile=${datadir}model${model}.eps
rayendsfile=${datadir}rayend.data
rayfile=${datadir}ray.data
raypsfile=${datadir}ray.eps
wavefile=${datadir}wave.data
wavepsfile=${datadir}wave.eps
wavefilet=${datadir}wavet.data

nangle=50 nxz=301 fangle=30 langle=45 nt=20

####   create  modelfile and modelpsfile
# ../Models/PSmodel${model}

#### Create Figure 5 in the GBmod paper
####   shoot the rays and get wavefronts
triray <$modelfile >$rayendsfile \
        rayfile=$rayfile wavefile=$wavefile \
        nangle=$nangle fangle=$fangle langle=$langle \
        nxz=$nxz nt=$nt xs=0.1 zs=0.1 \
        refseq=4,1

####   plot the rays 
psgraph <$rayfile >$raypsfile \
        nplot=`cat outpar` n=$nxz hbox=3.0 wbox=6.0  \
        x1beg=0 x1end=1 x2beg=0 x2end=2 \
        label1="Depth (km)" label2="Distance (km)" \
        style=seismic linegray=1

####    transpose the wavefile
transp <$wavefile >$wavefilet n1=$nt n2=$nangle nbpe=8

####    plot the wavefronts
psgraph <$wavefilet >$wavepsfile  \
	linewidth=0 mark=8 marksize=2 \
        nplot=$nt n=$nangle hbox=3.0 wbox=6.0 \
        x1beg=0 x1end=1 x2beg=0 x2end=2 \
        label1="Depth (km)" label2="Distance (km)" \
        style=seismic

####    display rays and waves
$CWPROOT/bin/psmerge in=$modelpsfile in=$raypsfile in=$wavepsfile | $PSPREVIEWER  &


exit 0
