      PROGRAM SUMMPLOTS
C
C PLOT HIGH RESOLUTION SUMMARY DATA 
C
C ORIGINAL SOURCE BY SANDY KRAMER - 11/09/93.
C MODIFIED FOR AUTOMATED DAILY PLOTS PLOTTING BY S. KRAMER - 01/03/93.
C MODIFIED TO PRODUCE 2 DAY PLOTS BY S. KRAMER
C Y2K COMPLIANCE ADDED 03/16/99 BY S. KRAMER
C Field range modified for early mission (0-10 nT)  SBK - 02/04/2009
C
C PROGRAMMING NOTES:
C
C  PLOT TIME AXIS LABELS IMPLEMENTED FROM MAIN PROGRAM USING
C  PSEUDO-DECIMAL YEAR FORMAT REAL(YYDDDHH) WHICH IS CONVERTED
C  INTERNALLY BY PROCEDURE "GRID" TO DECIMAL YEAR WITH 4 DIGIT
C  YEAR FOR Y2K COMPLIANCE.
C
C  DATA PLOTTED BY DECIMAL DAYS FOR PRECISION REASONS.
C
      PARAMETER (IARR=30000,NFIELD=11)
C
C PLOTTING PARAMETERS
C
      CHARACTER PLOTTIME*8,PLOTDATE*9,DSN*50,
     &          VALUE*10,XFORM*10,YFORM*10,CTIME1*10,CTIME2*10,
     &          XSTRING*50,YSTRING*50,TITLE*50,SUBTITLE*50
      INTEGER*4 NUM(IARR,NFIELD),NOUT(NFIELD),XLAB,YLAB,TYPE
      LOGICAL*1 DONE
      REAL*4 MAG(IARR,NFIELD),AVG(IARR,NFIELD),VAL,
     &       B(IARR),BX(IARR),BY(IARR),BZ(IARR),
     &       BDEL(IARR),BLAM(IARR),
     &       BRMS(IARR),BRMSX(IARR),BRMSY(IARR),BRMSZ(IARR),
     &       BS(IARR),BSX(IARR),BSY(IARR),BSZ(IARR),
     &       BSDEL(IARR),BSLAM(IARR),TDAY(IARR),XLOW,XHIGH,CONVERTD
      REAL*8 TIME(IARR),TOUT(IARR),T1,T2,T,CONVERT,DECYR,DECDY,TEMP
C
C SUMMARY FORMAT PARAMETERS
C
      CHARACTER RECTYPE*4,TELFMT*4,FLTID*4,TIMEFMT*4,
     &          COORD*2,MODE*9,RUNMONTH*4,RUNTIME*12
      INTEGER*2 DATE(6),MODCNT(3),DATAID(2)
      INTEGER*4 IB192(25),IB96(5),IB48,SIB48,
     &          RUNYEAR,RUNDAY
      LOGICAL*1 VOYSYS(32),MAGSYS(32),RECORD(10000)
      REAL*4 B192(25,3),DEL192(25),LAM192(25),FMOD192(25),
     &       FMAG192(25),RMS192(25,3),
     &       B96(3,5),DEL96(5),LAM96(5),FMOD96(5),
     &       FMAG96(5),RMS96(5,3),
     &       B48(3),DEL48,LAM48,FMOD48,
     &       FMAG48,RMS48(3),
     &       SRMS48(3),SX48,SY48,SZ48,SFMOD48,
     &       MTB(3,3),MTB5(3,3),MHG(3,3),POS(3),ANG(2),
     &       SCFLD(155),POSN(13),ATTD(27),
     &       HDR1(100),HDR(32),DATA(341),SUMOUT(568)
      REAL*8 DD,TD,TN,TP
C
      EQUIVALENCE (HDR(1),RECTYPE),    (HDR(2),TELFMT),
     &            (HDR(3),FLTID),      (HDR(4),DATE(1)),
     &            (HDR(7),DD),         (HDR(9),TD),
     &            (HDR(11),TIMEFMT),   (HDR(12),TIMEPD),
     &            (HDR(13),MODCNT(1)), (HDR(17),DATAID(1))
C
      EQUIVALENCE (HDR1(7),RUNYEAR),   (HDR1(9),RUNMONTH),
     &            (HDR1(10),RUNDAY),   (HDR1(69),VOYSYS(1)),
     &            (HDR1(77),MAGSYS(1))
C
      EQUIVALENCE (DATA(1),FMAG48),      (DATA(2),FMOD48),
     &            (DATA(3),DEL48),       (DATA(4),LAM48),     
     &            (DATA(5),B48(1)),      (DATA(8),RMS48(1)),
     &            (DATA(11),IB48),       (DATA(12),FMAG96(1)),
     &            (DATA(17),FMOD96(1)),  (DATA(22),DEL96(1)),  
     &            (DATA(27),LAM96(1)),   (DATA(32),B96(1,1)),  
     &            (DATA(47),RMS96(1,1)), (DATA(62),IB96(1)),
     &            (DATA(67),FMAG192(1)), (DATA(92),FMOD192(1)),
     &            (DATA(117),DEL192(1)), (DATA(142),LAM192(1)),
     &            (DATA(167),B192(1,1)), (DATA(242),RMS192(1,1)),
     &            (DATA(317),IB192(1))
C
      EQUIVALENCE (SCFLD(148),SRMS48(1)), (SCFLD(151),SIB48),
     &            (SCFLD(152),SX48),      (SCFLD(153),SY48),
     &            (SCFLD(154),SZ48),      (SCFLD(155),SFMOD48)
C
      EQUIVALENCE (POSN(1),TN),     (POSN(3),TP),
     &            (POSN(5),POS(1)), (POSN(11),RANGE),
     &            (POSN(12),ANG(1))
C
      EQUIVALENCE (ATTD(1),MTB(1,1)), (ATTD(10),MTB5(1,1)),
     &            (ATTD(19),MHG(1,1))
C
      EQUIVALENCE (SUMOUT(1),HDR(1)),     (SUMOUT(33),DATA(1)),
     &            (SUMOUT(342),SCFLD(1)), (SUMOUT(529),POSN(1)),
     &            (SUMOUT(542),ATTD(1)),  (SUMOUT(1),HDR1(1))
C
      EQUIVALENCE (RECORD(1),SUMOUT(1))
C
      WRITE(6,*)
      TYPE*,'ENTER LFM SUMMARY DSN'
      READ(5,'(A)') DSN
      OPEN(50,FILE=DSN,STATUS='OLD',FORM='FORMATTED',
     &     RECORDTYPE='VARIABLE',RECL=8191,READONLY)
C
C      WRITE(6,*)
C      WRITE(6,*) 'ENTER START TIME (YYDDDHH)'
C      READ(5,*) XLOW
C      WRITE(6,*)
C      WRITE(6,*) 'ENTER STOP TIME (YYDDDHH)'
C      READ(5,*) XHIGH
C      T1 = CONVERT(XLOW)
C      T2 = CONVERT(XHIGH)
C
      DONE = .FALSE.
      NPTS = 0
      TEMP = 0.0D0
    5 CONTINUE
C
C GET SYSTEM DATE AND TIME
C
      CALL GETDATE(PLOTDATE)
      CALL GETTIME(PLOTTIME)
C
C READ BINARY SUMMARY LFM DATASET 
C
   10 CONTINUE
       READ(50,'(Q,<LEN>A1)',END=100,ERR=10) LEN,(RECORD(I),I=1,LEN)
C
       IF (DATAID(1).EQ.8) THEN
C
C DETERMINE COORDINATE SYSTEM
C
        IF ( MAGSYS(9) ) THEN
         COORD = 'PL'
        ELSE IF ( MAGSYS(17) ) THEN
         COORD = 'S3'
        ELSE IF ( MAGSYS(18) ) THEN
         COORD = 'L1'
        ELSE IF ( MAGSYS(19) ) THEN
         COORD = 'U1'
        ELSE IF ( MAGSYS(20) ) THEN
         COORD = 'N1'
        ELSE
         COORD = 'HG'
        END IF
C
C READ MVS GENERATED SUMMARY DATA
C
        IF ( VOYSYS(8) ) COORD = 'HG'
        IF ( VOYSYS(9) ) COORD = 'L1'
        IF ( VOYSYS(10) ) COORD = 'S3'
        IF ( VOYSYS(11) ) COORD = 'U1'
        IF ( VOYSYS(12) ) COORD = 'N1'
C
C DETERMINE MAGNETOMETER MODE
C
        IF ( MAGSYS(28) ) THEN
         MODE = 'COMPOSITE'
        ELSE IF ( MAGSYS(21) .AND. .NOT.MAGSYS(22) ) THEN
         MODE = 'PRI'
        ELSE IF ( .NOT.MAGSYS(21) .AND. MAGSYS(22) ) THEN
         MODE = 'SEC'
        ELSE IF ( MAGSYS(21) .AND. MAGSYS(22) ) THEN
         MODE = 'PRI/SEC'
        ELSE IF ( MAGSYS(27) ) THEN
         MODE = 'SCF' ! DISTINGUISH AMBIENT FROM S/C FIELD DATA  SBK 10/20/94
        ELSE 
         MODE = 'DUAL'
        END IF
C
C EDR PROCESS DATE
C
        WRITE(RUNTIME,'(A4,''/'',I2.2,''/'',I4)') 
     &   RUNMONTH,RUNDAY,RUNYEAR
C
       END IF
C
       IF (DATAID(1).NE.1) GOTO 10
       T = DECYR(DATE) ! DECIMAL 4 DIGIT YEAR
C CHECK START TIME
C      IF (T.LT.T1) GOTO 10
C CHECK END TIME
C      IF (T.GE.T2) GOTO 100
C CHECK TIME ADVANCE FORWARD
       IF (T.LE.TEMP) GOTO 10
C CHECK TIME JUMP LESS THAN APPROXIMATELY 5 DAYS
       IF (NPTS.GT.1.AND.(T-TEMP).GT.0.137D0) GOTO 10
       TEMP = T
       NPTS = NPTS + 1
C GET STARTING INFORMATION
       IF (NPTS.EQ.1) THEN
C GET DATA START TIME FOR PLOT LABEL (YEAR & DAY ONLY).  REAL(YYDDDHH) FORMAT
        XLOW = FLOAT(DATE(1))*100000.0 + FLOAT(DATE(2))*100.0 
C SET WINDOW LENGTH TO TWO DAYS 
        XHIGH = XLOW + 200.0
C DOUBLE PRECISION DECIMAL YEAR
        T1 = CONVERT(XLOW)
        T2 = CONVERT(XHIGH)
C DOUBLE PRECISION DECIMAL DAY
        XL = CONVERTD(XLOW)
        XH = CONVERTD(XHIGH)
        IY1 = INT(T1)
        IY2 = INT(T2)
        IF ( IY2.GT.IY1 ) THEN
         IF ( MOD(IY1,4).EQ.0 ) THEN
          ROLLOVER = 366.0
         ELSE
          ROLLOVER = 365.0
         END IF
        ELSE
         ROLLOVER = 0.0
        END IF
C FLIGHT ID
        IF (FLTID.EQ.'FLT1') THEN
         TITLE = 'VOYAGER 1 '//COORD//' '//MODE//'\E'
        ELSE IF (FLTID.EQ.'FLT2') THEN 
         TITLE = 'VOYAGER 2 '//COORD//' '//MODE//'\E'
        ELSE
         WRITE(6,*) '***INVALID FLIGHT ID***'
         STOP
        END IF
C YEAR OF DATA
        IF ( DATE(1).LT.77 ) THEN
         ICENT = 2000
        ELSE
         ICENT = 1900
        END IF
        WRITE(SUBTITLE,'(I4,A2)') DATE(1) + ICENT,'\E'
C END STARTING INFORMATION
       END IF
C STORE DECIMAL YEAR
       TIME(NPTS) = T
       IY = INT(T)
C STORE DECIMAL DAY
       T = DECDY(DATE) 
       IF ( IY.GT.IY1 ) THEN
        TDAY(NPTS) = REAL(T) + ROLLOVER
       ELSE
        TDAY(NPTS) = REAL(T) 
       END IF
c       write(6,*) xl,xh,TDAY(NPTS)
C FIELD X COMPONENT
       MAG(NPTS,1) = DATA(5)
C FIELD Y COMPONENT
       MAG(NPTS,2) = DATA(6)
C FIELD Z COMPONENT
       MAG(NPTS,3) = DATA(7)
C FIELD MAGNITUDE
       MAG(NPTS,4) = DATA(1)
C RMS X FIELD COMPONENT
       MAG(NPTS,9) = DATA(8)
C RMS Y FIELD COMPONENT
       MAG(NPTS,10) = DATA(9)
C RMS Z FIELD COMPONENT
       MAG(NPTS,11) = DATA(10)
C S/C FIELD OR SECONDARY X COMPONENT
       MAG(NPTS,5) = SCFLD(152)
C S/C FIELD OR SECONDARY Y COMPONENT
       MAG(NPTS,6) = SCFLD(153)
C S/C FIELD OR SECONDARY Z COMPONENT
       MAG(NPTS,7) = SCFLD(154)
C S/C FIELD OR SECONDARY Z COMPONENT
       MAG(NPTS,8) = SCFLD(155)
C
C      WRITE(6,803) ID,DATE,(MAG(NPTS,K),K=1,4)
C
      IF (TIME(NPTS).LT.T2) GOTO 10
      GOTO 101
C
  100 CONTINUE
      DONE = .TRUE.
  101 CONTINUE
C
  803 FORMAT(1X,I1,6(1X,I3),4(1X,F7.3))
C
      WRITE(6,*)
      WRITE(6,'(1X,''START TIME: '',I7,3X,
     & ''END TIME: '',I7)') INT(XLOW),INT(XHIGH)
      WRITE(6,*)
      WRITE(6,'(1X,I5.5,'' RECORDS INPUTTED'')') NPTS
      IF (NPTS.LT.10.AND.DONE) THEN
       WRITE(6,*)
       WRITE(6,*) '***INSUFFICIENT NUMBER OF POINTS FOR PLOTTING***'
       WRITE(6,*)
       GOTO 999
      END IF
      IF (NPTS.GT.IARR) THEN
       WRITE(6,*)
       WRITE(6,*) 'REDIMENSION IARR GT ',NPTS
       WRITE(6,*)
       GOTO 999
      END IF
C
      FILL = 999.0
      NMAX = NPTS
      DO JJ = 1,NMAX
       IF (MAG(JJ,4).NE.FILL) THEN
        BX(JJ) = MAG(JJ,1)
        BY(JJ) = MAG(JJ,2)
        BZ(JJ) = MAG(JJ,3)
        BRMSX(JJ) = MAG(JJ,9)
        BRMSY(JJ) = MAG(JJ,10)
        BRMSZ(JJ) = MAG(JJ,11)
        B2 = SQRT(BX(JJ)**2+BY(JJ)**2+BZ(JJ)**2)
        B(JJ) = MAG(JJ,4)
        IF (B2.NE.0.0) THEN
         BDEL(JJ) = ASIN(BZ(JJ)/B2) * 180.0/3.14159
        ELSE
         BDEL(JJ) = FILL
        END IF
        IF (BX(JJ).NE.0.0.OR.BY(JJ).NE.0.0) THEN
         BLAM(JJ) = 180.0 - ATAN2(BY(JJ),-BX(JJ))  * 180.0/3.14159
        END IF
       ELSE
        BX(JJ) = FILL
        BY(JJ) = FILL
        BZ(JJ) = FILL
        BRMSX(JJ) = FILL
        BRMSY(JJ) = FILL
        BRMSZ(JJ) = FILL
        B(JJ) = FILL
        BDEL(JJ) = FILL
        BLAM(JJ) = FILL
       END IF
C
C S/C FIELD OR SECONDARY DATA
C
       IF (MAG(JJ,8).NE.FILL) THEN
        BSX(JJ) = MAG(JJ,5)
        BSY(JJ) = MAG(JJ,6)
        BSZ(JJ) = MAG(JJ,7)
        BS2 = SQRT(BSX(JJ)**2+BSY(JJ)**2+BSZ(JJ)**2)
        BS(JJ) = MAG(JJ,8)
        IF (BSZ(JJ).NE.0.0) THEN
         BSDEL(JJ) = ASIN(BSZ(JJ)/BS2) * 180.0/3.14159
        ELSE
         BSDEL(JJ) = FILL
        END IF
        IF (BSY(JJ).NE.0.0.OR.BSX(JJ).NE.0.0) THEN
         BSLAM(JJ) = 180.0 - ATAN2(BSY(JJ),-BSX(JJ))  * 180.0/3.14159
        ELSE
         BSLAM(JJ) = FILL
        END IF
       ELSE
        BSX(JJ) = FILL
        BSY(JJ) = FILL
        BSZ(JJ) = FILL
        BS(JJ) = FILL
        BSDEL(JJ) = FILL
        BSLAM(JJ) = FILL
       END IF
      END DO
      INDEX = NMAX
C
  900 CONTINUE
C
C MAG PLOT ROUTINE
C
C PLOT DISCONTINUITY INTERVAL (HOURS)
C
      DTIME = 1.0
      DTIME = DTIME/24.0
C
      CALL MGOINIT
      CALL MGOSETUP(-5)
      CALL MGOERASE
C
      XX1 = 400.
      XX2 = 2225.
C
C PLOT SUBTITLE
C
      YY1 = 3000.
      YY2 = 3100.
      CALL MGOSETEXPAND(1.0)
      CALL MGOSETLOC(XX1,YY1,XX2,YY2)
      CALL MGOSETLIM(0.,0.,100.,100.)
      CALL MGORELOCATE(50.,50.)
      CALL MGOPUTLABEL(50,SUBTITLE,5)
C
C FIELD MAGNITUDE
C
      YY1 = 2730.
      YY2 = 3000.
      XT1 = 100.
      XT2 = 200.
      YT1 = 100.
      YT2 = 200.
      TT1 = 3100.
      TT2 = 3200.
      XNUMH = .8
      YNUMH = .8
      XTXTH = 1.
      YTXTH = 1.
      TITLH = 1.3
      XLAB = 4
      YLAB = 3
      XTLEN = 50.
      YTLEN = 50.
      XSTRING = 'DAY/HOUR\E'
      YSTRING = 'B (nT)\E'
      TYPE = -2
C
C LARGE/SMALL TICK INTERVALS (HOURS)
C    
      BSPANX = 24.0
      SSPANX = 1.0
C CONVERT TICK SPACING FROM HOURS TO DAYS
      BSPANX = BSPANX/24.0
      SSPANX = SSPANX/24.0
C
C SET MIN AND MAX VALUE FOR Y AXIS
C
      YLOW = 0.
      YHIGH = 4000.
      SSPANY = 100.0
      BSPANY = 1000.0
      YFORM = 'F6.1'
C
C FLAG OUT OF BOUNDS POINTS AND REMOVE
C
      CALL MGOSETLOC(XX1,YY1,XX2,YY2)
      CALL MGOSETLIM(XL,YLOW,XH,YHIGH)
      write(6,*) xl,xh,TDAY(NPTS/2),B(NPTS/2)
      CALL MGOSETEXPAND(1.0)
      DO I = 1,NPTS
       IF (B(I).GT.YHIGH.AND.B(I).NE.FILL) THEN
        FOLD = MOD(B(I),YHIGH)
        CALL MGORELOCATE(TDAY(I),FOLD)
        CALL MGOPOINT(10,1)
        B(I) = FILL
       END IF
      END DO
      DO I = 1,NPTS
       IF (BS(I).GT.YHIGH.AND.BS(I).NE.FILL) THEN
        FOLD = MOD(BS(I),YHIGH)
        CALL MGORELOCATE(TDAY(I),FOLD)
        CALL MGOPOINT(10,1)
        BS(I) = FILL
       END IF
      END DO
C
C PLOT GRID AND LABELS
C
      CALL GRID(XX1,XX2,YY1,YY2,XT1,XT2,YT1,YT2,TT1,TT2,
     C          XNUMH,YNUMH,XTXTH,YTXTH,TITLH,XLAB,YLAB,
     C          XTLEN,YTLEN,XSTRING,YSTRING,TITLE,TYPE,
     C          XLOW,XHIGH,YLOW,YHIGH,SSPANX,BSPANX,
     C          SSPANY,BSPANY,XFORM,YFORM)
C
C PLOT DATA
C
      CALL MGOSETLOC(XX1,YY1,XX2,YY2)
C     CALL MGOPLOTID(F1,F2)
      CALL MGOSETEXPAND(0.5)
      CALL MGOGRELOCATE(XX2,YY2)
      CALL MGOPUTLABEL(20,PLOTDATE//'   '//PLOTTIME,7)
      CALL MGOGRELOCATE(XX1,YY2)
      CALL MGOPUTLABEL(31,'EDR PROCESS DATE:  '//RUNTIME,9)
      CALL MGOSETLIM(XL,YLOW,XH,YHIGH)
      CALL MGOSETEXPAND(0.01)
      IF ( COORD.EQ.'PL' .AND. MODE.EQ.'PRI/SEC' ) 
     & CALL MGOSETLWEIGHT(2)
      CALL MCONNECT(FILL,TDAY,B,NPTS,DTIME)
      CALL MGOSETLWEIGHT(1)
      IF ( COORD.EQ.'PL' .AND. MODE.EQ.'PRI/SEC' ) 
     & CALL MCONNECT(FILL,TDAY,BS,NPTS,DTIME)
      CALL MGOSETEXPAND(1.0)
C
C MAKE RMS PLOTS FOR PAYLOAD DATA
C
      IF ( COORD.EQ.'PL' ) THEN
C
C RMS BX
C
      YY1 = 2460.
      YY2 = 2730.
      XLAB = 4
      TITLH = 0.
      YSTRING = 'RMS\DX\E'
C
C SET MIN AND MAX VALUES FOR Y AXIS
C
      YLOW = 0.
      YHIGH = .04
      SSPANY = .0025
      BSPANY = .01
      YFORM = 'F5.2'
C
C PLOT GRID AND LABELS
C
      CALL GRID(XX1,XX2,YY1,YY2,XT1,XT2,YT1,YT2,TT1,TT2,
     C          XNUMH,YNUMH,XTXTH,YTXTH,TITLH,XLAB,YLAB,
     C          XTLEN,YTLEN,XSTRING,YSTRING,TITLE,TYPE,
     C          XLOW,XHIGH,YLOW,YHIGH,SSPANX,BSPANX,
     C          SSPANY,BSPANY,XFORM,YFORM)
C
C PLOT DATA
C
      CALL MGOSETLOC(XX1,YY1,XX2,YY2)
      CALL MGOSETLIM(XL,YLOW,XH,YHIGH)
      CALL MGOSETEXPAND(0.01)
C     CALL MGOSETLWEIGHT(2)
      CALL MCONNECT(FILL,TDAY,BRMSX,NPTS,DTIME)
      CALL MGOSETLWEIGHT(1)
      CALL MGOSETEXPAND(1.0)
C
C RMS BY
C
      YY1 = 2190.
      YY2 = 2460.
      XLAB = 4
      TITLH = 0.
      YSTRING = 'RMS\DY\E'
C
C PLOT GRID AND LABELS
C
      CALL GRID(XX1,XX2,YY1,YY2,XT1,XT2,YT1,YT2,TT1,TT2,
     C          XNUMH,YNUMH,XTXTH,YTXTH,TITLH,XLAB,YLAB,
     C          XTLEN,YTLEN,XSTRING,YSTRING,TITLE,TYPE,
     C          XLOW,XHIGH,YLOW,YHIGH,SSPANX,BSPANX,
     C          SSPANY,BSPANY,XFORM,YFORM)
C
C PLOT DATA
C
      CALL MGOSETLOC(XX1,YY1,XX2,YY2)
      CALL MGOSETLIM(XL,YLOW,XH,YHIGH)
      CALL MGOSETEXPAND(0.01)
C     CALL MGOSETLWEIGHT(2)
      CALL MCONNECT(FILL,TDAY,BRMSY,NPTS,DTIME)
      CALL MGOSETLWEIGHT(1)
      CALL MGOSETEXPAND(1.0)
C
C RMS BZ
C
      YY1 = 1920.
      YY2 = 2190.
      XLAB = 4
      TITLH = 0.
      YSTRING = 'RMS\DZ\E'
C
C PLOT GRID AND LABELS
C
      CALL GRID(XX1,XX2,YY1,YY2,XT1,XT2,YT1,YT2,TT1,TT2,
     C          XNUMH,YNUMH,XTXTH,YTXTH,TITLH,XLAB,YLAB,
     C          XTLEN,YTLEN,XSTRING,YSTRING,TITLE,TYPE,
     C          XLOW,XHIGH,YLOW,YHIGH,SSPANX,BSPANX,
     C          SSPANY,BSPANY,XFORM,YFORM)
C
C PLOT DATA
C
      CALL MGOSETLOC(XX1,YY1,XX2,YY2)
      CALL MGOSETLIM(XL,YLOW,XH,YHIGH)
      CALL MGOSETEXPAND(0.01)
C     CALL MGOSETLWEIGHT(2)
      CALL MCONNECT(FILL,TDAY,BRMSZ,NPTS,DTIME)
      CALL MGOSETLWEIGHT(1)
      CALL MGOSETEXPAND(1.0)
C
C MAKE DELTA AND LAMBDA PLOTS FOR ROTATED DATA
C
      ELSE IF ( COORD.NE.'PL' ) THEN
C
C LATITUDINAL FIELD ANGLE
C
      YY1 = 2460.
      YY2 = 2730.
      XLAB = 4
      TITLH = 0.
      YSTRING = '\Gd\DB\E'
C
C SET MIN AND MAX VALUES FOR Y AXIS
C
      YLOW = -90.
      YHIGH = 90.
      SSPANY = 15.
      BSPANY = 45.
      YFORM = 'F5.1'
C
C PLOT GRID AND LABELS
C
      CALL GRID(XX1,XX2,YY1,YY2,XT1,XT2,YT1,YT2,TT1,TT2,
     C          XNUMH,YNUMH,XTXTH,YTXTH,TITLH,XLAB,YLAB,
     C          XTLEN,YTLEN,XSTRING,YSTRING,TITLE,TYPE,
     C          XLOW,XHIGH,YLOW,YHIGH,SSPANX,BSPANX,
     C          SSPANY,BSPANY,XFORM,YFORM)
C
C PLOT DATA
C
      CALL MGOSETLOC(XX1,YY1,XX2,YY2)
      CALL MGOSETLIM(XL,YLOW,XH,YHIGH)
      CALL MGOSETEXPAND(0.01)
      IF ( COORD.EQ.'PL' .AND. MODE.EQ.'PRI/SEC' ) 
     & CALL MGOSETLWEIGHT(2)
      CALL MCONNECT(FILL,TDAY,BDEL,NPTS,DTIME)
      CALL MGOSETLWEIGHT(1)
      IF ( COORD.EQ.'PL' .AND. MODE.EQ.'PRI/SEC' ) 
     & CALL MCONNECT(FILL,TDAY,BSDEL,NPTS,DTIME)
      CALL MGOSETEXPAND(1.0)
C
C LONGITUDINAL FIELD ANGLE
C
      YY1 = 1920.
      YY2 = 2460.
      XLAB = 4
      TITLH = 0.
      YSTRING = '\Gl\DB\E'
C
C SET MIN AND MAX VALUES FOR Y AXIS
C
      YLOW = 0.
      YHIGH = 360.
      SSPANY = 15.
      BSPANY = 90.
C
C PLOT GRID AND LABELS
C
      CALL GRID(XX1,XX2,YY1,YY2,XT1,XT2,YT1,YT2,TT1,TT2,
     C          XNUMH,YNUMH,XTXTH,YTXTH,TITLH,XLAB,YLAB,
     C          XTLEN,YTLEN,XSTRING,YSTRING,TITLE,TYPE,
     C          XLOW,XHIGH,YLOW,YHIGH,SSPANX,BSPANX,
     C          SSPANY,BSPANY,XFORM,YFORM)
C
C PLOT DATA
C
      CALL MGOSETLOC(XX1,YY1,XX2,YY2)
      CALL MGOSETLIM(XL,YLOW,XH,YHIGH)
      CALL MGOSETEXPAND(0.01)
      IF ( COORD.EQ.'PL' .AND. MODE.EQ.'PRI/SEC' ) 
     & CALL MGOSETLWEIGHT(2)    
      CALL MCONNECT(FILL,TDAY,BLAM,NPTS,DTIME)
      CALL MGOSETLWEIGHT(1)
      IF ( COORD.EQ.'PL' .AND. MODE.EQ.'PRI/SEC' ) 
     & CALL MCONNECT(FILL,TDAY,BSLAM,NPTS,DTIME)
      CALL MGOSETEXPAND(1.0)
C
C END RMS/ANGLE PLOT DECISION BLOCK
C
      END IF
C
C X FIELD COMPONENT
C
      YY1 = 1380.
      YY2 = 1920.
      XLAB = 4
      TITLH = 0.
      IF ( COORD.EQ.'PL' ) THEN 
       YSTRING = 'B\DX (nT)\E'
      ELSE IF ( COORD.EQ.'HG' ) THEN
       YSTRING = 'B\DR (nT)\E'
      ELSE 
       YSTRING = 'B\DR (nT)\E'
      END IF
C
C SET MIN AND MAX VALUE FOR Y AXIS
C
      YLOW = -4000.
      YHIGH = 4000.
      SSPANY = 100.0
      BSPANY = 1000.0
      YFORM = 'F7.1'
C
C FLAG OUT OF BOUNDS POINTS AND REMOVE
C
      CALL MGOSETLOC(XX1,YY1,XX2,YY2)
      CALL MGOSETLIM(XL,YLOW,XH,YHIGH)
      CALL MGOSETEXPAND(1.0)
      DO I = 1,NPTS
       IF ((BX(I).LT.YLOW.OR.
     &      BX(I).GT.YHIGH).AND.
     &      BX(I).NE.FILL) THEN
        IF (BX(I).LT.YLOW) FOLD = MOD(BX(I),YLOW)
        IF (BX(I).GT.YHIGH) FOLD = MOD(BX(I),YHIGH)
        CALL MGORELOCATE(TDAY(I),FOLD)
        CALL MGOPOINT(10,1)
        BX(I) = FILL
       END IF
      END DO
      DO I = 1,NPTS
       IF ((BSX(I).LT.YLOW.OR.
     &      BSX(I).GT.YHIGH).AND.
     &      BSX(I).NE.FILL) THEN
        IF (BSX(I).LT.YLOW) FOLD = MOD(BSX(I),YLOW)
        IF (BSX(I).GT.YHIGH) FOLD = MOD(BSX(I),YHIGH)
        CALL MGORELOCATE(TDAY(I),FOLD)
        CALL MGOPOINT(10,1)
        BSX(I) = FILL
       END IF
      END DO
C
C PLOT GRID AND LABELS
C
      CALL GRID(XX1,XX2,YY1,YY2,XT1,XT2,YT1,YT2,TT1,TT2,
     C          XNUMH,YNUMH,XTXTH,YTXTH,TITLH,XLAB,YLAB,
     C          XTLEN,YTLEN,XSTRING,YSTRING,TITLE,TYPE,
     C          XLOW,XHIGH,YLOW,YHIGH,SSPANX,BSPANX,
     C          SSPANY,BSPANY,XFORM,YFORM)
C
C PLOT DATA
C
      CALL MGOSETLOC(XX1,YY1,XX2,YY2)
      CALL MGOSETLIM(XL,YLOW,XH,YHIGH)
      CALL MGOSETEXPAND(0.01)
      IF ( COORD.EQ.'PL'.AND. MODE.EQ.'PRI/SEC' ) 
     & CALL MGOSETLWEIGHT(2)
      CALL MCONNECT(FILL,TDAY,BX,NPTS,DTIME)
      CALL MGOSETLWEIGHT(1)
      IF ( COORD.EQ.'PL' .AND. MODE.EQ.'PRI/SEC' ) 
     & CALL MCONNECT(FILL,TDAY,BSX,NPTS,DTIME)
      CALL MGOSETEXPAND(1.0)
C
C Y FIELD COMPONENT
C
      YY1 = 840.
      YY2 = 1380.
      XLAB = 4
      TITLH = 0.
      IF ( COORD.EQ.'PL' ) THEN 
       YSTRING = 'B\DY (nT)\E'
      ELSE IF ( COORD.EQ.'HG' ) THEN
       YSTRING = 'B\DT (nT)\E'
      ELSE 
       YSTRING = 'B\DQ (nT)\E'
      END IF
C
C SET MIN AND MAX VALUE FOR Y AXIS
C
      YLOW = -4000.
      YHIGH = 4000.
      SSPANY = 100.0
      BSPANY = 1000.0
C
C FLAG OUT OF BOUNDS POINTS AND REMOVE
C
      CALL MGOSETLOC(XX1,YY1,XX2,YY2)
      CALL MGOSETLIM(XL,YLOW,XH,YHIGH)
      CALL MGOSETEXPAND(1.0)
      DO I = 1,NPTS
       IF ((BY(I).LT.YLOW.OR.
     &      BY(I).GT.YHIGH).AND.
     &      BY(I).NE.FILL) THEN
        IF (BY(I).LT.YLOW) FOLD = MOD(BY(I),YLOW)
        IF (BY(I).GT.YHIGH) FOLD = MOD(BY(I),YHIGH)
        CALL MGORELOCATE(TDAY(I),FOLD)
        CALL MGOPOINT(10,1)
        BY(I) = FILL
       END IF
      END DO
      DO I = 1,NPTS
       IF ((BSY(I).LT.YLOW.OR.
     &      BSY(I).GT.YHIGH).AND.
     &      BSY(I).NE.FILL) THEN
        IF (BSY(I).LT.YLOW) FOLD = MOD(BSY(I),YLOW)
        IF (BSY(I).GT.YHIGH) FOLD = MOD(BSY(I),YHIGH)
        CALL MGORELOCATE(TDAY(I),FOLD)
        CALL MGOPOINT(10,1)
        BSY(I) = FILL
       END IF
      END DO
C
C PLOT GRID AND LABELS
C
      CALL GRID(XX1,XX2,YY1,YY2,XT1,XT2,YT1,YT2,TT1,TT2,
     C          XNUMH,YNUMH,XTXTH,YTXTH,TITLH,XLAB,YLAB,
     C          XTLEN,YTLEN,XSTRING,YSTRING,TITLE,TYPE,
     C          XLOW,XHIGH,YLOW,YHIGH,SSPANX,BSPANX,
     C          SSPANY,BSPANY,XFORM,YFORM)
C
C PLOT DATA
C
      CALL MGOSETLOC(XX1,YY1,XX2,YY2)
      CALL MGOSETLIM(XL,YLOW,XH,YHIGH)
      CALL MGOSETEXPAND(0.01)
      IF ( COORD.EQ.'PL' .AND. MODE.EQ.'PRI/SEC' ) 
     & CALL MGOSETLWEIGHT(2)
      CALL MCONNECT(FILL,TDAY,BY,NPTS,DTIME)
      CALL MGOSETLWEIGHT(1)
      IF ( COORD.EQ.'PL' .AND. MODE.EQ.'PRI/SEC' ) 
     & CALL MCONNECT(FILL,TDAY,BSY,NPTS,DTIME)
      CALL MGOSETEXPAND(1.0)
C
C Z FIELD COMPONENT
C
      YY1 = 300.
      YY2 = 840.
      XLAB = 0
      TITLH = 0.
      IF ( COORD.EQ.'PL' ) THEN 
       YSTRING = 'B\DZ (nT)\E'
      ELSE IF ( COORD.EQ.'HG' ) THEN
       YSTRING = 'B\DN (nT)\E'
      ELSE 
       YSTRING = 'B\DF (nT)\E'
      END IF
C
C SET MIN AND MAX VALUE FOR Y AXIS
C
      YLOW = -4000.
      YHIGH = 4000.
      SSPANY = 100.0
      BSPANY = 1000.0
C
C FLAG OUT OF BOUNDS POINTS AND REMOVE
C
      CALL MGOSETLOC(XX1,YY1,XX2,YY2)
      CALL MGOSETLIM(XL,YLOW,XH,YHIGH)
      CALL MGOSETEXPAND(1.0)
      DO I = 1,NPTS
       IF ((BZ(I).LT.YLOW.OR.
     &      BZ(I).GT.YHIGH).AND.
     &      BZ(I).NE.FILL) THEN
        IF (BZ(I).LT.YLOW) FOLD = MOD(BZ(I),YLOW)
        IF (BZ(I).GT.YHIGH) FOLD = MOD(BZ(I),YHIGH)
        CALL MGORELOCATE(TDAY(I),FOLD)
        CALL MGOPOINT(10,1)
        BZ(I) = FILL
       END IF
      END DO
      DO I = 1,NPTS
       IF ((BSZ(I).LT.YLOW.OR.
     &      BSZ(I).GT.YHIGH).AND.
     &      BSZ(I).NE.FILL) THEN
        IF (BSZ(I).LT.YLOW) FOLD = MOD(BSZ(I),YLOW)
        IF (BSZ(I).GT.YHIGH) FOLD = MOD(BSZ(I),YHIGH)
        CALL MGORELOCATE(TDAY(I),FOLD)
        CALL MGOPOINT(10,1)
        BSZ(I) = FILL
       END IF
      END DO
C
C PLOT GRID AND LABELS
C
      CALL GRID(XX1,XX2,YY1,YY2,XT1,XT2,YT1,YT2,TT1,TT2,
     C          XNUMH,YNUMH,XTXTH,YTXTH,TITLH,XLAB,YLAB,
     C          XTLEN,YTLEN,XSTRING,YSTRING,TITLE,TYPE,
     C          XLOW,XHIGH,YLOW,YHIGH,SSPANX,BSPANX,
     C          SSPANY,BSPANY,XFORM,YFORM)
C
C PLOT DATA
C
      CALL MGOSETLOC(XX1,YY1,XX2,YY2)
      CALL MGOSETLIM(XL,YLOW,XH,YHIGH)
      CALL MGOSETEXPAND(0.01)
      IF ( COORD.EQ.'PL' .AND. MODE.EQ.'PRI/SEC' ) 
     & CALL MGOSETLWEIGHT(2)
      CALL MCONNECT(FILL,TDAY,BZ,NPTS,DTIME)
      CALL MGOSETLWEIGHT(1)
      IF ( COORD.EQ.'PL' .AND. MODE.EQ.'PRI/SEC' ) 
     & CALL MCONNECT(FILL,TDAY,BSZ,NPTS,DTIME)
      CALL MGOSETEXPAND(1.0)
C
  888 CONTINUE
C
      CALL MGOPRNTPLOT(NVEC)
      WRITE(6,*)
      WRITE(6,'(1X,''VECTORS PLOTTED: '',I6)') NVEC
      WRITE(6,*)
      DO I=1,1000000
         CALL GETTIME(PLOTTIME)
      END DO
C
      TIME(1) = TIME(NPTS)
      NPTS = 1
C INCREMENT WINDOW BY TWO DAYS
      XLOW = XHIGH
      XHIGH = XLOW + 200.0
C DOUBLE PRECISION DECIMAL YEAR
      T1 = CONVERT(XLOW)
      T2 = CONVERT(XHIGH)
C DOUBLE PRECISION DECIMAL DAY
      XL = CONVERTD(XLOW)
      XH = CONVERTD(XHIGH)
      IY1 = INT(T1)
      IY2 = INT(T2)
      IF ( IY2.GT.IY1 ) THEN
       IF ( MOD(IY1,4).EQ.0 ) THEN
        ROLLOVER = 366.0
       ELSE
        ROLLOVER = 365.0
       END IF
      ELSE
       ROLLOVER = 0.0
      END IF
C
      IF (.NOT.DONE) GOTO 5
C
  999 CONTINUE
      CLOSE(50)
C
      STOP
      END      
      SUBROUTINE GETDATE(SYSDATE)
C
C CALL DEC FORTRAN SUBROUTINE TO GET SYSTEM DATE 
C
      CHARACTER SYSDATE*9
C
      CALL DATE(SYSDATE)
C
      RETURN 
      END
      SUBROUTINE GETTIME(SYSTIME)
C
C CALL DEC FORTRAN SUBROUTINE TO GET SYSTEM TIME
C
      CHARACTER SYSTIME*8
C
      CALL TIME(SYSTIME)
C
      RETURN
      END
      SUBROUTINE HRAVG(T1,T2,FILL,BTIME,DTIME,T,DATA,NPTS,
     &                 TOUT,DOUT,PARR,NOUT)
C
C PRODUCE AVERAGES OF PERIOD DTIME FROM INPUTTED HOUR AVERAGE DATA.
C ROUTINE DOES NOT REMOVE TIME POINTS WITH NULL AVERAGES.
C
      INTEGER*4 PARR(1)
      REAL*4 DATA(1),DOUT(1)
      REAL*8 T(1),TIME,TIME1,TIME2,TOUT(1),T1,T2,TEND
C
C T1 START TIME
C T2 END TIME
C FILL  FILL DATA VALUE
C BTIME RESOLUTION OF INPUT DATA (HOURS)
C DTIME AVERAGING INTERVAL (HOURS)
C T     INPUT TIME ARRAY
C DATA  INPUT DATA ARRAY
C NPTS  NUMBER OF POINTS IN INPUT DATA
C TOUT  OUTPUT TIME ARRAY
C DOUT  OUTPUT ARRAY OF DTIME AVERAGES
C PARR  OUTPUT ARRAY - # OF AVGS/PT
C NOUT NUMBER OF POINTS IN OUTPUT DATA
C
      SUM = 0.0
      NUM = 0
      TEND = T2 + DBLE(BTIME/(100.0*24.0*DAYS(T2)))
C
      TIME1 = T1 - DBLE(BTIME/(100.0*24.0*DAYS(T1)))
      TIME2 = TIME1 + DBLE(DTIME/(24.0*DAYS(TIME1)))
      IARR = 0
      INDEX = 0
  10  CONTINUE
       IARR = IARR + 1
       IF (IARR.GT.NPTS) GOTO 20
       TIME = T(IARR)
       VAL = DATA(IARR)
C
C GET DATA FROM PROPER TIME INTERVAL
C
  15   CONTINUE
       IF (TIME.LT.TIME1) GOTO 10
       IF (TIME.GE.TIME2) GOTO 20
       IF (TIME.GT.TEND) GOTO 20
C
C ACCUMULATE DATA NOT INCLUDING FILL POINTS
C
       IF (VAL.NE.FILL) THEN
        SUM = SUM + VAL
        NUM = NUM + 1
       END IF
C
      GOTO 10
   20 CONTINUE
C
C AVERAGE DATA AND INCREMENT TIME INTERVAL
C
       AVG = FILL
       IF (NUM.NE.0) THEN
        AVG = SUM/FLOAT(NUM)
       END IF
       INDEX = INDEX + 1
       TOUT(INDEX) = TIME1 + DBLE(BTIME/(100.0*24.0*DAYS(TIME1)))
       DOUT(INDEX) = AVG
       PARR(INDEX) = NUM
       TIME1 = TIME2
       TIME2 = TIME1 + DBLE(DTIME/(24.0*DAYS(TIME1)))
       SUM = 0.0
       NUM = 0
      IF (IARR.LE.NPTS.AND.TIME.LE.TEND) GOTO 15
C
C RETURN FROM AVERAGING OVER DESIRED TIME INTERVAL
C
 100  CONTINUE
      WRITE(6,800) TOUT(1),TOUT(INDEX),INDEX
      CLOSE(10)
      NOUT = INDEX
C
      RETURN
 800  FORMAT(' START ',F12.8,'   STOP ',F12.8,'   POINTS ',I5)
      END
C
C MONGO PLOTTING ROUTINE WITH USER DEFINED SPECIFICS
C DEVELOPED BY SAUNDERS B. KRAMER JR., CODE 692
C
C XX1     PLOT LOWER X AXIS LIMIT IN DEVICE COORDINATES - R*4
C XX2     PLOT UPPER X AXIS LIMIT IN DEVICE COORDINATES - R*4
C YY1     PLOT LOWER Y AXIS LIMIT IN DEVICE COORDINATES - R*4
C YY2     PLOT UPPER Y AXIS LIMIT IN DEVICE COORDINATES - R*4
C XT1     Y LABEL LOWER X AXIS LIMIT IN DEVICE COORDINATES - R*4
C XT2     Y LABEL UPPER X AXIS LIMIT IN DEVICE COORDINATES - R*4
C YT1     X LABEL LOWER Y AXIS LIMIT IN DEVICE COORDINATES - R*4
C YT2     X LABEL UPPER Y AXIS LIMIT IN DEVICE COORDINATES - R*4
C TT1     TITLE LOWER Y AXIS LIMIT IN DEVICE COORDINATES - R*4
C TT2     TITLE UPPER Y AXIS LIMIT IN DEVICE COORDINATES - R*4
C XNUMH   SIZE OF X SCALE CHARACTERS - R*4
C YNUMH   SIZE OF Y SCALE CHARACTERS - R*4
C XTXTH   SIZE OF X LABEL CHARACTERS - R*4
C YTXTH   SIZE OF Y LABEL CHARACTERS - R*4
C TITLH   SIZE OF TITLE CHARACTERS - R*4
C XLAB    CONTROLS X AXIS SCALING AND LABELLING - I*4
C   0 PRINT ALL SCALES AND X LABEL
C   1 SKIP FIRST SCALE
C   2 SKIP LAST SCALE
C   3 SKIP FIRST AND LAST SCALE
C   4 NO SCALES 
C YLAB    PARALLELS XLAB ABOVE FOR Y AXIS - I*4
C XTLEN   SIZE OF X AXIS LARGE TICKS - R*4
C YTLEN   SIZE OF Y AXIS LARGE TICKS - R*4
C XSTRING X LABEL STRING - C*50
C YSTRING Y LABEL STRING - C*50
C TITLE TITLE STRING - C*50
C TYPE INDICATES INDEPENDENT AXIS FORMAT - I*4
C  -6 TIME FORMAT (YEAR//DAY) YEAR INTERVALS
C  -5 TIME FORMAT (YEAR)       *
C  -4 TIME FORMAT (HOUR)      DAY INTERVALS
C  -3 TIME FORMAT (DAY)        *  
C  -2 TIME FORMAT (DAY//HOUR)  *
C  -1 TIME FORMAT (YEAR//DAY)  *
C   0 STANDARD
C   1 LOGARITHMIC INDEPENDENT AXIS, STANDARD Y AXIS
C   2 LOGARITHMIC DEPENDENT AXIS, STANDARD X AXIS
C   3 LOGARITHMIC INDEPENDENT AND DEPENDENT AXIS
C XLOW    PLOT LOWER X AXIS LIMIT IN USER COORDINATES - R*4
C XHIGH   PLOT UPPER X AXIS LIMIT IN USER COORDINATES - R*4
C YLOW    PLOT LOWER Y AXIS LIMIT IN USER COORDINATES - R*4
C YHIGH   PLOT LOWER Y AXIS LIMIT IN USER COORDINATES - R*4
C SSPANX  SPACING OF INDEPENDENT AXIS SMALL TICKS - R*4
C BSPANX  SPACING OF INDEPENDENT AXIS LARGE TICKS - R*4
C SSPANY  SPACING OF DEPENDENT AXIS SMALL TICKS - R*4
C BSPANY  SPACING OF DEPENDENT AXIS LARGE TICKS - R*4
C XFORM   FORTRAN FORMAT OF INDEPENDENT AXIS SCALING - C*10
C YFORM   FORTRAN FORMAT OF DEPENDENT AXIS SCALING - C*10
C
      SUBROUTINE GRID(XX1,XX2,YY1,YY2,XT1,XT2,YT1,YT2,TT1,TT2,
     C                XNUMH,YNUMH,XTXTH,YTXTH,TITLH,XLAB,YLAB,
     C                XTLEN,YTLEN,XSTRING,YSTRING,TITLE,TYPE,
     C                XLOW,XHIGH,YLOW,YHIGH,SSPANX,BSPANX,
     C                SSPANY,BSPANY,XFORM,YFORM)
      CHARACTER*1 TCHAR
      CHARACTER*10 VALUE,XFORM,YFORM
      CHARACTER*50 XSTRING,YSTRING,TITLE
      INTEGER*4 XNUM,YNUM,XLAB,YLAB,TYPE
      LOGICAL TIME
      REAL*8 TLOW,THIGH,DEC,CONVERT
C
C DETERMINE IF INDEPENDENT AXIS HAS TIME SCALE
C
      TIME = .FALSE.
      IF (TYPE.LT.0) TIME = .TRUE.
      IF (.NOT.TIME) GOTO 30
C
C CONVERT 6 DIGIT CALENDAR YEAR (YYDDDHH) INTO TRUE DECIMAL YEAR 
C W/4 DIGITS TO LEFT OF DECIMAL POINT
C
       TLOW = CONVERT(XLOW)
       THIGH = CONVERT(XHIGH)
   30 CONTINUE
C
C DRAW BOX
C
      CALL MGOSETLOC(XX1,YY1,XX2,YY2)
      CALL MGOLINE(XX1,YY1,XX1,YY2)
      CALL MGOLINE(XX1,YY2,XX2,YY2)
      CALL MGOLINE(XX2,YY2,XX2,YY1)
      CALL MGOLINE(XX2,YY1,XX1,YY1)
C
C DRAW TIME (INDEPENDENT) AXIS TICK MARKS
C
      IF (.NOT.TIME) GOTO 40
       CALL MGOSETEXPAND(XNUMH)
       CALL BIGTICK(XX1,XX2,YY1,YY2,TLOW,THIGH,BSPANX,
     C              XLAB,XTLEN,TYPE)
       CALL SMALLTICK(XX1,XX2,YY1,YY2,TLOW,THIGH,SSPANX,
     C              XTLEN,TYPE)
      GOTO 45
   40 CONTINUE       
C
C DRAW INDEPENDENT AXIS LARGE TICK MARKS AND LABEL
C
      IF (TYPE.EQ.1.OR.TYPE.EQ.3) THEN
       ISCALE = 1
      ELSE 
       ISCALE = 0
      END IF
      CALL MGOSETEXPAND(XNUMH)
      CALL LRGTICK(YY1,YY2,XX1,XX2,XLOW,XHIGH,BSPANX,
     C             XFORM,XTLEN,XLAB,ISCALE)
C
C DRAW INDEPENDENT AXIS SMALL TICK MARKS
C
      IF (TYPE.EQ.1.OR.TYPE.EQ.3) THEN
       CALL LOGTICK(YY1,YY2,XX1,XX2,XLOW,XHIGH,SSPANX,XTLEN)
      ELSE
       CALL SMLTICK(YY1,YY2,XX1,XX2,XLOW,XHIGH,SSPANX,XTLEN)
      END IF
C
   45 CONTINUE
C 
C DRAW DEPENDENT AXIS LARGE TICK MARKS AND LABEL 
C
      IF (TYPE.EQ.2.OR.TYPE.EQ.3) THEN
       ISCALE = 1
      ELSE
       ISCALE = 0
      END IF
      CALL MGOSETEXPAND(YNUMH)
      XX1 = -XX1
      CALL LRGTICK(XX1,XX2,YY1,YY2,YLOW,YHIGH,
     C             BSPANY,YFORM,YTLEN,YLAB,ISCALE)
C
C DRAW DEPENDENT AXIS SMALL TICK MARKS 
C
      XX1 = -XX1
      IF (TYPE.EQ.2.OR.TYPE.EQ.3) THEN
       CALL LOGTICK(XX1,XX2,YY1,YY2,YLOW,YHIGH,SSPANY,YTLEN)
      ELSE
       CALL SMLTICK(XX1,XX2,YY1,YY2,YLOW,YHIGH,SSPANY,YTLEN)
      END IF
C
C LABEL TITLE
C
      CALL STRINGSET(XX1,TT1,XX2,TT2,TITLE,50,1,TITLH)
C
C LABEL INDEPENDENT AXIS
C
      CALL STRINGSET(XX1,YT1,XX2,YT2,XSTRING,50,1,XTXTH)
C
C LABEL DEPENDENT AXIS
C
      CALL MGOSETANGLE(90.)
      CALL STRINGSET(XT1,YY1,XT2,YY2,YSTRING,50,1,YTXTH)
      CALL MGOSETANGLE(0.)
C
C RESET LOCATION,SIZE AND WEIGHT
C
      CALL MGOSETLOC(XX1,YY1,XX2,YY2)
      CALL MGOSETEXPAND(1.0)
      CALL MGOSETLWEIGHT(1)
C
      RETURN
      END
      SUBROUTINE STRINGSET(XX1,YY1,XX2,YY2,STRING,LEN,WGHT,SIZE)
C
C PLACE LABEL CENTERED IN A BOX DEFINED BY DEVICE COORDINATES
C
      CHARACTER*(*) STRING
      INTEGER*4 WGHT
      CALL MGOSETLOC(XX1,YY1,XX2,YY2)
      XX = (XX1 + XX2) / 2.0
      YY = (YY1 + YY2) / 2.0
      CALL MGOGRELOCATE(XX,YY)
      CALL MGOSETEXPAND(SIZE)
      CALL MGOSETLWEIGHT(WGHT)
      IF (SIZE.NE.0.) CALL MGOPUTLABEL(LEN,STRING,5)
      CALL MGOSETEXPAND(1.0)
      CALL MGOSETLWEIGHT(1)
      RETURN
      END
      SUBROUTINE BIGTICK(XX1,XX2,YY1,YY2,XLOW,XHIGH,SPAN,
     C                   XLAB,XTLEN,TYPE)
C
C PRODUCE LARGE TICK MARKS AND LABELS ALONG INDEPENDENT AXIS (TIME)
C AT USER DEFINED INTERVALS IN DAYS.  RANGE OF AXIS VALUES AND USER 
C SPECIFIED SPACING DETERMINE NUMBER OF TICK MARKS.  DEVICE 
C COORDINATES USED FOR LOCATING MARKS.
C
C XX1   = LOWER INDEPENDENT AXIS BOUNDRY
C XX2   = UPPER INDEPENDENT AXIS BOUNDRY
C YY1   = LOWER DEPENDENT AXIS BOUNDRY
C YY2   = UPPER DEPENDENT AXIS BOUNDRY
C XLOW  = DECIMAL YEAR START
C XHIGH = DECIMAL YEAR STOP
C SPAN  = DISTANCE BETWEEN TICK MARKS
C XLAB  = BEGIN/END LABEL PRINT FLAG
C XTELN = INDEPENDENT AXIS TICK SIZE
C TYPE  = INDEPENDENT AXIS LABEL TYPE
C
      CHARACTER*10 XCHAR
      INTEGER*4 XLAB,TYPE
      REAL*8 XLOW,XHIGH,XVAL,VAL,YR1,YR2,FRAC1,FRAC2
C
      Y0 = YY1 - XTLEN/2.
      Y1 = YY1 + XTLEN
      Y2 = YY2 - XTLEN
      IF (TYPE.GT.-5) THEN
       XDIFF = REAL((XHIGH-XLOW)*365.25D0)
      ELSE
       XDIFF = REAL(XHIGH-XLOW)
      END IF
      YDIFF = YY2 - YY1
      ISTOP = NINT(XDIFF/SPAN) + 1
      XVAL = XLOW
      DO 50 I=1,ISTOP
       XLIN = XX1 + (XX2-XX1) * REAL((XVAL-XLOW)/(XHIGH-XLOW))
       ERR = XVAL/XHIGH
       IF (ERR.GT.1.0) GOTO 50
C       ERR = REAL(XVAL-XHIGH)/REAL(XHIGH-XLOW)
C       IF (ERR.GT.0.0001) GOTO 50
C       IF (XVAL.GT.XHIGH) GOTO 50
       CALL MGOLINE(XLIN,YY1,XLIN,Y1)
       CALL MGOLINE(XLIN,YY2,XLIN,Y2)
       IF ( XLAB.EQ.4 ) GOTO 30
       CALL CALENDAR(XVAL,IY,ID,IH)
C ONLY WANT LAST 2 DIGITS OF 4 DIGIT YEAR (1960-2059)
       IF (IY.GT.1999) THEN
        IY = IY - 2000
       ELSE
        IY = IY - 1900
       END IF
       IF (IY.GT.99) THEN
        WRITE(6,*) '***BIGTICK***STOP'
        STOP
       END IF
C BUILD SCALE STRING
       IF (TYPE.EQ.-1) WRITE(XCHAR,'(I2,I3.3)') IY,ID
       IF (TYPE.EQ.-2) WRITE(XCHAR,'(I3.3,I2.2)') ID,IH
       IF (TYPE.EQ.-3) WRITE(XCHAR,'(I3.3)') ID
       IF (TYPE.EQ.-4) WRITE(XCHAR,'(I2.2)') IH
       IF (TYPE.EQ.-5) WRITE(XCHAR,'(I2.2)') IY
       IF (TYPE.EQ.-6) WRITE(XCHAR,'(I2.2,I3.3)') IY,ID
       INUM = 0
       DO ILET=1,10
        LVAL = ICHAR(XCHAR(ILET:ILET))
        IF ( LVAL.NE.32 ) INUM = INUM + 1
       END DO
       IF ( XLAB.EQ.1.AND.I.EQ.1 ) GOTO 30
       IF ( XLAB.EQ.2.AND.I.EQ.ISTOP ) GOTO 30
       IF ( XLAB.EQ.3.AND.(I.EQ.1.OR.I.EQ.ISTOP) ) GOTO 30
       CALL MGOGRELOCATE(XLIN,Y0)
C       TYPE*,'X LABEL *'//XCHAR//'*'
       CALL MGOPUTLABEL(INUM,XCHAR,2)
   30  IF (I.EQ.ISTOP) GOTO 50
       IF (TYPE.GT.-5) THEN
C BALANCE CALENDAR ARITHMETIC ACROSS CHANGE OF YEAR
        YR1 = DBLE(INT(XVAL))
        VAL = XVAL + DBLE(SPAN) / DBLE(DAYS(XVAL))
        IF (INT(VAL).GT.INT(XVAL)) THEN
         YR2 = YR1 + 1.D0
         FRAC1 = (YR2 - XVAL) * DBLE(DAYS(XVAL))
         FRAC2 = DBLE(SPAN) - FRAC1 
         VAL = YR2 + FRAC2 / DBLE(DAYS(YR2))
        ENDIF
        XVAL = VAL
       ELSE
        XVAL = XVAL + SPAN
       END IF
   50 CONTINUE
C
      RETURN
      END
      SUBROUTINE SMALLTICK(XX1,XX2,YY1,YY2,XLOW,XHIGH,SPAN,
     C                     XTLEN,TYPE)
C
C PRODUCE SMALL TICK MARKS ALONG INDEPENDENT AXIS (TIME)
C AT USER DEFINED INTERVALS IN DAYS.
C
C XX1   = LOWER INDEPENDENT AXIS BOUNDRY
C XX2   = UPPER INDEPENDENT AXIS BOUNDRY
C YY1   = LOWER DEPENDENT AXIS BOUNDRY
C YY2   = UPPER DEPENDENT AXIS BOUNDRY
C XLOW  = DECIMAL YEAR START
C XHIGH = DECIMAL YEAR STOP
C SPAN  = DISTANCE BETWEEN TICK MARKS
C XLAB  = BEGIN/END LABEL PRINT FLAG
C XTELN = INDEPENDENT AXIS TICK SIZE
C TYPE  = INDEPENDENT AXIS LABEL TYPE
C
      INTEGER*4 TYPE
      REAL*8 XLOW,XHIGH,XVAL,VAL,YR1,YR2,FRAC1,FRAC2
C
      Y1 = YY1 + XTLEN/2.
      Y2 = YY2 - XTLEN/2.
      IF (TYPE.GT.-5) THEN
       XDIFF = REAL((XHIGH-XLOW)*365.25D0)
      ELSE
       XDIFF = REAL(XHIGH-XLOW)
      END IF
      YDIFF = YY2 - YY1
      ISTOP = NINT(XDIFF/SPAN) + 1
      XVAL = XLOW
      DO 50 I=1,ISTOP
       XLIN = XX1 + (XX2-XX1) * REAL((XVAL-XLOW)/(XHIGH-XLOW))
       IF (XVAL.GT.XHIGH) GOTO 50
       CALL MGOLINE(XLIN,YY1,XLIN,Y1)
       CALL MGOLINE(XLIN,YY2,XLIN,Y2)
       IF (TYPE.GT.-5) THEN
        YR1 = DBLE(INT(XVAL))
        VAL = XVAL + DBLE(SPAN) / DBLE(DAYS(XVAL))
        IF (INT(VAL).GT.INT(XVAL)) THEN
         YR2 = YR1 + 1.D0
         FRAC1 = (YR2 - XVAL) * DBLE(DAYS(XVAL))
         FRAC2 = DBLE(SPAN) - FRAC1 
         VAL = YR2 + FRAC2 / DBLE(DAYS(YR2))
        ENDIF
        XVAL = VAL
       ELSE
        XVAL = XVAL + SPAN
       END IF
   50 CONTINUE
C
      RETURN
      END
      SUBROUTINE LRGTICK(VV1,VV2,WW1,WW2,VLOW,VHIGH,
     C                   SPAN,FORM,TLEN,LAB,ISCALE)
      CHARACTER*1 TCHAR
      CHARACTER*5 EXPONENT
      CHARACTER*10 FORM,NCHAR
      REAL*4 NUMH
C
C DRAW LARGE TICK MARKS AND LABEL
C
      IF (VV1.LT.0) THEN
       VV1 = ABS(VV1)
       IAX = 2
      ELSE 
       IAX = 1
      END IF
      DIFF = VHIGH - VLOW
      ISTOP = NINT(DIFF/SPAN) + 1
      V0 = VV1 - TLEN/2.
      V1 = VV1 + TLEN
      V2 = VV2 - TLEN
      VAL = VLOW
      DO I=1,ISTOP 
       VLIN = WW1 + (WW2-WW1)*(VAL-VLOW)/DIFF
       ERR = (VLIN-WW2)/(WW2-WW1)
       IF (ERR.GT.0.0001) GOTO 100
       IF (IAX.EQ.1) THEN
C
C DRAW INDEPENDENT AXIS LARGE TICKS
C
        CALL MGOLINE(VLIN,VV1,VLIN,V1)
        CALL MGOLINE(VLIN,VV2,VLIN,V2)
       ELSE IF (IAX.EQ.2) THEN
C
C DRAW DEPENDENT AXIS LARGE TICKS
C
        CALL MGOLINE(VV1,VLIN,V1,VLIN)
        CALL MGOLINE(VV2,VLIN,V2,VLIN)
       END IF
C
       IF ( LAB.EQ.4 ) GOTO 60
C
C LINEAR SCALING
C
       IF (ISCALE.EQ.0) THEN
        TCHAR = FORM(1:1)
        IF (TCHAR.EQ.'I') THEN
         WRITE(NCHAR,'('//FORM//')') INT(VAL)
        ELSE
         WRITE(NCHAR,'('//FORM//')') VAL
        END IF
       END IF
C
C LOGARITHMIC SCALING
C
       IF (ISCALE.EQ.1) THEN
        WRITE(EXPONENT,'(I5)') INT(VAL)
        WRITE(NCHAR,'(A10)') '10\\U'//EXPONENT
       END IF
C
C SHIFT LEFT CHARACTER STRING
C
       NUM = 0
       DO 55 ILET=1,10
        TCHAR = NCHAR(ILET:ILET)
        LVAL = ICHAR(TCHAR)
        IF ( LVAL.LT.32 .OR. LVAL.EQ.127 ) GOTO 55
        IF ( LVAL.NE.32 ) THEN
         NUM = NUM + 1
         NCHAR(NUM:NUM) = TCHAR
         IF (NUM.NE.ILET) NCHAR(ILET:ILET) = CHAR(32)
        END IF
   55  CONTINUE
C
       IF ( LAB.EQ.1.AND.I.EQ.1 ) GOTO 60
       IF ( LAB.EQ.2.AND.I.EQ.ISTOP ) GOTO 60
       IF ( LAB.EQ.3.AND.(I.EQ.1.OR.I.EQ.ISTOP) ) GOTO 60
       IF (IAX.EQ.1) THEN
C        TYPE*,'X LABEL *'//NCHAR//'*'
        CALL MGOGRELOCATE(VLIN,V0)
        CALL MGOPUTLABEL(NUM,NCHAR,2)
       END IF
       IF (IAX.EQ.2) THEN
C        TYPE*,'Y LABEL *'//NCHAR//'*'
        CALL MGOGRELOCATE(V0,VLIN)
        CALL MGOPUTLABEL(NUM,NCHAR,4)
       END IF
   60  VAL = VAL + SPAN
      END DO
  100 CONTINUE
      RETURN 
      END
      SUBROUTINE SMLTICK(VV1,VV2,WW1,WW2,VLOW,VHIGH,SPAN,TLEN)
C
C DRAW SMALL TICK MARKS
C
      IF (VV1.LT.0) THEN
       VV1 = ABS(VV1)
       IAX = 2
      ELSE 
       IAX = 1
      END IF
      DIFF = VHIGH - VLOW
      ISTOP = NINT(DIFF/SPAN) + 1
      V1 = VV1 + TLEN/2.
      V2 = VV2 - TLEN/2.
      VAL = VLOW
      DO I=1,ISTOP
       VAL = VAL + SPAN
       VLIN = WW1 + (WW2-WW1)*(VAL-VLOW)/DIFF
       IF (VLIN.GT.WW2) GOTO 45
       IF (IAX.EQ.1) THEN
C
C DRAW INDEPENDENT AXIS SMALL TICKS
C
        CALL MGOLINE(VLIN,VV1,VLIN,V1)
        CALL MGOLINE(VLIN,VV2,VLIN,V2)
       ELSE IF (IAX.EQ.2) THEN
C
C DRAW DEPENDENT AXIS SMALL TICKS
C
        CALL MGOLINE(VV1,VLIN,V1,VLIN)
        CALL MGOLINE(VV2,VLIN,V2,VLIN)
       END IF
C
      END DO
   45 CONTINUE
      RETURN
      END
      SUBROUTINE LOGTICK(VV1,VV2,WW1,WW2,VLOW,VHIGH,SPAN,TLEN)
C
C DRAW SMALL TICK MARKS
C
      IF (VV1.LT.0) THEN
       VV1 = ABS(VV1)
       IAX = 2
      ELSE 
       IAX = 1
      END IF
      DIFF = VHIGH - VLOW
      ISTOP = NINT(DIFF/SPAN) + 1
      V1 = VV1 + TLEN/2.
      V2 = VV2 - TLEN/2.
      DO DECADE = VLOW,VHIGH,1.0
       CURDEC = DECADE
       DO I = 1,10
        VAL = ALOG10(REAL(I)*10.0**CURDEC)
        VLIN = WW1 + (WW2-WW1)*(VAL-VLOW)/DIFF
        IF (VLIN.GT.WW2) GOTO 45
        IF (IAX.EQ.1) THEN
C
C DRAW INDEPENDENT AXIS SMALL TICKS
C
         CALL MGOLINE(VLIN,VV1,VLIN,V1)
         CALL MGOLINE(VLIN,VV2,VLIN,V2)
        ELSE IF (IAX.EQ.2) THEN
C
C DRAW DEPENDENT AXIS SMALL TICKS
C
         CALL MGOLINE(VV1,VLIN,V1,VLIN)
         CALL MGOLINE(VV2,VLIN,V2,VLIN)
        END IF
C
       END DO
      END DO
   45 CONTINUE
      RETURN
      END
      REAL*8 FUNCTION DECYR(DATE)
C
C CONVERT HIGH RESOLUTION CALENDAR TIME INTO DECIMAL YEAR TIME
C
C MODIFIED TO RETURN 4 DIGIT WHOLE YEAR  SBK 06/12/98
C
      INTEGER*2 DATE(6)
      REAL*8 YEAR,DAY,HOUR,MIN,SEC,MSEC,LEAP
      LEAP = 365.0D0
      IF (MOD(DATE(1),4).EQ.0) LEAP = 366.0D0
      MSEC = DFLOAT(DATE(6))
      SEC = DFLOAT(DATE(5)) + MSEC/1000.0D0
      MIN = DFLOAT(DATE(4)) + SEC/60.0D0
      HOUR = DFLOAT(DATE(3)) + MIN/60.0D0
      DAY = DFLOAT(DATE(2)-1) + HOUR/24.0D0
      YEAR = DFLOAT(DATE(1)) + DAY/LEAP
C
      IF ( DATE(1).LT.60 ) THEN
       YEAR = YEAR + 2000.0D0
      ELSE
       YEAR = YEAR + 1900.0D0
      END IF
C
      DECYR = YEAR
C
      RETURN
      END
      REAL*8 FUNCTION DECDY(DATE)
C
C CONVERT HIGH RESOLUTION CALENDAR DATA TO DECIMAL DAY (JAN. 1 = 1)
C
      INTEGER*2 DATE(6)
      REAL*8 DAY,HOUR,MIN,SEC,MSEC
      MSEC = DFLOAT(DATE(6))
      SEC = DFLOAT(DATE(5)) + MSEC/1000.0D0
      MIN = DFLOAT(DATE(4)) + SEC/60.0D0
      HOUR = DFLOAT(DATE(3)) + MIN/60.0D0
      DAY = DFLOAT(DATE(2)) + HOUR/24.0D0
      DECDY = DAY
      RETURN
      END
      SUBROUTINE CALENDAR(TIME,IYEAR,IDAY,IHOUR)
C 
C CONVERT DOUBLE PRECISION DECIMAL YEAR INTO INTEGER YEAR, DAY AND HOUR 
C

      REAL*8 TIME,DYEAR,LEAP
C
      IYEAR = INT(TIME)
      DYEAR = DBLE(IYEAR)
      LEAP = DBLE(DAYS(TIME))
      DDAY = REAL((TIME-DYEAR)*LEAP) + 1.
      IDAY = INT(DDAY)
      DHOUR = (DDAY-REAL(IDAY))*24.
      IHOUR = INT(DHOUR)
      DMIN = DHOUR - REAL(IHOUR)
      IF (DMIN.GT.0.8) IHOUR = IHOUR + 1
C
      IF (IHOUR.EQ.24) THEN
       IDAY = IDAY + 1
       IHOUR = 0
      END IF
C
      ILEAP = INT(LEAP)
      IF (IDAY.GT.ILEAP) THEN
       IYEAR = IYEAR + 1
       IDAY = IDAY - ILEAP
      END IF
C
      RETURN
      END
      REAL*8 FUNCTION CONVERT(XVAL)
      CHARACTER*10 VALUE
      REAL*8 DEC
C
C CONVERT CONCATENATED CALENDAR DATE STORED AS REAL (YYDDDHH) INTO DECIMAL YEAR
C
C 1960-1999 = 60-99   
C 2000-2059 = 100-159
C
      IVAL = INT(XVAL)
      WRITE(VALUE,'(I7)') IVAL
      READ(VALUE,'(I2,I3,I2)') IY,ID,IH
      IF ( IY.LT.0 ) THEN
       WRITE(6,*)
     &  '***CONVERT***INVALID NEGATIVE YEAR.  STOPPING EXECUTION.' 
       STOP
      END IF
      CONVERT = DEC(IY,ID,IH)
      RETURN
      END
      REAL*4 FUNCTION CONVERTD(XVAL)
      CHARACTER*10 VALUE
C
C CONVERT CONCATENATED CALENDAR DATE STORED AS REAL(YYDDDHH) INTO DECIMAL DAY
C
C 1960-1999 = 60-99   
C 2000-2059 = 100-159
C
      IVAL = INT(XVAL)
      WRITE(VALUE,'(I7)') IVAL
      READ(VALUE,'(I2,I3,I2)') IY,ID,IH
      IF ( IY.LT.0 ) THEN
       WRITE(6,*)
     &  '***CONVERT***INVALID NEGATIVE YEAR.  STOPPING EXECUTION.' 
       STOP
      END IF
C
      CONVERTD = FLOAT(ID) + FLOAT(IH)/24.0
C    
      RETURN
      END
      REAL*8 FUNCTION DEC(NYR,NDY,NHR)
C
C CONVERT CALENDAR TIME INTO DECIMAL YEAR,
C REQUIRES TO DIGIT YEAR (00-99) AS INPUT
C OUTPUT DECIMAL YEAR W/4 DIGIT YEAR.
C
      REAL*8 YEAR,DAY,HOUR
C
      IF ( MOD(NYR,4).EQ.0 ) THEN
       LEAP = 366
      ELSE
       LEAP = 365
      END IF
C
      IF ( NDY.GT.LEAP ) THEN
       NYR = NYR + 1
       NDY = NDY - LEAP
      END IF
C
      YEAR = DFLOAT(NYR)
      DAY = DFLOAT(NDY)
      HOUR= DFLOAT(NHR)
      DEC = YEAR + (DAY-1.D0 + HOUR/24.D0) / DBLE(DAYS(YEAR))
C
      IF (NYR.LT.60) THEN
       DEC = DEC + 2000.0D0
      ELSE
       DEC = DEC + 1900.0D0
      END IF
C
      RETURN
      END
      REAL*4 FUNCTION DAYS(TIME)
C
C COMPUTE LENGTH OF YEAR IN DAYS
C
      REAL*8 TIME
      NYR = INT(TIME)
      DAYS = 365.
      IF (MOD(NYR,4).EQ.0) DAYS = 366.
C
      RETURN
      END
      SUBROUTINE MCONNECT(FILL,XVAL,YVAL,PTS,DELTA)
C
C ROUTINE TO CONNECT ADJACENT ARRAY POINTS WITH INDEPENDENT DELTAS
C LESS THAN <DELTA> WHILE REMOVING FILL VALUES.
C
      INTEGER*4 PTS
      LOGICAL*1 PRIME
      REAL*4 X(2),Y(2),XVAL(1),YVAL(1)
C
C FILL  FILL DATA VALUE
C XVAL  INDEPENDENT AXIS ARRAY
C YVAL  DEPENDENT AXIS ARRAY
C PTS   NUMBER OF ELEMENTS IN PARALLEL ARRAYS XVAL & VAL
C DELTA MAXIMUM GAP BETWEEN TWO ADJACENT POINTS TO BE CONNECTED
C
      PRIME = .FALSE.
      IPT = 0
      N = PTS - 1
      DO I = 1,N
       IF ( YVAL(I).NE.FILL ) THEN
        X(1) = XVAL(I)
        Y(1) = YVAL(I)
        PRIME = .TRUE.
       END IF
       IF ( YVAL(I+1).NE.FILL .AND. PRIME ) THEN
        X(2) = XVAL(I+1)
        Y(2) = YVAL(I+1)
        DIFF = X(2) - X(1)
        IF ( DIFF.LT.DELTA ) THEN
         CALL MGOCONNECT(X,Y,2)
        ELSE
         CALL MGORELOCATE(XVAL(I+1),YVAL(I+1))
         CALL MGOPOINT(10,3)
        END IF
        IPT = IPT + 1
        X(1) = X(2)
        Y(1) = Y(2)
       END IF
      END DO
C
      RETURN
      END
      REAL*8 FUNCTION REALTIME(TIME)
C
C CONVERT INTEGER CALENDAR TIME INTO DECIMAL YEAR REAL TIME.
C
      INTEGER*2 TIME(6)
      REAL*8 DAYS
C
      DAYS = 365.0D0
      IF (MOD(TIME(1),4).EQ.0) DAYS = 366.0D0
      REALTIME = DBLE(TIME(1)) + 
     &           DBLE(TIME(2)-1)/DAYS +
     &           DBLE(TIME(3))/24.0D0/DAYS +
     &           DBLE(TIME(4))/60.0D0/24.0D0/DAYS +
     &           DBLE(TIME(5))/60.0D0/60.0D0/24.0D0/DAYS +
     &           DBLE(TIME(6))/1000.0D0/60.0D0/60.0D0/24.0D0/DAYS
C
C ASSUME 2 DIGIT YEAR.  ANY YEAR BEFORE VOYAGER LAUNCH (77) IS INTO NEXT
C CENTURY.
C
      IF (TIME(1).LT.77) REALTIME = REALTIME + 100.0D0
C
      RETURN                                                            
      END
