#! /bin/sh
WIDTH=420
HEIGHT=400
WIDTHOFF1=0
WIDTHOFF2=430
WIDTHOFF3=860
HEIGHTOFF1=50
HEIGHTOFF2=500

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 
suxwigb <modeldata title="Data: 64ms reverbs" \
	windowtitle="Data" \
	label1="Time (sec)" label2="Trace" \
	wbox=$WIDTH hbox=$HEIGHT xbox=$WIDTHOFF1 ybox=$HEIGHTOFF1 &
 
# Plot the autocorrelation with gain
sugain <modeldata tpow=$tpow |
suacor ntout=31 |
suxwigb title="Autocorrelation" \
	windowtitle="AutoCorr" \
	label1="Time (sec)" label2="Trace" \
	wbox=$WIDTH hbox=$HEIGHT xbox=$WIDTHOFF2 ybox=$HEIGHTOFF1 &

# 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 |
suxwigb label1="Time"  label2="Trace" title="Spiking Decon: 4,40ms" \
	windowtitle="Spike" \
	wbox=$WIDTH hbox=$HEIGHT xbox=$WIDTHOFF3 ybox=$HEIGHTOFF1 &

# Plot the autocorrelation after spike
supef <modeldata maxlag=.04 |
sugain tpow=$tpow |
suacor ntout=31 |
suxwigb title="Autocor after spike" \
	windowtitle="AutoCorr/Spike" \
	label1="Time (sec)" label2="Trace" \
	wbox=$WIDTH hbox=$HEIGHT xbox=$WIDTHOFF1 ybox=$HEIGHTOFF2 &

# Second, go after the reverberations:  note reflector!
supef <modeldata maxlag=.04 |
supef minlag=.05 maxlag=.16 |
suxwigb label1="Time"  label2="Trace" \
	title="PEF: 50,160ms (Note dipping reflector!)" \
	windowtitle="PEF" \
	wbox=$WIDTH hbox=$HEIGHT xbox=$WIDTHOFF2 ybox=$HEIGHTOFF2 &

# Bandpass
supef <modeldata maxlag=.04 |
supef minlag=.05 maxlag=.16 |
sufilter f=5,15,80,100 |
suxwigb label1="Time"  label2="Trace" title="Spike, PEF, BandPass" \
	windowtitle="Final" \
	wbox=$WIDTH hbox=$HEIGHT xbox=$WIDTHOFF3 ybox=$HEIGHTOFF2 &

exit 0
