Makefile 9.78 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
	  | tidofi -type bin -Type bin -v -o -cs $@ $<
113
cpsdbgmod1.bin: cpsdbgraw.bin
114
	printf "fac  $(CPSD_AMPBG)\nlp1  2.\navg  0\nhpb  100.,4\nend\n" \
115
	  | tidofi -type bin -Type bin -v -o -cs $@ $<
116
117
118
cpsdbgmod2.bin: cpsdbgmod1.bin
	printf "dif  0.\nend\n" \
	  | tidofi -type bin -Type bin -v -o -cs $@ $<
119
120

# add common noise signal to channel specific noise floor
121
122
123
124
125
126
127
cpsd1add.bin: cpsd1mod.bin cpsdbgmod1.bin
	teseco -v -o -a -type bin -Type bin $@ $^
cpsd2add.bin: cpsd2mod.bin cpsdbgmod1.bin
	teseco -v -o -a -type bin -Type bin $@ $^
cpsd3add.bin: cpsd3mod.bin cpsdbgmod1.bin
	teseco -v -o -a -type bin -Type bin $@ $^
cpsd4add.bin: cpsd4mod.bin cpsdbgmod2.bin
128
129
	teseco -v -o -a -type bin -Type bin $@ $^
	
130
131
132
133
134
135
136
# ----------------------------------------------------------------------
# adjust time series headers
# --------------------------
#
#  recorded data is expected to be taken with slightly different sampling;
#  croposp is expected to be able to handle this
#
137
138
139
140
141
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
# 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) $@ $<

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

187
188
# define string summarizing essential croposp options
# this will be included in the title of diagrams
189
190
191
192
CPSD_CROPOSP_OPT=--nsegments=$(CPSD_NSEGMENTS) \
		     --divisor=$(CPSD_DIVISOR) \
		     --padfactor=$(CPSD_PADFACTOR) \
		     --overlap=$(CPSD_OVERLAP) $(CPSD_LOG)
193

194
195
196
197
198
199
200
201
202
203
204
# 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
205
206
207
208
209
210
211
212
213
214
215
216
$(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 \
	  $^

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

226
227
# compute PSD with foutra, which is used as a reference to confirm correct
# calibration of algorithms
228
229
230
231
232
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 $^

233
# convert foutra output to the file format expected by croposplot.py
234
235
236
237
psd%foutra.xxx: psd%.001.asc
	printf "# #1: foutra (reference): %s\n" $< > $@
	cat $< >> $@

238
239
240
241
242
243
244
# ======================================================================
# display results
# ===============
#
# location of plot tool
CROPOSPLOT=../croposplot.py

245
# common options
246
247
CROPOSPLOT_MATCH=.
CROPOSPLOT_OPT=--usemarkers --match='$(CROPOSPLOT_MATCH)'
248

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

261
262
# display PSD of raw data (apparent instrumental noise and seismic
# background noise separately)
263
264
PSDraw.pdf: psdraw.xxx $(CPSD_RAWREFFILES)
	$(CROPOSPLOT) --title="$@ \n $(CPSD_CROPOSP_OPT)" \
265
	  --ylim=1.e-14:100. \
266
	  $(CROPOSPLOT_OPT) \
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
	$(CROPOSPLOT) --title="$@ \n $(CPSD_CROPOSP_OPT)" \
273
	  $(CROPOSPLOT_OPT) \
274
	  --titfontsize='small' -v --grid -o $@ $^
275

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

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

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

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