Commit 13c902bc authored by thomas.forbriger's avatar thomas.forbriger
Browse files

sleeman [FEATURE]: provide control script

tricrosp.sh is a control script to run a complete cross spectral coherency
analyses. The final result of the run is a graphical display of the signal
values.
parent 129fba9e
#!/bin/sh
# this is <tricrosp.sh>
# ----------------------------------------------------------------------------
#
# Copyright (c) 2016 by Thomas Forbriger (BFO Schiltach)
#
# run a full cross spectral analysis
#
# method published by:
# Reinoud Sleeman, Arie van Wettum, and Jeannot Trampert, 2006. Three-Channel
# Correlation Analysis: A New Technique to Measure Instrumental Noise of
# Digitizers and Seismic Sensors. Bull. Seism. Soc. Am. vol. 96(1), 258–271,
# doi: 10.1785/0120050032
#
# use programs implemented by E. Wielandt and provided for Linux at
# https://git.scc.kit.edu/Seitosh/software-for-seismometry-linux
#
# REVISIONS and CHANGES
# 14/09/2016 V1.0 Thomas Forbriger
#
# ============================================================================
#
VERSION="14/09/2016 V1.0"
# ============================================================================
# shell functions
# ===============
# usage function
function usage {
cat << HERE
This is <tricrosp.sh>
version: $VERSION
usage: tricrosp.sh [-v] [-d] [-t type] [-o name]
[-n nband] [-g gain] [-G gain]
T0 h file1 file2 file3
or: tricrosp.sh -h
Run a full cross spectral analysis.
HERE
}
# ----------------------------------------------------------------------------
# usage function
function longusage {
cat << HERE
-h print detailed usage information
-v be verbose
-d print debug output
-t type input data files are of file type \"type\"
-o name use \"name\" as a basename for all results
-n nband use \"nband\" frequency bands per decade
-g gain gain of seismometer is \"gain\" (in Vs/m)
-G gain gain of digitizer is \"gain\" (in uV/count)
T0 eigenperiod of reference seismometer (in seconds)
h damping of reference seismometer as a fraction of critical damping
file1 seismic time series of reference seismometer
file2 seismic time series of second seismometer
file3 seismic time series of third seismometer
The method is published by:
Reinoud Sleeman, Arie van Wettum, and Jeannot Trampert, 2006. Three-Channel
Correlation Analysis: A New Technique to Measure Instrumental Noise of
Digitizers and Seismic Sensors. Bull. Seism. Soc. Am. vol. 96(1), 258–271,
doi: 10.1785/0120050032
The script uses programs implemented by E. Wielandt and provided for Linux
at https://git.scc.kit.edu/Seitosh/software-for-seismometry-linux
HERE
}
# ============================================================================
# read command line
# =================
verbose=0
VFLAG=
VFLAGL=
debug=0
TYPE=seife
BASE=tcresult
NBAND=6
SGAIN=1.
DGAIN=1.
# fetch commandline arguments
while getopts 'vdht:o:n:g:G:' Option;
do
case ${Option} in
h) hflag=1
usage; longusage; exit 0;;
v) verbose=1; VFLAG=-v; VFLAGL=--verbose;;
d) debug=1;;
t) TYPE="$OPTARG";;
o) BASE="$OPTARG";;
n) NBAND="$OPTARG";;
g) SGAIN="$OPTARG";;
G) DGAIN="$OPTARG";;
*) echo -e "ERROR: Unimplemented option chosen.\n"
usage; exit 2;;
esac
done
shift $((${OPTIND} - 1))
MINPARAMS=5
# check number of commandline arguments
if [ $# -lt ${MINPARAMS} ]
then
usage
exit 2
fi
T0=$1
h=$2
FILE1=$3
FILE2=$4
FILE3=$5
# ============================================================================
# start action
# ============
# report
if test $verbose -gt 0
then
cat <<HERE
-----------------------------------------------------------------------
$0
Version: $VERSION
Run a full cross spectral analysis.
Reference signal: $FILE1
Signal of 2nd seimometer: $FILE2
Signal of 3rd seimometer: $FILE3
Eigenperiod of reference seismometer: $T0 s
Damping of reference seismometer: $h of critical damping
Seismometer's gain: $SGAIN Vs/m
Digitizer's gain: $DGAIN uV/count
Input file type: $TYPE
Output basename: $BASE
-----------------------------------------------------------------------
HERE
fi
DATA1=${BASE}.seis1.seife
DATA2=${BASE}.seis2.seife
DATA3=${BASE}.seis3.seife
any2any $VFLAGL --itype=$TYPE --otype=seife --overwrite $DATA1 $FILE1
any2any $VFLAGL --itype=$TYPE --otype=seife --overwrite $DATA2 $FILE2
any2any $VFLAGL --itype=$TYPE --otype=seife --overwrite $DATA3 $FILE3
bandnois $DATA1 $DGAIN $T0 $h $SGAIN $NBAND
tricrosp $DATA1 $DATA2 $DATA3 $NBAND
tri2db
for d in bandnois.out tri2db.out tricrosp.out
do
/bin/mv $VFLAG $d ${BASE}.$d
done
bandnois2table.awk ${BASE}.bandnois.out > ${BASE}.bandnois.dat
tricrosp2table.awk ${BASE}.tricrosp.out > ${BASE}.tricrosp.dat
tridb2table.awk ${BASE}.tri2db.out > ${BASE}.tri2db.dat
cat > ${BASE}.gpt << HERE
set terminal postscript color solid enhanced
set output '${BASE}.ps'
set logscale x
set grid
set style line 10 lt 1 lw 3
set style line 1 lt 1
set style line 2 lt 3
set style line 3 lt 4
set style line 4 lt 5
set style line 5 lt 7
set ylabel 'PSD acceleration (rel. to 1 m^{2} s^{-4} Hz^{-1}) / dB'
set xlabel 'frequency / Hz'
tomps(x)=x*1.e-6*$DGAIN/$SGAIN
om0=2*pi/$T0
om(x)=2*pi*x
Trec(x)=sqrt((om0**2-om(x)**2)**2+(2*$h*om0*om(x))**2)
tompss(f,x)=x*Trec(f)/om(f)
RBW=(exp(log(10.)/$NBAND)-1.)/exp(log(10.)/(2.*$NBAND))
cntdb(f,x)=10.*log10(((tompss(f,tomps(x)))**2)/(2.*f*RBW))
set title 'Power spectral density'
plot '${BASE}.tri2db.dat' u 1:2 t 'NLNM' w l ls 10, \
'${BASE}.bandnois.dat' u 1:5 t 'reference signal' w l ls 2, \
'${BASE}.tri2db.dat' u 1:3 t 'incoherent part of ref. seismometer' w l ls 3, \
'${BASE}.tri2db.dat' u 1:5 t 'incoherent part of 2nd seismometer' w l ls 4, \
'${BASE}.tri2db.dat' u 1:7 t 'incoherent part of 3rd seismometer' w l ls 5
plot '${BASE}.tri2db.dat' u 1:2 t 'NLNM' w l ls 10, \
'${BASE}.bandnois.dat' u 1:5 t 'reference signal' w l ls 2, \
'${BASE}.tricrosp.dat' u (\$1):(cntdb(\$1,\$2)) \
t 'coherent part of signals' w l ls 3
set title 'Power spectral density relative to NLNM'
plot '${BASE}.tri2db.dat' u 1:4 t 'incoherent part of ref. seismometer' \
w l ls 3, \
'${BASE}.tri2db.dat' u 1:6 t 'incoherent part of 2nd seismometer' w l ls 4, \
'${BASE}.tri2db.dat' u 1:8 t 'incoherent part of 3rd seismometer' w l ls 5
set ylabel 'coherency'
set yrange [*:1.01]
set title 'Coherency'
plot '${BASE}.tricrosp.dat' u 1:(sqrt(\$8)) t '2nd seism. to ref.' w l ls 1, \
'${BASE}.tricrosp.dat' u 1:(sqrt(\$9)) t '3rd seism. to ref.' w l ls 2
set yrange[*:*]
set title 'Transfer function'
set ylabel 'amplitude factor'
plot '${BASE}.tricrosp.dat' u 1:3 t '2nd seism. to ref.' w l ls 1, \
'${BASE}.tricrosp.dat' u 1:4 t '3rd seism. to ref.' w l ls 2
HERE
gnuplot ${BASE}.gpt
gv ${BASE}.ps
# ----- END OF tricrosp.sh -----
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment