subroutine sum2cnss (card, ounit, sunit, eunit, date,cmag,amag, 1 zmag, lmag, wmag, iyear, imonth, iday, ihour, imin, sec, ncodv, 2 ncodz, nampl, nampx, y2k, ierr, typmag) c c reformat NCSN summary lines into CNSS unified format c 2008/07/02 DSN Add dc eventid field (copy of local eventid) c implicit none real amag ! (output) Eaton amplitude magnitude real cmag ! (output) Eaton coda magnitude real zmag ! (output) Hirshorn Z coda mag real lmag ! (output) USGS computed WA mag real wamag ! UCB external WA mag real wmag ! (output) Moment mag integer nampx ! (input) number of Samp readings used in Eaton amplitude magnitude integer nampl ! (input) number of WA readings used in local (ML) magnitude integer ncodv ! (input) number of coda readings used in high-gain magnitude integer ncodz ! (input) number of coda readings used in low-gain magnitude logical y2k ! y2k format found for hypoinverse file character*2 auxrmk integer azgap character*(*) card character*8 date real depth real dmin character*1 east real e(3,3) real errh character*5 errmag real errz integer eventid integer eunit integer i integer iaz(3) integer iday integer idip(3) integer ierr integer ihour integer imin integer imonth integer iq integer iyear integer j integer k real latdeg real latd real latmin character*3 locrmk character*2 loctyp real londeg real lond real lonmin real mag(5) ! magnitudes corresponding to magtype character*1 magtyp(5) ! 1=Xmag,2=Fmag,3=Ext_mag (eg.,UCB ML),4=Alt_amp (eg.,NCSN ML),5=Alt_coda (e.g,Zmag) real magerr(5) real magobs(5) character*1 magpref character*3 model character*1 lstauth character*1 most_amp character*1 most_dur character*1 most_p integer npfm integer nps_all integer nps_sol integer ns character*4 obsmag integer ounit real pmag real pmagobs character*1 pmagtyp c character*1 qual real qd real qs real rdeg real rms real rqual real sec real serr(3) character*1 south character*3 src integer sunit real t1 character*2 typmag character*1 versn1 character*1 versn2 integer ios c rdeg = 180./(4.*atan(1.0)) loctyp = 'H ' c do 20 i = 1, 5 mag(i) = 0. magtyp(i) = ' ' magerr(i) = 0. magobs(i) = 0 20 continue if (.not. y2k) then read (card, 100, iostat = ios) 1 iyear, imonth, iday, ihour, imin, sec, 2 latdeg, south, latmin, londeg, east, 3 lonmin, depth, mag(1), nps_sol, azgap, 4 dmin, rms, (iaz(i), idip(i), serr(i), i = 1, 2), 5 mag(2), locrmk, 6 serr(3), auxrmk, ns, errh, errz, 7 npfm, magobs(1), magobs(2), magerr(1), magerr(2), 8 model, lstauth, most_p, most_dur, most_amp, 9 magtyp(2), nps_all, magtyp(1), magtyp(3), mag(3), 1 magobs(3), magtyp(4), mag(4), magobs(4), eventid, 2 pmagtyp, pmag, pmagobs, magtyp(5), mag(5), magobs(5) 100 format( 1 5i2, f4.2, 2 f2.0, a1, f4.2, f3.0, a1, 3 f4.2, f5.2, f2.1, i3, i3, 4 f3.0, f4.2, 2(i3, i2, f4.2), 5 f2.1, a3, 6 f4.2, a2, i2, f4.2, f4.2, 7 i2, f3.1, f3.1, f3.2, f3.2, 8 a3, a1, a1, a1, a1, 9 a1, i3, a1, a1, f3.2, 1 f3.1, a1, f3.2, f3.1, i10, 2 a1, f3.2, f3.1,t146, a1, f3.2, f3.1) else read (card, 110, iostat = ios) 1 iyear, imonth, iday, ihour, imin, sec, 2 latdeg, south, latmin, londeg, east, 3 lonmin, depth, mag(1), nps_sol, azgap, 4 dmin, rms, (iaz(i), idip(i), serr(i), i = 1, 2), 5 mag(2), locrmk, 6 serr(3), auxrmk, ns, errh, errz, 7 npfm, magobs(1), magobs(2), magerr(1), magerr(2), 8 model, lstauth, most_p, most_dur, most_amp, 9 magtyp(2), nps_all, magtyp(1), magtyp(3), mag(3), 1 magobs(3), magtyp(4), mag(4), magobs(4), eventid, 2 pmagtyp, pmag, pmagobs, magtyp(5), mag(5), magobs(5), 3 versn1, versn2 110 format ( 1 i4, 4i2, f4.2, 2 f2.0, a1, f4.2, f3.0, a1, 3 f4.2, f5.2, f3.2, i3, i3, 4 f3.0, f4.2, 2(i3, i2, f4.2), 5 f3.2, a3, 6 f4.2, a2, i3, 2f4.2, 7 i3, 2f4.1, 2f3.2, 8 a3, 4a1, 9 a1, i3, a1, a1, f3.2, 1 f3.1, a1, f3.2, f3.1, i10, 2 2(a1, f3.2, f4.1), 3 2a1) endif if (ios .ne. 0) then ierr = 1 return endif if (.not. y2k) then if (iyear .gt. 0 .and. iyear .lt. 67) then iyear = 2000 + iyear else iyear = 1900 + iyear endif endif call datenm (iyear, imonth, iday, ihour, imin, sec) if (south .eq. 'S' .or. south .eq. 's') latdeg = -latdeg if (east .eq. 'W' .or. east .eq. ' ' .or. east .eq. 'w') 1 londeg = -londeg lond = sign(abs(londeg) + abs(lonmin)/60., londeg) latd = sign(abs(latdeg) + abs(latmin)/60., latdeg) c c paste in an "L" for local earthquake if not something else. Ok to overwrite '*' and '-' c if (auxrmk(1:1) .eq. 'Q' .or. auxrmk(2:2) .eq. 'Q') then auxrmk(1:2) = 'Q ' elseif (auxrmk(1:1) .eq. 'B' .or. auxrmk(2:2) .eq. 'B') then auxrmk(1:2) = 'B ' elseif (auxrmk(1:1) .eq. 'N' .or. auxrmk(2:2) .eq. 'N') then auxrmk(1:2) = 'N ' elseif (auxrmk(1:1) .eq. 'R' .or. auxrmk(2:2) .eq. 'R') then auxrmk(1:2) = 'R ' else auxrmk(1:2) = 'L ' endif c complete error ellipse. orientation of 3rd std error is cross product c of first two do i=1,2 e(i,3)=sin(idip(i)/rdeg) t1=sqrt(1.-e(i,3)**2) e(i,1)=cos(iaz(i)/rdeg)*t1 e(i,2)=sin(iaz(i)/rdeg)*t1 end do do i=1,3 j=i+1 k=i+2 if (j.gt.3) j=j-3 if (k.gt.3) k=k-3 e(3,i)=e(1,j)*e(2,k)-e(1,k)*e(2,j) end do c get the orientation of serr(3) from the unit vectors if (e(3,2).eq.0. .and. e(3,1).eq.0.) then iaz(3)=0. else iaz(3)=rdeg* atan2 (e(3,2), e(3,1)) end if if (iaz(3).lt.0) iaz(3)=iaz(3)+360 idip(3)=rdeg* abs( atan2( e(3,3), sqrt( e(3,1)**2 +e(3,2)**2))) c compute quality if (rms .lt. 0.15 .and. errh .le. 1.0 .and. errz .le. 2.0) then qs = 1 elseif (rms .lt. 0.30 .and. errh .le. 2.5 .and. errz .le. 5.0) 1 then qs = 2 elseif (rms .lt. 0.50 .and. errh .le. 5.0) then qs = 3 else qs = 4 endif if (nps_sol .ge. 6 .and. azgap .le. 90.0 .and. 1 (dmin .le. depth .or. dmin .le. 5.0)) then qd = 1 elseif (nps_sol .ge. 6 .and. azgap .le. 135.0 .and. 1 (dmin .le. 2.*depth .or. dmin .le. 10.0)) then qd = 2 elseif (nps_sol .ge. 6 .and. azgap .le. 180.0 .and. 1 (dmin .le. 50.0)) then qd = 3 else qd = 4 endif iq = int((qs + qd)/2.) if (iq .eq. 1) then c qual = 'A' rqual = 1. elseif (iq .eq. 2) then c qual = 'B' rqual = .75 elseif (iq .eq. 3) then c qual = 'C' rqual = .50 else c qual = 'D' rqual = .25 endif if (most_p .eq. 'B') then src = 'BK ' else src = 'NC ' endif c c write out $loc line c write (ounit, 200) iyear, imonth, iday, ihour, imin, sec, 1 latd, lond, depth, loctyp, src, nps_sol, azgap, dmin, rms, 2 errh, errz, auxrmk, rqual, eventid, date, eventid write (sunit, 200) iyear, imonth, iday, ihour, imin, sec, 1 latd, lond, depth, loctyp, src, nps_sol, azgap, dmin, rms, 2 errh, errz, auxrmk, rqual, eventid, date, eventid 200 format ('$loc ', i4, 4i2.2, f7.4, 1 f9.5, f10.5, f8.4, a2, a3, i4, i3, f10.4, f7.4, 2 t88, 2f7.4, a2, f7.4, i12, a8, i12) c write out $add line write (ounit, 300) nps_all, ns, npfm, 1 (iaz(i), idip(i), serr(i), i = 3, 1, -1) write (sunit, 300) nps_all, ns, npfm, 1 (iaz(i), idip(i), serr(i), i = 3, 1, -1) 300 format ('$add$loc', 3i4, 3(i3, i2, f10.4)) c PREFERRED MAGNITUDE SELECTION ORDER c c Order label label mag sum_mag_wts mag_MAD min # min type c code column columns columns columns readings mag c ---------------------------------------------------------------------- c 1. L 115 116-118 119-121 N/A 4 4.0 "external" c 2. Z 146 147-149 150-152 N/A 6 4.0 "alternate" c 3. D 110 68- 69 94- 96 100-102 1 0.0 "primary" c 4. X or A 114 35- 36 91- 93 97- 99 1 0.0 "primary" c 5. L 122 123-125 126-128 N/A 4 4.0 "alternate" c 6. L 115 116-118 119-121 N/A 0 0.0 "external" c 7. L 122 123-125 126-128 N/A 0 0.0 "alternate" c magtype = (1=Xmag, 2=Fmag, 3=Ext_mag (eg., UCB ML or 1966 mags), 4=Alt_amp (eg., NCSN ML), 5=Alt_coda (e.g, Zmag)) c c c No info on number or error of ML, and no error info on Zmag c c output $mag line. Only output info if mag is non-zero c Don't output value for magerr, magobs if unknown c do i=1,5 if (mag(i) .gt. 0.) then c c Is this the preferred magnitude? If so, sub in the preferred magnitude since it has more precision than D or X mag c if (pmagtyp .eq. magtyp(i)) then magobs(i) = pmagobs magpref = 'P' mag(i) = pmag else magpref = ' ' endif c c Set the CNSS magtyp; return the mag value c if (magtyp(i) .eq. 'X') then typmag = 'a ' amag = mag(i) elseif (magtyp(i) .eq. 'D') then typmag = 'd ' cmag = mag(i) elseif (magtyp(i) .eq. 'L' .and. i .eq. 3) then typmag = 'l ' wamag = mag(i) elseif (magtyp(i) .eq. 'L' .and. i .eq. 4) then typmag = 'l3' lmag = mag(i) elseif (magtyp(i) .eq. 'W') then typmag = 'w ' wmag = mag(i) elseif (magtyp(i) .eq. 'Z') then typmag = 'z ' zmag = mag(i) elseif (magtyp(i) .eq. 'A') then typmag = 'a ' amag = mag(i) elseif (magtyp(i) .eq. 'G') then typmag = 'a ' amag = mag(i) elseif (pmagtyp .ne. 'H') then write (eunit, '(a)') 'Unknown magtyp' return endif c plug in number of observations instead of weighted # if available if (magobs(i) .eq. 0) then obsmag = " " else if (typmag .eq. 'd ') then magobs(i) = ncodv else if (typmag .eq. 'z ') then magobs(i) = ncodz else if (typmag .eq. 'a ') then magobs(i) = nampx else if (typmag .eq. 'l3' .and. i .eq. 4) then c if nampl = 0, use the value on the summary card (see hypoinverse command XCH) if (nampl .ne. 0) magobs(i) = nampl endif write (obsmag, '(i4)') nint(magobs(i)) endif c c write magnitude error if non-zero c if (magerr(i) .eq. 0) then errmag = " " else write (errmag, '(f5.2)') magerr(i) endif c c make external UCB mags have BK src c if ((typmag .eq. 'l ' .and. i .eq. 3) .or. 1 typmag .eq. 'w ') then src = 'BK ' else src = 'NC ' endif c c output the $mag line c write (ounit, 400) magpref, mag(i), typmag, src, 1 obsmag, errmag, date write (sunit, 400) magpref, mag(i), typmag, src, 1 obsmag, errmag, date 400 format ('$mag', a1, f5.2, a2, a3, 1 a4, a5, t45, a8) endif end do c c if there is only a "Hand" magnitude in the preferred column, output it c if (pmagtyp .eq. 'H') then magpref = 'P' typmag = 'un' write (obsmag, '(i4)') nint(pmagobs) errmag = " " write (ounit, 400) magpref, pmag, typmag, src, 1 obsmag, errmag, date write (sunit, 400) magpref, pmag, typmag, src, 1 obsmag, errmag, date endif return end