****************************************************************************
c
c      SHAREWARE NOTICE              Copyright (C) D.I. Hoyer, 1988.
c     ==================             GRAPHLIB.for  v1.0
c
c     If you find these Fortran subroutines useful, you are requested to 
c     register by sending US$20-00 to the address below :-
c
c                    David I Hoyer
c                    31 Rossian Place
c                    Cherrybrook  NSW 2120
c                    AUSTRALIA
c
c     This subroutine library is a shareware product. Copies of the original
c     unmodified programs and manual on disk may be made and distributed as 
c     required, as long as they are not charged for. You may not modify the 
c     source code or manual except for your personal use. 
c
c     Requirements  :  ANSI Fortran 77, full language.
c     
****************************************************************************
c
c     No responsibility is accepted for any errors in this software, 
c     or for any loss or damage resulting from using it.  
c
****************************************************************************

      SUBROUTINE PRTGRF(IOFF,MMAXX,MMAXY,LOHI,IGRAPH)
*
* Print graph IGRAPH on the dot matrix printer.
*
* Set the following character constants after referring to the printer manual:
*
* LSP7 and LSP12 = Set line spacing to 7/72" and 12/72" respectively.
* LSP1  = Set line spacing to 1/216"  [or  1/144"]   (for hi-res plotting)
* LSP20 = Set line spacing to 20/216" [or 13/144"]   ( "    "       "    )
* GFXON = Set graphics mode on, with number of dot columns across page
*             Epson, Star NX-10 : <ESC> 'K'...  =  60 dpi (lo-res)
*                                 <ESC> 'L'...  = 120 dpi (hi-res)
*
* LOHI : Two options are available, viz low and high density.
* 1= Low density = 72 dots per inch vertical and 60 horizontal. Each
*                  row of integers is printed as two lines of 7 rows
*                  each on the printer. 8"x10" = 41k. IGRAPH(480,43)
* 2= High   "    = 144 dpi vert and 120 horiz. This is printed in two
*                  interleaved rows, with 1/216" or 1/144" line feed.
*                  8"x10" requires about 164 kbytes : IGRAPH(960,86)
*
* If the printer cannot do line feeds of 1/216" or 1/144" then only lo-res
* graphs can be printed.
*
      INTEGER*2 IGRAPH
      CHARACTER*5 LSP7,LSP12
      CHARACTER*4 GFXON
      CHARACTER*3 LSP1,LSP20
      CHARACTER BLANK*100,LINE1*1200,LINE2*1200
      DIMENSION IGRAPH(MMAXX,MMAXY)
      DATA BLANK/'
     $                                              '/
*
* This subroutine is set up for IBM, Epson, Star type dot matrix printers
* which have a smallest line feed of 1/216". Some printers have a smallest 
* line feed of 1/144", in which case you should re-define LSP20 to :
*      LSP20 = CHAR(27)//'3'//CHAR(13)
*
* If you have a printer which uses different graphics commands, or if the
* graph is not printing out correctly, consult the printer manual and change
* the appropriate parameters for LSP7, LSP12, LSP1, LSP20, GFXON :
*
      LSP7 = CHAR(27)//'A'//CHAR(7)//CHAR(27)//'2'
      LSP12 = CHAR(27)//'A'//CHAR(12)//CHAR(27)//'2'
      LSP1 = CHAR(27)//'3'//CHAR(1)
      LSP20 = CHAR(27)//'3'//CHAR(20)
*
      IOFF = MAX(1,IOFF)
      I2 = MMAXX/256
      I1 = MMAXX - 256*I2
      IF(LOHI.EQ.1) THEN  ! Set 60 dots per inch across the page (lo-res)
        GFXON = CHAR(27)//'K'//CHAR(I1)//CHAR(I2)
      ELSE  ! Set 120 dots per inch across the page (hi-res)
        GFXON = CHAR(27)//'L'//CHAR(I1)//CHAR(I2)
      ENDIF
      WRITE(6,303) LSP7
 303  FORMAT(1X,A5)
      DO 10 IROW=MMAXY, 1, -1
        DO 30 ICOL=1, MMAXX
          ICH1 = IGRAPH(ICOL,IROW)/128
          ICH2 = 2*(IGRAPH(ICOL,IROW) - 128*ICH1)
          ICH1 = 2*ICH1
          IF(ICH1.EQ.26) ICH1=18
          IF(ICH2.EQ.26) ICH2=18
          LINE1(ICOL:ICOL) = CHAR(ICH1)
          LINE2(ICOL:ICOL) = CHAR(ICH2)
  30    CONTINUE
        IF(LOHI.EQ.1) THEN
          WRITE(6,101) BLANK(1:IOFF),GFXON,LINE1(1:MMAXX),BLANK(1:IOFF),
     $     GFXON,LINE2(1:MMAXX)
 101      FORMAT(1X,A,A5,A/1X,A,A5,A)
        ELSE
          WRITE(6,404) BLANK(1:IOFF),GFXON,LINE1(1:MMAXX),LSP1
          WRITE(6,404) BLANK(1:IOFF),GFXON,LINE2(1:MMAXX),LSP20
 404      FORMAT(1X,A,A5,A,A3)
        ENDIF 
  10  CONTINUE
      WRITE(6,202) LSP12
 202  FORMAT(1X,A5)
      RETURN
      END

*-----------------------------------------------------------------------

      SUBROUTINE SETBIT(I,J)
*
* Set bit j in integer i. This method is marginally faster than using
* the IOR function which is available with some Fortran compilers.
*
      INTEGER*2 I
      DIMENSION MASK(0:14)
      DATA MASK/16384,8192,4096,2048,1024,512,256,128,64,32,16,8,4,2,1/
      IF(MOD(I,MASK(J-1)).LT.MASK(J)) I = I + MASK(J)
      RETURN
      END

*-----------------------------------------------------------------------

      SUBROUTINE PREP(XL,YL,NDIVX,NDIVY,XMINA,XMAXA,YMINA,YMAXA,
     $ XMIN,YMIN,XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,IVH)
*
* Prepare the input data for plotting.
* Calc graph sizes and limits, log conversions etc.
* This subroutine must be called before any plotting is done.
*
* If log scales are specified, then take logs of axis details
      IF(NDIVX.EQ.0) THEN
        XMIN = ALOG10(XMIN)
        XMAX = ALOG10(XMAX)
      ENDIF
      IF(NDIVY.EQ.0) THEN
        YMIN = ALOG10(YMIN)
        YMAX = ALOG10(YMAX)
      ENDIF
* Allow extra space above, below, left & right, for labels and headings.
* The numbers refers to cm here.
      XMINA = XMIN - 1./XL*(XMAX-XMIN)
      XMAXA = XMAX + 1./XL*(XMAX-XMIN)
      YMINA = YMIN - 1./YL*(YMAX-YMIN)
      YMAXA = YMAX + 1./YL*(YMAX-YMIN)
      XL = XL + 2.
      YL = YL + 2.
      IF(IVH.EQ.1) THEN
*------ Vertical format (portrait, or normal)
        WIDTH = XL
        HEIGHT = YL
      ELSE
*------ If horizontal (landscape) format selected, transform width/height
        WIDTH = YL
        HEIGHT = XL
      ENDIF
* Set the number of dots vert. and horiz. for required graph size
      IF(LOHI.EQ.1) THEN
        MAXX = INT(WIDTH*23.62+1.5)
        MAXY = INT(HEIGHT*28.35+1.5)
      ELSE
        MAXX = INT(WIDTH*47.25+1.5)
        MAXY = INT(HEIGHT*56.7+1.5)
      ENDIF
      MMAXX = MAXX
      MMAXY = MAXY/14+1
      IF(IVH.EQ.2) THEN
        I = MAXX
        MAXX = MAXY
        MAXY = I
      ENDIF
      RETURN
      END

*-----------------------------------------------------------------------

      SUBROUTINE AXES(IGRID,NDIVX,NDIVY,NSDIVX,NSDIVY,NDPX,NDPY,
     $ GLABEL,XLABEL,YLABEL,XMINA,XMAXA,YMINA,YMAXA,XMIN,YMIN,XMAX,YMAX,
     $ MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,IVH)
*
* Draw a set of axes on linear or logarithmic scales. 
* IGRID = 0 for no grid lines, 
*       > 0 to draw grid lines, in line type specified by IGRID
*         (eg. IGRID = 1 for solid line, 2 for dotted .....)
* NDIVX, NDIVY = No of divisions along axes 
*                (>0 for linear scale, 0 for log.)
*
      INTEGER*2 IGRAPH
      CHARACTER*80 GLABEL, XLABEL, YLABEL, AXVAL, FMT
      DIMENSION IGRAPH(MMAXX,MMAXY)
      CALL LINE(1,XMIN,YMIN,XMAX,YMIN,XMINA,XMAXA,YMINA,YMAXA,XMIN,YMIN,
     $ XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,IVH,NPP,0)
      CALL LINE(1,XMIN,YMIN,XMIN,YMAX,XMINA,XMAXA,YMINA,YMAXA,XMIN,YMIN,
     $ XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,IVH,NPP,0)
      CALL LINE(1,XMIN,YMAX,XMAX,YMAX,XMINA,XMAXA,YMINA,YMAXA,XMIN,YMIN,
     $ XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,IVH,NPP,0)
      CALL LINE(1,XMAX,YMAX,XMAX,YMIN,XMINA,XMAXA,YMINA,YMAXA,XMIN,YMIN,
     $ XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,IVH,NPP,0)
* Plot graph heading 
      ITXSIZ = LOHI+1
      Y = YMAX + FLOAT(LOHI*15)/FLOAT(MAXY)*(YMAXA-YMINA)
      X = (XMAX+XMIN)/2. - FLOAT((LENG(GLABEL)-1)*(LOHI*2-1)*3)/
     $     FLOAT(MAXX)*(XMAXA-XMINA)
      CALL TEXT(X,Y,1,ITXSIZ,GLABEL,XMINA,XMAXA,YMINA,YMAXA,XMINA,YMINA,
     $ XMAXA,YMAXA,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,IVH)
      DO 30 IAX=0, 1
* Plot x-axis or y-axis label
        IF(IAX.EQ.0) THEN
          NDIVS = NDIVX
          NSDIVS = NSDIVX
          AXMIN = XMIN
          AXMAX = XMAX
          Y = YMIN - FLOAT(LOHI*18)/FLOAT(MAXY)*(YMAXA-YMINA)
          X = (XMAX+XMIN)/2. - FLOAT((LENG(XLABEL)-1)*LOHI*3)/
     $         FLOAT(MAXX)*(XMAXA-XMINA)
          ITXSIZ = LOHI
          CALL TEXT(X,Y,1,ITXSIZ,XLABEL,XMINA,XMAXA,YMINA,YMAXA,XMINA,
     $     YMINA,XMAXA,YMAXA,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,IVH)
        ELSE
          NDIVS = NDIVY
          NSDIVS = NSDIVY
          AXMIN = YMIN
          AXMAX = YMAX
          X = XMIN - FLOAT(LOHI*20)/FLOAT(MAXX)*(XMAXA-XMINA)
          Y = (YMAX+YMIN)/2. - FLOAT((LENG(YLABEL)-1)*LOHI*3)/
     $         FLOAT(MAXY)*(YMAXA-YMINA)
          ITXSIZ = LOHI
          CALL TEXT(X,Y,2,ITXSIZ,YLABEL,XMINA,XMAXA,YMINA,YMAXA,XMINA,
     $     YMINA,XMAXA,YMAXA,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,IVH)
        ENDIF
* Set size of tick marks along axis, if no grid required
        IF(IGRID.EQ.0) THEN
          IF(IAX.EQ.0) AXLL = (YMAXA-YMINA)*3./FLOAT(MAXY-1)
          IF(IAX.EQ.1) AXLL = (XMAXA-XMINA)*3./FLOAT(MAXX-1)
          IF(LOHI.EQ.2) AXLL = 2.*AXLL
        ENDIF
* Calc number of major/minor tick marks on axis
        JDIV = NSDIVS
        IF(JDIV.LT.1) JDIV = 1
        IF(NDIVS.GT.0) THEN  ! Linear axis scale
          I1 = 0
          I2 = NDIVS
        ELSE  ! Log scale
          IF(JDIV.NE.1) JDIV = 9
          I1 = INT(AXMIN) - 1
          I2 = INT(AXMAX) + 1
        ENDIF
* Plot tick mark(s) on the axis
        DO 10 I=I1, I2
          DO 20 J = 1, JDIV
            AXL = AXLL
            IF(J.EQ.1) AXL = AXLL*1.5
            IF(NDIVS.GT.0) THEN
              XY = AXMIN+(AXMAX-AXMIN)/FLOAT(NDIVS)*(I+(J-1)/
     $             FLOAT(NSDIVS))
            ELSE
              XY = ALOG10(10.**FLOAT(I)*FLOAT(J))
            ENDIF
            IF(IGRID.EQ.0) THEN
              IF(IAX.EQ.0) THEN
                CALL LINE(1,XY,YMAX,XY,YMAX-AXL,XMINA,XMAXA,YMINA,YMAXA,
     $           XMIN,YMIN,XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,
     $           IVH,NPP,0)
                CALL LINE(1,XY,YMIN,XY,YMIN+AXL,XMINA,XMAXA,YMINA,YMAXA,
     $           XMIN,YMIN,XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,
     $           IVH,NPP,0)
              ELSE
                CALL LINE(1,XMAX,XY,XMAX-AXL,XY,XMINA,XMAXA,YMINA,YMAXA,
     $           XMIN,YMIN,XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,
     $           IVH,NPP,0)
                CALL LINE(1,XMIN,XY,XMIN+AXL,XY,XMINA,XMAXA,YMINA,YMAXA,
     $           XMIN,YMIN,XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,
     $           IVH,NPP,0)
              ENDIF
            ELSE
              IF(IAX.EQ.0) CALL LINE(IGRID,XY,YMIN,XY,YMAX,XMINA,XMAXA,
     $         YMINA,YMAXA,XMIN,YMIN,XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,
     $         LOHI,IGRAPH,IVH,NPP,0)
              IF(IAX.EQ.1) CALL LINE(IGRID,XMIN,XY,XMAX,XY,XMINA,XMAXA,
     $         YMINA,YMAXA,XMIN,YMIN,XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,
     $         LOHI,IGRAPH,IVH,NPP,0)
            ENDIF
* Plot the axis value at this position
            IF(J.EQ.1.AND.IAX.EQ.0.AND.XY.GE.XMIN.AND.XY.LE.XMAX) THEN
              IF(NDIVS.EQ.0) NDPX = MAX0(0,INT(0.1-XY))
              WRITE(FMT,"('(F20.',I1,')')") NDPX
              WRITE(AXVAL,FMT) XY
              IF(NDIVS.EQ.0) WRITE(AXVAL,FMT) 10.**(XY)
              IB = 0
  40          IB = IB + 1
              IF(AXVAL(IB:IB).EQ.' '.OR.AXVAL(IB:IB).EQ.'0') GOTO 40
              IF(IB.GT.1) AXVAL = AXVAL(IB:LENG(AXVAL))
              IF(LENG(AXVAL).EQ.1.AND.AXVAL(1:1).EQ.'.') AXVAL = '0'
              IB = LENG(AXVAL)
              IF(AXVAL(IB:IB).EQ.'.') AXVAL = AXVAL(1:IB-1)
              Y = YMIN - FLOAT(LOHI*8)/FLOAT(MAXY)*(YMAXA-YMINA)
              X = XY - FLOAT((LENG(AXVAL)-1)*LOHI*3)/
     $             FLOAT(MAXX)*(XMAXA-XMINA)
              CALL TEXT(X,Y,1,ITXSIZ,AXVAL,XMINA,XMAXA,YMINA,YMAXA,XMINA
     $         ,YMINA,XMAXA,YMAXA,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,IVH)
            ELSE IF(J.EQ.1.AND.IAX.EQ.1.AND.XY.GE.YMIN.AND.XY.LE.YMAX) 
     $      THEN
              IF(NDIVS.EQ.0) NDPY = MAX0(0,INT(0.1-XY))
              WRITE(FMT,"('(F20.',I1,')')") NDPY
              WRITE(AXVAL,FMT) XY
              IF(NDIVS.EQ.0) WRITE(AXVAL,FMT) 10.**(XY)
              IB = 0
  50          IB = IB + 1
              IF(AXVAL(IB:IB).EQ.' '.OR.AXVAL(IB:IB).EQ.'0') GOTO 50
              IF(IB.GT.0) AXVAL = AXVAL(IB:LENG(AXVAL))
              IF(LENG(AXVAL).EQ.1.AND.AXVAL(1:1).EQ.'.') AXVAL = '0'
              IB = LENG(AXVAL)
              IF(AXVAL(IB:IB).EQ.'.') AXVAL = AXVAL(1:IB-1)
              X = XMIN - FLOAT(LOHI*8)/FLOAT(MAXX)*(XMAXA-XMINA)
              Y = XY - FLOAT((LENG(AXVAL)-1)*LOHI*3)/
     $             FLOAT(MAXY)*(YMAXA-YMINA)
              CALL TEXT(X,Y,2,ITXSIZ,AXVAL,XMINA,XMAXA,YMINA,YMAXA,XMINA
     $         ,YMINA,XMAXA,YMAXA,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,IVH)
            ENDIF
  20      CONTINUE
  10    CONTINUE
  30  CONTINUE
      RETURN
      END
   
*-----------------------------------------------------------------------
 
      SUBROUTINE PLOTD(NPTS,LTYP,MARK,MSIZE,LEGEND,IORI,ITXSIZ,
     $  XX,YY,XMINA,XMAXA,YMINA,YMAXA,XMIN,YMIN,XMAX,YMAX,MAXX,MAXY,
     $  MMAXX,MMAXY,LOHI,IGRAPH,IVH,STRNG,IFN,X1,X2,P,X,Y,ILEGND,
     $  LEGPOS,NDIVX,NDIVY)
*
* Plot the data for this curve, or this data set.
*
* First plot the legend for this data set, if required
*
      CHARACTER*80 LEGEND
      DIMENSION X(*), Y(*), P(*)
      IF(NPTS.GE.0.AND.(LENG(LEGEND).GT.1.OR.LEGEND(1:1).NE.' ')) THEN
        ILEGND = ILEGND + 1
        YY = YMAX - FLOAT(ILEGND*LOHI*15)/FLOAT(MAXY)*(YMAX-YMIN)
        IF(LEGPOS.LT.1.OR.LEGPOS.GT.8) LEGPOS = 1
        LEGP = LEGPOS
        IF(LEGPOS.GT.4) LEGP = LEGPOS - 4
        IF(LEGPOS.GT.4) YY = YMIN + FLOAT(ILEGND*LOHI*15)/FLOAT(MAXY)
     $   *(YMAX-YMIN)
        DX = (XMAX-XMIN)/FLOAT(MAXX)*FLOAT(LOHI)
        XX = XMIN + 20.*DX + FLOAT(LEGP-1)*(XMAX-XMIN)/4.
        CALL MARKPT(MARK,MSIZE,XX,YY,XMINA,XMAXA,YMINA,YMAXA,XMIN,
     $   YMIN,XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,IVH)
        XX1 = XX - 15.*DX
        XX2 = XX + 15.*DX
        IF(LTYP.NE.0) CALL LINE(IABS(LTYP),XX1,YY,XX2,YY,XMINA,XMAXA,
     $   YMINA,YMAXA,XMIN,YMIN,XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,
     $   IGRAPH,IVH,NPP,0)
        XX = XX + 20.*DX
        CALL TEXT(XX,YY,1,LOHI,LEGEND,XMINA,XMAXA,YMINA,YMAXA,XMIN,
     $   YMIN,XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,IVH)
      ENDIF
      IF(NPTS.LT.0) THEN
* If NPTS<0 then text is to be plotted.
        IF(NDIVX.EQ.0) XX = ALOG10(XX)
        IF(NDIVY.EQ.0) YY = ALOG10(YY)
        CALL TEXT(XX,YY,IORI,ITXSIZ,STRNG,XMINA,XMAXA,YMINA,YMAXA,
     $   XMIN,YMIN,XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,IVH)
      ELSEIF(NPTS.EQ.0) THEN  
* If NPTS=0 then a function plot is required.
        IF(NDIVX.EQ.0) THEN
*-------- Log scale on x-axis, so take logs of X1, X2
          X1 = ALOG10(X1)
          X2 = ALOG10(X2)
        ENDIF
        CALL FUNCT(IFN,P,NDIVX,NDIVY,LTYP,X1,X2,XMINA,XMAXA,YMINA,
     $   YMAXA,XMIN,YMIN,XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,
     $   IVH)
      ELSE
* Plot the set of points
        IF(NDIVX.EQ.0) THEN
*-------- Log scale for x-axis, so take logs of x values
          DO 20 J=1, NPTS
            X(J) = ALOG10(X(J))
  20      CONTINUE
        ENDIF
        IF(NDIVY.EQ.0) THEN
*-------- Log scale for y-axis...
          DO 30 J=1, NPTS
            Y(J) = ALOG10(Y(J))
  30      CONTINUE
        ENDIF
*------ Plot the points
        CALL POINTS(MARK,MSIZE,LTYP,NPTS,X,Y,XMINA,XMAXA,YMINA,YMAXA,
     $   XMIN,YMIN,XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,IVH) 
      ENDIF
      RETURN
      END

*-----------------------------------------------------------------------

      SUBROUTINE POINT(X,Y,XMINA,XMAXA,YMINA,YMAXA,XMIN,YMIN,
     $ XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,IVH)
*
* Plot a single point at the coordinate (X,Y) on IGRAPH
* (XMIN,YMIN) and (XMAX,YMAX) are the lower left and upper right points
* of the graph being plotted.
*   MAXX  = no. of dots along x-axis
*   MAXY  = no. of dots along y-axis
*   MMAXX = no. of columns to be printed out;
*   MMAXY = (no. of graph rows)/14. 
*          Each line printed prints 7 rows of the graph.
*
* Each integer in the array IGRAPH has 16 bits, and each bit which 
* equals 1 is to be printed as a dot.  However, only 14 are used 
* for various technical reasons. In other words the integer at
* IGRAPH(5,2) would contain 14 dot positions corresponding to column
* 5 and rows 17 to 32 of the graph to be printed on the dot matrix 
* printer.  MASKHI converts for Hi-res 144 dpi vertical printing, in 
* which the row is printed as two interleaved rows 1/144" apart.
*
      INTEGER*2 IGRAPH
      DIMENSION IGRAPH(MMAXX,MMAXY), MASKHI(14)
      DATA MASKHI/1,8,2,9,3,10,4,11,5,12,6,13,7,14/
      X = AMAX1(XMIN,AMIN1(XMAX,X))
      Y = AMAX1(YMIN,AMIN1(YMAX,Y))
      IF(IVH.EQ.1) THEN
*------ Vertical format (normal)
        ICOL = INT((X-XMINA)/(XMAXA-XMINA)*FLOAT(MAXX-1)+1.5)
        IY = INT((Y-YMINA)/(YMAXA-YMINA)*FLOAT(MAXY-1)+0.5)
      ELSE
*------ If horizontal format selected, then transform values
        ICOL = MAXY - INT((Y-YMINA)/(YMAXA-YMINA)*FLOAT(MAXY-1)+0.5)
        IY = INT((X-XMINA)/(XMAXA-XMINA)*FLOAT(MAXX-1)+0.5)
      ENDIF
      IROW = IY/14+1
      IBIT = IROW*14-IY
      IF(LOHI.EQ.2) IBIT = MASKHI(IBIT)
      CALL SETBIT(IGRAPH(ICOL,IROW),IBIT)
      RETURN
      END

*-----------------------------------------------------------------------

      SUBROUTINE MARKPT(IPT,ISIZE,X,Y,XMINA,XMAXA,YMINA,YMAXA,
     $ XMIN,YMIN,XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,IVH)
*
* Plot a mark (symbol) of size ISIZE, centred at the point (X,Y) on IGRAPH.
* The symbol shapes are stored as characters, and accessed from TEXT,
* so changes to TEXT may alter these definitions.
*
*   IPT  =  1 : point           
*           2 : open octagon    
*           3 : filled  "       
*           4 : open square     
*           5 : filled  "       
*           6 : open triangle 
*           7 : filled  "
*           8 : cross
*           9 : plus
*          10 : star
*          11 : open diamond
*          12 : filled  "         
*      13-31  :  undefined? Might be used later. Blank for now.
*      32-126 :  plot the corresponding ASCII character
*                (Orientation = 1,  ISIZE = Text size)
*
      INTEGER*2 IGRAPH
      CHARACTER*80 CH, BLNK
      DIMENSION IGRAPH(MMAXX,MMAXY)
      DATA BLNK /' '/
      CH = BLNK
      CALL POINT(X,Y,XMINA,XMAXA,YMINA,YMAXA,XMIN,YMIN,XMAX,YMAX,MAXX,
     $ MAXY,MMAXX,MMAXY,LOHI,IGRAPH,IVH)
      IORI = 1
      IF(IPT.GT.1.AND.IPT.LE.126) THEN 
        CH(1:1) = CHAR(IPT)
        CALL TEXT(X,Y,IORI,ISIZE,CH,XMINA,XMAXA,YMINA,YMAXA,XMIN,YMIN,
     $   XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,IVH)
      ENDIF
      IF(IPT.EQ.10) THEN  
        CH(1:1) = CHAR(8)
        CALL TEXT(X,Y,IORI,ISIZE,CH,XMINA,XMAXA,YMINA,YMAXA,XMIN,YMIN,
     $   XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,IVH)
      ENDIF
      RETURN
      END        

*-----------------------------------------------------------------------

      SUBROUTINE TEXT(X,Y,IORI,ITXSIZ,STRNG,XMINA,XMAXA,YMINA,
     $ YMAXA,XMIN,YMIN,XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,IVH)
*
* Plot a string of text. Also used to plot graph symbols, using character
* positions 1 to 31 in the ASCII sequence to define these symbols.
* (X,Y) is the centre of the capital letters (or symbols)
* IORI = orientation.  1 = Normal (vertical)
*                      2 = Rotated 90 deg anti-clockwise
*                      3 = Upside down
*                      4 = Rotated 90 deg clockwise
*
* IORIX = psuedo-orientation, corrected for IVH
* ITXSIZ= Text size, 1 to n. (Steps of two for lo-res : 2, 4, 6, ...)
*
      INTEGER*2 IGRAPH, ICH, NB, NI, CHARS
      CHARACTER*80 STRNG
      CHARACTER CH
      DIMENSION IGRAPH(MMAXX,MMAXY), CHARS(3,126), MASK2(16), 
     $ IDOTS(10,0:6)
      DATA MASK2/1,2,4,8,16,32,64,128,256,512,1024,2048,4096,8192,16384,
     $           32768/
      DATA CHARS/    0,    0,    0, 4208, 1090, 7185,28784, 1987, 7199, 
     $            4344, 1090,15889,
     $           28920, 1987,15903, 8240, 1153, 3082,24624, 1921, 3086,
     $            8328,  257, 8714,16416, 1984, 2052,16416, 1984, 2052,
     $            8224, 1089, 2058,24608, 1985, 2062,    0,    0,    0,
     $               0,    0,    0,    0,    0,    0,    0,    0,    0,
     $               0,    0,    0,    0,    0,    0,    0,    0,    0,
     $               0,    0,    0,    0,    0,    0,    0,    0,    0,
     $               0,    0,    0,    0,    0,    0,    0,    0,    0,
     $               0,    0,    0,    0,    0,    0,    0,    0,    0,
     $               0,    0,    0,    0,    0,    0,    0,    0,    0,
     $               0,    0,    0,    0, 4000,    0,    0,    7,   56,
     $           30800,17031, 5183,20552, 4066, 9237, 8584,16646, 8969,
     $           18648, 2724, 1297,    0, 3077,    0,    0,18368,   32,
     $            2048, 1988,    0,24648, 2016, 4614,16416, 1984, 2052,
     $            6656,  112,    0,16416,  256, 2052, 6144,   96,    0,
     $            8200,  256, 8200, 2296,18468,15904, 2048,20450,    0,
     $            6276,18596,12580,18692,19236,17972, 8240,17537, 1087,
     $            2504,18981,20008,18552,18722, 1572,  256, 2532,24616,
     $           18648,18724,13860,18624, 2340,15397,22528,  865,    0,
     $           23040,  881,    0, 8224,17473,   32, 8272,  641, 5130,
     $            2048, 1092, 2058,  128, 2212,12324, 2296,19364,15658,
     $            8316, 2178, 7954,30980,18727,13860, 2296,18468, 8736,
     $           30980,18471,15904,18940,18724,16676,16892, 2308,16420,
     $            2296,18468,20004,16892,  256,32516, 2048,20452,   32,
     $            2056, 2080,16447,16892,  640,16657, 2556,16416,  256,
     $             508,  770,32528,  508,  514,32516, 2296,18468,15904,
     $           16892, 2308,12324, 2296,18596,16161,16892, 2436,12581,
     $           18632,18724, 9764,  256, 4068,16416, 2552,16416,32256,
     $            4592,   32,31745, 2552,16832,32256, 8588,  257,25354,
     $             384,  481,24584,10508,18724,24872,30720,18471,   32,
     $             128,  257,  514, 2048,18468,   63,   64, 2050, 4112,
     $             513, 4104,   64,    0, 3072,   40,10248,17057, 3850,
     $            2556,16929, 3592, 2104,16929, 4360, 2104,16929,32520,
     $           10296,17057, 3338,30752, 2307,   36, 2617,21033, 8072,
     $             508,  513, 3848, 2048,19425,    0,  513,25096,   47,
     $            8700,  384,  265, 2048,20452,    0,   60,  481, 3848,
     $             124,  513, 3848, 2104,16929, 3592, 2175,16929, 3592,
     $            2104,16929, 8136,16508,  512, 4104,10276,17057, 4618,
     $           28736,16935,  520, 2168,16416, 7936, 4208,   32, 7169,
     $            2168,16576, 7680,20548,  128, 4357, 2680,20520, 8064,
     $            6212,17057, 4364,16384,18112,16672,    0, 3808,    0,
     $            2308, 1732,    4,  128, 1028, 8200/
      ISIZE = MAX0(1,(ITXSIZ*LOHI+1)/2)
      IF(IORI.LT.1.OR.IORI.GT.4) IORI = 1
* Calc centre of first character, by analogy with POINT.
* First calc integer limits of plotting area
      IF(IVH.EQ.1) THEN
        IMINX = INT((XMIN-XMINA)/(XMAXA-XMINA)*FLOAT(MAXX-1)+1.5)
        IMAXX = INT((XMAX-XMINA)/(XMAXA-XMINA)*FLOAT(MAXX-1)+1.5)
        IMINY = INT((YMIN-YMINA)/(YMAXA-YMINA)*FLOAT(MAXY-1)+0.5)
        IMAXY = INT((YMAX-YMINA)/(YMAXA-YMINA)*FLOAT(MAXY-1)+0.5)
        IORIX = IORI
        IX = INT((X-XMINA)/(XMAXA-XMINA)*FLOAT(MAXX-1)+1.5)
        IY = INT((Y-YMINA)/(YMAXA-YMINA)*FLOAT(MAXY-1)+0.5)
      ELSE
        IMINX = INT((XMIN-XMINA)/(XMAXA-XMINA)*FLOAT(MAXY-1)+0.5)
        IMAXX = INT((XMAX-XMINA)/(XMAXA-XMINA)*FLOAT(MAXY-1)+0.5)
        IMINY = INT((YMIN-YMINA)/(YMAXA-YMINA)*FLOAT(MAXX-1)+1.5)
        IMAXY = INT((YMAX-YMINA)/(YMAXA-YMINA)*FLOAT(MAXX-1)+1.5)
        IORIX = IORI + 1
        IF(IORIX.GT.4) IORIX = 1
        IX = MAXY - INT((Y-YMINA)/(YMAXA-YMINA)*FLOAT(MAXY-1)+0.5)
        IY = INT((X-XMINA)/(XMAXA-XMINA)*FLOAT(MAXX-1)+0.5)
      ENDIF
      DO 60 I=1, LENG(STRNG)
* Calc bottom left of this character
        GOTO(10,20,30,40) IORIX
  10    IC = IX + (I*6-9)*ISIZE
        IR = IY - 6*ISIZE
        GOTO 50
  20    IC = IX + 6*ISIZE
        IR = IY + (I*6-9)*ISIZE
        GOTO 50
  30    IC = IX + (9-I*6)*ISIZE
        IR = IY + 6*ISIZE
        GOTO 50
  40    IC = IX - 6*ISIZE
        IR = IY + (9-I*6)*ISIZE
* Plot the character
  50    IF(ISIZE.GT.1) THEN
          DO 200 J=1, 10
            DO 210 K=0, 6
              IDOTS(J,K) = 0
 210        CONTINUE
 200      CONTINUE
        ENDIF  
        CH = STRNG(I:I)
        ICH = ICHAR(CH)
        DO 70 J=1, 5
          DO 80 K=1, 9
            NI = (J*9+K+5)/15
            NB = J*9+K+6-15*NI
            IF(MOD(CHARS(NI,ICH),MASK2(NB+1)).GE.MASK2(NB)) THEN
              CALL CHRDOT(IORIX,IC,IR,J*ISIZE,K*ISIZE,LOHI,
     $         IMINX,IMAXX,IMINY,IMAXY,MMAXX,MMAXY,IGRAPH)
              IDOTS(K,J) = 1
            ENDIF
  80      CONTINUE
  70    CONTINUE
        IF(ISIZE.GT.1) THEN
          DO 320 ISZ=1, ISIZE-1
            DO 300 J=1, 5
              DO 310 K=1, 9
                JJ = ISIZE*J
                KK = ISIZE*K
                IICH = 0
                IF(ICH.EQ.3.OR.ICH.EQ.5.OR.ICH.EQ.7.OR.ICH.EQ.12) IICH=1
                IF (IDOTS(K,J).EQ.1) THEN
                  IF(IDOTS(K,J+1).EQ.1) CALL CHRDOT(IORIX,IC,IR,JJ+ISZ,
     $             KK,LOHI,IMINX,IMAXX,IMINY,IMAXY,MMAXX,MMAXY,IGRAPH)
                  IF(IDOTS(K+1,J).EQ.1) THEN
                    CALL CHRDOT(IORIX,IC,IR,JJ,KK+ISZ,LOHI,
     $               IMINX,IMAXX,IMINY,IMAXY,MMAXX,MMAXY,IGRAPH)
                  ELSE
                    IF(IDOTS(K,J+1).NE.1.AND.IDOTS(K+1,J+1).EQ.1) CALL
     $               CHRDOT(IORIX,IC,IR,JJ+ISZ,KK+ISZ,LOHI,
     $               IMINX,IMAXX,IMINY,IMAXY,MMAXX,MMAXY,IGRAPH)
                    IF(IDOTS(K,J-1).NE.1.AND.IDOTS(K+1,J-1).EQ.1) CALL
     $               CHRDOT(IORIX,IC,IR,JJ-ISZ,KK+ISZ,LOHI,
     $               IMINX,IMAXX,IMINY,IMAXY,MMAXX,MMAXY,IGRAPH)
                  ENDIF  
                  IF(IICH.EQ.1) THEN
                    IF(IDOTS(K+1,J+1).EQ.1) CALL CHRDOT(IORIX,IC,IR,
     $   JJ+ISZ,KK+ISZ,LOHI,IMINX,IMAXX,IMINY,IMAXY,MMAXX,MMAXY,IGRAPH)
                    IF(IDOTS(K+1,J-1).EQ.1) CALL CHRDOT(IORIX,IC,IR,
     $   JJ-ISZ,KK+ISZ,LOHI,IMINX,IMAXX,IMINY,IMAXY,MMAXX,MMAXY,IGRAPH)
                  ENDIF
                ENDIF
 310          CONTINUE
 300        CONTINUE
 320      CONTINUE
        ENDIF  
  60  CONTINUE
      RETURN
      END

*-----------------------------------------------------------------------

      SUBROUTINE LINE(LTYP,X1,Y1,X2,Y2,XMINA,XMAXA,YMINA,YMAXA,
     $ XMIN,YMIN,XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,IVH,NPP,
     $ IFPBL)
*
* Plot a line from (X1,Y1) to (X2,Y2) on IGRAPH.
* Calls POINT to plot the individual points.
*
* LTYP  = Line type : 1 =  continuous line
*                     2 =  .................
*                     3 =  . . . . . . . . .
*                     4 =  - - - - - - - - -
*                     5 =  -- . -- . -- . --
*
* NPP is the number of points plotted.
* IFPBL = 0 to print both end-points of the line, 
*         else blank first end-point (for plotting chains of lines)
* DOT(ijk) = Array of patterns for dotted lines
*            i=0 or 1 for dot or no dot, j=LTYP, k=LOHI
*
      INTEGER*2 IGRAPH
      DIMENSION IGRAPH(MMAXX,MMAXY), DOT(12,5,2)
      DATA DOT/1,1,1,1,1,1,1,1,1,1,1,1,
     $         1,0,1,0,1,0,1,0,1,0,1,0,
     $         1,0,0,1,0,0,1,0,0,1,0,0,
     $         1,1,1,1,0,0,1,1,1,1,0,0,
     $         1,1,1,0,1,0,1,1,1,0,1,0,
     $         1,1,1,1,1,1,1,1,1,1,1,1,
     $         1,0,0,1,0,0,1,0,0,1,0,0,
     $         1,0,0,0,0,0,1,0,0,0,0,0,
     $         1,1,1,1,1,1,1,1,0,0,0,0,
     $         1,1,1,1,1,1,1,0,0,1,0,0/
      IF(NPP.LT.1) NPP = 1
      IF(LTYP.LT.0.OR.LTYP.GT.5) LTYP = 1
* Calc min. no. of points for an unbroken line
      PNTS = ABS((X2-X1)/(XMAXA-XMINA)*(MAXX-1))
      PNTS = AMAX1(1.,PNTS,ABS((Y2-Y1)/(YMAXA-YMINA)*(MAXY-1)))
      IPNTS = INT(PNTS+0.5)
      DX = (X2-X1)/PNTS
      DY = (Y2-Y1)/PNTS
      ISTART = IFPBL
      IF(ISTART.NE.0) ISTART = 1
      DO 10 I=ISTART, IPNTS
        IF(NPP.GT.12) NPP = 1
        IF(DOT(NPP,LTYP,LOHI).EQ.1) THEN
          X = X1 + I*DX
          Y = Y1 + I*DY
          CALL POINT(X,Y,XMINA,XMAXA,YMINA,YMAXA,XMIN,YMIN,XMAX,YMAX,
     $     MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,IVH)
        ENDIF
        NPP = NPP + 1
  10  CONTINUE
      RETURN
      END

*-----------------------------------------------------------------------

      SUBROUTINE SORT(NPTS,X,Y)
*
* Sort the data pairs into ascending order
*
      DIMENSION X(*),Y(*)
      M = NPTS
      L = M/2 + 1
  10  IF(L.GT.1) THEN
        L = L - 1
        XX = X(L)
        YY = Y(L)
      ELSE
        XX = X(M)
        YY = Y(M)
        X(M) = X(1)
        Y(M) = Y(1)
        M = M - 1
        IF(M.EQ.1) GOTO 30
      ENDIF
      I = L
      J = 2*L
      IF(J.LE.M) THEN
  20    IF(J.LT.M.AND.X(J).LT.X(J+1)) J = J+1
        IF(XX.LT.X(J)) THEN
          X(I) = X(J)
          Y(I) = Y(J)
          I = J
          J = 2*J
        ELSE
          J = M + 1
        ENDIF
        IF(J.LE.M) GO TO 20
      ENDIF
      X(I) = XX
      Y(I) = YY
      GO TO 10
  30  X(1) = XX
      Y(1) = YY
      RETURN
      END

*-----------------------------------------------------------------------

      SUBROUTINE SPLINE(NPIN,NPOUT,X,Y,XX,YY)
*
* Fit a cubic spline to the npin (x,y) pairs, and return npout evenly 
* spaced fitted (xx,yy) data pairs. Also returns sorted (x,y)           
*
      DIMENSION X(*), Y(*), XX(*), YY(*), S(500), C(500)
      CALL SORT(NPIN,X,Y)
      AI = X(2)-X(1) 
      PI = (Y(2)-Y(1))/AI 
      S(1) = -1 
      C(1) = 0 
      DI = -AI 
      CI = 0
      DO 10 I=2,NPIN-1
        A1 = X(I+1)-X(I) 
        Z = 2*(A1+AI)-DI 
        P2 = (Y(I+1)-Y(I))/A1
        C(I) = (6*(P2-PI)-CI)/Z 
        CI = C(I)*A1 
        PI = P2
        S(I) = A1/Z 
        DI = S(I)*A1 
        AI = A1
  10  CONTINUE
      S(NPIN) = C(NPIN-1)/(1+S(NPIN-1)) 
      J = NPIN
      DO 20 I=1,NPIN-1
        J = J-1 
        S(J) = C(J)-S(J)*S(J+1) 
  20  CONTINUE
      IF(NPOUT.GT.500) NPOUT = 500
      DX = (X(NPIN)-X(1))/(NPOUT-1) 
      XX(1) = X(1) 
      YY(1) = Y(1) 
      J = 1
      DO 30 I=1,NPIN-1 
        IF(XX(J)+DX.LT.X(I+1)) THEN
  40      J = J+1 
          XX(J) = XX(J-1)+DX
          XV = XX(J)-X(I) 
          T = 2*S(I)+XV*(S(I)-S(I+1))/(X(I)-X(I+1))+S(I+1)
          YY(J)=Y(I)+XV*((Y(I)-Y(I+1))/(X(I)-X(I+1))+(XX(J)-X(I+1))*T/6)
          IF(XX(J)+DX.LT.X(I+1)) GOTO 40
        ENDIF
  30  CONTINUE
      XX(NPOUT) = X(NPIN) 
      YY(NPOUT) = Y(NPIN)
      RETURN
      END

*-----------------------------------------------------------------------

      SUBROUTINE POINTS(IPT,ISIZE,LTYP,NPTS,X,Y,XMINA,XMAXA,
     $ YMINA,YMAXA,XMIN,YMIN,XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,
     $ IGRAPH,IVH)
*
* Plot an array of points [X(i),Y(i)] on IGRAPH
*
* IPT   = Type of mark to be plotted at each point (see MARKPT)
* ISIZE = Size of mark 
*
* LTYP  = Line type.  0 = plot points only
*                    >0 = fit data with straight line segments
*                    <0 = fit data with a cubic spline
*        
* NPTS  = No of points to be plotted
*
      INTEGER*2 IGRAPH
      DIMENSION X(*), Y(*), IGRAPH(MMAXX,MMAXY), XX(500), YY(500)
      DO 10 I=1, NPTS
        CALL MARKPT(IPT,ISIZE,X(I),Y(I),XMINA,XMAXA,YMINA,YMAXA,XMIN,
     $   YMIN,XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,IVH)
  10  CONTINUE 
      NPP = 1
      IF(LTYP.GT.0) THEN
        DO 20 I=2,NPTS
          CALL LINE(LTYP,X(I-1),Y(I-1),X(I),Y(I),XMINA,XMAXA,YMINA,
     $     YMAXA,XMIN,YMIN,XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,
     $     IVH,NPP,1)
  20    CONTINUE
      ELSEIF(LTYP.LT.0) THEN   ! Spline fit
*------ First sort into ascending x values, then calc and plot the spline
        NP = INT(ABS((X(NPTS)-X(1))/(XMAXA-XMINA)*(MAXX-1)/5.))
        CALL SPLINE(NPTS,NP,X,Y,XX,YY)
        DO 30 I=2, NP
          CALL LINE(-LTYP,XX(I-1),YY(I-1),XX(I),YY(I),XMINA,XMAXA,
     $     YMINA,YMAXA,XMIN,YMIN,XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,
     $     LOHI,IGRAPH,IVH,NPP,1)
  30    CONTINUE
      ENDIF         
      RETURN
      END

*-----------------------------------------------------------------------

      SUBROUTINE FUNCT(IFN,P,NDIVX,NDIVY,LTYP,X1,X2,XMINA,
     $ XMAXA,YMINA,YMAXA,XMIN,YMIN,XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,
     $ LOHI,IGRAPH,IVH) 
*
* Plot a continuous function from X1 to X2 on IGRAPH
*
* NDIVX, NDIVY = 0 for log. scale on X or Y axis, >0 for linear scale.
* IFN   = Function number
* LTYP  = Line type (solid, dotted etc)
*
      INTEGER*2 IGRAPH
      DIMENSION IGRAPH(MMAXX,MMAXY), P(*)
      IPNTS = INT(ABS((X2-X1)/(XMAXA-XMINA)*(MAXX-1)/4.))
      IF(IPNTS.LT.1) IPNTS = 1
      XB = X1
* Take logs as necessary if log scales specified
      IF(NDIVX.GT.0) YB = FGRAPH(IFN,P,X1)
      IF(NDIVX.EQ.0) YB = FGRAPH(IFN,P,10.**X1)
      IF(NDIVY.EQ.0) YB = ALOG10(YB)
      DX = (X2-X1)/FLOAT(IPNTS)
      NPP = 1
      DO 10 I=1, IPNTS
        XA = XB
        YA = YB
        XB = X1 + I*DX
        IF(NDIVX.GT.0) YB = FGRAPH(IFN,P,XB)
        IF(NDIVX.EQ.0) YB = FGRAPH(IFN,P,10.**XB)
        IF(NDIVY.EQ.0) YB = ALOG10(YB)
        CALL LINE(LTYP,XA,YA,XB,YB,XMINA,XMAXA,YMINA,YMAXA,XMIN,YMIN,
     $   XMAX,YMAX,MAXX,MAXY,MMAXX,MMAXY,LOHI,IGRAPH,IVH,NPP,1)
  10  CONTINUE
      RETURN
      END

*-----------------------------------------------------------------------

      SUBROUTINE CLRGRF(MMAXX,MMAXY,IGRAPH)
*
* Set IGRAPH to zero - ie clear the graph
*
      INTEGER*2 IGRAPH
      DIMENSION IGRAPH(MMAXX,MMAXY)
      DO 10 ICOL=1, MMAXX
        DO 20 IROW=1, MMAXY
          IGRAPH(ICOL,IROW) = 0
  20    CONTINUE
  10  CONTINUE
      RETURN
      END

*-----------------------------------------------------------------------

      SUBROUTINE CHRDOT(IORIX,IC,IR,J,K,LOHI,IMINX,IMAXX,IMINY,
     $ IMAXY,MMAXX,MMAXY,IGRAPH)
      INTEGER*2 IGRAPH
      DIMENSION IGRAPH(MMAXX,MMAXY), MASKHI(14)
      DATA MASKHI/1,8,2,9,3,10,4,11,5,12,6,13,7,14/
*
* Set this dot to on (for plotting characters). Called by TEXT
* IC,IR : Col, row for bottom left of character
* J,K   : Offset in dots to current dot position
*
      GOTO(110,120,130,140) IORIX
 110  ICOL = IC + J
      IROW = IR + K
      GOTO 150
 120  ICOL = IC - K
      IROW = IR + J
      GOTO 150
 130  ICOL = IC - J
      IROW = IR - K
      GOTO 150
 140  ICOL = IC + K
      IROW = IR - J
 150  IF(ICOL.GT.IMINX.AND.ICOL.LT.IMAXX.AND.IROW.GT.IMINY.AND.IROW.
     $ LT.IMAXY) THEN
        IRR = IROW/14+1
        IBIT = IRR*14-IROW
        IF(LOHI.EQ.2) IBIT = MASKHI(IBIT)
        CALL SETBIT(IGRAPH(ICOL,IRR),IBIT)
      ENDIF
      RETURN
      END

*-----------------------------------------------------------------------

      FUNCTION LENG(STRING)
*
* Returns the length of STRING, excluding trailing blanks
*
      CHARACTER*80 STRING
      I = LEN(STRING)
      IF(I.LE.1) GOTO 20
  10  I = I - 1
      IF(STRING(I:I).EQ.' '.AND.I.GT.1) GOTO 10
  20  LENG = I
      RETURN
      END

*-----------------------------------------------------------------------