#! /bin/sh
demo=decon1

tpow=1		# gain power before examining autocorrelations
		# tpow=2 might be better for field data

###PURPOSE###
# We make some simple data: 8 traces with 2-way reverberations and
# a weak tilting reflector hidden underneath the reverberations.  The
# spike data is convolved with the minimum phase wavelet specified below.
# Then a tiny bit of noise is added--tiny because: (1) the definition
# of signal to noise is stringent--based on the biggest amplitude
# on the trace and (2) these are single traces, usually you'd have
# a cmp gather and the ultimate goal of stacking to combat noise.

# In this demo, we do spiking decon first and then do gapped decon
# (prediction error filtering).  You are to assume that the decon
# parameters are estimated from the autocorrelation shown in frame two.

# Construction of the minimum phase wavelet with Mathematica:
# (2-z)(3+z)(3-2z)(4+z)(4-2z)(4+3z)//(CoefficientList[#,z])&
#     {1152, -384, -904, 288, 174, -34, -12}


# First make the synthetic data for the deconvolution demo.
######START OF MODEL######
I=${CWPROOT}/include
L=${CWPROOT}/lib

make

./traces |
suaddhead ns=512 |
sushw key=dt a=4000 |
suaddnoise sn=200 |
suconv filter=1152,-384,-904,288,174,-34,-12 >modeldata

rm traces
######END OF MODEL######


# Plot the model
supswigp <modeldata title="Model Data: 64ms reverbs" \
	label1="Time (sec)" label2="Trace" >${demo}.1.eps

# Plot the autocorrelation with gain
sugain <modeldata tpow=$tpow |
suacor ntout=31 |
supswigp title="Autocorrelation (with gain)" \
	label1="Time (sec)" label2="Trace" >${demo}.2.eps

# The autocorrelation shows that:
#       wavelet extends to about 30ms
#	first repetition ends at about 160ms
# this accords well with the actual 64ms reverberations and wavelet of 24ms


######DECON EXAMPLES######
# First, spike it up
supef <modeldata maxlag=.04 |
supswigp label1="Time" label2="Trace" \
	 title="Spiking Decon: 4,40ms" >${demo}.3.eps

# Plot the autocorrelation after spike
supef <modeldata maxlag=.04 |
sugain tpow=$tpow |
suacor ntout=31 |
supswigp title="Autocorrelation after spike (with gain)" \
	label1="Time (sec)" label2="Trace" >${demo}.4.eps

# Second, go after the reverberations
supef <modeldata maxlag=.04 |
supef minlag=.05 maxlag=.16 |
supswigp label1="Time"  label2="Trace" \
title="PEF: 50,160ms (Note dipping reflector)" >${demo}.5.eps


# Bandpass
supef <modeldata maxlag=.04 |
supef minlag=.05 maxlag=.16 |
sufilter f=5,15,80,100 |
supswigp label1="Time"  label2="Trace" \
title="Spike, PEF, BandPass" >${demo}.6.eps

echo 
echo "Output PostScript files:"

echo -n "     "
for i in ${demo}.*.eps
do
	 echo -n "$i "
done
echo
echo "You may view these with your PostScript Previewer"

exit 0
