Modified recfunk_pick_through.f

click title to see Google Drive directory JParkCodes

More options for specifying the start time of the data record are recognised.

c program recfunk_pick_through
c version of recfunk_pick that reads seismic record start times from
c the sac header, and writes to in_recpick the filenames that have been picked
c using EITHER the A-marker or the T1, T2, or T3 markers in the SAC header, placed by running plotpk
c within SAC. Set by pressing "T" and 1,2 or 3 in plotpk sequentially
c 10/13/05 JJP
c input file: in_recpick_raw
c output file: in_recpick
c
c xf77 -o /park/backup/bin/recfunk_pick_through recfunk_pick_through.f /park/backup/Ritz/jlib.a
c xf77 -o /Users/jjpark/bin/recfunk_pick_through recfunk_pick_through.f /Users/jjpark/Ritz/jlib.a
c
......implicit real*4(a-h,o-z)
......real*8 ar,ai,el(12),aaa
......complex*8 zc,zero
......complex*8 afft,rf,rfs
......character*80 name,subtext,string
......character*28 title(2)
......character*10 cmps(3)
......character comp(3)
......common/nnaamm/iwflag
......common/npiprol/anpi
......common/datastuff/a(99000,3),pdat(16400),tdat(16400),drf(4100,2)
......common/header/ahead(158)
......dimension iah(158),chead(158),fr(3),tt(3),ttt(3)
......equivalence (ahead,iah),(ahead,chead)
......data comp/'r','t','z'/
......print *,'This code reads start of interval from sac header'
......print *,'it writes filenames of those files that have been picked'
......open(10,file='in_recpick_raw',form='formatted')
......open(11,file='in_recpick',form='formatted')
10 print *,'enter filename e.g. 1998.361.bhz'
......read(10,101) name
101 format(a)
......if(name(1:4).eq.'stop') go to 111
......subtext=name
c we assume that the last character is either 'r','t',or 'z'
......do i=80,1,-1
......if(name(i:i).ne.' ') then
......namlen=i
......go to 98
......endif
......end do
......print *,name
......print *,'no filename?'
......stop
98 continue
......do kcmp=3,3
......name(namlen:namlen)=comp(kcmp)
......print *,name
......call sacread(name,a(1,kcmp),ierr)
......if(ierr.eq.1) stop
......dt=ahead(1)
......tt(1)=0.
......tt(2)=dt
......tt(3)=dt
......nscan=iah(80)
......if(nscan.gt.99000) then
......print *,'Careful! data length greater than 99000:', nscan
......pause
......endif
......end do
c end kluge
c if(baz.lt.70) then
c call plotit_axes(0.,0.,0.,0.)
c call manyplot(99000,1,nscan,3,tt,a,1.1,name,'time(sec)',' ')
c endif
c plot kluge
c KLUGE: CHANGES FOR READING TIME INTERVAL FROM SAC HEADER
c ahead(9) is SAC's A parameter, for arrival time. Set by pressing A in plotpk
c ahead(12,13,14) are SAC's T1,T2,T3 parameters, for arrival time. Set by pressing "T" and "N" in plotpk
......if(ahead(9).gt.1.0.or.ahead(12).gt.1.0.or.ahead(13).gt.1.0.
......x or.ahead(14).gt.1.0) write(11,101) subtext
......pstart=ahead(9)-ahead(6)-3.
......print *,'sach params b and a, pstart',ahead(6),ahead(9),pstart,
......x (ahead(n),n=12,14)
......go to 10
111 write(11,101) name
......close(10)
......close(11)
......stop
......end