Initial commit

1 files changed, 1313 insertions, 0 deletions

diff --git a/sys/gio/ncarutil/conrec.f b/sys/gio/ncarutil/conrec.f
new file mode 100644
index 00000000..b3e246c1
--- /dev/null
+++ b/sys/gio/ncarutil/conrec.f
@@ -0,0 +1,1313 @@
+      SUBROUTINE CONREC (Z,L,M,N,FLO,HI,FINC,NSET,NHI,NDOT)
+C
+C                                                                               
+C +-----------------------------------------------------------------+           
+C |                                                                 |           
+C |                Copyright (C) 1986 by UCAR                       |           
+C |        University Corporation for Atmospheric Research          |           
+C |                    All Rights Reserved                          |           
+C |                                                                 |           
+C |                 NCARGRAPHICS  Version 1.00                      |           
+C |                                                                 |           
+C +-----------------------------------------------------------------+           
+C                                                                               
+C                                                                               
+C
+C
+C
+C     DIMENSION OF           Z(L,N)
+C     ARGUMENTS
+C
+C     LATEST REVISION        JUNE 1984
+C
+C     PURPOSE                CONREC DRAWS A CONTOUR MAP FROM DATA STORED
+C                            IN A RECTANGULAR ARRAY, LABELING THE LINES.
+C
+C     USAGE                  IF THE FOLLOWING ASSUMPTIONS ARE MET, USE
+C
+C                              CALL EZCNTR (Z,M,N)
+C
+C                              ASSUMPTIONS:
+C                                  --ALL OF THE ARRAY IS TO BE CONTOURED.
+C                                  --CONTOUR LEVELS ARE PICKED
+C                                    INTERNALLY.
+C                                  --CONTOURING ROUTINE PICKS SCALE
+C                                    FACTORS.
+C                                  --HIGHS AND LOWS ARE MARKED.
+C                                  --NEGATIVE LINES ARE DRAWN WITH A
+C                                    DASHED LINE PATTERN.
+C                                  --EZCNTR CALLS FRAME AFTER DRAWING THE
+C                                    CONTOUR MAP.
+C
+C                            IF THESE ASSUMPTIONS ARE NOT MET, USE
+C
+C                              CALL CONREC (Z,L,M,N,FLO,HI,FINC,NSET,
+C                                           NHI,NDOT)
+C
+C     ARGUMENTS
+C
+C     ON INPUT               Z
+C     FOR EZCNTR               M BY N ARRAY TO BE CONTOURED.
+C
+C                            M
+C                              FIRST DIMENSION OF Z.
+C
+C                            N
+C                              SECOND DIMENSION OF Z.
+C
+C     ON OUTPUT              ALL ARGUMENTS ARE UNCHANGED.
+C     FOR EZCNTR
+C
+C     ON INPUT               Z
+C     FOR CONREC               THE (ORIGIN OF THE) ARRAY TO BE
+C                              CONTOURED.  Z IS DIMENSIONED L BY N.
+C
+C                            L
+C                              THE FIRST DIMENSION OF Z IN THE CALLING
+C                              PROGRAM.
+C
+C                            M
+C                              THE NUMBER OF DATA VALUES TO BE CONTOURED
+C                              IN THE X-DIRECTION (THE FIRST SUBSCRIPT
+C                              DIRECTION).  WHEN PLOTTING AN ENTIRE
+C                              ARRAY, L = M.
+C
+C                            N
+C                              THE NUMBER OF DATA VALUES TO BE CONTOURED
+C                              IN THE Y-DIRECTION (THE SECOND SUBSCRIPT
+C                              DIRECTION).
+C
+C                            FLO
+C                              THE VALUE OF THE LOWEST CONTOUR LEVEL.
+C                              IF FLO = HI = 0., A VALUE ROUNDED UP FROM
+C                              THE MINIMUM Z IS GENERATED BY CONREC.
+C
+C                            HI
+C                              THE VALUE OF THE HIGHEST CONTOUR LEVEL.
+C                              IF HI = FLO = 0., A VALUE ROUNDED DOWN
+C                              FROM THE MAXIMUM Z IS GENERATED BY
+C                              CONREC.
+C
+C                            FINC
+C                              > 0  INCREMENT BETWEEN CONTOUR LEVELS.
+C                              = 0  A VALUE, WHICH PRODUCES BETWEEN 10
+C                                   AND 30 CONTOUR LEVELS AT NICE VALUES,
+C                                   IS GENERATED BY CONREC.
+C                              < 0  THE NUMBER OF LEVELS GENERATED BY
+C                                   CONREC IS ABS(FINC).
+C
+C                            NSET
+C                              FLAG TO CONTROL SCALING.
+C                              = 0  CONREC AUTOMATICALLY SETS THE
+C                                   WINDOW AND VIEWPORT TO PROPERLY
+C                                   SCALE THE FRAME TO THE STANDARD
+C                                   CONFIGURATION.
+C                                   THE GRIDAL ENTRY PERIM IS
+C                                   CALLED AND TICK MARKS ARE PLACED
+C                                   CORRESPONDING TO THE DATA POINTS.
+C                              > 0  CONREC ASSUMES THAT THE USER
+C                                   HAS SET THE WINDOW AND VIEWPORT
+C                                   IN SUCH A WAY AS TO PROPERLY
+C                                   SCALE THE PLOTTING
+C                                   INSTRUCTIONS GENERATED BY CONREC.
+C                                   PERIM IS NOT CALLED.
+C                              < 0  CONREC GENERATES COORDINATES SO AS
+C                                   TO PLACE THE (UNTRANSFORMED) CONTOUR
+C                                   PLOT WITHIN THE LIMITS OF THE
+C                                   USER'S CURRENT WINDOW AND
+C                                   VIEWPORT.  PERIM IS NOT CALLED.
+C
+C                            NHI
+C                              FLAG TO CONTROL EXTRA INFORMATION ON THE
+C                              CONTOUR PLOT.
+C                              = 0  HIGHS AND LOWS ARE MARKED WITH AN H
+C                                   OR L AS APPROPRIATE, AND THE VALUE
+C                                   OF THE HIGH OR LOW IS PLOTTED UNDER
+C                                   THE SYMBOL.
+C                              > 0  THE DATA VALUES ARE PLOTTED AT
+C                                   EACH Z POINT, WITH THE CENTER OF
+C                                   THE STRING INDICATING THE DATA
+C                                   POINT LOCATION.
+C                              < 0  NEITHER OF THE ABOVE.
+C
+C                            NDOT
+C                              A 10-BIT CONSTANT DESIGNATING THE DESIRED
+C                              DASHED LINE PATTERN.
+C                              IF ABS(NDOT) = 0, 1, OR 1023, SOLID LINES
+C                              ARE DRAWN.
+C                              > 0 NDOT PATTERN IS USED FOR ALL LINES.
+C                              < 0 ABS(NDOT) PATTERN IS USED FOR NEGA-
+C                                TIVE-VALUED CONTOUR LINES, AND SOLID IS
+C                                USED FOR POSITIVE-VALUED CONTOURS.
+C                                CONREC CONVERTS NDOT
+C                                TO A 16-BIT PATTERN AND DASHDB IS USED.
+C                                SEE DASHDB COMMENTS IN THE DASHLINE
+C                                DOCUMENTATION FOR DETAILS.
+C
+C
+C
+C     ON OUTPUT              ALL ARGUMENTS ARE UNCHANGED.
+C     FOR CONREC
+C
+C
+C     ENTRY POINTS           CONREC, CLGEN, REORD, STLINE, DRLINE,
+C                            MINMAX, PNTVAL, CALCNT, EZCNTR, CONBD
+C
+C     COMMON BLOCKS          INTPR, RECINT, CONRE1, CONRE2, CONRE3,
+C                            CONRE4,CONRE5
+C
+C     REQUIRED LIBRARY       STANDARD VERSION: DASHCHAR, WHICH AT
+C     ROUTINES               NCAR ISLOADED BY DEFAULT.
+C                            SMOOTH VERSION: DASHSMTH WHICH MUST BE
+C                            REQUESTED AT NCAR.
+C                            BOTH VERSIONS REQUIRE  GRIDAL, THE
+C                            ERPRT77 PACKAGE, AND THE SPPS.
+C
+C     I/O                    PLOTS CONTOUR MAP.
+C
+C     PRECISION              SINGLE
+C
+C     LANGUAGE               FORTRAN 77
+C
+C     HISTORY                REPLACES OLD CONTOURING PACKAGE CALLED
+C                            CALCNT AT NCAR.
+C
+C     ALGORITHM              EACH LINE IS FOLLOWED TO COMPLETION.  POINTS
+C                            ALONG A LINE ARE FOUND ON BOUNDARIES OF THE
+C                            (RECTANGULAR) CELLS. THESE POINTS ARE
+C                            CONNECTED BY LINE SEGMENTS USING THE
+C                            SOFTWARE DASHED LINE PACKAGE, DASHCHAR.
+C                            DASHCHAR IS ALSO USED TO LABEL THE
+C                            LINES.
+C
+C     NOTE                   TO DRAW NON-UNIFORM CONTOUR LEVELS, SEE
+C                            THE COMMENTS IN CLGEN.  TO MAKE SPECIAL
+C                            MODIFICATIONS FOR SPECIFIC NEEDS SEE THE
+C                            EXPLANATION OF THE INTERNAL PARAMETERS
+C                            BELOW.
+C
+C     TIMING                 VARIES WIDELY WITH SIZE AND SMOOTHNESS OF
+C                            Z.
+C
+C     INTERNAL PARAMETERS    NAME   DEFAULT         FUNCTION
+C                            ----   -------         --------
+C
+C                            ISIZEL   1      SIZE OF LINE LABELS,
+C                                            AS PER THE SIZE DEFINITIONS
+C                                            GIVEN IN THE SPPS
+C                                            DOCUMENTATION FOR WTSTR.
+C
+C                            ISIZEM   2      SIZE OF LABELS FOR MINIMUMS
+C                                            AND MAXIMUMS,
+C                                            AS PER THE SIZE DEFINITIONS
+C                                            GIVEN IN THE SPPS
+C                                            DOCUMENTATION FOR WTSTR.
+C
+C                            ISIZEP   0      SIZE OF LABELS FOR DATA
+C                                            POINT VALUES AS PER THE SIZE
+C                                            DEFINITIONS GIVEN IN THE SPPS
+C                                            DOCUMENTATION FOR WTSTR.
+C
+C                            NLA      16     APPROXIMATE NUMBER OF
+C                                            CONTOUR LEVELS WHEN
+C                                            INTERNALLY GENERATED.
+C
+C                            NLM      40     MAXIMUM NUMBER OF CONTOUR
+C                                            LEVELS.  IF THIS IS TO BE
+C                                            INCREASED, THE DIMENSIONS
+C                                            OF CL AND RWORK IN CONREC
+C                                            MUST BE INCREASED BY THE
+C                                            SAME AMOUNT.
+C
+C                            XLT      .05    LEFT HAND EDGE OF THE PLOT
+C                                            (0.0 IS THE LEFT EDGE OF
+C                                            THE FRAME AND 1.0 IS THE
+C                                            RIGHT EDGE OF THE FRAME.)
+C
+C                            YBT      .05    BOTTOM EDGE OF THE PLOT
+C                                            (0.0 IS THE BOTTOM OF THE
+C                                            FRAME AND 1.0 IS THE TOP
+C                                            OF THE FRAME.)
+C
+C                            SIDE     0.9    LENGTH OF LONGER EDGE OF
+C                                            PLOT (SEE ALSO EXT).
+C
+C                            NREP     6      NUMBER OF REPETITIONS OF
+C                                            THE DASH PATTERN BETWEEN
+C                                            LINE LABELS.
+C
+C                            NCRT     2      NUMBER OF CRT UNITS PER
+C                                            ELEMENT (BIT) IN THE DASH
+C                                            PATTERN.
+C +NOAO - Value of ncrt changed from 4 to 2 in conbd.
+C -NOAO
+C
+C                            ILAB     1      FLAG TO CONTROL THE DRAWING
+C                                            OF LINE LABELS.
+C                                            . ILAB NON-ZERO MEANS LABEL
+C                                              THE LINES.
+C                                            . ILAB = 0 MEANS DO NOT
+C                                              LABEL THE LINES.
+C
+C                            NULBLL   3      NUMBER OF UNLABELED LINES
+C                                            BETWEEN LABELED LINES.  FOR
+C                                            EXAMPLE, WHEN NULBLL = 3,
+C                                            EVERY FOURTH LEVEL IS
+C                                            LABELED.
+C
+C                            IOFFD    0      FLAG TO CONTROL
+C                                            NORMALIZATION OF LABEL
+C                                            NUMBERS.
+C                                            . IOFFD = 0 MEANS INCLUDE
+C                                              DECIMAL POINT WHEN
+C                                              POSSIBLE (DO NOT
+C                                              NORMALIZE UNLESS
+C                                              REQUIRED).
+C                                            . IOFFD NON-ZERO MEANS
+C                                              NORMALIZE ALL LABEL
+C                                              NUMBERS AND OUTPUT A
+C                                              SCALE FACTOR IN THE
+C                                              MESSAGE BELOW THE GRAPH.
+C
+C                            EXT      .0625  LENGTHS OF THE SIDES OF THE
+C                                            PLOT ARE PROPORTIONAL TO M
+C                                            AND N (WHEN CONREC SETS
+C                                            THE WINDOW AND VIEWPORT).
+C                                            IN EXTREME CASES, WHEN
+C                                            MIN(M,N)/MAX(M,N) IS LESS
+C                                            THAN EXT, CONREC
+C                                            PRODUCES A SQUARE PLOT.
+C
+C                            IOFFP      0    FLAG TO CONTROL SPECIAL
+C                                            VALUE FEATURE.
+C                                            . IOFFP = 0 MEANS SPECIAL
+C                                              VALUE FEATURE NOT IN USE.
+C                                            . IOFFP NON-ZERO MEANS
+C                                              SPECIAL VALUE FEATURE IN
+C                                              USE.  (SPVAL IS SET TO THE
+C                                              SPECIAL VALUE.)  CONTOUR
+C                                              LINES WILL THEN BE
+C                                              OMITTED FROM ANY CELL
+C                                              WITH ANY CORNER EQUAL TO
+C                                              THE SPECIAL VALUE.
+C
+C                            SPVAL    0.     CONTAINS THE SPECIAL VALUE
+C                                            WHEN IOFFP IS NON-ZERO.
+C
+C                            IOFFM    0      FLAG TO CONTROL THE MESSAGE
+C                                            BELOW THE PLOT.
+C                                            . IOFFM = 0  IF THE MESSAGE
+C                                              IS TO BE PLOTTED.
+C                                            . IOFFM NON-ZERO IF THE
+C                                              MESSAGE IS TO BE OMITTED.
+C
+C                            ISOLID   1023   DASH PATTERN FOR
+C                                            NON-NEGATIVE CONTOUR LINES.
+C
+C
+C +NOAO - Block data conbd rewritten as run time initialization.
+C     EXTERNAL        CONBD
+C -NOAO
+C
+      SAVE
+      CHARACTER*1     IGAP       ,ISOL       ,RCHAR
+      CHARACTER       ENCSCR*22  ,IWORK*126
+C +NOAO - Character variable added for improved label processing.
+      character*25    string(5)
+C -NOAO
+      DIMENSION       LNGTHS(5)  ,HOLD(5)    ,WNDW(4)    ,VWPRT(4)
+      DIMENSION       Z(L,N)     ,CL(40)     ,RWORK(40)  ,LASF(13)
+      COMMON /INTPR/ PAD1, FPART, PAD(8)
+      COMMON /CONRE1/ IOFFP      ,SPVAL
+      COMMON /CONRE3/ IXBITS     ,IYBITS
+      COMMON /CONRE4/ ISIZEL     ,ISIZEM     ,ISIZEP     ,NREP       ,
+     1                NCRT       ,ILAB       ,NULBLL     ,IOFFD      ,
+     2                EXT        ,IOFFM      ,ISOLID     ,NLA        ,
+     3                NLM        ,XLT        ,YBT        ,SIDE
+      COMMON /CONRE5/ SCLY
+      COMMON /RECINT/ IRECMJ     ,IRECMN     ,IRECTX
+C +NOAO - Value of LNGTHS have been changed from original defaults.  Additional
+C         common block noaolb added for communication with calling routine.
+C
+      common /noaolb/ hold
+      DATA  LNGTHS(1),LNGTHS(2),LNGTHS(3),LNGTHS(4),LNGTHS(5)
+     1      /  13,       4,        21,       10,        19       /
+      DATA  ISOL, IGAP /'$', ''''/
+C
+C -NOAO
+C
+C ISOL AND IGAP (DOLLAR-SIGN AND APOSTROPHE) ARE USED TO CONSTRUCT PAT-
+C TERNS PASSED TO ROUTINE DASHDC IN THE SOFTWARE DASHED-LINE PACKAGE.
+C
+C
+C
+C +NOAO - Blockdata conbd called as run time initialization subroutine
+      call conbd
+C -NOAO
+C
+C THE FOLLOWING CALL IS FOR GATHERING STATISTICS ON LIBRARY USE AT NCAR
+C
+      CALL Q8QST4 ('GRAPHX','CONREC','CONREC','VERSION 01')
+C
+C NONSMOOTHING VERSION
+C
+C
+C
+C  CALL RESET FOR COMPATIBILITY WITH ALL DASH ROUTINES(EXCEPT DASHLINE)
+C
+      CALL RESET
+C
+C  GET NUMBER OF BITS IN INTEGER ARITHMETIC
+C
+      IARTH = I1MACH(8)
+      IXBITS = 0
+      DO 101 I=1,IARTH
+         IF (M .LE. (2**I-1)) GO TO 102
+         IXBITS = I+1
+  101 CONTINUE
+  102 IYBITS = 0
+      DO 103 I=1,IARTH
+         IF (N .LE. (2**I-1)) GO TO 104
+         IYBITS = I+1
+  103 CONTINUE
+  104 IF ((IXBITS*IYBITS).GT.0 .AND. (IXBITS+IYBITS).LE.24) GO TO 105
+C
+C REPORT ERROR NUMBER ONE
+C
+      IWORK =    'CONREC  - DIMENSION ERROR - M*N .GT. (2**IARTH)    M =
+     +            N = '
+C +NOAO
+C
+C     WRITE (IWORK(56:62),'(I6)') M
+      call encode (6, '(i6)', iwork(56:62), m)
+C     WRITE (IWORK(73:79),'(I6)') N
+      call encode (6, '(i6)', iwork(73:79), n)
+C -NOAO
+C
+      CALL SETER( IWORK, 1, 1 )
+      RETURN
+  105 CONTINUE
+C
+C INQUIRE CURRENT TEXT AND LINE COLOR INDEX
+C
+      CALL GQTXCI ( IERR, ITXCI )
+      CALL GQPLCI ( IERR, IPLCI )
+C
+C SET LINE AND TEXT ASF TO INDIVIDUAL
+C
+      CALL GQASF ( IERR, LASF )
+      LSV3  = LASF(3)
+      LSV10 = LASF(10)
+      LASF(3)  = 1
+      LASF(10) = 1
+      CALL GSASF ( LASF )
+C
+      GL = FLO
+      HA = HI
+      GP = FINC
+      MX = L
+      NX = M
+      NY = N
+      IDASH = NDOT
+      NEGPOS = ISIGN(1,IDASH)
+      IDASH = IABS(IDASH)
+      IF (IDASH.EQ.0 .OR. IDASH.EQ.1) IDASH = ISOLID
+C
+C SET CONTOUR LEVELS.
+C
+      CALL CLGEN (Z,MX,NX,NY,GL,HA,GP,NLA,NLM,CL,NCL,ICNST)
+C
+C FIND MAJOR AND MINOR LINES
+C
+      IF (ILAB .NE. 0) CALL REORD (CL,NCL,RWORK,NML,NULBLL+1)
+      IF (ILAB .EQ. 0) NML = 0
+C
+C SAVE CURRENT NORMALIZATION TRANS NUMBER NTORIG AND LOG SCALING FLAG
+C
+      CALL GQCNTN ( IERR, NTORIG )
+      CALL GETUSV ('LS',IOLLS)
+C
+C SET UP SCALING
+C
+      CALL GETUSV ( 'YF' , IYVAL )
+      SCLY = 1.0 / ISHIFT ( 1, 15 - IYVAL )
+C
+      IF (NSET) 106,107,111
+  106 CALL GQNT ( NTORIG,IERR,WNDW,VWPRT )
+      X1 = VWPRT(1)
+      X2 = VWPRT(2)
+      Y1 = VWPRT(3)
+      Y2 = VWPRT(4)
+C
+C SAVE NORMALIZATION TRANS 1
+C
+      CALL GQNT (1,IERR,WNDW,VWPRT)
+C
+C DEFINE NORMALIZATION TRANS AND LOG SCALING
+C
+      CALL SET(X1, X2, Y1, Y2, 1.0, FLOAT(NX), 1.0, FLOAT(NY), 1)
+      GO TO 111
+  107 CONTINUE
+      X1 = XLT
+      X2 = XLT+SIDE
+      Y1 = YBT
+      Y2 = YBT+SIDE
+      X3 = NX
+      Y3 = NY
+      IF (AMIN1(X3,Y3)/AMAX1(X3,Y3) .LT. EXT) GO TO 110
+      IF (NX-NY) 108,110,109
+  108 X2 = SIDE*X3/Y3+XLT
+      GO TO 110
+  109 Y2 = SIDE*Y3/X3+YBT
+C
+C SAVE NORMALIZATION TRANS 1
+C
+  110 CALL GQNT ( 1, IERR, WNDW, VWPRT )
+C
+C DEFINE NORMALIZATION TRANS 1 AND LOG SCALING
+C
+      CALL SET(X1,X2,Y1,Y2,1.0,X3,1.0,Y3,1)
+C
+C DRAW PERIMETER
+C
+      CALL PERIM (NX-1,1,NY-1,1)
+  111 IF (ICNST .NE. 0) GO TO 124
+C
+C SET UP LABEL SCALING
+C
+      IOFFDT = IOFFD
+      IF (GL.NE.0.0 .AND. (ABS(GL).LT.0.1 .OR. ABS(GL).GE.1.E5))
+     1    IOFFDT = 1
+      IF (HA.NE.0.0 .AND. (ABS(HA).LT.0.1 .OR. ABS(HA).GE.1.E5))
+     1    IOFFDT = 1
+      ASH = 10.**(3-IFIX(ALOG10(AMAX1(ABS(GL),ABS(HA),ABS(GP)))-5000.)-
+     1                                                             5000)
+      IF (IOFFDT .EQ. 0) ASH = 1.
+      HOLD(1) = GL
+      HOLD(2) = HA
+      HOLD(3) = GP
+      HOLD(4) = Z(3,3)
+      HOLD(5) = ASH
+      NCHAR = 0
+      IF (IOFFM .NE. 0) GO TO 115
+C +NOAO -  This label generation has been reworked to eliminate the large
+C          spaces in between fields of the label.
+C     IWORK = 'CONTOUR FROM              TO              CONTOUR INTERVAL
+C    1 OF              PT(3,3)=              LABELS SCALED BY'
+      string(1)(1:13) = 'CONTOUR FROM '
+      string(2)(1:4)  = ' TO '
+      string(3)(1:21) = '; CONTOUR INTERVAL = '
+      string(4)(1:11) = '; PT(3,3)= '
+      string(5)(1:19) = '; LABELS SCALED BY '
+C
+      DO 114 I=1,5
+C (NOAO) WRITE ( ENCSCR, '(G13.5)' ) HOLD(I)
+         call encd (hold(i), ash, encscr, nc, ioffd)
+         do 1113 k = 1, lngths(i)
+            nchar = nchar + 1
+ 1113    iwork(nchar:nchar) = string(i)(k:k)
+C
+C (NOAO) NCHAR = NCHAR+LNGTHS(I)
+C (NOAO) DO 113 J=1,13
+         do 113 j = 1, nc
+            NCHAR = NCHAR+1
+            IWORK(NCHAR:NCHAR) = ENCSCR(J:J)
+  113    CONTINUE
+  114 CONTINUE
+C
+C +NOAO IF (ASH .EQ. 1.) NCHAR = NCHAR-13-LNGTHS(5)
+         if (ash .eq. 1.) nchar = nchar - nc - lngths(5)
+C -NOAO
+C
+C SET TEXT INTENSITY TO LOW, AND WRITE TITLE USING NORMALIZATION
+C TRANS NUMBER 0
+C
+      CALL GSTXCI (IRECTX)
+      CALL GETUSV('LS',LSO)
+      CALL SETUSV('LS',1)
+      CALL GSELNT (0)
+C +NOAO - following text output centered on current viewport
+C     CALL WTSTR ( 0.5, 0.015625, IWORK(1:NCHAR), 0, 0, 0 )
+      CALL WTSTR ( ((x1+x2)/2.0), y1 - 0.03, IWORK(1:NCHAR), 0, 0, 0 )
+C -NOAO
+      CALL SETUSV('LS',LSO)
+      CALL GSELNT (1)
+C
+C
+C
+C                       * * * * * * * * * *
+C                            * * * * * * * * * *
+C
+C
+C PROCESS EACH LEVEL
+C
+  115 FPART = .5
+C
+      DO 123 I=1,NCL
+         CONTR = CL(I)
+         NDASH = IDASH
+         IF (NEGPOS.LT.0 .AND. CONTR.GE.0.) NDASH = ISOLID
+C
+C CHANGE 10 BIT PATTERN TO 10 CHARACTER PATTERN.
+C
+         DO 116 J=1,10
+            IBIT = IAND(ISHIFT(NDASH,(J-10)),1)
+            RCHAR = IGAP
+            IF (IBIT .NE. 0) RCHAR = ISOL
+            IWORK(J:J) = RCHAR
+  116    CONTINUE
+         IF (I .GT. NML) GO TO 121
+C
+C SET UP MAJOR LINE (LABELED)
+C
+C SET LINE INTENSITY TO HIGH
+C
+         CALL GSPLCI ( IRECMJ )
+C
+C NREP REPITITIONS OF PATTERN PER LABEL.
+C
+         NCHAR = 10
+         IF (NREP .LT. 2) GO TO 119
+         DO 118 J=1,10
+            NCHAR = J
+            RCHAR = IWORK(J:J)
+            DO 117 K=2,NREP
+               NCHAR = NCHAR+10
+               IWORK(NCHAR:NCHAR) = RCHAR
+  117       CONTINUE
+  118    CONTINUE
+  119    CONTINUE
+C
+C PUT IN LABEL.
+C
+         CALL ENCD (CONTR,ASH,ENCSCR,NCUSED,IOFFDT)
+         DO 120 J=1,NCUSED
+            NCHAR = NCHAR+1
+            IWORK(NCHAR:NCHAR) = ENCSCR(J:J)
+  120    CONTINUE
+         GO TO 122
+C
+C SET UP MINOR LINE (UNLABELED).
+C
+  121    CONTINUE
+C
+C SET LINE INTENSITY TO LOW
+C
+         CALL GSPLCI ( IRECMN )
+         NCHAR = 10
+  122    CALL DASHDC ( IWORK(1:NCHAR),NCRT, ISIZEL )
+C
+C
+C DRAW ALL LINES AT THIS LEVEL.
+C
+         CALL STLINE (Z,MX,NX,NY,CONTR)
+C
+C
+  123 CONTINUE
+C
+C FIND RELATIVE MINIMUMS AND MAXIMUMS IF WANTED, AND MARK VALUES IF
+C WANTED.
+C
+      IF (NHI .EQ. 0) CALL MINMAX (Z,MX,NX,NY,ISIZEM,ASH,IOFFDT)
+      IF (NHI .GT. 0) CALL MINMAX (Z,MX,NX,NY,ISIZEP,-ASH,IOFFDT)
+      FPART = 1.
+      GO TO 127
+  124 CONTINUE
+         IWORK = 'CONSTANT FIELD'
+C +NOAO 
+C     WRITE( ENCSCR, '(G22.14)' ) GL
+      i = gl
+      call encode (22, '(g22.14)', encscr, i)
+C -NOAO
+      DO 126 I=1,22
+         IWORK(I+14:I+14) = ENCSCR(I:I)
+  126 CONTINUE
+C
+C WRITE TITLE USING NORMALIZATION TRNS 0
+C
+      CALL GETUSV('LS',LSO)
+      CALL SETUSV('LS',1)
+      CALL GSELNT (0)
+C +NOAO
+C     CALL WTSTR ( 0.09765, 0.48825, IWORK(1:36), 3, 0, -1 )
+      CALL WTSTR ( x1+0.03, (y1+y2)/2.0, IWORK(1:36), 3, 0, -1 )
+C -NOAO
+C
+C RESTORE NORMALIZATION TRANS 1, LINE AND TEXT INTENSITY TO ORIGINAL
+C
+  127 IF (NSET.LE.0) THEN
+          CALL SET(VWPRT(1),VWPRT(2),VWPRT(3),VWPRT(4),
+     -             WNDW(1),WNDW(2),WNDW(3),WNDW(4),IOLLS)
+      END IF
+      CALL GSPLCI ( IPLCI )
+      CALL GSTXCI ( ITXCI )
+C
+C SELECT ORIGINAL NORMALIZATION TRANS NUMBER NTORIG, AND RESTORE ASF
+C
+      CALL GSELNT ( NTORIG )
+      LASF(3)  = LSV3
+      LASF(10) = LSV10
+      CALL GSASF ( LASF )
+C
+      RETURN
+C
+C
+      END
+      SUBROUTINE CLGEN (Z,MX,NX,NNY,CCLO,CHI,CINC,NLA,NLM,CL,NCL,ICNST)
+      SAVE
+      DIMENSION       CL(NLM)    ,Z(MX,NNY)
+      COMMON /CONRE1/ IOFFP      ,SPVAL
+C
+C CLGEN PUTS THE VALUES OF THE CONTOUR LEVELS IN CL.
+C VARIABLE NAMES MATCH THOSE IN CONREC, WITH THE FOLLOWING ADDITIONS.
+C         NCL     -NUMBER OF CONTOUR LEVELS PUT IN CL.
+C         ICNST   -FLAG TO TELL CONREC IF A CONSTANT FIELD WAS DETECTED.
+C                 .ICNST=0 MEANS NON-CONSTANT FIELD.
+C                 .ICNST NON-ZERO MEANS CONSTANT FIELD.
+C
+C TO PRODUCE NON-UNIFORM CONTOUR LEVEL SPACING, REPLACE THE CODE IN THIS
+C SUBROUTINE WITH CODE TO PRODUCE WHATEVER SPACING IS DESIRED.
+C
+      ICNST = 0
+      NY = NNY
+      CLO = CCLO
+      GLO = CLO
+      HA = CHI
+      FANC = CINC
+      CRAT = NLA
+      IF (HA-GLO) 101,102,111
+  101 GLO = HA
+      HA = CLO
+      GO TO 111
+  102 IF (GLO .NE. 0.) GO TO 120
+      GLO = Z(1,1)
+      HA = Z(1,1)
+      IF (IOFFP .EQ. 0) GO TO 107
+      DO 106 J=1,NY
+         DO 105 I=1,NX
+            IF (Z(I,J) .EQ. SPVAL) GO TO 105
+            GLO = Z(I,J)
+            HA = Z(I,J)
+            DO 104 JJ=J,NY
+               DO 103 II=1,NX
+                  IF (Z(II,JJ) .EQ. SPVAL) GO TO 103
+                  GLO = AMIN1(Z(II,JJ),GLO)
+                  HA = AMAX1(Z(II,JJ),HA)
+  103          CONTINUE
+  104       CONTINUE
+            GO TO 110
+  105    CONTINUE
+  106 CONTINUE
+      GO TO 110
+  107 DO 109 J=1,NY
+         DO 108 I=1,NX
+            GLO = AMIN1(Z(I,J),GLO)
+            HA = AMAX1(Z(I,J),HA)
+  108    CONTINUE
+  109 CONTINUE
+  110 IF (GLO .GE. HA) GO TO 119
+  111 IF (FANC) 112,113,114
+  112 CRAT = AMAX1(1.,-FANC)
+  113 FANC = (HA-GLO)/CRAT
+      P = 10.**(IFIX(ALOG10(FANC)+5000.)-5000)
+      FANC = AINT(FANC/P)*P
+  114 IF (CHI-CLO) 116,115,116
+  115 GLO = AINT(GLO/FANC)*FANC
+      HA = AINT(HA/FANC)*FANC*(1.+SIGN(1.E-6,HA))
+  116 DO 117 K=1,NLM
+         CC = GLO+FLOAT(K-1)*FANC
+         IF (CC .GT. HA) GO TO 118
+         KK = K
+         CL(K) = CC
+  117 CONTINUE
+  118 NCL = KK
+      CCLO = CL(1)
+      CHI = CL(NCL)
+      CINC = FANC
+      RETURN
+  119 ICNST = 1
+      NCL = 1
+      CCLO = GLO
+      RETURN
+  120 CL(1) = GLO
+      NCL = 1
+      RETURN
+      END
+      SUBROUTINE DRLINE (Z,L,MM,NN)
+      SAVE
+      DIMENSION       Z(L,NN)
+C
+C THIS ROUTINE TRACES A CONTOUR LINE WHEN GIVEN THE BEGINNING BY STLINE.
+C TRANSFORMATIONS CAN BE ADDED BY DELETING THE STATEMENT FUNCTIONS FOR
+C FX AND FY IN DRLINE AND MINMAX AND ADDING EXTERNAL FUNCTIONS.
+C X=1. AT Z(1,J), X=FLOAT(M) AT Z(M,J). X TAKES ON NON-INTEGER VALUES.
+C Y=1. AT Z(I,1), Y=FLOAT(N) AT Z(I,N). Y TAKES ON NON-INTEGER VALUES.
+C
+      COMMON /CONRE2/ IX         ,IY         ,IDX        ,IDY        ,
+     1                IS         ,ISS        ,NP         ,CV         ,
+     2                INX(8)     ,INY(8)     ,IR(80000)  ,NR
+c + noao:  dimension of ir array in conre2 changed from 500 to 20000 6March87
+c + noao:  dimension of ir array in conre2 changed from 20000 to 80000 6-93
+      COMMON /CONRE1/ IOFFP      ,SPVAL
+      COMMON /CONRE3/ IXBITS     ,IYBITS
+      LOGICAL         IPEN       ,IPENO
+      DATA IPEN,IPENO/.TRUE.,.TRUE./
+C
+      FX(X,Y) = X
+      FY(X,Y) = Y
+      IXYPAK(IXX,IYY) = ISHIFT(IXX,IYBITS)+IYY
+      C(P1,P2) = (P1-CV)/(P1-P2)
+C
+      M = MM
+      N = NN
+      IF (IOFFP .EQ. 0) GO TO 101
+      ASSIGN 110 TO JUMP1
+      ASSIGN 115 TO JUMP2
+      GO TO 102
+  101 ASSIGN 112 TO JUMP1
+      ASSIGN 117 TO JUMP2
+  102 IX0 = IX
+      IY0 = IY
+      IS0 = IS
+      IF (IOFFP .EQ. 0) GO TO 103
+      IX2 = IX+INX(IS)
+      IY2 = IY+INY(IS)
+      IPEN = Z(IX,IY).NE.SPVAL .AND. Z(IX2,IY2).NE.SPVAL
+      IPENO = IPEN
+  103 IF (IDX .EQ. 0) GO TO 104
+      Y = IY
+      ISUB = IX+IDX
+      X = C(Z(IX,IY),Z(ISUB,IY))*FLOAT(IDX)+FLOAT(IX)
+      GO TO 105
+  104 X = IX
+      ISUB = IY+IDY
+      Y = C(Z(IX,IY),Z(IX,ISUB))*FLOAT(IDY)+FLOAT(IY)
+  105 CALL FRSTD (FX(X,Y),FY(X,Y))
+  106 IS = IS+1
+      IF (IS .GT. 8) IS = IS-8
+      IDX = INX(IS)
+      IDY = INY(IS)
+      IX2 = IX+IDX
+      IY2 = IY+IDY
+      IF (ISS .NE. 0) GO TO 107
+      IF (IX2.GT.M .OR. IY2.GT.N .OR. IX2.LT.1 .OR. IY2.LT.1) GO TO 120
+  107 IF (CV-Z(IX2,IY2)) 108,108,109
+  108 IS = IS+4
+      IX = IX2
+      IY = IY2
+      GO TO 106
+  109 IF (IS/2*2 .EQ. IS) GO TO 106
+      GO TO JUMP1,(110,112)
+  110 ISBIG = IS+(8-IS)/6*8
+      IX3 = IX+INX(ISBIG-1)
+      IY3 = IY+INY(ISBIG-1)
+      IX4 = IX+INX(ISBIG-2)
+      IY4 = IY+INY(ISBIG-2)
+      IPENO = IPEN
+      IF (ISS .NE. 0) GO TO 111
+      IF (IX3.GT.M .OR. IY3.GT.N .OR. IX3.LT.1 .OR. IY3.LT.1) GO TO 120
+      IF (IX4.GT.M .OR. IY4.GT.N .OR. IX4.LT.1 .OR. IY4.LT.1) GO TO 120
+  111 IPEN = Z(IX,IY).NE.SPVAL .AND. Z(IX2,IY2).NE.SPVAL .AND.
+     1       Z(IX3,IY3).NE.SPVAL .AND. Z(IX4,IY4).NE.SPVAL
+  112 IF (IDX .EQ. 0) GO TO 113
+      Y = IY
+      ISUB = IX+IDX
+      X = C(Z(IX,IY),Z(ISUB,IY))*FLOAT(IDX)+FLOAT(IX)
+      GO TO 114
+  113 X = IX
+      ISUB = IY+IDY
+      Y = C(Z(IX,IY),Z(IX,ISUB))*FLOAT(IDY)+FLOAT(IY)
+  114 GO TO JUMP2,(115,117)
+  115 IF (.NOT.IPEN) GO TO 118
+      IF (IPENO) GO TO 116
+C
+C END OF LINE SEGMENT
+C
+      CALL LASTD
+      CALL FRSTD (FX(XOLD,YOLD),FY(XOLD,YOLD))
+C
+C CONTINUE LINE SEGMENT
+C
+  116 CONTINUE
+  117 CALL VECTD (FX(X,Y),FY(X,Y))
+  118 XOLD = X
+      YOLD = Y
+      IF (IS .NE. 1) GO TO 119
+      NP = NP+1
+      IF (NP .GT. NR) GO TO 120
+      IR(NP) = IXYPAK(IX,IY)
+  119 IF (ISS .EQ. 0) GO TO 106
+      IF (IX.NE.IX0 .OR. IY.NE.IY0 .OR. IS.NE.IS0) GO TO 106
+C
+C END OF LINE
+C
+  120 CALL LASTD
+      RETURN
+      END
+      SUBROUTINE MINMAX (Z,L,MM,NN,ISSIZM,AASH,JOFFDT)
+C
+C THIS ROUTINE FINDS RELATIVE MINIMUMS AND MAXIMUMS.  A RELATIVE MINIMUM
+C (OR MAXIMUM) IS DEFINED TO BE THE LOWEST (OR HIGHEST) POINT WITHIN
+C A CERTAIN NEIGHBORHOOD OF THE POINT.  THE NEIGHBORHOOD USED HERE
+C IS + OR - MN IN THE X DIRECTION AND + OR - NM IN THE Y DIRECTION.
+C
+C ORIGINATOR       DAVID KENNISON
+C
+      SAVE
+      CHARACTER*6     IA
+      DIMENSION       Z(L,NN)
+C
+C
+C
+      COMMON /CONRE1/ IOFFP      ,SPVAL
+      COMMON /CONRE5/ SCLY
+C
+      FX(X,Y) = X
+      FY(X,Y) = Y
+C
+      M = MM
+      N = NN
+C
+C SET UP SCALING FOR LABELS
+C
+      SIZEM = (ISSIZM + 1)*256*SCLY
+      ISIZEM = ISSIZM
+C
+      ASH = ABS(AASH)
+      IOFFDT = JOFFDT
+C
+      IF (AASH .LT. 0.0) GO TO 128
+C
+      MN = MIN0(15,MAX0(2,IFIX(FLOAT(M)/8.)))
+      NM = MIN0(15,MAX0(2,IFIX(FLOAT(N)/8.)))
+      NM1 = N-1
+      MM1 = M-1
+C
+C LINE LOOP FOLLOWS - THE COMPLETE TWO-DIMENSIONAL TEST FOR A MINIMUM OR
+C MAXIMUM OF THE FIELD IS ONLY PERFORMED FOR POINTS WHICH ARE MINIMA OR
+C MAXIMA ALONG SOME LINE - FINDING THESE CANDIDATES IS MADE EFFICIENT BY
+C USING A COUNT OF CONSECUTIVE INCREASES OR DECREASES OF THE FUNCTION
+C ALONG THE LINE
+C
+      DO 127 JP=2,NM1
+C
+         IM = MN-1
+         IP = -1
+         GO TO 126
+C
+C CONTROL RETURNS TO STATEMENT 10 AS LONG AS THE FUNCTION IS INCREASING
+C ALONG THE LINE - WE SEEK A POSSIBLE MAXIMUM
+C
+  101    IP = IP+1
+         AA = AN
+         IF (IP .EQ. MM1) GO TO 104
+         AN = Z(IP+1,JP)
+         IF (IOFFP.NE.0 .AND. AN.EQ.SPVAL) GO TO 125
+         IF (AA-AN) 102,103,104
+  102    IM = IM+1
+         GO TO 101
+  103    IM = 0
+         GO TO 101
+C
+C FUNCTION DECREASED - TEST FOR MAXIMUM ON LINE
+C
+  104    IF (IM .GE. MN) GO TO 106
+         IS = MAX0(1,IP-MN)
+         IT = IP-IM-1
+         IF (IS .GT. IT) GO TO 106
+         DO 105 II=IS,IT
+            IF (AA .LE. Z(II,JP)) GO TO 112
+  105    CONTINUE
+  106    IS = IP+2
+         IT = MIN0(M,IP+MN)
+         IF (IS .GT. IT) GO TO 109
+         DO 108 II=IS,IT
+            IF (IOFFP.EQ.0 .OR. Z(II,JP).NE.SPVAL) GO TO 107
+            IP = II-1
+            GO TO 125
+  107       IF (AA .LE. Z(II,JP)) GO TO 112
+  108    CONTINUE
+C
+C WE HAVE MAXIMUM ON LINE - DO TWO-DIMENSIONAL TEST FOR MAXIMUM OF FIELD
+C
+  109    JS = MAX0(1,JP-NM)
+         JT = MIN0(N,JP+NM)
+         IS = MAX0(1,IP-MN)
+         IT = MIN0(M,IP+MN)
+         DO 111 JK=JS,JT
+            IF (JK .EQ. JP) GO TO 111
+            DO 110 IK=IS,IT
+               IF (Z(IK,JK).GE.AA .OR.
+     1             (IOFFP.NE.0 .AND. Z(IK,JK).EQ.SPVAL)) GO TO 112
+  110       CONTINUE
+  111    CONTINUE
+C
+         X = FLOAT(IP)
+         Y = FLOAT(JP)
+         CALL WTSTR ( FX(X,Y),FY(X,Y),'H',ISIZEM,0,0 )
+         CALL FL2INT ( FX(X,Y),FY(X,Y),IFX,IFY )
+C
+C SCALE TO USER SET RESOLUTION
+C
+         IFY = IFY*SCLY
+         CALL ENCD (AA,ASH,IA,NC,IOFFDT)
+         MY = IFY - SIZEM
+         TMY = CPUY ( MY )
+         CALL WTSTR ( FX(X,Y),TMY,IA(1:NC),ISIZEM,0,0 )
+  112    IM = 1
+         IF (IP-MM1) 113,127,127
+C
+C CONTROL RETURNS TO STATEMENT 20 AS LONG AS THE FUNCTION IS DECREASING
+C ALONG THE LINE - WE SEEK A POSSIBLE MINIMUM
+C
+  113    IP = IP+1
+         AA = AN
+         IF (IP .EQ. MM1) GO TO 116
+         AN = Z(IP+1,JP)
+         IF (IOFFP.NE.0 .AND. AN.EQ.SPVAL) GO TO 125
+         IF (AA-AN) 116,115,114
+  114    IM = IM+1
+         GO TO 113
+  115    IM = 0
+         GO TO 113
+C
+C FUNCTION INCREASED - TEST FOR MINIMUM ON LINE
+C
+  116    IF (IM .GE. MN) GO TO 118
+         IS = MAX0(1,IP-MN)
+         IT = IP-IM-1
+         IF (IS .GT. IT) GO TO 118
+         DO 117 II=IS,IT
+            IF (AA .GE. Z(II,JP)) GO TO 124
+  117    CONTINUE
+  118    IS = IP+2
+         IT = MIN0(M,IP+MN)
+         IF (IS .GT. IT) GO TO 121
+         DO 120 II=IS,IT
+            IF (IOFFP.EQ.0 .OR. Z(II,JP).NE.SPVAL) GO TO 119
+            IP = II-1
+            GO TO 125
+  119       IF (AA .GE. Z(II,JP)) GO TO 124
+  120    CONTINUE
+C
+C WE HAVE MINIMUM ON LINE - DO TWO-DIMENSIONAL TEST FOR MINIMUM OF FIELD
+C
+  121    JS = MAX0(1,JP-NM)
+         JT = MIN0(N,JP+NM)
+         IS = MAX0(1,IP-MN)
+         IT = MIN0(M,IP+MN)
+         DO 123 JK=JS,JT
+            IF (JK .EQ. JP) GO TO 123
+            DO 122 IK=IS,IT
+               IF (Z(IK,JK).LE.AA .OR.
+     1             (IOFFP.NE.0 .AND. Z(IK,JK).EQ.SPVAL)) GO TO 124
+  122       CONTINUE
+  123    CONTINUE
+C
+         X = FLOAT(IP)
+         Y = FLOAT(JP)
+         CALL WTSTR ( FX(X,Y),FY(X,Y),'L',ISIZEM,0,0 )
+         CALL FL2INT( FX(X,Y),FY(X,Y),IFX,IFY )
+         IFY = SCLY*IFY
+         CALL ENCD (AA,ASH,IA,NC,IOFFDT)
+         MY = IFY - SIZEM
+         TMY = CPUY ( MY )
+         CALL WTSTR ( FX(X,Y),TMY,IA(1:NC),ISIZEM,0,0 )
+  124    IM = 1
+         IF (IP-MM1) 101,127,127
+C
+C SKIP SPECIAL VALUES ON LINE
+C
+  125    IM = 0
+  126    IP = IP+1
+         IF (IP .GE. MM1) GO TO 127
+         IF (IOFFP.NE.0 .AND. Z(IP+1,JP).EQ.SPVAL) GO TO 125
+         IM = IM+1
+         IF (IM .LE. MN) GO TO 126
+         IM = 1
+         AN = Z(IP+1,JP)
+         IF (Z(IP,JP)-AN) 101,103,113
+C
+  127 CONTINUE
+C
+      RETURN
+C
+C ****************************** ENTRY PNTVAL **************************
+C     ENTRY PNTVAL (Z,L,MM,NN,ISSIZM,AASH,JOFFDT)
+C
+  128 CONTINUE
+      II = (M-1+24)/24
+      JJ = (N-1+48)/48
+      NIQ = 1
+      NJQ = 1
+      DO 130 J=NJQ,N,JJ
+         Y = J
+         DO 129 I=NIQ,M,II
+            X = I
+            ZZ = Z(I,J)
+            IF (IOFFP.NE.0 .AND. ZZ.EQ.SPVAL) GO TO 129
+            CALL ENCD (ZZ,ASH,IA,NC,IOFFDT)
+            CALL WTSTR (FX(X,Y),FY(X,Y),IA(1:NC),ISIZEM,0,0 )
+  129    CONTINUE
+  130 CONTINUE
+      RETURN
+      END
+      SUBROUTINE REORD (CL,NCL,C1,MARK,NMG)
+      SAVE
+      DIMENSION       CL(NCL)    ,C1(NCL)
+C
+C THIS ROUTINE PUTS THE MAJOR (LABELED) LEVELS IN THE BEGINNING OF CL
+C AND THE MINOR (UNLABELED) LEVELS IN END OF CL.  THE NUMBER OF MAJOR
+C LEVELS IS RETURNED IN MARK.  C1 IS USED AS A WORK SPACE.  NMG IS THE
+C NUMBER OF MINOR GAPS (ONE MORE THAN THE NUMBER OF MINOR LEVELS BETWEEN
+C MAJOR LEVELS).
+C
+      NL = NCL
+      IF (NL.LE.4 .OR. NMG.LE.1) GO TO 113
+      NML = NMG-1
+      IF (NL .LE. 10) NML = 1
+C
+C CHECK FOR ZERO OR OTHER NICE NUMBER FOR A MAJOR LINE
+C
+      NMLP1 = NML+1
+      DO 101 I=1,NL
+         ISAVE = I
+         IF (CL(I) .EQ. 0.) GO TO 104
+  101 CONTINUE
+      L = NL/2
+      L = ALOG10(ABS(CL(L)))+1.
+      Q = 10.**L
+      DO 103 J=1,3
+         Q = Q/10.
+         DO 102 I=1,NL
+            ISAVE = I
+            IF (AMOD(ABS(CL(I)+1.E-9*CL(I))/Q,FLOAT(NMLP1)) .LE. .0001)
+     1          GO TO 104
+  102    CONTINUE
+  103 CONTINUE
+      ISAVE = NL/2
+C
+C PUT MAJOR LEVELS IN C1
+C
+  104 ISTART = MOD(ISAVE,NMLP1)
+      IF (ISTART .EQ. 0) ISTART = NMLP1
+      NMAJL = 0
+      DO 105 I=ISTART,NL,NMLP1
+         NMAJL = NMAJL+1
+         C1(NMAJL) = CL(I)
+  105 CONTINUE
+      MARK = NMAJL
+      L = NMAJL
+C
+C PUT MINOR LEVELS IN C1
+C
+      IF (ISTART .EQ. 1) GO TO 107
+      DO 106 I=2,ISTART
+         ISUB = L+I-1
+         C1(ISUB) = CL(I-1)
+  106 CONTINUE
+  107 L = NMAJL+ISTART-1
+      DO 109 I=2,NMAJL
+         DO 108 J=1,NML
+            L = L+1
+            ISUB = ISTART+(I-2)*NMLP1+J
+            C1(L) = CL(ISUB)
+  108    CONTINUE
+  109 CONTINUE
+      NLML = NL-L
+      IF (L .EQ. NL) GO TO 111
+      DO 110 I=1,NLML
+         L = L+1
+         C1(L) = CL(L)
+  110 CONTINUE
+C
+C PUT REORDERED ARRAY BACK IN ORIGINAL PLACE
+C
+  111 DO 112 I=1,NL
+         CL(I) = C1(I)
+  112 CONTINUE
+      RETURN
+  113 MARK = NL
+      RETURN
+      END
+      SUBROUTINE STLINE (Z,LL,MM,NN,CONV)
+      SAVE
+      DIMENSION       Z(LL,NN)
+C
+C THIS ROUTINE FINDS THE BEGINNINGS OF ALL CONTOUR LINES AT LEVEL CONV.
+C FIRST THE EDGES ARE SEARCHED FOR LINES INTERSECTING THE EDGE (OPEN
+C LINES) THEN THE INTERIOR IS SEARCHED FOR LINES WHICH DO NOT INTERSECT
+C THE EDGE (CLOSED LINES).  BEGINNINGS ARE STORED IN IR TO PREVENT RE-
+C TRACING OF LINES.  IF IR IS FILLED, THE SEARCH IS STOPPED FOR THIS
+C CONV.
+C
+      COMMON /CONRE2/ IX         ,IY         ,IDX        ,IDY        ,
+     1                IS         ,ISS        ,NP         ,CV         ,
+     2                INX(8)     ,INY(8)     ,IR(80000)  ,NR
+c + noao:  dimension of ir array in conre2 changed from 500 to 20000 6March87
+c + noao:  dimension of ir array in conre2 changed from 20000 to 80000 6-93
+      COMMON /CONRE3/ IXBITS     ,IYBITS
+C
+C
+C
+C
+C
+C
+      IXYPAK(IXX,IYY) = ISHIFT(IXX,IYBITS)+IYY
+C
+      L = LL
+      M = MM
+      N = NN
+      CV = CONV
+      NP = 0
+      ISS = 0
+      DO 102 IP1=2,M
+         I = IP1-1
+         IF (Z(I,1).GE.CV .OR. Z(IP1,1).LT.CV) GO TO 101
+         IX = IP1
+         IY = 1
+         IDX = -1
+         IDY = 0
+         IS = 1
+         CALL DRLINE (Z,L,M,N)
+  101    IF (Z(IP1,N).GE.CV .OR. Z(I,N).LT.CV) GO TO 102
+         IX = I
+         IY = N
+         IDX = 1
+         IDY = 0
+         IS = 5
+         CALL DRLINE (Z,L,M,N)
+  102 CONTINUE
+      DO 104 JP1=2,N
+         J = JP1-1
+         IF (Z(M,J).GE.CV .OR. Z(M,JP1).LT.CV) GO TO 103
+         IX = M
+         IY = JP1
+         IDX = 0
+         IDY = -1
+         IS = 7
+         CALL DRLINE (Z,L,M,N)
+  103    IF (Z(1,JP1).GE.CV .OR. Z(1,J).LT.CV) GO TO 104
+         IX = 1
+         IY = J
+         IDX = 0
+         IDY = 1
+         IS = 3
+         CALL DRLINE (Z,L,M,N)
+  104 CONTINUE
+      ISS = 1
+      DO 108 JP1=3,N
+         J = JP1-1
+         DO 107 IP1=2,M
+            I = IP1-1
+            IF (Z(I,J).GE.CV .OR. Z(IP1,J).LT.CV) GO TO 107
+            IXY = IXYPAK(IP1,J)
+            IF (NP .EQ. 0) GO TO 106
+            DO 105 K=1,NP
+               IF (IR(K) .EQ. IXY) GO TO 107
+  105       CONTINUE
+  106       NP = NP+1
+            IF (NP .GT. NR) THEN
+C
+C THIS PRINTS AN ERROR MESSAGE IF THE LOCAL ARRAY IR IN SUBROUTINE
+C STLINE HAS AN OVERFLOW
+C THIS MESSAGE IS WRITTEN BOTH ON THE FRAME AND ON THE STANDARD ERROR
+C UNIT
+C
+C +NOAO - Message is written only to stderr, not to the plotting frame.
+C         Error is written with uliber, not FTN write statement.
+C
+      call uliber (1, 'STLINE (CONREC) - WORK ARRAY OVERFLOW', 80)
+      call uliber (1,'STLINE - ***WARNING -- PICTURE INCOMPLETE***',80)
+C     IUNIT = I1MACH(4)
+C     WRITE(IUNIT,1000)
+C1000 FORMAT(
+C    1' WARNING FROM ROUTINE STLINE IN CONREC--WORK ARRAY OVERFLOW')
+C     CALL GETSET(VXA,VXB,VYA,VYB,XA,XB,YA,YB,LTYPE)
+C     Y = (YB - YA) / 2.
+C     X = (XB - XA) / 2.
+C     CALL PWRIT(X,Y,
+C    1'**WARNING--PICTURE INCOMPLETE**',
+C    2 31,3,0,0)
+C     Y = Y * .7
+C     CALL PWRIT(X,Y,
+C    1'WORK ARRAY OVERFLOW IN STLINE',
+C    2 29,3,0,0)
+C -NOAO
+        RETURN
+        ENDIF
+            IR(NP) = IXY
+            IX = IP1
+            IY = J
+            IDX = -1
+            IDY = 0
+            IS = 1
+            CALL DRLINE (Z,L,M,N)
+  107    CONTINUE
+  108 CONTINUE
+      RETURN
+      END
+      SUBROUTINE CALCNT (Z,M,N,A1,A2,A3,I1,I2,I3)
+C
+C THIS ENTRY POINT IS FOR USERS WHO ARE TOO LAZY TO SWITCH OLD DECKS
+C TO THE NEW CALLING SEQUENCE.
+C
+      DIMENSION       Z(M,N)
+      SAVE
+C
+C THE FOLLOWING CALL IS FOR GATHERING STATISTICS ON LIBRARY USE AT NCAR
+C
+      CALL Q8QST4 ('GRAPHX','CONREC','CALCNT','VERSION 01')
+C
+      CALL CONREC (Z,M,M,N,A1,A2,A3,I1,I2,-IABS(I3))
+      RETURN
+      END
+      SUBROUTINE EZCNTR (Z,M,N)
+C
+C CONTOURING VIA SHORTEST POSSIBLE ARGUMENT LIST
+C ASSUMPTIONS --
+C     ALL OF THE ARRAY IS TO BE CONTOURED,
+C     CONTOUR LEVELS ARE PICKED INTERNALLY,
+C     CONTOURING ROUTINE PICKS SCALE FACTORS,
+C     HIGHS AND LOWS ARE MARKED,
+C     NEGATIVE LINES ARE DRAWN WITH A DASHED LINE PATTERN,
+C     EZCNTR CALLS FRAME AFTER DRAWING THE CONTOUR MAP.
+C IF THESE ASSUMPTIONS ARE NOT MET, USE CONREC.
+C
+C ARGUMENTS
+C     Z   ARRAY TO BE CONTOURED
+C     M   FIRST DIMENSION OF Z
+C     N   SECOND DIMENSION OF Z
+C
+      SAVE
+      DIMENSION       Z(M,N)
+      DATA NSET,NHI,NDASH/0,0,682/
+C
+C                             682=1252B
+C THE FOLLOWING CALL IS FOR GATHERING STATISTICS ON LIBRARY USE AT NCAR
+C
+      CALL Q8QST4 ('GRAPHX','CONREC','EZCNTR','VERSION 01')
+C
+      CALL CONREC (Z,M,M,N,0.,0.,0.,NSET,NHI,-NDASH)
+C +NOAO - EZCNTR no longer calls frame.
+C     CALL FRAME
+C -NOAO
+      RETURN
+      END
+C
+C REVISION HISTORY---
+C
+C JANUARY 1980     ADDED REVISION HISTORY AND CHANGED LIBRARY NAME
+C                  FROM CRAYLIB TO PORTLIB FOR MOVE TO PORTLIB
+C
+C MAY 1980         ARRAYS IWORK AND ENCSCR, PREVIOUSLY TOO SHORT FOR
+C                  SHORT-WORD-LENGTH MACHINES, LENGTHENED.  SOME
+C                  DOCUMENTATION CLARIFIED AND CORRECTED.
+C
+C JUNE 1984        CONVERTED TO FORTRAN 77 AND TO GKS
+C
+C JUNE 1985        ERROR HANDLING LINES ADDED; IF OVERFLOW HAPPENS TO
+C                  WORK ARRAY IN STLINE, A WARNING MESSAGE IS WRITTEN
+C                  BOTH ON PLOT FRAME AND ON STANDARD ERROR MESSAGE.
+C-------------------------------------------------------------------
+C