Waveform filter programs
susei.f
Go to the documentation of this file.
1 c this is <susei.f>
2 c------------------------------------------------------------------------------
3 c
4 c Copyright 1998, 2010 by Thomas Forbriger (IfG Stuttgart)
5 c
6 c SUm of SEIsmograms
7 c
8 c ----
9 c This program is free software; you can redistribute it and/or modify
10 c it under the terms of the GNU General Public License as published by
11 c the Free Software Foundation; either version 2 of the License, or
12 c (at your option) any later version.
13 c
14 c This program is distributed in the hope that it will be useful,
15 c but WITHOUT ANY WARRANTY; without even the implied warranty of
16 c MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 c GNU General Public License for more details.
18 c
19 c You should have received a copy of the GNU General Public License
20 c along with this program; if not, write to the Free Software
21 c Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
22 c ----
23 c
24 c
25 c uses code from coro.f
26 c
27 c REVISIONS and CHANGES
28 c 16/06/98 V1.0 Thomas Forbriger
29 c 27/02/99 V1.1 add more than two
30 c 22/10/01 V1.2 calculate arithmetic mean of traces
31 c 27/01/04 V1.3 switch for normalization
32 c 28/10/13 V1.4 use libfapidxx to support additional file formats
33 c
34 c==============================================================================
35 c
36  program susei
37 c
38  character*79 version
39  parameter(version='SUSEI V1.4 SUm of SEIsmograms')
40 c
41 c dimensions
42  integer maxsamps, maxfree
43  parameter(maxsamps=150000, maxfree=20)
44 c
45  integer i, nsumfiles,ifile
46 c input data
47  real xdata(maxsamps), ydata(maxsamps)
48  integer ixdata(maxsamps), iydata(maxsamps)
49  equivalence(xdata,ixdata)
50  equivalence(ydata,iydata)
51  character*80 inxname, inyname, inputformat
52  integer lui
53  parameter(lui=11)
54 c output data
55  real data(maxsamps)
56  integer idata(maxsamps)
57  equivalence(data,idata)
58  character*80 outname, outputformat
59  character*132 wid2line
60  integer luo
61  parameter(luo=10)
62 c sff extras
63  character*80 free(maxfree)
64  integer nfree
65  real libversion
66  character timestamp*13, code*10, type*20, cs*1, date*6, time*10
67  character code2*10, cs2*1
68  real c1,c2,c3,c12,c22,c32
69  integer ierr,nsamp,nsamp2,nstack,nstack2
70  logical last,last2
71  real tanf,dt
72  character*132 wid2line2
73 c commandline
74  integer maxopt, lastarg, iargc
75  character*80 argument
76  parameter(maxopt=5)
77  character*3 optid(maxopt)
78  character*40 optarg(maxopt)
79  logical optset(maxopt), opthasarg(maxopt)
80 c
81  logical nonormalize
82 c debugging
83  logical debug
84 c here are the keys to our commandline options
85  data optid/'-d','-nn','-ty','-it','-ot'/
86  data opthasarg/2*.false.,3*.true./
87  data optarg/2*'-','sff','sff','sff'/
88 c
89 c------------------------------------------------------------------------------
90 c basic information
91 c
92  print *,version
93  print *,'Usage: susei [-nn] [-ty f] [-it f] [-ot f]'
94  print *,' file1 file2 ... outfile'
95  print *,' or: susei -help'
96 c
97  if (iargc().lt.1) stop 'ERROR: missing arguments'
98  call getarg(1, argument)
99  if (argument(1:5).eq.'-help') then
100  print *,' '
101  print *,'SUm of SEIsmograms'
102  print *,' '
103  print *,'file1 first input data file'
104  print *,'file2 second input data file to be stacked to file1'
105  print *,'... further input data files to be stacked'
106  print *,'outfile output data file'
107  print *,' '
108  print *,'-nn do not normalize'
109  print *,' The program is inteded to be used for stacking.'
110  print *,' In that case the amplitude should be normalized'
111  print *,' by the number of traces stacked. This option'
112  print *,' turns off amplitude normalization.'
113  print *,'-ty f choose file format ''f'' instead of SFF'
114  print *,' for input and output files'
115  print *,'-it f choose input file format ''f'' instead of SFF'
116  print *,'-ot f choose output file format ''f'' instead of SFF'
117  print *,' '
118  print *,'The data files are expected to have the same number of'
119  print *,'traces each. The number of traces is determined'
120  print *,'from the file with the least number of traces.'
121  print *,'All traces of equal index within each file are'
122  print *,'expected to have the same number of samples. The'
123  print *,'programs aborts, if this is not the case. No further'
124  print *,'consistency checks are applied.'
125  print *,' '
126  call sff_help_formats
127  stop
128  endif
129  if (iargc().lt.3) stop 'ERROR: too few arguments'
130 c
131 c
132 c------------------------------------------------------------------------------
133 c read command line arguments
134 c
135  call tf_cmdline(1, lastarg, maxopt, optid,
136  & optarg, optset, opthasarg)
137  debug=optset(1)
138  nonormalize=optset(2)
139  inputformat=optarg(3)
140  outputformat=optarg(3)
141  if (optset(4)) inputformat=optarg(4)
142  if (optset(5)) outputformat=optarg(5)
143 c
144  nsumfiles=iargc()-lastarg-2
145 c
146  call getarg(lastarg+1, inxname)
147  call getarg(lastarg+nsumfiles+2, outname)
148 c
149 c
150  call sff_select_input_format(inputformat, ierr)
151  if (ierr.ne.0) stop 'ERROR: selected input format is not supported'
152  call sff_select_output_format(outputformat, ierr)
153  if (ierr.ne.0) stop 'ERROR: selected output format is not supported'
154 c
155 c------------------------------------------------------------------------------
156 c
157 c
158  free(1)=version
159  free(2)='first data from '//inxname(1:index(inxname, ' '))
160  free(3)='second data from '//inyname(1:index(inyname, ' '))
161  free(4)='all information is taken from first data file'
162  free(5)='no plausibility checks were performed'
163  nfree=5
164 c
165  print *,'open ',inxname(1:index(inxname,' ')),'for input'
166  call sff_ropens(lui, inxname,
167  & libversion, timestamp, code,
168  & type, cs, c1, c2, c3, date, time, ierr)
169  if (ierr.ne.0) stop 'ERROR: opening x-component'
170  call sff_new(luo, outname, ierr)
171  print *,'open ',outname(1:index(outname,' ')),'for output'
172  call sff_wopenfs(luo, outname,
173  & free, nfree,
174  & type, cs, c1, c2, c3, date, time, ierr)
175  if (ierr.ne.0) stop 'ERROR: opening output file'
176  do ifile=1,nsumfiles
177  call getarg(lastarg+ifile+1, inyname)
178  print *,'open ',inyname(1:index(inyname,' ')),'for input'
179  call sff_ropens(lui+ifile, inyname,
180  & libversion, timestamp, code,
181  & type, cs, c1, c2, c3, date, time, ierr)
182  if (ierr.ne.0) stop 'ERROR: opening y-component'
183  enddo
184 c call sff_ROpenS(luz, inzname,
185 c & libversion, timestamp, code,
186 c & type, cs, c1, c2, c3, date, time, ierr)
187 c if (ierr.ne.0) stop 'ERROR: opening z-component'
188 c
189  last=.false.
190  do while (.not.(last))
191  nsamp=maxsamps
192  call sff_rtracei(lui, tanf, dt, wid2line, nsamp, xdata, ixdata,
193  & code, last, cs, c1, c2, c3, nstack, ierr)
194  if (ierr.ne.0) stop 'ERROR: reading x-data'
195  do i=1,nsamp
196  data(i)=xdata(i)
197  enddo
198  do ifile=1,nsumfiles
199  nsamp2=maxsamps
200  call sff_rtracei(lui+ifile, tanf, dt, wid2line2,
201  & nsamp2, ydata, iydata,
202  & code2, last2, cs2, c12, c22, c32, nstack2, ierr)
203  if (ierr.ne.0) then
204  print *,'ERROR: reading y-data no. ',ifile
205  stop
206  endif
207  if (nsamp.ne.nsamp2) then
208  print *,'ERROR: wrong number of samples in y-data no. ',ifile
209  stop
210  endif
211  last=(last.or.last2)
212  do i=1,nsamp
213  data(i)=(data(i)+ydata(i))
214  enddo
215  enddo
216  if (.not.nonormalize) then
217  do i=1,nsamp
218  data(i)=data(i)/float(nsumfiles+1)
219  enddo
220  endif
221  wid2line(36:38)='NSP'
222 c
223  call sff_wtracei(luo, wid2line, nsamp, data, idata, last,
224  & cs, c1, c2, c3, nstack, ierr)
225  if (ierr.ne.0) stop 'ERROR: writing trace'
226  enddo
227 c
228  stop
229  end
230 c
231 c ----- END OF susei.f -----
program susei
Definition: susei.f:36