Makefile 9.69 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
# this is <Makefile>
# ----------------------------------------------------------------------------
# 
# Copyright (c) 2019 by Thomas Forbriger (BFO Schiltach) 
# 
# compile croposp
# 
# ----
# This program is free software; you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation; either version 2 of the License, or
# (at your option) any later version. 
# 
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU General Public License for more details.
# 
# You should have received a copy of the GNU General Public License
# along with this program. If not, see <http://www.gnu.org/licenses/>.
# ----
#
# REVISIONS and CHANGES 
#    06/02/2019   V1.0   Thomas Forbriger
# 
# ============================================================================
#

alldata: cpsd1.bin cpsd2.bin cpsd3.bin cpsd4.bin

31
32
33
34
35
36
.PHONY: clean
clean: ; 
	-find . -name \*.bak | xargs --no-run-if-empty /bin/rm -v
	-/bin/rm -vf flist *.bin *.o *.xxx *.ps *.pdf *~ 
	-/bin/rm -vf *.asc

37
38
39
40
41
42
# ======================================================================
# create time series
# ==================
#
# set time series parameters
# --------------------------
43
CPSD_DT=0.01
44
CPSD_DUR=8640.
45
#CPSD_DUR=864.
46

47
# ----------------------------------------------------------------------
48
# create raw time series
49
50
51
52
53
54
# ----------------------
# Several time series of normally distributed random samples with zero mean
# are computed with siggenx. To make them incoherent, each run requires its
# own specific seed value for the random number generator. We use the
# nanosecond field of the system clock.
#
55
CPSD_SIGGENOPT=-v -o -d $(CPSD_DT) -T $(CPSD_DUR) -ot bin
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78

cpsd1raw.bin:
	/bin/bash -c 'export GSL_RNG_SEED=$$(date +"%N"); \
	  siggenx 14 $@ $(CPSD_SIGGENOPT); \
	  echo seed: $$GSL_RNG_SEED'
cpsd2raw.bin:
	/bin/bash -c 'export GSL_RNG_SEED=$$(date +"%N"); \
	  siggenx 14 $@ $(CPSD_SIGGENOPT); \
	  echo seed: $$GSL_RNG_SEED'
cpsd3raw.bin:
	/bin/bash -c 'export GSL_RNG_SEED=$$(date +"%N"); \
	  siggenx 14 $@ $(CPSD_SIGGENOPT); \
	  echo seed: $$GSL_RNG_SEED'
cpsd4raw.bin:
	/bin/bash -c 'export GSL_RNG_SEED=$$(date +"%N"); \
	  siggenx 14 $@ $(CPSD_SIGGENOPT); \
	  echo seed: $$GSL_RNG_SEED'
# create time series for common noise signal
cpsdbgraw.bin:
	/bin/bash -c 'export GSL_RNG_SEED=$$(date +"%N"); \
	  siggenx 14 $@ $(CPSD_SIGGENOPT); \
	  echo seed: $$GSL_RNG_SEED'

79
# provide a time series plot to confirm incoherency
80
81
82
rawsig.ps: cpsd1raw.bin cpsd2raw.bin cpsd3raw.bin cpsd4raw.bin cpsdbgraw.bin
	stuplox -d $@/cps -s x -i -ty bin $^

83
84
85
86
87
88
89
90
91
92
93
# ----------------------------------------------------------------------
# scale and filter time series
# ----------------------------
#
#  The raw random sequences are adjusted to represent four sets of incoherent
#  noise and one random sequence ('bg') which will be common in all final data
#  files. We understand the first four as a simulation of instrumental noise
#  and the bg-sequence as a representation of background ground noise to be
#  common in all reacordings.
#
# set amplitudes
94
95
96
97
98
CPSD_AMP1=0.2
CPSD_AMP2=4.
CPSD_AMP3=100.
CPSD_AMP4=1.e-2
CPSD_AMPBG=10.
99

100
101
# modify white noise to channel specific frequency characteristic
cpsd1mod.bin: cpsd1raw.bin
102
	printf "fac  $(CPSD_AMP1)\nhpb  1.,2\nend\n" \
103
104
	  | tidofi -type bin -Type bin -v -o -cs $@ $<
cpsd2mod.bin: cpsd2raw.bin
105
	printf "fac  $(CPSD_AMP2)\nhpb 4.,4\nlpb  2.,4\nend\n" \
106
107
	  | tidofi -type bin -Type bin -v -o -cs $@ $<
cpsd3mod.bin: cpsd3raw.bin
108
	printf "fac  $(CPSD_AMP3)\nlp2  100.,3\nend\n" \
109
110
	  | tidofi -type bin -Type bin -v -o -cs $@ $<
cpsd4mod.bin: cpsd4raw.bin
111
	printf "fac  $(CPSD_AMP4)\nend\n" \
112
113
	  | tidofi -type bin -Type bin -v -o -cs $@ $<
cpsdbgmod.bin: cpsdbgraw.bin
114
	printf "fac  $(CPSD_AMPBG)\nlp1  2.\navg  0\nhpb  100.,4\nend\n" \
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
	  | tidofi -type bin -Type bin -v -o -cs $@ $<

# add common noise signal to channel specific noise floor
cpsd%addraw.bin: cpsd%mod.bin cpsdbgmod.bin
	teseco -v -o -a -type bin -Type bin $@ $^

# modify sensor response
cpsd1add.bin: cpsd1addraw.bin
	printf "fac  1.\nend\n" \
	  | tidofi -type bin -Type bin -v -o -cs $@ $<
cpsd2add.bin: cpsd2addraw.bin
	printf "fac  1.\nend\n" \
	  | tidofi -type bin -Type bin -v -o -cs $@ $<
cpsd3add.bin: cpsd3addraw.bin
	printf "fac  1.\nend\n" \
	  | tidofi -type bin -Type bin -v -o -cs $@ $<
cpsd4add.bin: cpsd4addraw.bin
	printf "lpb  10.,2\nend\n" \
	  | tidofi -type bin -Type bin -v -o -cs $@ $<
	
135
136
137
138
139
140
141
# ----------------------------------------------------------------------
# adjust time series headers
# --------------------------
#
#  recorded data is expected to be taken with slightly different sampling;
#  croposp is expected to be able to handle this
#
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
# set header fields to reasonable time an channel identifiers
CPSD_DATE=2018/12/16_14:24:16
CPSD_DATE1=$(CPSD_DATE).010980
CPSD_DATE2=$(CPSD_DATE).011050
CPSD_DATE3=$(CPSD_DATE).010760
CPSD_DATE4=$(CPSD_DATE).010760

CPSD_HDOPTS=-verbose -overwrite -itype bin -otype bin -ss TST -sa XX -sc BHZ
cpsd1hd.bin: cpsd1add.bin
	sehefixx $< $@ $(CPSD_HDOPTS) -st $(CPSD_DATE1) -si SYN1
cpsd2hd.bin: cpsd2add.bin
	sehefixx $< $@ $(CPSD_HDOPTS) -st $(CPSD_DATE2) -si SYN2
cpsd3hd.bin: cpsd3add.bin
	sehefixx $< $@ $(CPSD_HDOPTS) -st $(CPSD_DATE3) -si SYN3
cpsd4hd.bin: cpsd4add.bin
	sehefixx $< $@ $(CPSD_HDOPTS) -st $(CPSD_DATE4) -si SYN4

# extract inconsistent time windows
CPSD_TOPT1=-samplesf 3 -samplesl 7
CPSD_TOPT2=-samplesf 25
CPSD_TOPT3=-samplesl 81
CPSD_TOPT4=-samplesl 1
cpsd1.bin: cpsd1hd.bin
	tijerasxx -verbose -overwrite -iformat bin -oformat bin \
	  $(CPSD_TOPT1) $@ $<
cpsd2.bin: cpsd2hd.bin
	tijerasxx -verbose -overwrite -iformat bin -oformat bin \
	  $(CPSD_TOPT2) $@ $<
cpsd3.bin: cpsd3hd.bin
	tijerasxx -verbose -overwrite -iformat bin -oformat bin \
	  $(CPSD_TOPT3) $@ $<
cpsd4.bin: cpsd4hd.bin
	tijerasxx -verbose -overwrite -iformat bin -oformat bin \
	  $(CPSD_TOPT4) $@ $<

177
178
179
180
181
# ======================================================================
# analyze time series
# ===================
#
# set croposp options
182
CPSD_NSEGMENTS=50
183
184
185
186
CPSD_OVERLAP=0.5
CPSD_DIVISOR=100
CPSD_PADFACTOR=1
CPSD_NLOG=20
187
CPSD_LOG=--log $(CPSD_NLOG) --avgout
188
189
190
191
CPSD_OPT=--DEBUG
CPSD_OPT=
OUTPUTS=psd npsd coherency transfer phase

192
193
# define string summarizing essential croposp options
# this will be included in the title of diagrams
194
195
196
197
CPSD_CROPOSP_OPT=--nsegments=$(CPSD_NSEGMENTS) \
		     --divisor=$(CPSD_DIVISOR) \
		     --padfactor=$(CPSD_PADFACTOR) \
		     --overlap=$(CPSD_OVERLAP) $(CPSD_LOG)
198

199
200
201
202
203
204
205
206
207
208
209
# select sets to be analyzed
# four are available, three must be used at least
CPSD_SET=1 2 3
CPSD_SET=1 2 3 4
# compile list of file names as input to croposp
CPSD_BINFILES=$(addprefix cpsd,$(addsuffix .bin,$(CPSD_SET)))
# compile list of raw files (apparent instrumental noise and seismic
# background noise separately)
CPSD_RAWFILES=$(addprefix cpsd,$(addsuffix mod.bin,$(CPSD_SET) bg)) \

# run complete analysis on the apparent seismometer output
210
211
212
213
214
215
216
217
218
219
220
221
$(addsuffix _out.xxx,$(OUTPUTS)): $(CPSD_BINFILES)
	croposp $(CPSD_OPT) --verbose --overwrite --itype bin --trim \
	  --pattern="%I" \
	  $(CPSD_CROPOSP_OPT) \
	  --datetolerance 0.1 \
	  --psd psd_out.xxx \
	  --npsd npsd_out.xxx \
	  --coherency coherency_out.xxx \
	  --transfer transfer_out.xxx \
	  --phase phase_out.xxx \
	  $^

222
223
# compute PSD for raw files (apparent instrumental noise and seismic
# background noise separately)
224
225
226
227
228
229
230
psdraw.xxx: $(CPSD_RAWFILES)
	croposp $(CPSD_OPT) --verbose --overwrite --itype bin --trim \
	  --pattern="%F" \
	  $(CPSD_CROPOSP_OPT) \
	  --datetolerance 0.1 \
	  --psd $@ $^

231
232
# compute PSD with foutra, which is used as a reference to confirm correct
# calibration of algorithms
233
234
235
236
237
psd%.001.asc: cpsd%.bin
	foutra -v -o -type bin -ASCII=$(patsubst cpsd%.bin,psd%,$<) \
	  -power -rbw=0.1 -demean -detrend -logascii=0.1 -avgascii \
	  -nsegments=20 junk.xxx $^

238
# convert foutra output to the file format expected by croposplot.py
239
240
241
242
psd%foutra.xxx: psd%.001.asc
	printf "# #1: foutra (reference): %s\n" $< > $@
	cat $< >> $@

243
244
245
246
247
248
249
250
251
252
253
254
255
256
# ======================================================================
# display results
# ===============
#
# location of plot tool
CROPOSPLOT=../croposplot.py

# ----------------------------------------------------------------------
# reference plots
# ---------------
#
#  Display PSD together with results from foutra. 
#
# compile list of file names as analyzed with foutra
257
CPSD_REFFILES=$(addprefix psd,$(addsuffix foutra.xxx,$(CPSD_SET)))
258
259
# compile list of raw files (apparent instrumental noise and seismic
# background noise separately) analyzed with foutra
260
261
CPSD_RAWREFFILES=$(addprefix psd,$(addsuffix modfoutra.xxx,$(CPSD_SET) bg))

262
263
# display PSD of raw data (apparent instrumental noise and seismic
# background noise separately)
264
265
PSDraw.pdf: psdraw.xxx $(CPSD_RAWREFFILES)
	$(CROPOSPLOT) --title="$@ \n $(CPSD_CROPOSP_OPT)" \
266
	  --ylim=1.e-14:100. \
267
	  --titfontsize='small' -v --grid -o $@ $^
268

269
270
# display PSD of apparent seismometer output in comparison with foutra
# analysis
271
PSD.pdf: psd_out.xxx $(CPSD_REFFILES)
272
273
	$(CROPOSPLOT) --title="$@ \n $(CPSD_CROPOSP_OPT)" \
	  --titfontsize='small' -v --grid -o $@ $^
274

275
276
277
# ----------------------------------------------------------------------
# pattern rules for quick plots
# -----------------------------
278
%.pdf: %_out.xxx
279
280
	$(CROPOSPLOT) --title='$@ \n $(CPSD_CROPOSP_OPT)' \
	  --titfontsize='small' -v --grid -o $@ $<
281
282

%_lin.pdf: %_out.xxx
283
284
	$(CROPOSPLOT) --title="$@ \n $(CPSD_CROPOSP_OPT)" \
	  --titfontsize='small' -v --grid --nology -o \
285
286
	  $@ $<

287
288
289
# ----------------------------------------------------------------------
#  previewer commands
#  ------------------
290
291
292
293
%.pdp: %.pdf; evince $<; /bin/rm -fv $<
%.psp: %.ps; gv $<; /bin/rm -fv $<

# ----- END OF Makefile -----