c program sac_wavelet <br> 
c JJP 11/12/97, spruced up for distribution 9/26/00 <br> 
c f77 -o sac_wavelet sac_wavelet.f -fast -native -O5   <br> 
c <br> 
c code reads SAC-format data, filename.bh? <-- routine subs e,n,z for "?" <br> 
c <br> 
c generates wavelets with eigenvector decomposition + interpolation <br> 
c calculate padded FFT of data for convolution with wavelet FFTs <br> 
c general use - calculates set of wavelets for a specified number of points <br> 
c interpolates and fft other lengths on the fly <br> 
c max=12 real-wavelets max = 6 complex-valued wavelets <br> 
c <br> 
c code is based on the algorithms described in <br> 
c <br> 
c Lilly, J., and J. Park, Multiwavelet spectral and polarization analysis <br> 
c of seismic records, Geophysical J. International, v122, 1001-1021, 1995. <br> 
c <br> 
c Park, J. and M. E. Mann, Interannual temperature events and shifts in <br> 
c global temperature: A multiwavelet correlation approach, Earth Interactions, <br> 
c v. 4, online paper 4-001, 2000. <br> 
c <br> 
C SEE ALSO: <br> 
c Bear, L. K., and G. L. Pavlis, Estimation of slowness vectors and <br> 
c their uncertainties using multi-wavelet seismic array processing, <br> 
c Bull. Seismol. Soc. Am., v87, 755--769, 1997. <br> 
c <br> 
c this code generates 2-D arrays of spectra, polarization and statistical stuff <br> 
c it writes GMT scripts to display these arrays in dazzling color <br> 
c <br> 
c for an input filename sacfile.bh? the GMT scripts are written to <br> 
c GMT_sacfile.bh, GM[123]_sacfile.bh <br> 
c <br> 
c the tarfile for distribution should have binary SACfiles nilcut.sh? <br> 
c -- a P wave coda <br> 
c you can test the program on this data <br> 
c <br> 
c the ASCII files to plot are named <br> 
c cpt_nilcut.sh lam_nilcut.sh1 spw_nilcut.shz <br> 
c dat_nilcut.sh lam_nilcut.sh2 spw_nilcut.she <br> 
c ell_nilcut.sht lam_nilcut.sh3 spw_nilcut.shn <br> 
c ell_nilcut.shz spw_nilcut.shs <br> 
c cp2_nilcut.sh lam_nilcut.sh0 spw_nilcut.shv <br> 
c <br> 
c dat_* is data to plot <br> 
c ell_* contains polarization ellipses in Z/R (z) and T/R (t) planes <br> 
c spw_* are wavelet spectra fot *z *r *t <br> 
c spw_*s contains splitting ellipses and spw_*v has rectilinearity estimate <br> 
c cpt_* and cp2_* are colorfiles for GMT <br> 
c lam_* are singular values of the wavelet polariation matrices <br> 
c <br> 
c you can fill a lot of disk space if you are not tidy with old files! <br> 
c <br> 
c user must make scripts executable ("chmod +x GM*_sacfile.bh"), <br> 
c execute the scripts and display the postscript files <br> 
c <br> 
c <br> 
<br> 
EXAMPLE OF GMT SCRIPT WRITTEN BY sac_wavelet.f <br> 
<br> 
colorfil=cpt_nilcut.sh <br> 
infil=spw_nilcut.shz <br> 
/bin/rm aaplot <br> 
xyz2grd $infil -Gblobg -I0.15/0.3800 -F -R0.00/30.15/-4.2489/1.8311 -V <br> 
grdimage blobg -C$colorfil -B5.00:'time (sec)':/1.0enSW -JX5.5/2.3 -X1.0 -Y1.0 -K -P > aaplot <br> 
pstext -JX -R -W -O -P -K <<END>>aaplot <br> 
1.500 -2.749 15 0 15 1 Vertical <br> 
END <br> 
infil=spw_nilcut.shn <br> 
xyz2grd $infil -Gblobg -I0.15/0.3800 -F -R0.00/30.15/-4.2489/1.8311 -V <br> 
grdimage blobg -C$colorfil -B5.00/1.0:'period (2@+y@+ sec)':ensW -JX5.5/2.3 -X0.0 -Y2.3 -O -K -P >> aaplot <br> 
pstext -JX -R -W -O -P -K <<END>>aaplot <br> 
1.500 -2.749 15 0 15 1 Transverse <br> 
END <br> 
psscale -C$colorfil -D5.8/1.15/3/0.25 -L -O -K -P >> aaplot <br> 
infil=spw_nilcut.she <br> 
xyz2grd $infil -Gblobg -I0.15/0.3800 -F -R0.00/30.15/-4.2489/1.8311 -V <br> 
grdimage blobg -C$colorfil -B5.00/1.0ensW -JX5.5/2.3 -X0.0 -Y2.3 -O -K -P >> aaplot <br> 
pstext -JX -R -W -O -P -K <<END>>aaplot <br> 
1.500 -2.749 15 0 15 1 Radial <br> 
END <br> 
psxy dat_nilcut.sh -R/0.000/30.150/-64700.434/188000.172 -B5.00:'.nilcut':/50000.000ensW -M -JX5.5/2.0 -X0.0 -Y2.3 -O -K -P >> aaplot <br> 
pstext -JX -R -O -P <<END>>aaplot <br> 
1.500 6250.000 10 0 15 1 Vertical <br> 
1.500 114759.695 10 0 15 1 Transverse <br> 
1.500 156296.094 10 0 15 1 Radial <br> 
END <br>