stuplo.f 46.7 KB
Newer Older
1
2
c this is <stuplo.f> by Thomas Forbriger 1996
c
thomas.forbriger's avatar
thomas.forbriger committed
3
c $Id$
thomas.forbriger's avatar
thomas.forbriger committed
4
c
thomas.forbriger's avatar
thomas.forbriger committed
5
c This is a simple plotting tool for seismic time series in
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
31
32
33
34
35
36
37
38
39
c SFF format
c 
c (sff is the Stuttgart File Format and is defined in sff.doc)
c 
c BE SHURE YOUR COMPUTER USES 4 BYTES OF MEMORY TO STORE
c INTEGER AND REAL VARIABLES 
c (look at the definition of idata and data in the subroutine checkmem)
c
c REVISIONS and CHANGES
c  V1.0   14/11/96   first running version
c  V1.1   15/11/96   still under major construction
c  V1.2   19/11/96   first complete release
c  V1.3   20/11/96   correct labelling of time axis
c  V1.4   26/11/96   added array traceinfile
c  V1.5   03/12/96   changed control output when reading traces
c                    and added keep option
c  V1.6   13/12/96   set correct file status on open command
c  V1.7   07/01/97   corrected timeax-calculation and timeax min/max
c  V1.8   21/01/97   correct nsamp parameter to libstuff V1.07 reading rouintes
c  V1.9   23/01/97   fixed problem with traces of zero amplitude
c  V1.10  03/02/97   * inlcuded chaining for overlay plots
c                    * included scaling factor for y- and x-size
c                    * included y-centering
c  V1.11  09/07/98   introduced time marker plotting
c  V1.12  15/07/98   * do not provide usage information in case of correct
c                      commandline parameters
c                    * introduced verbose flag
c  V1.13  10/02/99   * allow non-black colour for first set
c                    * removed old tflib calls (changed to tf_ calls)
c  V1.14  23/04/99   * linewidth and character height may be set
c                      from command line
c                    * write y-axis numeric labels perpendicular to axis
c                      if desired
c  V1.15  05/11/99   * use correct TBOX option when using full time x-labels
thomas.forbriger's avatar
thomas.forbriger committed
40
c  V1.16  19/07/01   * descriptive help text for time marker
thomas.forbriger's avatar
thomas.forbriger committed
41
c  V1.17  09/10/01   * allow long postscript file names
thomas.forbriger's avatar
thomas.forbriger committed
42
43
44
45
c  V1.18  04/04/02   new options:
c                    -T set title
c                    -X set x label
c                    -S shift time axis
thomas.forbriger's avatar
thomas.forbriger committed
46
c  V1.19             new options:
thomas.forbriger's avatar
thomas.forbriger committed
47
c                    -E force pgenv and pglab
thomas.forbriger's avatar
thomas.forbriger committed
48
49
c                    -x select time interval
c                    -y select amplitude interval
thomas.forbriger's avatar
thomas.forbriger committed
50
c  V1.19a 05/04/02   -E character scaling was incomplete
thomas.forbriger's avatar
thomas.forbriger committed
51
c  V1.19b 17/06/02   longer annotation string
thomas.forbriger's avatar
thomas.forbriger committed
52
c  V1.20  26/06/02   new option -W
53
c  V1.21  28/12/03   new skipdata subroutine
thomas.forbriger's avatar
thomas.forbriger committed
54
c  V1.21a 13/01/04   improved time scale
thomas.forbriger's avatar
thomas.forbriger committed
55
c  V1.21b 27/01/04   increased title length
thomas.forbriger's avatar
thomas.forbriger committed
56
c  V1.22  04/03/04   new option -c1
thomas.forbriger's avatar
thomas.forbriger committed
57
c  V1.23  26/09/05   new option -oS (stops line style cycling)
58
c  V1.24  29/11/05   provide new option for clean captions
thomas.forbriger's avatar
thomas.forbriger committed
59
c  V1.25  29/06/07   provide station location in caption
thomas.forbriger's avatar
thomas.forbriger committed
60
c  V1.26  02/07/07   provide option to set tick intervals
thomas.forbriger's avatar
thomas.forbriger committed
61
62
c  V1.27  18/01/08   do not use write statement for annotations
c                    just assign to character variable line
thomas.forbriger's avatar
thomas.forbriger committed
63
c  V1.28  04/12/09   use correct DIN notation for units
64
65
66
67
68
69
70
71
72
c 
c======================================================================
      program stuplo
c 
c variable declaration
c --------------------
c 
c version
      character*77 version, creator
73
      parameter(version=
thomas.forbriger's avatar
thomas.forbriger committed
74
     &  'STUPLO  V1.28  plot seismic time series (SFF format)')
75
76
77
      parameter(creator='1996 by Thomas Forbriger (IfG Stuttgart)')
c parameter definitions
      integer maxsamples, maxselect, lu, maxtraces, maxchain, maxstyle
thomas.forbriger's avatar
ohj    
thomas.forbriger committed
78
      parameter(maxsamples=6 000 000)
thomas.forbriger's avatar
thomas.forbriger committed
79
      parameter(maxtraces=200)
80
81
82
83
84
85
86
87
88
89
90
91
92
      parameter(maxselect=200)
      parameter(maxchain=10)
      parameter(maxstyle=5)
      parameter(lu=10)
      real capheight, labheight, numheight
      parameter(capheight=1.5, labheight=1.5, numheight=2.)
c chaining of traces
      integer nchain, ichain, npanels, ipanel, istyle
      integer firstic(maxchain), tracesic(maxchain)
c counters
      integer trace, filep, ftrace
c variables
      logical debug, verbose, goahead
thomas.forbriger's avatar
thomas.forbriger committed
93
      character line*80
94
95
96
97
98
99
100
101
102
103
104
      integer ntrim, ntraces, i, n
c here is the one big array to contain all seismic traces
c and this big array will hold real and integer data together
c the way we do this both real and integer type variables
c must allocate 4 bytes of memory each! 
      real data(maxsamples)
      integer idata(maxsamples)
      equivalence (idata, data)
c arras to hold information on data
      integer nsamples(maxtraces), firstsample(maxtraces)
      integer traceinfile(maxtraces)
thomas.forbriger's avatar
thomas.forbriger committed
105
      character*80 filename(maxtraces), firstfree(maxtraces)
106
107
108
109
110
111
      character*10 date(maxtraces)
      character*12 time(maxtraces)
      character*5 station(maxtraces)
      character*3 channel(maxtraces)
      character*4 auxid(maxtraces)
      character*6 instype(maxtraces)
thomas.forbriger's avatar
thomas.forbriger committed
112
113
      character*1 loccs(maxtraces)
      real locc1(maxtraces), locc2(maxtraces)
114
115
116
      real sectime(maxtraces), dt(maxtraces) 
      real maxval(maxtraces), average(maxtraces), minval(maxtraces)
c file reading variables
thomas.forbriger's avatar
thomas.forbriger committed
117
      character*200 infile
118
119
120
121
122
123
      logical moretraces
c using selections
      logical useselect
      logical selection(maxselect)
c plot caption
      character captionsel*20, captionstr*250, unitslabel*40
thomas.forbriger's avatar
thomas.forbriger committed
124
      character annotation*100
125
126
127
128
129
130
131
      integer caplen
c marker
      logical optmarker
      real xmark
c plot style
      logical optgrid, optscalex, optscaley, optinline, optabstime
      logical opttbox, ftoplab, fbotlab, optfixed, optcolor, optcenter
thomas.forbriger's avatar
thomas.forbriger committed
132
      logical optblack, optvertlab, opttitle, optxlabel, optpgenv
thomas.forbriger's avatar
thomas.forbriger committed
133
      logical optnoxlabels, optsetxrange, optsetyrange, optwhitepaper
134
      logical optcapfromfirst, optcyclestyle, optcleancaprect
135
136
137
138
      real lwf_std,lwf_txt,lwf_cur
      real chf_std,chf_xlab,chf_ylab,chf_cap
      integer lw_std,lw_txt,lw_cur
      real ch_std,ch_xlab,ch_ylab,ch_cap
thomas.forbriger's avatar
thomas.forbriger committed
139
      character*200 partitle, parxlabel
thomas.forbriger's avatar
thomas.forbriger committed
140
      real partimeoff, parxmin, parxmax, parymin, parymax
thomas.forbriger's avatar
thomas.forbriger committed
141
142
      real majorxticks,minorxticks
      integer nminorxticks,nmajorxticks
143
144
145
146
147
148
149
150
151
152
c cursor action
      logical optkeep
      integer curresult, pgcurs
      character curchar*2
      real cursx, cursy
c plot range
      real xmaxp, xminp, ymaxp, yminp, xmaxt, xmint, ymaxt, ymint
      real nlabels, Yfac, yamp, ymid, xlfac, xrfac, xamp
      real labelheight,fractop, fracbot, vptop, vpbot, boxratio
      real yvalheight
153
154
155
156
c box coordinates
      real boxxcoor(4), boxycoor(4)
c variables to do some intermediate calculations
      real jx1,jx2,jx3,jx4,jx5
157
158
159
160
161
c variables related to pgplot
      character*80 device, botlabel
      real timeax(maxsamples)
c commandline
      integer maxopt, lastarg, iargc
thomas.forbriger's avatar
thomas.forbriger committed
162
      parameter(maxopt=36)
thomas.forbriger's avatar
thomas.forbriger committed
163
      character*3 optid(maxopt)
thomas.forbriger's avatar
thomas.forbriger committed
164
      character*200 optarg(maxopt)
165
166
167
      logical optset(maxopt), opthasarg(maxopt)
c here are the keys to our commandline options
      data optid/2h-d,2h-D,2h-g,2h-c,2h-u,2h-s,2h-i,2h-a,2h-t,2h-f,2h-A,2h-k,
thomas.forbriger's avatar
thomas.forbriger committed
168
     &  2h-C,2h-Y,2h-R,2h-L,2h-z,2h-m,2h-v,2h-N,2h-l,2h-h,2h-V,2h-X,
thomas.forbriger's avatar
thomas.forbriger committed
169
     &  2h-T,2h-S,2h-E,3h-nT,2h-x,2h-y,2h-W,3h-Fc,3h-oS,3h-wc,3h-xt,3h-nx/
170
      data opthasarg/.TRUE.,2*.FALSE.,3*.TRUE.,4*.FALSE.,.TRUE.,2*.FALSE.,
thomas.forbriger's avatar
thomas.forbriger committed
171
     &  3*.TRUE.,.FALSE.,.TRUE.,2*.FALSE.,2*.TRUE.,.FALSE.,3*.TRUE.,
thomas.forbriger's avatar
thomas.forbriger committed
172
     &  2*.FALSE.,2*.TRUE.,4*.FALSE.,2*.TRUE./
173
      data optarg/3hx11,2*1h-,3hfTt,1h ,1hy,4*1h-,1h*,2*1h-,2h1.,2*2h0.,1h-,
thomas.forbriger's avatar
thomas.forbriger committed
174
     &  2h0.,2*1h-,8h1.,1.,1.,11h1.,1.,1.,1.,1h-,2*1h-,2h0.,2*1h-,
thomas.forbriger's avatar
thomas.forbriger committed
175
     &  2*5h0.,1.,4*1h-,4h0.,0,4h0,0./
176
177
178
179
180
181
182
183
184
185
186
187
188
189
c----------------------------------------------------------------------
c give basic information and help
      line=' '
      if (iargc().eq.1) call getarg(1, line)
c----------------------------------------------------------------------
c give help information
      if ((line(1:5).eq.'-help').or.(iargc().lt.1)) then
        print 80,version
        print *,creator
        print 80,'Usage: stuplo [-d dev] [-g] [-a] [-i] [-t] [-k] [-N]'
        print 80,'              [-Y fac] [-R fac] [-L fac] [-z] [-C]'
        print 80,'              [-c options] [-s x|y|xy] [-u units]'
        print 80,'              [-f] [-A comment] [-m time] [-v] [-V]'
        print 80,'              [-l std,lab,cur] [-h std,xlab,ylab,cap]'
thomas.forbriger's avatar
thomas.forbriger committed
190
        print 80,'              [-X label] [-T title] [-S sec] [-E]'
thomas.forbriger's avatar
thomas.forbriger committed
191
        print 80,'              [-nT] [-x f,t] [-y f,t] [-W] [-Fc] [-oS]'
thomas.forbriger's avatar
thomas.forbriger committed
192
        print 80,'              [-wc] [-xt i,n] [-nx n,i]'
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
        print 80,'              file [list] [nc:] ...'
        print 80,'or:    stuplo -help'
        if (iargc().lt.1) stop 'ERROR: missing parameters\n'
        print *,' '
        print *,'STUPLO plots seismic time series that are stored'
        print *,'in SFF format (Stuttgart File Format).'
        print *,' '
        print *,'-d dev       This option selects the pgplot output'
        print *,'             device. Information on available ouput'
        print *,'             devices is given below.'
        print *,'-v           be verbose'
        print *,'-k           keep program and window active'
        print *,'             This is useful for opening more than one'
        print *,'             pgplot window. Stop action by pressing x'
        print *,'             or the right mouse button.'
        print *,' '
        print *,'markers'
        print *,'-------'
        print *,' '
        print *,'-g           plot grid'
        print *,'-m time      plot a marker at specified time'
thomas.forbriger's avatar
thomas.forbriger committed
214
215
216
        print *,'             ''time'' is fed directly into the PGPLOT'
        print *,'             coordinate system. It has to be given in'
        print *,'             seconds on the abscissa-scale.'
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
        print *,' '
        print *,'labels'
        print *,'------'
        print *,' '
        print *,'-t           use full time labels for x-axis (tbox style)'
        print *,'             rather than counting seconds'
        print *,'-u units     define label for y-axis'
        print *,'-V           write numeric labels along left or right'
        print *,'             viewport boundary perpendicular to axis'
        print *,'             (vertically orientated) rather than parallel'
        print *,'             to the axis.'
        print *,'-c options   define caption string'
        print *,'             the options string may consist of the'
        print *,'             following letters:'
        print *,'             f   filename'
        print *,'             t   time'
        print *,'             d   date'
        print *,'             n   nsamples'
        print *,'             i   sampling interval'
        print *,'             s   station'
        print *,'             I   instrument'
        print *,'             c   channel'
        print *,'             a   GSE2.0 auxid'
        print *,'             m   minimum and maximum value'
        print *,'             F   first line of FREE block'
        print *,'             T   number of trace'
        print *,'             A   annotation string given with -A option'
thomas.forbriger's avatar
thomas.forbriger committed
244
        print *,'             L   location of station'
245
        print *,'             default is: fTt'
246
        print *,'-wc          print caption on white box, if inside panel'
thomas.forbriger's avatar
thomas.forbriger committed
247
        print *,'-Fc          create captions only from first trace in panel'
248
249
        print *,'-A comment   annotate comment to caption string'
        print *,'-i           plot caption inside box'
thomas.forbriger's avatar
thomas.forbriger committed
250
251
252
        print *,'-T title     use this string instead of the automatically'
        print *,'             generated title'
        print *,'-X label     use this label for the time axis'
thomas.forbriger's avatar
thomas.forbriger committed
253
        print *,'-nT          use this to ommit any time scale'
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
        print *,' '
        print *,'scales and ranges'
        print *,'-----------------'
        print *,' '
        print *,'-a           use absolute time scale'
        print *,'-s x|y|xy    use same scale for all boxes'
        print *,'             x   all x-axis will have the same range'
        print *,'             y   all y-axis will have the same range'
        print *,'             xy  all x-axis and y-axis will have'
        print *,'                 the same range'
        print *,'-Y fac       The calculated maximum amplitude (rule is'
        print *,'             given by -s option) will fill the fraction'
        print *,'             fac of the plot panel (default if fac=1.).'
        print *,'-R fac       The fraction of the seismogram to cut off'
        print *,'             at the end (default fac=0.).'
        print *,'-L fac       The fraction of the seismogram to cut off'
        print *,'             at the beginning (default fac=0.).'
        print *,'-z           Set zero-line always to the center of the'
        print *,'             plot panel.'
thomas.forbriger's avatar
thomas.forbriger committed
273
        print *,'-S sec       add sec (in seconds) to the time values'
thomas.forbriger's avatar
thomas.forbriger committed
274
275
        print *,'-x f,t       set time scale from f to t (s)'
        print *,'-y f,t       set amplitude scale from f to t (s)'
thomas.forbriger's avatar
thomas.forbriger committed
276
277
278
        print *,'-xt i,n      set major tickmark interval of time axis'
        print *,'             to ''i'' and the number minor tick marks'
        print *,'             between major tick marks to ''n'' '
thomas.forbriger's avatar
thomas.forbriger committed
279
280
281
282
        print *,'-nx n,i      use ''n'' major tick marks with a minor'
        print *,'             tick marks with a minor tick mark'
        print *,'             interval of ''i'' '
        print *,'             This options overrides option -xt'
283
284
285
286
        print *,' '
        print *,'styles'
        print *,'------'
        print *,' '
thomas.forbriger's avatar
thomas.forbriger committed
287
        print *,'-W           plot black on white (also on screen)'
288
289
        print *,'-C           use different colours instead of different'
        print *,'             linestyles for marking different chains.'
thomas.forbriger's avatar
thomas.forbriger committed
290
291
        print *,'-oS          use only one style for all curves (no'
        print *,'             cycling of line style or colour).'
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
        print *,'-N           use non-black colour (or linestyle) for'
        print *,'             first dataset.'
        print *,'-l std,txt,cur'
        print *,'             The default linewidth for the following'
        print *,'             elements is scaled by the given factor:'
        print *,'             txt   text strings'
        print *,'             cur   time series curve'
        print *,'             std   all others'
        print *,'-h std,xlab,ylab,cap'
        print *,'             The default character height for the following'
        print *,'             elements is scaled by the given factor:'
        print *,'             xlab  y-axis labels'
        print *,'             ylab  x-axis labels'
        print *,'             cap   caption string'
        print *,'             std   all others'
        print *,'-f           used fixed character size'
thomas.forbriger's avatar
thomas.forbriger committed
308
309
310
311
        print *,'-E           This is a workaround for difficult cases in'
        print *,'             which the -h option does not lead to the'
        print *,'             desired result. This option forces stuplo'
        print *,'             to use a single panel and pgenv and pglab.'
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
        print *,' '
        print *,' '
        print *,'Each datafile name may be followed by a list of'
        print *,'traces. This list selects a range of traces in'
        print *,'the file which will be processed. This list may'
        print *,'contain no blank (which is the separator to the'
        print *,'next filname). The traces will always be processed'
        print *,'in the order they appear in the data file.'
        print *,' '
        print *,'Examples:'
        print *,'  t:2           will select only trace 2'
        print *,'  t:4-6,2,4     will select traces 2, 4, 5 and 6'
        print *,'  t:9,8,10,14   will select traces 8, 9, 10 and 14'
        print *,' '
        print *,'For comparison purposes it is usefull to do overlay plots.'
        print *,'This may be achieved be building severl chains of seismic'
        print *,'traces. To start a new chain in the list of files use'
        print *,'''nc:'' as a separator. Using this feature the first'
        print *,'panel will contain all first files of each chain, the'
        print *,'second one all second files, etc. Each trace will be plotted'
        print *,'with a different linestyle within each panel. Using the'
        print *,'option -C there will be a different color for each trace.'
        print *,' '
        print *,'Example for building two chains of different length with'
        print *,'the datafiles ''data1'' and ''data2'':'
        print *,' '
        print *,'    stuplo data1 t:1,2 data2 t:5 nc: data1 t:4-7'
        print *,' '
        print *,'  This is what the plotting panels will contain:'
        print *,' '
        print *,'    panel  file.trace'
        print *,'        1  data1.1 data1.4'
        print *,'        2  data1.2 data1.5'
        print *,'        3  data2.5 data1.6'
        print *,'        4  data1.7'
        print *,' '
        call pgp_showdevices
        print *,' '
        print *,'Compiled array dimensions are:'
        print *,'                 total number of samples:',maxsamples
        print *,'                 (all traces together)'
        print *,'      maximum number of traces to select:',maxselect
        print *,'        maximum number of traces to hold:',maxtraces
        print *,'       maximum number of chains to build:',maxchain
        print *,' maximum number of different plot styles:',maxstyle
        print *,' '
        stop
      endif
c----------------------------------------------------------------------
c set standard options
      device='x11'
      captionsel='fTt'
      unitslabel=' '
      optinline=.false.
      optscalex=.false.
      optscaley=.false.
      optabstime=.false.
      opttbox=.false.
c first read the commandline
      call tf_cmdline(1, lastarg,
     &     maxopt, optid, optarg, optset, opthasarg)
      if (optset(1)) device=optarg(1)
      debug=optset(2)
      optgrid=optset(3)
      if (optset(4)) captionsel=optarg(4)
      if (optset(5)) unitslabel=optarg(5)
      if (optset(6)) then
        if (index(optarg(6),'x').gt.0) optscalex=.true.
        if (index(optarg(6),'y').gt.0) optscaley=.true.
      endif
      optinline=optset(7)
      optabstime=optset(8)
      opttbox=optset(9)
      optfixed=optset(10)
      annotation=optarg(11)
      optkeep=optset(12)
      optcolor=optset(13)
thomas.forbriger's avatar
thomas.forbriger committed
389
      read(optarg(14), *, end=96, err=97) Yfac
390
391
392
393
      if (Yfac.lt.0.01) then
        print *,'WARNING: Y-fac is set to 1. (was too small: ',Yfac,')'
        Yfac=1.
      endif
thomas.forbriger's avatar
thomas.forbriger committed
394
395
      read(optarg(15), *, end=96, err=97) xrfac
      read(optarg(16), *, end=96, err=97) xlfac
396
397
398
399
400
401
      optcenter=optset(17)
      if ((xrfac+xlfac).ge.1.) 
     &  stop 'ERROR: check setting of -R and -L (nothing will be left)'
      if (debug) print *,'DEBUG: debug messages are switched on'
      if (debug) call checkmem(idata, data, maxsamples)
      optmarker=optset(18)
thomas.forbriger's avatar
thomas.forbriger committed
402
      if (optmarker) read(optarg(18), *, end=96, err=97) xmark
403
404
      verbose=optset(19)
      optblack=optset(20)
thomas.forbriger's avatar
thomas.forbriger committed
405
406
407
      read(optarg(21), *, end=96, err=97) lwf_std,lwf_txt,lwf_cur
      read(optarg(22), *, end=96, err=97) 
     &  chf_std,chf_xlab,chf_ylab,chf_cap
408
      optvertlab=optset(23)
thomas.forbriger's avatar
thomas.forbriger committed
409
410
411
412
413
      optxlabel=optset(24)
      parxlabel=optarg(24)
      opttitle=optset(25)
      partitle=optarg(25)
      read(optarg(26), *, end=96, err=97) partimeoff
thomas.forbriger's avatar
thomas.forbriger committed
414
415
      optpgenv=optset(27)
      optnoxlabels=optset(28)
thomas.forbriger's avatar
thomas.forbriger committed
416
417
418
419
      optsetxrange=optset(29)
      read(optarg(29), *, end=96, err=97) parxmin,parxmax
      optsetyrange=optset(30)
      read(optarg(30), *, end=96, err=97) parymin,parymax
thomas.forbriger's avatar
thomas.forbriger committed
420
      optwhitepaper=optset(31)
thomas.forbriger's avatar
thomas.forbriger committed
421
      optcapfromfirst=optset(32)
thomas.forbriger's avatar
thomas.forbriger committed
422
      optcyclestyle=(.not.optset(33))
423
      optcleancaprect=optset(34)
thomas.forbriger's avatar
thomas.forbriger committed
424
425
      read(optarg(35), *, end=96, err=97) majorxticks, nminorxticks
      read(optarg(36), *, end=96, err=97) nmajorxticks, minorxticks
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
c      device='/krm3'
c----------------------------------------------------------------------
c initialize
      do ichain=1,maxchain
        firstic(ichain)=-1
        tracesic(ichain)=0
      enddo
c----------------------------------------------------------------------
c 
c enter main program body
c -----------------------
c
c first of all read in the time series
c
      trace=0
      nchain=1
      moretraces=.false.
      filep=lastarg+1
      goahead=.true.
      if (iargc().lt.filep) stop 'ERROR: missing data file\n'
c start file loop
    1 continue
c open new file 
        call getarg(filep, infile)
c is there a new chain to start
        if (infile(1:3).eq.'nc:') then
          nchain=nchain+1
          if (nchain.gt.maxchain) stop 'ERROR: opened too many chains'
          goto 4
        endif
        call sff_TrimLen(infile,ntrim)
        if (verbose) print 81,'  opening data file ',infile(1:ntrim)
        open(lu, file=infile, err=99, status='old')
        call sffopen(lu, filep, useselect, selection, maxselect, debug)
        ftrace=0
c start trace loop
    2   continue
          ftrace=ftrace+1
          if (((.not.(useselect)).or.(selection(ftrace))).and.(goahead)) then
            if (trace.lt.maxtraces) then
              trace=trace+1
              if (firstic(nchain).lt.0) then
c start a new chain with this trace
                firstic(nchain)=trace
              endif
              tracesic(nchain)=tracesic(nchain)+1
              traceinfile(trace)=ftrace
              if (verbose) print 82,'    reading trace no.',ftrace
              call sffread(lu, infile, trace, moretraces, debug,
     &                     maxsamples, maxtraces, filename, dt, timeax,
     &                     firstsample, nsamples, firstfree, date,
     &                     time, sectime, station, channel, auxid,
     &                     instype, data, idata, maxval, average, minval,
thomas.forbriger's avatar
thomas.forbriger committed
479
480
     &                     optabstime, verbose, partimeoff,
     &                     loccs, locc1, locc2)
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
              if (debug) then
                print *,'DEBUG: dataset parameters:'
                write(6,'(" DEBUG:",6(2h >,a,1h<))')
     &            date(trace),time(trace),station(trace),
     &            channel(trace),auxid(trace),instype(trace)
                print *,'DEBUG: dt:',dt(trace),' sectime:',sectime(trace)
                print *,'DEBUG: firstfree: ',firstfree(trace)
                print *,'DEBUG: file name: ',filename(trace)
                print *,'DEBUG: first sample:',firstsample(trace),
     &                  ' last sample:',(firstsample(trace)-1+
     &                                   nsamples(trace)),
     &                  ' samples:',nsamples(trace)
                print *,'DEBUG: maxval:',maxval(trace),
     &                  ' minval:',minval(trace),
     &                  ' average:',average(trace)
              endif
            else
              print *,'WARNING: reached limit of',maxtraces,' traces'
              goahead=.false.
            endif
          else
            if (verbose) print 82,'    skipping trace no.',ftrace
            call skipdata(lu, moretraces)
          endif
          if (moretraces) goto 2
c close data file
        if (verbose) print 81,'  closing data file ',infile(1:ntrim)
        close(lu, err=98)
    4   filep=filep+1
        if (filep.le.iargc()) goto 1
      ntraces=trace
c----------------------------------------------------------------------
c maximum number of plot panels will be:
      npanels=0
      do ichain=1,nchain
        npanels=max(npanels,tracesic(ichain))
      enddo
c----------------------------------------------------------------------
c 
c all traces are read
c
c prepare data for plotting now
c
c calculate axis range
c
      xminp=timeax(firstsample(1))
      xmaxp=timeax(firstsample(1)-1+nsamples(1))
      yminp=minval(1)
      ymaxp=maxval(1)
      do trace=1,ntraces
        xminp=min(xminp,timeax(firstsample(trace)))
        xmaxp=max(xmaxp,timeax(firstsample(trace)-1+nsamples(trace)))
        yminp=min(yminp,minval(trace))
        ymaxp=max(ymaxp,maxval(trace))
      enddo
c 
c determine fractioning of surface
c 
c the way of labeling:
c
c   sx: optscalex
c   il: optinline
c   cp: caption
c   nu: numeric labels
c   la: bottom label
c
c   sx    il       cp        nu        la
c
c    x     x        -      last      last
c    -     x        -       all       all
c    x     -      all       all      last
c    -     -      all       all       all
c
c we have to check headers and scale labels
      if ((optscalex).and.(optinline)) then
c no header and only one scale label line at the bottom
        nlabels=numheight+labheight
      elseif (optinline) then
c each trace has a scale label and the last one has a label at the bottom
        nlabels=npanels*(numheight+labheight)
      elseif (optscalex) then
c each trace has a caption and numeric labels and
c the last one has a label at the bottom
        nlabels=npanels*(capheight+numheight)+labheight
      else
c each trace has caption, numeric labels and bottom label
        nlabels=npanels*(capheight+numheight+labheight)
      endif
      boxratio=40./npanels
      labelheight=1./(nlabels+npanels*boxratio)
      yvalheight=labelheight*20.
      if (optfixed) then
        labelheight=1./40.
        boxratio=((1/labelheight)-nlabels)/npanels
      endif
c no viewport up to now --> set last viewport value to top
      vpbot=0.99
c 
c----------------------------------------------------------------------
c plot style
c   linewidth
      lw_txt=lwf_txt
      lw_cur=lwf_cur
      lw_std=lwf_std
c   character height
thomas.forbriger's avatar
thomas.forbriger committed
586
587
588
589
590
591
592
593
594
595
596
      if (optpgenv) then
        ch_std=chf_std
        ch_xlab=chf_xlab
        ch_ylab=chf_ylab
        ch_cap=chf_cap
      else
        ch_std=chf_std*labelheight*30.
        ch_xlab=chf_xlab*labelheight*30.
        ch_ylab=chf_ylab*yvalheight
        ch_cap=chf_cap*labelheight*30.
      endif
597
598
599
600
601
602
603
604
605
606
607
c
c----------------------------------------------------------------------
c 
c do plot
c
      call pgp_setdevice(device,1,1)
      call pgask(.false.)
      call pgscr(3, 0.,0.,1.)
      call pgscr(4, 0.,1.,0.)
      call pgslw(lw_std)
      call pgsch(ch_std)
thomas.forbriger's avatar
thomas.forbriger committed
608
609
610
611
612
613
614
615
616
c 
c prepare background
c ------------------
c
      if (optwhitepaper) then
        call pgscr(0, 1.,1.,1.)
        call pgscr(1, 0.,0.,0.)
      endif
c
617
618
619
620
621
622
623
624
625
      do ipanel=1,npanels
c
c build caption string
c
        caplen=1
        call sff_TrimLen(captionsel,n)
        do ichain=1,nchain
c is there a trace to plot in this chain?
          if (ipanel.gt.tracesic(ichain)) goto 5
thomas.forbriger's avatar
thomas.forbriger committed
626
c or skip if only first is selcted
thomas.forbriger's avatar
thomas.forbriger committed
627
          if ((ichain.ne.1).and.(optcapfromfirst)) goto 5
628
629
630
            trace=firstic(ichain)-1+ipanel
            do i=1,n
              if (captionsel(i:i).eq.'f') then
thomas.forbriger's avatar
thomas.forbriger committed
631
                line=filename(trace)
632
              elseif(captionsel(i:i).eq.'t') then
thomas.forbriger's avatar
thomas.forbriger committed
633
                line=time(trace)
634
635
                if (debug) print *,'DEBUG: time ',time(trace)
              elseif(captionsel(i:i).eq.'d') then
thomas.forbriger's avatar
thomas.forbriger committed
636
                line=date(trace)
637
638
639
640
641
642
              elseif(captionsel(i:i).eq.'T') then
                write(line,'("(",i2,")")') traceinfile(trace)
              elseif(captionsel(i:i).eq.'n') then
                write(line,'("n=",i8)') nsamples(trace)
                if (debug) print *,'DEBUG: nsamples ',nsamples(trace)
              elseif(captionsel(i:i).eq.'s') then
thomas.forbriger's avatar
thomas.forbriger committed
643
                line=station(trace)
644
645
646
              elseif(captionsel(i:i).eq.'i') then
                write(line,'("dt=",e12.3)') dt(trace)
              elseif(captionsel(i:i).eq.'I') then
thomas.forbriger's avatar
thomas.forbriger committed
647
                line=instype(trace)
648
              elseif(captionsel(i:i).eq.'c') then
thomas.forbriger's avatar
thomas.forbriger committed
649
                line=channel(trace)
650
              elseif(captionsel(i:i).eq.'a') then
thomas.forbriger's avatar
thomas.forbriger committed
651
                line=auxid(trace)
652
              elseif(captionsel(i:i).eq.'A') then
thomas.forbriger's avatar
thomas.forbriger committed
653
                line=annotation
thomas.forbriger's avatar
thomas.forbriger committed
654
655
656
657
658
659
660
661
              elseif(captionsel(i:i).eq.'L') then
                if (loccs(trace).eq.'S') then
                  write(line,'(f8.3,"N; ",f8.3,"E")') 
     &              locc1(trace), locc2(trace)
                elseif (loccs(trace).eq.'C') then
                  write(line,'("x="f8.3,"m; y=",f8.3,"m")') 
     &              locc1(trace), locc2(trace)
                else
thomas.forbriger's avatar
thomas.forbriger committed
662
                  line='no location'
thomas.forbriger's avatar
thomas.forbriger committed
663
                endif
664
665
666
667
              elseif(captionsel(i:i).eq.'m') then
                write(line,'("min/max:",e12.3,"/",e12.3)') 
     &          minval(trace),maxval(trace)
              elseif(captionsel(i:i).eq.'F') then
thomas.forbriger's avatar
thomas.forbriger committed
668
                line=firstfree(trace)
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
              else
                print *,'WARNING: unknown caption character: ',captionsel(i:i)
                line='*'
              endif
              call sff_TrimLen(line, ntrim)
              if (debug) print *,'DEBUG: ',captionsel(i:i)
     &          ,'>',line(1:ntrim),'<', line
              if ((caplen+ntrim).lt.len(captionstr)) then
                captionstr(caplen:ntrim+caplen-1)=line(1:ntrim)
                caplen=caplen+ntrim
                if ((caplen+ntrim+3).lt.len(captionstr)) then
                  captionstr(caplen:caplen+2)='   '
                  caplen=caplen+3
                endif
              else
                print *,'WARNING: caption is truncated'
              endif
c enddo i
            enddo
    5     continue
c enddo ichain
        enddo
thomas.forbriger's avatar
thomas.forbriger committed
691
692
693
694
695
696
        if (opttitle) then
          captionstr=partitle
        else
          do i=caplen,len(captionstr)
            captionstr(i:i)=' '
          enddo
thomas.forbriger's avatar
thomas.forbriger committed
697
        endif
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
c 
c determine range of box and open plot space
c
        ichain=0
    6   continue
          ichain=ichain+1
          if (ipanel.gt.tracesic(ichain)) goto 6
        trace=firstic(ichain)-1+ipanel
        xmint=timeax(firstsample(trace))
        xmaxt=timeax(firstsample(trace)-1+nsamples(trace))
        ymint=minval(trace)
        ymaxt=maxval(trace)
        do ichain=1,nchain
c is there a trace to plot in this chain?
          if (ipanel.le.tracesic(ichain)) then
            trace=firstic(ichain)-1+ipanel
            xmint=min(xmint,timeax(firstsample(trace)))
            xmaxt=max(xmaxt,timeax(firstsample(trace)-1+nsamples(trace)))
            ymint=min(ymint,minval(trace))
            ymaxt=max(ymaxt,maxval(trace))
          endif
c enddo ichain
        enddo
c may be we have to use global scales
        if (optscalex) then
          xmint=xminp
          xmaxt=xmaxp
        endif
        if (optscaley) then
          ymint=yminp
          ymaxt=ymaxp
        endif
c modify scales by selected factors
        if (optcenter) then
          ymaxt=max(abs(ymaxt),abs(ymint))
          ymint=-max(abs(ymaxt),abs(ymint))
        endif
        yamp=(ymaxt-ymint)/Yfac
        ymid=(ymaxt+ymint)/2.
        ymint=ymid-0.5*yamp
        ymaxt=ymid+0.5*yamp
        xamp=(xmaxt-xmint)
        xmint=xmint+xamp*xlfac
        xmaxt=xmaxt-xamp*xrfac
thomas.forbriger's avatar
thomas.forbriger committed
742
743
744
745
        if (optsetxrange) then
          xmint=parxmin
          xmaxt=parxmax
        endif
746
747
748
749
750
751
752
753
754
755
756
c be aware of zero amplitudes
        if (ymaxt.eq.ymint) then
          if(ymaxt.eq.0.) then
            ymaxt=.5e-30
            ymint=-.5e-30
            print *,'NOTICE: trace ',trace,' has amplitude 0.'
          else
            ymaxt=ymint+abs(.5*ymint)
            ymint=ymint-abs(.5*ymint)
          endif
        endif
thomas.forbriger's avatar
thomas.forbriger committed
757
758
759
760
        if (optsetyrange) then
          ymint=parymin
          ymaxt=parymax
        endif
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
        if (debug) print *,'DEBUG: trace, ymaxt, ymint ',trace,ymaxt,ymint
        call pgswin(xmint, xmaxt, ymint, ymaxt)
c 
c calc values for new vp
c 
c start at botom of last box
        vptop=vpbot
c is there a caption
        fractop=vptop-(capheight*labelheight)
        if (optinline) fractop=vptop
c get bottom of box
        fracbot=fractop-boxratio*labelheight
c is there a numeric label
        vpbot=fracbot
        if (.not.((optscalex).and.(optinline))) then
          vpbot=vpbot-(numheight*labelheight)
        elseif (trace.eq.ntraces) then
          vpbot=vpbot-(numheight*labelheight)
        endif
        if (.not.(optscalex)) then
          vpbot=vpbot-(labheight*labelheight)
        elseif (trace.eq.ntraces) then
          vpbot=vpbot-(labheight*labelheight)
        endif
        call pgsch(ch_std)
thomas.forbriger's avatar
thomas.forbriger committed
786
787
788
789
790
        if (optpgenv) then
          call pgvstd
        else
          call pgsvp(chf_ylab*labelheight*2.5,0.99,fracbot,fractop)
        endif
791
792
793
794
795
796
797
798
799
800
801
802
803
c
c determine labelling style
c
        fbotlab=.true.
        ftoplab=(.not.(optinline))
        if ((optinline).and.(optscalex)) then
          fbotlab=.false.
          ftoplab=.false.
        endif
        if (trace.eq.ntraces) fbotlab=.true.
        botlabel=' '
        if ((trace.eq.ntraces).or.(.not.(optscalex))) then
          if (optabstime) then
thomas.forbriger's avatar
thomas.forbriger committed
804
            botlabel='time since midnight / sec'
805
806
            if (opttbox) botlabel='time since midnight'
          else
thomas.forbriger's avatar
thomas.forbriger committed
807
            botlabel='time since first sample / sec'
808
809
            if (opttbox) botlabel='time since first sample'
          endif
thomas.forbriger's avatar
thomas.forbriger committed
810
          if (optxlabel) botlabel=parxlabel
811
812
813
814
        endif
c 
c determine box style and plot box
c
thomas.forbriger's avatar
thomas.forbriger committed
815
816
817
818
819
820
        if (nmajorxticks.ne.0) then
          if (minorxticks.lt.((xmaxt-xmint)/1.e3))
     &      stop 'ERROR: minor tick mark interval too small!'
          nminorxticks=int((xmaxt-xmint)/(nmajorxticks*minorxticks))
          majorxticks=nminorxticks*minorxticks
        endif
821
        if (opttbox) then
thomas.forbriger's avatar
thomas.forbriger committed
822
          if (fbotlab.and.(.not.optnoxlabels)) then
thomas.forbriger's avatar
thomas.forbriger committed
823
            call pgtbox('ZABCTSHO',majorxticks,nminorxticks,'BCTS',0.0,0)
824
825
            call pgslw(lw_txt)
            call pgsch(ch_xlab)
thomas.forbriger's avatar
thomas.forbriger committed
826
            call pgtbox('ZNHO',majorxticks,nminorxticks,' ',0.0,0)
827
828
829
            call pgslw(lw_std)
            call pgsch(ch_std)
          else
thomas.forbriger's avatar
thomas.forbriger committed
830
            call pgtbox('ZABCTSHO',majorxticks,nminorxticks,'BCTS',0.0,0)
831
832
          endif
        else
thomas.forbriger's avatar
thomas.forbriger committed
833
          if (fbotlab.and.(.not.optnoxlabels)) then
thomas.forbriger's avatar
thomas.forbriger committed
834
            call pgtbox('ABCTSYHO',majorxticks,nminorxticks,'BCTS',0.0,0)
835
836
            call pgslw(lw_txt)
            call pgsch(ch_xlab)
thomas.forbriger's avatar
thomas.forbriger committed
837
            call pgtbox('NY',majorxticks,nminorxticks,' ',0.0,0)
838
839
840
            call pgslw(lw_std)
            call pgsch(ch_std)
          else
thomas.forbriger's avatar
thomas.forbriger committed
841
            call pgtbox('ABCTSYHO',majorxticks,nminorxticks,'BCTS',0.0,0)
842
843
844
845
846
847
          endif
        endif
        call pgsch(yvalheight)
        call pgslw(lw_txt)
        call pgsch(ch_ylab)
        if (optvertlab) then
thomas.forbriger's avatar
thomas.forbriger committed
848
          call pgtbox(' ',majorxticks,nminorxticks,'NV',0.0,0)
849
        else
thomas.forbriger's avatar
thomas.forbriger committed
850
          call pgtbox(' ',majorxticks,nminorxticks,'N',0.0,0)
851
852
853
854
855
856
857
        endif
        call pgslw(lw_std)
        call pgsch(ch_std)
        call pgsch(labelheight*30.)
        if (optgrid) then
          call pgsls(4)
          if (opttbox) then
thomas.forbriger's avatar
thomas.forbriger committed
858
            call pgtbox('ZG',majorxticks,nminorxticks,'G',0.0,0)
859
          else
thomas.forbriger's avatar
thomas.forbriger committed
860
            call pgtbox('G',majorxticks,nminorxticks,'G',0.0,0)
861
862
863
864
865
866
867
868
          endif
          call pgsls(1)
        endif
c 
c plot a marker at specified time if desired
c 
      if (optmarker) then
        call pgsave
thomas.forbriger's avatar
thomas.forbriger committed
869
        call pgsch(ch_std)
870
871
872
873
874
875
876
877
878
879
880
881
        call pgpt1(xmark, ymaxt, -31)
        call pgpt1(xmark, ymint, -31)
        call pgslw(3)
        call pgsls(4)
        call pgmove(xmark,ymint)
        call pgdraw(xmark,ymaxt)
        call pgunsa
      endif
c 
c plot labels
c 
        call pgslw(lw_txt)
thomas.forbriger's avatar
thomas.forbriger committed
882
883
884
885
886
887
888
889
890
        if (optpgenv) then
          call pgsch(ch_xlab)
          call pglab(botlabel, ' ', ' ')
          call pgsch(ch_ylab)
          call pglab(' ', unitslabel, ' ')
          call pgsch(ch_cap)
          if (ftoplab) then
            call pglab(' ', ' ', captionstr)
          else
891
892
893
            if (.not.optcleancaprect) then
              call pgmtxt('LV',-1.,0.9,0.0,captionstr)
            endif
thomas.forbriger's avatar
thomas.forbriger committed
894
          endif
895
        else
thomas.forbriger's avatar
thomas.forbriger committed
896
897
898
899
900
901
902
903
          call pgsch(ch_xlab)
          call pgmtxt('B',(2.9-(chf_xlab-1.)),0.5,0.5,botlabel)
          call pgsch(ch_ylab)
          call pgmtxt('L',(3.2-(chf_ylab-1.)),0.5,0.5,unitslabel)
          call pgsch(ch_cap)
          if (ftoplab) then
            call pgmtxt('T',0.7,0.5,0.5,captionstr)
          else
904
905
906
            if (.not.optcleancaprect) then
              call pgmtxt('LV',-1.,(0.9-(chf_cap-1.)*0.025),0.0,captionstr)
            endif
thomas.forbriger's avatar
thomas.forbriger committed
907
          endif
908
909
910
911
912
913
914
915
916
917
918
919
        endif
        call pgslw(lw_std)
        call pgsch(ch_std)
c 
c plot time series
c
        call pgslw(lw_cur)
        istyle=0
        if (optblack) istyle=1
        do ichain=1,nchain
          istyle=istyle+1
          if (istyle.gt.maxstyle) istyle=1
thomas.forbriger's avatar
thomas.forbriger committed
920
921
922
923
924
925
          if (optcyclestyle) then
            if (optcolor) then
              call pgsci(istyle)
            else
              call pgsls(istyle)
            endif
926
          else
thomas.forbriger's avatar
thomas.forbriger committed
927
928
929
            istyle=1
            if (optblack) istyle=2
            call pgsci(istyle)
930
931
932
933
934
935
936
937
938
939
940
          endif
c is there a trace to plot in this chain?
          if (ipanel.le.tracesic(ichain)) then
            trace=firstic(ichain)-1+ipanel
            call pgline(nsamples(trace), 
     &        timeax(firstsample(trace)),
     &        data(firstsample(trace)))
            call pgupdt
          endif
c enddo ichain 
        enddo
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
c if optcleancaprect print caption here
        if (optcleancaprect) then
c   world coordinate ranges
          jx1=xmaxt-xmint
          jx2=ymaxt-ymint
c   world coordinate offsets
          jx3=0.02*jx1
          jx4=0.02*jx2
c   first read size of box
          call pgqtxt(xmint+jx3,ymint+jx4,0.,0.,captionstr,
     &      boxxcoor, boxycoor)
          jx5=boxycoor(2)-boxycoor(1)
          call pgqtxt(xmint+jx3,ymaxt-jx4-jx5,0.,0.,captionstr,
     &      boxxcoor, boxycoor)
          call pgsave
          call pgsfs(1)
          call pgsci(0)
c   clear background of cation's box
          call pgpoly(4,boxxcoor,boxycoor)
          call pgunsa
c   write caption
          call pgptxt(xmint+jx3,ymaxt-jx4-jx5,0.,0.,captionstr)
        endif
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
        call pgsls(1)
        call pgsci(1)
        call pgslw(lw_std)
c enddo ipanel
      enddo
c 
c keep alive
c
      if (optkeep) then
        cursx=0.
        cursy=0.
    3   curresult=pgcurs(cursx, cursy, curchar)  
        if ((curresult.eq.1).and.(curchar.ne.'x').and.(curchar.ne.'X'))
     &    goto 3
      endif
      call pgend
      stop 
c----------------------------------------------------------------------
c formats
   80 format(a)
   81 format(2a)
thomas.forbriger's avatar
thomas.forbriger committed
985
   82 format(a,i4)
986
c stop sequences
thomas.forbriger's avatar
thomas.forbriger committed
987
988
989
990
   99 stop 'ERROR: opening data file!'
   98 stop 'ERROR: closing data file!'
   97 stop 'ERROR: in command line option arguments!'
   96 stop 'ERROR: too few command line option arguments!'
991
992
993
994
995
996
997
998
999
1000
      end
      
c======================================================================
c 
c some subroutines
c
c----------------------------------------------------------------------
c
c checkmem
c
For faster browsing, not all history is shown. View entire blame