w3fp05.f

Go to the documentation of this file.

C> @file
C> @brief Printer contour subroutine.
C> @author Ralph Jones @date 1989-10-13
 
C> Prints a two-dimensional grid of any shape, with
C> contouring, if desired. grid values are scaled according to
C> to constants specified by the programer, rounded, and printed
C> as 4,3, or 2 digit integers with sign, the sign marking the
C> grid position of the printed number. if contouring is requested,
C> bessel's interpolation formula is used to optain the contour lines.
C> contours are indicated by alphabetic characters ranging from a to
C> h or numeric characters from 0 to 9. contour origin and interval
C> are specified by the programmer in terms of printed values.
C>
C> Program history log:
C> - Ralph Jones 1989-10-13
C> - Ralph Jones 1992-05-02 Add save
C>
C> @param[in] RDATA  Real array of grid data to be printed.
C> @param[in] KTBL Integer array with shape of array.
C> @param[in] CNST Real array of four elements, used in
C> scaling for printing and contouring.
C> @param[in] TITLE  Is a array of 132 characters or less of
C> hollerith data, 1st char. must be blank.
C> printed at bottom of the map.
C> @param[in] KRECT 1 if grid is rectangular, 0 otherwise.
C> @param[in] KCONTR 1 for contouring , 0 otherwise.
C> @param[in] LINEV 0 is for 6 lines per vertical inch,
C> non-zero 8 lines per vertical inch.
C> @param[in] IWIDTH Number of characters in print line,
C> 132 is standard printer.
C>
C> @note Normal subroutine return, unless number of rows is greater than 200,
C> prints error message and exits.
C>
C> @author Ralph Jones @date 1989-10-13
      SUBROUTINE w3fp05(RDATA,KTBL,CNST,TITLE,KRECT,KCONTR,LINEV,IWIDTH)
C
       REAL            CNST(4)
       REAL            RDATA(1)
       REAL            RWA(28)
       REAL            RWB(28)
       REAL            RWC(28)
       REAL            RWD(28)
       REAL            VDJA(29)
       REAL            VDJB(28)
       REAL            VDJC(28)
C
       INTEGER         KALFA(16)
       INTEGER         KALPH(20)
       INTEGER         KHTBL(10)
       INTEGER         KLINE(126)
       INTEGER         KLINES(132)
       INTEGER         KNUMB(20)
       INTEGER         KRLOC(200)
       INTEGER         KTBL(407)
       INTEGER         OUTPUT
       INTEGER         PAGNL
       INTEGER         PAGNR
       INTEGER         PAGN3
       INTEGER         PCCNT
       INTEGER         PCFST
       INTEGER         PGCNT
       INTEGER         PGCNTA
       INTEGER         PGFST
       INTEGER         PGFSTA
       INTEGER         PGMAX
C
       LOGICAL         DONE
       LOGICAL         LCNTR
       LOGICAL         RECT
C
       CHARACTER*1     TITLE(*)
C
       equivalence(crmx,vdja(29))
       equivalence(kline(1),klines(8))
       equivalence(vdjc(1),rwa(1))
C
C      ... THE VAULUE CRMX IS MACHINE DEPENDENT, IT SHOULD BE
C      ... SET TO A VALUE A LITTLE LESS THAN THE LARGEST POSITIVE
C      ... FLOATING  POINT NUMBER FOR THE COMPUTER.
C
       SAVE
C
       DATA  crmx  /10.e70/
       DATA  kalfa/
     a 1ha,1h ,1hb,1h ,1hc,1h ,1hd,1h ,1he,1h ,1hf,
     b 1h ,1hg,1h ,1hh,1h /
       DATA  khastr/1h*/
       DATA  khblnk/1h /
       DATA  khdolr/1h$/
       DATA  khmns /1h-/
       DATA  khplus/1h+/
       DATA  khrstr/1h1/
       DATA  khtbl /1h0,1h1,1h2,1h3,1h4,1h5,1h6,1h7,1h8,1h9/
C
C     ... LIMNRW IS LIMIT ON NUMBER OF ROWS ALLOWED
C     ...   AND IS DIMENSION OF KRLOC ...
C
       DATA  limnrw/200/
       DATA  knumb /1h0,1h ,1h1,1h ,1h2,1h ,1h3,1h ,1h4,1h ,
     1              1h5,1h ,1h6,1h ,1h7,1h ,1h8,1h ,1h9,1h /
       DATA  output/6/
       DATA  r5    /.2/
       DATA  r50   /.02/
C
 8000  FORMAT (1h0,10x,44herror from w3fp05 ... number of rows in your,
     1         9h array = ,i4,24h which exceeds limit of ,i4)
 8100  FORMAT (1ht)
 8200  FORMAT (1hs)
 8300  FORMAT (1h /1h /1h )
 8400  FORMAT (1h /1h )
 8500  FORMAT (132a1)
 8600  FORMAT (132a1)
C
C           COMPUTE VALUES FOR PRINTER WIDTH
C
         IF (iwidth.GE.132.OR.iwidth.LE.0) pgmax = 25
         IF (iwidth.GE.1.AND.iwidth.LE.22) pgmax = 3
         IF (iwidth.GT.22.AND.iwidth.LT.132) pgmax = (iwidth-7)/5
         pagn3 = pgmax + 3
         lw = pgmax * 5 + 7
         vdja(pagn3 + 1) = crmx
         mxpg = pgmax * 5 + 7
C
         IF (linev .EQ. 0) GO TO  100
C     ...OTHERWISE LINEV IS NON-ZERO, SO 8 LINES/INCH IS DESIRED...
           linate = 1
           r4 = 0.250
           r32 = 0.03125
           con2 = 10.0
           nbtwn = 3
           GO TO  200
C
 100     CONTINUE
           linate = 2
           r4 = 0.33333333
           r32 = 1.0/18.0
           con2 = 6.0
           nbtwn = 2
C
 200     CONTINUE
         pgcnta = 0
         pgfsta = 0
         rect  =  .false.
         done = .false.
         kz = 0
         kza = 1000
         a = cnst(1)
         kca = 2*(1-krect)
C           TO SET NO. OF DIGITS TO BE PRINTED
C           WHICH IS A FUNCTION OF THE TENS POSITION IN KCONTR
         nodig = iabs(kcontr/10)
         nodig = 3 - nodig
C           WHERE C(NODIG) + 1 IS NO. OF DIGITS TO BE PRINTED
         IF (nodig.LT.1 .OR. nodig.GT.3) nodig = 3
C           ANY OUT-OF-RANGE WILL GET 4 DIGITS
         lcntr = .false.
         nconq = iabs(mod(kcontr,10))
         IF (nconq .EQ. 0) GO TO 400
         IF (nconq .LE. 2) GO TO 300
C           OTHERWISE RESET NCONQ
         nconq = 0
         GO TO 400
 300     CONTINUE
         lcntr = .true.
C           WITH NCONQ=1 FOR LETTERS,AND =2 FOR NUMBERS IN CONTOUR BANDS
 400     CONTINUE
         IF (nconq .EQ. 2) GO TO 600
C        OTHERWISE SET AS LETTERS
C
         kcow = 16
         DO 500  j = 1,kcow
           kalph(j) = kalfa(j)
 500     CONTINUE
         GO TO  800
C
 600     CONTINUE
         kcow = 20
         DO 700  j = 1,kcow
           kalph(j) = knumb(j)
 700     CONTINUE
C
800      CONTINUE
         radj = 4 * kcow
         kd=1
C *** SET UP TABLE OF INDICES CORRESPONDING TO FIRST ITEM IN EACH ROW
C ***     THIS IS KRLOC
C *** PICK OUT SIZE AND ROW NUMBER OF LARGEST ROW (KCMX AND KCLMX)
C *** KZA LEFT-JUSTIFIES MAP IF ALL ROWS HAVE COMMON MINIMAL OFFSET
         IF (ktbl(1 ).EQ.(-1)) GO TO  1100
C *** ONE-DIMENSIONAL FORM
         ktf=3
         kza=0
         imin = ktbl(2)
         jmax = ktbl(3)+ktbl(1)-1
         nrws = ktbl(1)
         IF (nrws .GT. limnrw) GO TO 1200
         kc = kca * (nrws-1) + 1
C
         DO 1000 j = 1,nrws
           k = nrws-j+1
           krloc(k) = kd
           IF (ktbl(kc+4)+ktbl(kc+3).LE.kz ) GO TO 900
             kclmx = k
             imax = ktbl(kc+4)+ktbl(kc+3)
             kz = imax
             kcmx = krloc(k)+ktbl(kc+4)
 900       CONTINUE
           kd = kd+ktbl(kc+4)
           kc = kc-kca
 1000    CONTINUE
         GO TO  1600
C *** TWO-DIMENSIONAL FORM
C *** THE TWO-DIMENSIONAL FORM IS COMPILER-DEPENDENT
C *** IT DEPENDS ON THE TWO-DIMENSIONAL ARRAY BEING STORED COLUMN-WISE
C *** THAT IS, WITH THE FIRST INDEX VARYING THE FASTEST
 1100    CONTINUE
         imin = ktbl(6)
         jmin = ktbl(7)
         nrws = ktbl(5)
         IF (nrws .LE. limnrw) GO TO  1300
C     ... ELSE, NRWS EXCEEDS LIMIT ALLOWED ...
 1200    CONTINUE
         WRITE (output,8000) nrws,limnrw
         GO TO  7400
C
 1300    CONTINUE
         jmax  = ktbl(7) +ktbl(5)-1
         kc = 1
         DO 1500 j = 1,nrws
           krloc(j) = ktbl(2)*(ktbl(4)-j)+ktbl(kc+7)+1
           IF (ktbl(kc+7)+ktbl(kc+8).LE.kz) GO TO 1400
             imax = ktbl(kc+7)+ktbl(kc+8)
             kz = imax
             kcmx = krloc(j)+ktbl(kc+8)
             kclmx = j
 1400      CONTINUE
           IF (ktbl(kc+7).LT.kza) kza = ktbl(kc+7)
           kc = kc + kca
 1500    CONTINUE
         imax = imax-kza
         ktf = 7
 1600    CONTINUE
         pagnl = 0
         pagnr = pgmax
         IF (.NOT.lcntr) GO TO 1700
           adc = (cnst(1)-cnst(4))/cnst(3)+radj
           bc = cnst(2)/cnst(3)
C *** PRINT I-LABELS ACROSS TOP OF MAP
 1700    CONTINUE
C ***    WHICH PREPARES CDC512 PRINTER FOR 8 LINES PER INCH
         IF (linate.EQ.1) WRITE (output,8100)
C     ...WHICH PREPARES PRINTER FOR 6 LINES PER INCH
         IF (linate.EQ.2) WRITE (output,8200)
         klines(1) = khrstr
         assign 1800 to kbr
         GO TO 6900
C
 1800    CONTINUE
         IF (.NOT.lcntr) GO TO 2000
C *** INITIALIZE CONTOUR WORKING AREA
           DO 1900 j=1,pagn3
             rwc(j)=crmx
             rwd(j)=crmx
 1900      CONTINUE
C *** SET UP CONTOUR DATA AND PAGE LIMITERS FOR FIRST TWO ROWS
C
 2000    CONTINUE
         kra  =  1
         kc = ktf+1
         assign 2100 to kbr
         GO TO 5900
C
 2100    CONTINUE
         kra = 2
         kc = kc+kca
         assign 2200 to kbr
         GO TO 5900
C
 2200    CONTINUE
         kr = 0
C *** TEST IF THIS IS LAST PAGE
         IF (imax.GT.pgmax-1) GO TO 2300
           lmr = imax*5 + 2
           done = .true.
C *** DO LEFT J-LABELS
 2300    CONTINUE
         jcurr = jmax
C
 2400    CONTINUE
         kr = kr + 1
         kra = kr+2
         kc = kc+kca
         kta = mod(jcurr,10)
         ktb = mod(jcurr,100)/10
         ktc  =  mod(jcurr,1000)/100
         IF (kr .EQ. 1 .OR. (.NOT. lcntr)) GO TO 2500
           GO TO  2600
 2500      CONTINUE
           IF (linate.EQ.1) WRITE (output,8300)
           IF (linate.EQ.2) WRITE (output,8400)
 2600    CONTINUE
         klines(2) = khplus
         klines(1) = khblnk
         IF (jcurr.LT.0) klines(2)=khmns
         kta=iabs(kta)
         ktb=iabs(ktb)
         ktc = iabs(ktc)
         IF (ktc .EQ. 0) GO TO 2700
           klines(3) = khtbl(ktc+1)
           klines(4) = khtbl(ktb+1)
           klines(5) = khtbl(kta+1)
           GO TO 2800
C
 2700      CONTINUE
           klines(3) = khtbl(ktb+1)
           klines(4) = khtbl(kta+1)
           klines(5) = khblnk
C
 2800      CONTINUE
         DO 2900  j = 6,mxpg
           klines(j) = khblnk
 2900    CONTINUE
         IF (.NOT.done) GO TO 3000
C *** DO RIGHT J-LABELS IF LAST PAGE OF MAP
           kline(lmr) = klines(2)
           kline(lmr+1) = klines(3)
           kline(lmr+2) = klines(4)
           kline(lmr+3) = klines(5)
C *** FETCH AND CONVERT GRID VALUES TO A1 FORMAT FOR WHOLE LINE
 3000    CONTINUE
         krx = krloc(kr)
         klx = 5*pgfst+1
         IF (pgcnt.EQ.0) GO TO 4000
           DO 3800 kk = 1,pgcnt
           temp = rdata(krx)*cnst(2)+a
           ktemp = abs(temp)+.5
           kline(klx) = khplus
           IF (temp.LT.0.0) kline(klx) = khmns
           GO TO (3300,3200,3100),nodig
 3100        CONTINUE
             kta = mod(ktemp,10000)/1000
C
 3200        CONTINUE
             ktb = mod(ktemp,1000)/100
C
 3300      CONTINUE
           ktc = mod(ktemp,100)/10
           ktd = mod(ktemp,10)
           GO TO (3400,3500,3600),nodig
 3400        CONTINUE
             kline(klx+1) = khtbl(ktc+1)
             kline(klx+2) = khtbl(ktd+1)
             GO TO  3700
 3500        CONTINUE
             kline(klx+1) = khtbl(ktb+1)
             kline(klx+2) = khtbl(ktc+1)
             kline(klx+3) = khtbl(ktd+1)
             GO TO  3700
 3600        CONTINUE
             kline(klx+1) = khtbl(kta+1)
             kline(klx+2) = khtbl(ktb+1)
             kline(klx+3) = khtbl(ktc+1)
             kline(klx+4) = khtbl(ktd+1)
 3700    CONTINUE
         klx = klx + 5
         krx = krx+1
 3800    CONTINUE
C *** FOLLOWING CHECKS FOR POLE POINT AND INSERTS PROPER CHARACTER.
         IF (jcurr.NE.0) GO TO 4000
         IF (imin.LT.(-25).OR.imin.GT.0) GO TO 4000
         kx = -imin
         IF (kx.LT.pgfst.AND.kx.GT.pgcnt+pgfst) GO TO 4000
         kx = 5*kx
         IF (kline(kx+1).EQ.khmns) GO TO 3900
           kline(kx) = khdolr
           GO TO 4000
 3900    CONTINUE
         kline(kx+1) = khastr
C *** PRINT LINE OF MAP DATA
 4000    CONTINUE
           WRITE (output,8500) (klines(ii),ii=1,mxpg)
           krloc(kr) = krx
           jcurr = jcurr - 1
C *** TEST BOTTOM OF MAP
         IF (kr.EQ.nrws) GO TO 5700
C *** SET UP CONTOUR DATA AND PAGE LIMITERS FOR NEXT ROW
         assign 4100 to kbr
         GO TO 5900
C
 4100    CONTINUE
         IF (.NOT.lcntr) GO TO 2400
C *** DO CONTOURING
         DO 4200 jj=1,mxpg
             klines(jj)=khblnk
 4200    CONTINUE
C *** VERTICAL INTERPOLATIONS
         DO 4700 kk = 1,pagn3
           IF (rwb(kk).LT.crmx.AND.rwc(kk).LT.crmx) GO TO  4300
             vdjb(kk) = crmx
             vdjc(kk) = crmx
             GO TO 4600
 4300      CONTINUE
           IF (rwa(kk).LT.crmx.AND.rwd(kk).LT.crmx) GO TO  4400
             vdjc(kk) = 0.
             GO TO  4500
 4400    CONTINUE
         vdjc(kk) = r32*(rwa(kk)+rwd(kk)-rwb(kk)-rwc(kk))
 4500    CONTINUE
         vdjb(kk) = r4*(rwc(kk)-rwb(kk)-con2*vdjc(kk))
 4600    CONTINUE
         vdja(kk)=rwb(kk)
 4700    CONTINUE
C     ...DO 2 OR 3 ROWS OF CONTOURING BETWEEN GRID ROWS...
         DO 5600 ll = 1,nbtwn
           DO 4800 kk = 1,pagn3
             vdjb(kk) = vdjc(kk) + vdjb(kk)
             vdja(kk) = vdjb(kk) + vdja(kk)
 4800      CONTINUE
C     ...WHERE VDJA HAS THE INTERPOLATED VALUE FOR THIS INTER-ROW
C *** HORIZONTAL INTERPOLATIONS
           hdc = 0.0
           IF (vdja(1).GE.crmx) GO TO 4900
           hdc = r50*(vdja(4)+vdja(1)-vdja(2)-vdja(3))
 4900    CONTINUE
         kxb = 0
         DO 5200 kk = 1,pgmax
           IF (vdja(kk+1).GE.crmx) GO TO 5100
           hda = vdja(kk+1)
           IF (vdja(kk+2).GE.crmx) GO TO 5500
           IF (vdja(kk+3).GE.crmx) hdc = 0.
        hdb = r5*(vdja(kk+2)-vdja(kk+1)-15.*hdc)
C *** COMPUTE AND STORE CONTOUR CHARACTERS, 5 PER POINT
           khda=hda
           kdb = iabs(mod(khda,kcow))
           kline(kxb+1) = kalph(kdb+1)
           DO 5000 jj=2,5
           hdb = hdb+hdc
           hda = hda+hdb
           khda = hda
           kdb = iabs(mod(khda,kcow))
           kxa = kxb+jj
           kline(kxa) = kalph(kdb+1)
 5000    CONTINUE
         hdc = r50*(vdja(kk+4)+vdja(kk+1)-vdja(kk+2)-vdja(kk+3))
         IF (vdja(kk+4).GE.crmx) hdc = 0.
 5100    CONTINUE
         kxb = kxb+5
 5200    CONTINUE
 5300    CONTINUE
         WRITE (output,8500) (klines(ii),ii=1,mxpg)
         DO 5400 kk = 1,mxpg
           klines(kk) = khblnk
 5400    CONTINUE
         GO TO 5600
C
 5500    CONTINUE
         khda = hda
         kdb = iabs(mod(khda,kcow))
         kline(kxb+1) = kalph(kdb+1)
         GO TO 5300
 5600    CONTINUE
         GO TO 2400
C
 5700    CONTINUE
         IF (linate.EQ.1)  WRITE (output,8300)
         IF (linate.EQ.2)  WRITE (output,8400)
         klines(1) = khblnk
C *** PRINT I-LABELS ACROSS BOTTOM OF PAGE
         assign 5800 to kbr
         GO TO 6900
C
 5800    CONTINUE
         IF (linate.EQ.1)  WRITE (output,8300)
         IF (linate.EQ.2)  WRITE (output,8400)
C *** PRINT TITLE
         WRITE (output,8600) (title(ii),ii=1,lw)
C *** TEST END OF MAP
           IF (krloc(kclmx).EQ.kcmx) RETURN
C *** ADJUST PAGE LINE BOUNDARIES
C
         IF (imax.GT.pgmax)imax = imax-pgmax
         imin = ka
         pagnl = pagnl + pgmax
         pagnr = pagnr + pgmax
         GO TO   1700
C *** ROUTINE TO PRE-STORE ROWS FOR CONTOURING AND COMPUTE LINE LIMITERS
C
 5900    CONTINUE
         pgfst = pgfsta
         pgcnt = pgcnta
         IF (kra.GT.nrws) GO TO 6800
           krfst = ktbl(kc)-kza
           krcnt = ktbl(kc+1)
           kfx = krloc(kra)
           IF (rect) GO TO  6100
           IF (krfst-pagnl.LE.(-1)) GO TO 6400
             pcfst = krfst-pagnl+1
             IF (pcfst.GE.pagn3) GO TO 6700
               pgfsta = pcfst-1
               pccnt = min(pagnr-krfst+2,krcnt)
               IF (pgfsta.EQ.0) GO TO 6600
               pgcnta = min(pagnr-krfst,krcnt)
               IF (pgcnta.GT.0) GO TO 6000
                 pgcnta = 0
                 GO TO 6100
 6000          CONTINUE
               rect = krect.EQ.1.AND.pgcnta.LE.krcnt
 6100    CONTINUE
         IF (.NOT.lcntr) GO TO kbr,(1800,2100,2200,4100,5800)
         DO 6200 kk = 1,pagn3
           rwa(kk) = rwb(kk)
           rwb(kk) = rwc(kk)
           rwc(kk) = rwd(kk)
           rwd(kk) = crmx
 6200    CONTINUE
C
         IF (pccnt.EQ.0) GO TO kbr,(1800,2100,2200,4100,5800)
           kpc = pcfst+1
           kpd = pccnt
         DO 6300 kk = 1,pccnt
           rwd(kpc) = rdata(kfx)*bc+adc
           kfx = kfx+1
           kpc = kpc + 1
 6300    CONTINUE
         GO TO kbr,(1800,2100,2200,4100,5800)
C
 6400      CONTINUE
           pcfst = 0
           pgfsta = 0
           kfx = kfx-1
           pccnt = krfst+krcnt-pagnl+1
           IF (pccnt.LT.pagn3) GO TO 6500
             pccnt = pagn3
             pgcnta = pgmax
             GO TO 6100
 6500      CONTINUE
           IF (pccnt.GT.0) GO TO 6600
           pgcnta = 0
           pccnt = 0
           GO TO 6100
C
 6600    CONTINUE
         pgcnta = min(pgmax,krcnt+krfst-pagnl)
         GO TO 6100
C
 6700        CONTINUE
             pgcnta = 0
 6800    CONTINUE
         pccnt = 0
         GO TO 6100
C
C *** ROUTINE TO PRINT I-LABELS
C
 6900    CONTINUE
         DO 7000 kk = 2,mxpg
           klines(kk) =  khblnk
 7000    CONTINUE
C
C
         kk = 1
         ka = imin
         lbl = min(imax,pgmax)
C
         DO 7300 jj = 1,lbl
           kline(kk) = khplus
           IF (ka.LT.0) kline(kk) = khmns
           kta = iabs(mod(ka,100))/10
           ktb = iabs(mod(ka,10))
           ktc = iabs(mod(ka,1000))/100
           IF (ktc .EQ. 0) GO TO 7100
             kline(kk+1) = khtbl(ktc+1)
             kline(kk+2) = khtbl(kta+1)
             kline(kk+3) = khtbl(ktb+1)
             GO TO  7200
C
 7100      CONTINUE
           kline(kk+1) = khtbl(kta+1)
           kline(kk+2) = khtbl(ktb+1)
C
 7200    CONTINUE
         kk = kk + 5
         ka = ka+1
 7300    CONTINUE
C
         WRITE (output,8500) (klines(ii),ii=1,mxpg)
C
         GO TO kbr,(1800,2100,2200,4100,5800)
C
 7400    RETURN
C
       END