program dops(demo,input,output);
(* dops: pascal graphics library and preprocessor for PostScript

  Thomas D. Schneider, Ph.D.
  toms@alum.mit.edu
  toms@alum.mit.edu (permanent)
  https://alum.mit.edu/www/toms (permanent)

module libraries required: delman, prgmods *)

label 1; (* end of program *)

const
(* begin module version *)
   version = 2.80; (* of dops.p 2010 Nov 03
2010 Nov 03, 2.80: fix doaxis to give highest tic! remove duplicated funcs.
2010 Nov 02, 2.79: fix graphstring
2010 Oct 26, 2.78: allow font control by external code in startpic
2008 Aug 05, 2.77: protect postscript from () by using \ before them
2007 Aug 30, 2.76: doaxis: subtic marks beyond last tic mark! fix
2007 Aug 30, 2.75: doaxis: subtic marks beyond last tic mark!
2005 Aug  6, 2.74: cleanup
2003 Aug 28, 2.73: important upgrade to stringreal
2000 Apr 19, 2.72: alter axes routines to remove line
origin 1988 jan 28 from dosun *)
(* end module version *)

(* begin module describe.dops *)
(*
name
   dops: pascal graphics library and preprocessor for postscript

synopsis
   dops(demo: in, input: in, output: out)

files
   demo: a file for demonstration of the program.
      Start dops interactively.  Start a picture with
          .PS 81 2 2
      then type
          demo
      Graphics instructions will be read from the file 'demo', and the
      corresponding postscript will appear on the output.
      You can try instructions by hand.  Then type
          .PE
          ^d (control-d)
      to conclude.
   input: Graphics instructions.  Portions surrounded by .PS (with
      the appropriate parameters) and .PE (.PS =picture start and .PE =
      picture end) are searched for function names.  When a function
      name is found, the parameters on the same line are read.
   output: the functions detected within .PS to .PE
      are translated into PostScript graphics

description
   Dops converts the graphical instructions made by modules from
   domod.p and produces graphics in the language PostScript.

examples
   To demonstrate the 3-D graphics, use
     .PS 81 2 2
     test3d
     .PE
     (control-d to leave the program)
   A complete test file is called 'demo', which should be run
   non-interactively.

documentation

@article{Schneider1982,
author = "T. D. Schneider
 and G. D. Stormo
 and J. S. Haemer
 and L. Gold",
title = "A design for computer nucleic-acid sequence storage, retrieval and
manipulation",
year = "1982",
journal = "Nucleic Acids Research",
volume = "10",
pages = "3013-3024"}

@article{Schneider1984,
author = "T. D. Schneider
 and G. D. Stormo
 and M. A. Yarus
 and L. Gold",
title = "Delila system tools",
year = "1984",
journal = "Nucleic Acids Research",
volume = "12",
pages = "129-140"}

see also
   doodle.p, domod.p, dosun.p
   PostScript Language Tutorial and Cookbook,
   PostScript Language Language Reference Manual
   both from Addison Wesley, 1985
   dops.demo - file that demonstrates all functions
   xyplo.p - uses colors in plots.

author
  Dr. Thomas D. Schneider
  toms@alum.mit.edu
  permanent email: toms@alum.mit.edu
  https://alum.mit.edu/www/toms/

bugs
   none known

technical notes

   NONSTANDARD is a comment that means that this portion of the code is
   dependent on non-standard pascal (or graphics) for its function.

   See routine setcolor for a discussion of how colors are converted to a
   spectrum.

*)
(* end module describe.dops *)

(* begin module interact.const *)
(* begin module string.const *)
      maxstring = 2000; (* the maximum string *)
(* end module string.const version = 5.27; (@ of prgmod.p 2005 Aug 06 *)
(* end module interact.const version = 5.27; (@ of prgmod.p 2005 Aug 06 *)

(* begin module dops.filler.const *)
      fillermax = 20; (* the size of the filler array for a string *)
(* end module dops.filler.const *)

(* begin module pic.const *)
   pi = 3.14159265354; (* circumference divided by diameter of circle *)
   picwidth = 8; (* width of numbers printed to the file *)
   picdecim = 5; (* number of decimal places for numbers *)
   charwidth = 0.15875; (* the width of characters in cm (ie, cm/char)
                        this allows centering of strings. *)
   (* note: for the Times-Roman  font, 0.0625 is a good value.
            for the Courier-Bold font, 0.08   is a good value. *)
   dotfactor = 0.015875; (* the size of dots *)
{
   defscale = 72; (* default scale factor.  coordinate units per in *)
}
   defscale = 28.35; (* default scale factor.  coordinate units per cm *)
{ making this change would be a big shock to all the programs that
use it, unfortunately.  A major user is xyplo. }
(* end module pic.const version = 9.22; (@ of xyplo.p 2010 Nov 02 *)

type
(* begin module pic.3d.type *)
(* these types are used by the three dimensional graphics routines *)

threevector = array[1..3] of real; (* a point in 3 space *)

tbtarray = array[1..3,1..3] of real; (* a three by three array *)

screen = record (* define a screen for viewing a 3d object *)
   a: threevector; (* center of screen *)
   b: threevector; (* screen x coordinate direction *)
   c: threevector; (* screen y coordinate direction *)
   v: threevector; (* the position of the viewer *)
   g: threevector; (* gaze: viewing direction *)
   smag: real; (* the magnification factor for the screen *)
   range: real; (* 1/smag; the half width of the screen *)
end;
(* end module pic.3d.type *)

(* begin module interact.type *)
(* begin module string.type *)
   stringptr = ^string; (* pointer to a string *)
   string = record (* a string of characters *)
      letters: array[1..maxstring] of char; (* the letters in the string *)
      length: integer; (* the number of characters in the string *)
      current: integer; (* the letter we are working on *)
      next: stringptr; (* the next string in a series *)
   end;
(* end module string.type version = 5.27; (@ of prgmod.p 2005 Aug 06 *)
(* end module interact.type version = 9.22; (@ of xyplo.p 2010 Nov 02 *)

(* begin module trigger.type *)
      trigger = record (* an object to be searched for *)
         seek: string; (* the characters looked for *)
         state: integer; (* how close to triggering we are *)
         skip: boolean; (* trigger not found- skip the line *)
         found: boolean (* the trigger was found *)
      end;
(* end module trigger.type version = 5.27; (@ of prgmod.p 2005 Aug 06 *)

(* begin module filler.type *)
      (* the following is an array used to fill a string.
      it is convenient to have it much shorter than the maxstring, so that
      it is easy to fill the string using procedure fillstring.
      the user must declare the value of constant fillermax. *)
      filler = packed array[1..fillermax] of char;
(* end module filler.type version = 5.27; (@ of prgmod.p 2005 Aug 06 *)

var
(* begin module doodle.var *)
   demo: text; (* demonstration file *)
(* end module doodle.var *)

(* begin module pic.var *)
   inpicture: boolean; (* true if we are drawing the picture,
                          ie, startpic has been called *)
   picxglobal, picyglobal: real; (* absolute location in the graph *)
   pictolerance: real; (* 10 raised to the picwidth,
      to detect values close to zero *)

   scale: real; (* scale factor.  graphic coordinate units per cm *)

(* NONSTANDARD for efficient use of postscript, keep track of
whether there is a current path *)
   inpath: boolean;

(* NONSTANDARD keep track of number of segments drawn so that
they can be stroked.  This (probably) solves the problem of the
Apple printer dying because it can't handle the data. *)
   segments: integer;

   xsideold, ysideold: real; (* current size of a rectangle.  see rectsize *)
(* end module pic.var version = 9.22; (@ of xyplo.p 2010 Nov 02 *)

(* begin module halt *)
procedure halt;
(* stop the program.  the procedure performs a goto to the end of the
   program.  you must have a label:
      label 1;
   declared, and also the end of the program must have this label:
      1: end.
   examples are in the module libraries.
   this is the only goto in the delila system. *)
begin
      writeln(output,' program halt.');
      goto 1
end;
(* end module halt version = 9.22; (@ of xyplo.p 2010 Nov 02 *)

(* begin module interact.clearstring *)
(* begin module clearstring *)
(* These modules clear strings in various ways *)
(* ---- *)
procedure emptystring(var ribbon: string);
(* empty the contents of the string but do NOT
remove the pointer.  This is useful for clearing
one string within a linked list of them. *)
var   index: integer; (* to the ribbon *)
begin (* clearstring *)
   with ribbon do begin
      for index := 1 to maxstring do letters[index] := ' ';
      length := 0;
      current := 0;
   end
end; (* emptystring *)
(* ---- *)
procedure clearstring(var ribbon: string);
(* empty the string and remove the pointer *)
begin (* clearstring *)
   with ribbon do begin
      emptystring(ribbon);
      next := nil;
   end
end; (* clearstring *)
(* ---- *)
procedure initializestring(var ribbon: string);
(* start the string with a nil pointer.   This routine
should be called before doing linked list work.  This allows
the standard string routines to clear the string without
killing the pointer.  This is now deprecated, do not use it
since clearstring still clears the next pointer. *)
begin (* initializestring *)
   writeln(output,'remove initializestring routine!');
   writeln(output,'replace it with clearstring routine!');
   halt; (* to force deprecation *)
   clearstring(ribbon);
   ribbon.next := nil;
end; (* initializestring *)
(* end module clearstring version = 5.27; (@ of prgmod.p 2005 Aug 06 *)
(* end module interact.clearstring version = 9.22; (@ of xyplo.p 2010 Nov 02 *)

(* begin module interact.writestring *)
(* begin module writestring *)
procedure writestring(var tofile: text; var s: string);
(* write the string s to file tofile, no writeln *)
var   i: integer; (* index to s *)
begin (* writestring *)
      with s do for i := 1 to length do write(tofile, letters[i])
end; (* writestring *)
(* end module writestring version = 5.27; (@ of prgmod.p 2005 Aug 06 *)
(* end module interact.writestring version = 9.22; (@ of xyplo.p 2010 Nov 02 *)

(* begin module trigger.proc *)
(* this module allows one to scan a series of characters, as from
an array or a file, and to "trigger" or detect a simple string
in the series.  the advantage of the trigger is that several triggers
can "observe" a stream of characters at once, each looking for a
different thing.
some other modules required: interact.const, interact.type *)

procedure resettrigger(var t: trigger);
(* reset the trigger to ground state *)
begin (* resettrigger *)
      with t do begin
         state := 0;
         skip := false;
         found := false
      end
end; (* resettrigger *)

procedure testfortrigger(ch: char; var t: trigger);
(* look at the character ch.
   if it is part of the trigger (at the current trigger state),
       then the trigger state goes higher.
   if it is not part of the trigger then the trigger state is reset,
      skip is true and one should skip onward to find the trigger.
   if the trigger is found, found is true.
1996 Sep 12: Bug found!  In the case of a trigger "ab", the program
used to miss it for situations like "aab".  This was because at the
first a it would step up.  Then it would see the second a and recognize
that was not part of ab.  It would fail to realize that it could be
the start of a new one.  The code now accounts for that possibility. *)
begin (* testfortrigger *)
      with t do begin
         state := succ(state);
{
writestring(list,seek);
writeln(list,'testfortrigger seek.letters[',state:1,']:',
seek.letters[state],' ch:',ch);
}
         if seek.letters[state] = ch
         then begin
            skip := false;
            if state = seek.length then found := true
                                   else found := false
         end
         else begin
            (* it failed.  But wait!  It could be the beginning
            of a NEW trigger string! *)
            if seek.letters[1] = ch then begin
               state := 1;
               skip := false;
               found := false
            end
            else begin (* reset trigger *)
               state := 0;
               skip := true;
               found := false
            end
         end
      end
end; (* testfortrigger *)
(* end module trigger.proc version = 5.27; (@ of prgmod.p 2005 Aug 06 *)

(* begin module filler.fillstring *)
procedure fillstring(var s: string; a: filler);
(* this procedure makes it reasonably easy to fill the string s with
characters.  one calls the procedure as: *)
(*                           1         2         3         4         5 *)
(*                  12345678901234567890123456789012345678901234567890 *)
(*   fillstring(s, 'this-is-the-string                                ');
the two comments make it easy to line the characters up. also, for this
example, it was assumed that the length of filler as defined by the
constant fillermax was 50. *)
var
      length: integer; (* of the string without trailing blanks *)
      index: integer; (* of s *)
begin (* fillstring *)
      clearstring(s);
      length := fillermax;
      while (length > 1) and (a[length] = ' ') do length := pred(length);
      if (length = 1) and (a[length] = ' ') then begin
         writeln(output, 'fillstring: the string is empty');
         halt
      end;

      for index := 1 to length do s.letters[index] := a[index];
      s.length := length;
      s.current := 1
end; (* fillstring *)
(* end module filler.fillstring version = 5.27; (@ of prgmod.p 2005 Aug 06 *)

(* begin module filler.filltrigger *)
procedure filltrigger(var t: trigger;
                      a: filler);
(* fill the trigger t *)
begin (* filltrigger *)
         fillstring(t.seek,a)
end; (* fillstring *)
(* end module filler.filltrigger version = 5.27; (@ of prgmod.p 2005 Aug 06 *)

(* begin module copyaline *)
procedure copyaline(var fin, fout: text);
(* copy a line from file fin to file fout *)
begin (* copyaline *)
      while not eoln(fin) do begin
         fout^ := fin^;
         put(fout);
         get(fin)
      end;
      readln(fin);
      writeln(fout);
end; (* copyaline *)
(* end module copyaline version = 9.22; (@ of xyplo.p 2010 Nov 02 *)

(* ********************************************************************** *)
(* ********************************************************************** *)
(* ********************************************************************** *)

(* begin module pic.functions *)
(* ********************************************************************** *)
(* begin module pic.await *)
(* end module pic.await *)

(* begin module pic.startpic *)
(* end module pic.startpic *)

(* begin module pic.stoppic *)
(* end module pic.stoppic *)

(* begin module pic.drawr *)
(* end module pic.drawr *)

(* begin module pic.mover *)
(* end module pic.mover *)

(* begin module pic.liner *)
(* end module pic.liner *)


(* begin module pic.drawa *)
(* end module pic.drawa *)

(* begin module pic.movea *)
(* end module pic.movea *)

(* begin module pic.linea *)
(* end module pic.linea *)

(* begin module pic.graphstring *)
(* end module pic.graphstring *)

(* begin module pic.stringinteger *)
(* end module pic.stringinteger *)

(* begin module pic.stringreal *)
(* end module pic.stringreal *)

(* begin module pic.picnumber *)
(* end module pic.picnumber *)

(* begin module pic.xtic *)
(* end module pic.xtic *)

(* begin module pic.ytic *)
(* end module pic.ytic *)

(* begin module pic.doaxis *)
(* end module pic.doaxis *)

(* begin module pic.xaxis *)
(* end module pic.xaxis *)

(* begin module pic.yaxis *)
(* end module pic.yaxis *)
(* ********************************************************************** *)
(* end module pic.functions *)

(* ********************************************************************** *)
(* ********************************************************************** *)
(* ********************************************************************** *)

(* begin module pic.3d.package *)
(* ********************************************************************** *)
(* begin module pic.3d.determinant *)
(* end module pic.3d.determinant *)

(* begin module pic.3d.d32 *)
(* end module pic.3d.d32 *)

(* begin module pic.3d.view *)
(* end module pic.3d.view *)

(* begin module pic.3d.makescreen *)
(* end module pic.3d.makescreen *)

(* begin module pic.3d.project3d *)
(* end module pic.3d.project3d *)
(* ********************************************************************** *)
(* end module pic.3d.package *)

(* ********************************************************************** *)
(* ********************************************************************** *)
(* ********************************************************************** *)

(* begin module pic.startpic *)
procedure startpic(var afile:text; setscale,x,y: real; thefont: char);
(* open the graphics field, with the given scale, and at (x,y)
in that scale.  scale is in device coordinates per cm.
The font is chosen with thefont; t = Times-Roman, c = Courier-Bold *)
(* start pic output to file afile, set the globals *)
(* NONSTANDARD *)
(* this is the actual "world" coordinates used: *)
(*                        xmin,  xmax,  ymin, ymax *)
(*    ns; if (setwindow(-5.0/scale, +5.0/scale, -5.0/scale, +5.0/scale)*)
begin
   writeln(afile,'gsave'); (* save the current graphics state *)

{2005 Aug 6: get rid of these finally
   writeln(afile,'% initgraphics'); (* make sure the printer
           is ready to print, without this, sometimes an Apple laserwriter
           will print the graph upside down, tiny and backwards! *)
   writeln(afile,'% clear erasepage'); (* clean residue from before *)
}

   scale := setscale; (* set the global scale *)

   case thefont of
      'c': begin
           writeln(afile,'/Courier-Bold findfont'); (* locate the font *)
           writeln(afile,10:1,' scalefont'); (* set the font size in points*)
           writeln(afile,'setfont'); (* put the font into the current font *)
      end;
      't': begin
           writeln(afile,'/Times-Roman findfont'); (* locate the font *)
           writeln(afile,12:1,' scalefont'); (* set the font size in points*)
           writeln(afile,'setfont'); (* put the font into the current font *)
      end;
      'n': begin
           writeln(afile,'% /Times-Roman findfont'); (* locate the font *)
           writeln(afile,'%', 12:1,' scalefont'); (* set the font size in points*)
           writeln(afile,'% setfont'); (* put the font into the current font *)
      end;
   end;


(* If the following statement is done then it will work on the sun, but will
kill the applewriter!!!!

Sun's non-standard PostScript extension, setlinewidth has default 0,
as stated in the Read This First and the NeWS Manual.
This draws very quickly with 1 bit wide lines.
If '1 setlinequality' is not done, then one cannot set the width of lines.
So to use PostScript on the screen, I must first do
'1 setlinequality'.  However, if I send this code to the Applewriter,
it kills PostScript on the Applewriter and I get no output whatsoever!
(It took me several hours to figure this out, since once PostScript is killed
on the Applewriter, the NEXT output is also smashed and I had to figure that
out also...)

So a standard PostScript program will not work correctly with the
default.  "Correcting" the PostScript program so that it works on the
Sun means that it BOMBS on the Applewriter.

The default for setlinequality should be '1 setlinequality' so that
the same PostScript code can be used both on the Sun and with other devices.
If you want speed, use the nonstandard form.
An alternative is to redefine setlinequality so that '0 setlinequality'
does give correct results with standard PostScript.

Please review this, Randy.  I think that Sun should fix it.
   writeln(afile,'1 setlinequality');  makes lines at least 1 bit wide *)

   (* set the scale to cm
   writeln(afile,
                 scale:picwidth:picdecim,' ',
                 scale:picwidth:picdecim,' scale'); *)

   (* define some things in postscript *)
   (* doline allows less stuff to be put in the output file.
      it takes two numbers off the stack, copies them, draws a line
      to them as coordinates. *)
(* replaced by 'currentpoint translate'
   writeln(afile,'/doline { 2 copy lineto } def');
*)
   (* define a function that makes cm out of a number *)
(* do this all internally here, it's faster
   writeln(afile,'/i { ',scale:picwidth:picdecim,' mul} def');
*)

   (* move to the start point on the page *)
   writeln(afile,
                     (x*scale):picwidth:picdecim,
                 ' ',(y*scale):picwidth:picdecim,
                 ' translate');

    writeln(afile);
    writeln(afile,'% Define functions so the text produced is smaller');
    writeln(afile,'/a {stroke newpath 0 0} def     % special for arc');
    writeln(afile,'/c {stroke 0 0 moveto} def      % current point');
    writeln(afile,'/f {findfont 10 scalefont setfont} def');
    writeln(afile,'   % to set fonts simply use the f function.  Example:');
    writeln(afile,'   %/Symbol f (\142) /Courier-Bold f (-galactosidase');
    writeln(afile,'/l {lineto} def');
    writeln(afile,'/m {moveto} def');
    writeln(afile,'/n {stroke newpath 0 0 moveto} def'); (* new segment *)
    writeln(afile,'/rl {rlineto} def');
    writeln(afile,'/rm {rmoveto} def');
    writeln(afile,'/s {newpath 0 0 moveto} def     % Start path ');
    writeln(afile,'/t {currentpoint translate} def % translate ');
    writeln(afile,'/x {show} def                   % show teXt ');
    writeln(afile);

   (* start out the pathway *)
   inpath := false;
   (* start the number of segments written: *)
   segments := 0;

   (* now for the normal pic stuff: *)
   inpicture := true;
   picxglobal := 0.0;
   picyglobal := 0.0;
   pictolerance := trunc(exp(picwidth*ln(10))+0.5)
(*;writeln(output,'pictolerance = ',pictolerance:picwidth:picdecim);*)
end;
(* end module pic.startpic version = 2.63; (@ of dops.p 1991 November 2 *)

(* begin module pic.await *)
procedure await;
(* Wait for user to type a carriage return.  the routine assumes that
there is a global file called input. *)
begin
(* the infinite way:
    writeln(output);
    writeln(output,'*********************************');
    writeln(output,'* Use control-c to kill program *');
    writeln(output,'*********************************');
    while true do begin end;*)
end;
(* end module pic.await version = 2.63; (@ of dops.p 1991 November 2 *)

(* begin module pic.stoppic *)
procedure stoppic(var afile:text);
(* stop pic output to file afile *)
(* NONSTANDARD *)
begin
    if inpath then begin
       writeln(afile,'stroke');
       inpath := false
    end;

    writeln(afile,'showpage');
    writeln(afile,'grestore'); (* restore the current graphics state
                                  to what it was before the startpic *)

    await;

    inpicture := false;
end;
(* end module pic.stoppic version = 2.63; (@ of dops.p 1991 November 2 *)

(* begin module pic.drawr *)
procedure drawr(var afile: text; dx,dy: real; visibility: char;
               spacing: real);
(* make a line to file afile by relative draw of dx,dy with visibility
  i invisible
  - dashed
  . dotted
  l line
with the dashes or dots separated by the spacing given
(this has no effect with invisible and line). *)
(* NONSTANDARD *)
var
   ddx,ddy: real; (* changes in dx and dy for dots and dashes *)
   dr: real; (* the hypotenuse, the distance actually drawn *)
   on: boolean; (* draw linesegment if true *)
   y: real; (* the variable for tracking dots and dashes *)
   r: integer; (* number of times to cycle for dots and dashes *)
   ss: real; (* precalculated value to make things a bit faster *)
   theta: real; (* angle of the line *)
procedure checkseg(var afile: text);
(* NONSTANDARD checks how many segments have been written, if
more than 'buffer', stroke them to the postscript page *)
const buffer = 10;
begin
   if segments >= buffer
   then begin
   (* New segment: writeln(afile,'stroke newpath 0 0 moveto'); *)
      writeln(afile,'n');
      segments := 0
   end
   else segments := segments + 1;
end;
begin (* drawr *)
   if not inpath then begin
      (* starts from current coordinates *)
      (* Start path: writeln(afile,'newpath 0 0 moveto'); *)
      writeln(afile,'s');
      inpath := true
   end
   else checkseg(afile);
(* checks
if not (visibility in ['l','i','.','-'])
then writeln(afile,'%YELLLLLL!!!',visibility,'!');
writeln(afile,'% ',visibility,' line');*)

   (* put these on the stack, they will always be used *)
   write(afile,    (dx*scale):picwidth:picdecim,
               ' ',(dy*scale):picwidth:picdecim);
   case visibility of
  'l','i': begin
              case visibility of
                 'i': write(afile,' m');
                 'l': write(afile,' l');
              end
           end;
  '.','-': begin (* make up our own dots and dashes *)
              writeln(afile); (* move away from the (dx,dy) on the stack *)
              if spacing <= 0.0 then begin
                 writeln(output,'drawr: spacing zero with . or - line');
                 halt
              end;
              if dx = 0.0
              then begin
                 ddx := 0.0; (* avoid division by zero *)
                 ddy := scale*spacing;
                 if dy < 0 then ddy := - ddy; (* this makes sure that
                    we draw lines straight down if that was the request *)
              end
              else begin
                 (* find out the angle of the slope, intentionally
                    lose the sign *)
                 theta := arctan(abs(dy/dx));

                 ddx := scale*spacing*cos(theta);
                 ddy := scale*spacing*sin(theta);

                 (* return the sign to the little buggers *)
                 if dx < 0 then ddx := -ddx;
                 if dy < 0 then ddy := -ddy;
              end;
              y := 0;
              case visibility of
                 '.': ss := scale*dotfactor;
                 '-': on := true;
              end;

              dr := sqrt(dx*dx+dy*dy);
              for r := 1 to round(dr/spacing) do begin
                 case visibility of
                 '-': begin
                         write(afile,    (ddx):picwidth:picdecim,
                                     ' ',(ddy):picwidth:picdecim);
                         if on
                         then writeln(afile,' rl')
                         else writeln(afile,' rm');
                         on := not on
                      end;
                 '.': begin
                         (* put out a dot like in dotr *)
                         write(afile,
                                    +ss:picwidth:picdecim,' 0 rl');
                         write(afile,' ',
                                    -ss:picwidth:picdecim,' 0 rl');
                         write(afile,' ',(ddx):picwidth:picdecim,
                                     ' ',(ddy):picwidth:picdecim);
                         writeln(afile,' rm');
                      end;
                 end
              end;
              (* let's make really sure we got there!! *)
              writeln(afile,' m'); (* pulled from the stack *)
           end;
   end;

   (* an elegant way to make postscript keep a global record is
   to translate the coordinates! *)
   (* writeln(afile,' currentpoint translate'); *)
   writeln(afile,' t');

   picxglobal := picxglobal + dx;
   picyglobal := picyglobal + dy;
end;
(* end module pic.drawr version = 2.63; (@ of dops.p 1991 November 2 *)

(* begin module pic.mover *)
procedure mover(var afile: text; dx,dy: real);
(* move relative the amount (dx, dy). *)
begin
   drawr(afile,dx,dy,'i',0.0);
end;
(* end module pic.mover version = 2.63; (@ of dops.p 1991 November 2 *)

(* begin module pic.liner *)
procedure liner(var afile: text; dx,dy: real);
(* draw a line the relative amount (dx, dy). *)
begin
   drawr(afile,dx,dy,'l',0.0);
end;
(* end module pic.liner version = 2.63; (@ of dops.p 1991 November 2 *)

(* begin module pic.drawa *)
procedure drawa(var afile: text; x,y: real; visibility: char;
               spacing: real);
(* make a line to file afile to absolute coordinate x,y with visibility
  i invisible
  - dashed
  . dotted
  l line
with the dashes or dots separated by the spacing given
(this has no effect with invisible and line). *)
var
   dx, dy: real; (* differences between current and desired
                    locations *)
begin
   dx := x - picxglobal;
   dy := y - picyglobal;

   drawr(afile,dx,dy,visibility,spacing)
end;
(* end module pic.drawa version = 2.63; (@ of dops.p 1991 November 2 *)

(* begin module pic.movea *)
procedure movea(var afile: text; x,y: real);
(* move to absolute x and y *)
begin
   drawa(afile,x,y,'i',0.0);
end;
(* end module pic.movea version = 2.63; (@ of dops.p 1991 November 2 *)

(* begin module pic.linea *)
procedure linea(var afile: text; x,y: real);
(* draw a line from current position to absolute x and y *)
begin
   drawa(afile,x,y,'l',0.0);
end;
(* end module pic.linea version = 2.63; (@ of dops.p 1991 November 2 *)

(* begin module pic.graphstring *)
procedure graphstring(var tofile: text; var s: string; justification: char);
(* graph the string s.  If it is recognized as a quoted string (surrounded
by double quotes), graph it without the quotes and center it.
Otherwise justify it based on the justification character:
'l' left, 'c' centered, 'r' right.
For right and centered justification, the drawing point is the same
as before the string was done.  For left justification it is at the
right of the string to allow more to be added on there.
If not in picture (global variable inpicture), there is no output.
2008 Aug 05.  '(' and ')' need to be protected by putting '\' in
front of them.
2010 Nov 02 ONLY IF THERE IS NO PROCEEDING '\'!
*)
(* NONSTANDARD: PostScript dependent code.  Since different fonts
have different sized characters, one must rely on the PostScript
to handle the justification of the string. *)
var i: integer; (* index to s, and temporary storage *)
    quoted: boolean; (* true if the string is quoted *)
    skipping: boolean; (* true if skipping leading blanks *)
procedure postscriptprotect(var afile: text; c: char);
(* 2008 Aug 05: protect parenthesis *)
begin
{
2010 Nov 02 DO NOT DO THIS!!!  It messes up beta
   if (c = '(') or
      (c = ')')
   then write(afile,'!\');
Instead, protect only if \ is given.
}
end;
begin
   if (inpicture and (s.length > 0))
   then with s do begin
      if length > 2
      then if (letters[1]='"') and (letters[length]='"')
           then quoted := true
           else quoted := false
      else quoted := false;

      (* override so quoted strings are always centered *)
      if quoted then justification := 'c';

      (* do the non-standard postscript: *)
      if justification <> 'l' then write(tofile,'gsave ');

      (* do postscript to complete previous path *)
      (* set current point: writeln(tofile,'stroke 0 0 moveto'); *)
      writeln(tofile,'c');

      if justification = 'c' then begin
         (* when centering, skip leading blanks *)
         if letters[1] = ' ' then skipping := true
                             else skipping := false;
      end
      else skipping := false;

      write(tofile,'('); (* begin postscript literal *)
      if quoted (* take it literally *)
      then for i := 2 to length-1 do begin
               if letters[i-1] <> '\'
               then postscriptprotect(tofile,letters[i]);
               write(tofile,letters[i]);
           end
      else for i := 1 to length
           do if skipping then begin (* skip leading blanks *)
                 if letters[i] <> ' '
                 then begin
                      skipping := false;
                      if i = 1 then postscriptprotect(tofile,letters[i])
                               else if letters[i-1] <> '\'
                                    then postscriptprotect(tofile,letters[i]);
                      write(tofile,letters[i]);
                 end
              (* else skip the blank by not writing it *)
              end
              else begin
                  if i = 1 then postscriptprotect(tofile,letters[i])
                           else if letters[i-1] <> '\'
                                then postscriptprotect(tofile,letters[i]);
                  write(tofile,letters[i]);
              end;
      write(tofile,')'); (* end postscript literal *)

      if justification = 'c' (* center the string *)
      then write(tofile,' dup stringwidth pop neg 2 div 0 rmoveto')
      else if justification = 'r' (* rigth justify the string *)
      then write(tofile,' dup stringwidth pop neg 0 rmoveto');

      writeln(tofile,' x'); (* show the literal *)
      inpath := false; (* force new path from here *)

      if justification <> 'l' then write(tofile,'grestore ');
   end
(* There is no output if not in picture
   else begin
      writestring(tofile,s);
      writeln(tofile)
   end
 *)
end;
(* end module pic.graphstring *)

(* begin module pic.stringinteger *)
procedure stringinteger(number: integer; var name: string;
                        width: integer; leadingzeros: boolean);
(* make the string from the number, start putting characters in
after the current length point. use width characters.
if leadingzeros is true, trail zeros before the number. *)
var
   bigdigit: integer; (* the location of the biggest digit *)
   dig: integer; (* number of digits in the number *)
   place: integer; (* place to write the next digit of the number *)
   sign: integer; (* the sign of the number *)
begin
   with name do begin
      if number < 0
      then begin
         sign := -1;
         length := length + 1; (* provide room for the sign!! *)
         number := -number;
         if leadingzeros then begin
            writeln(output,'WARNING: stringinteger: the sign of a negative',
                           ' number with leading zeros is lost');
         end
      end
      else sign := +1;

      (* log 10 of the number plus 1 is the number of digits in the number.
      On this sun computer ln(1000)/ln(10) is 2.9999, which when
      truncated gives 2, rather than the desired 3.  To avoid this
      kind of problem, 0.1 is added. *)
      if number > 9
      then dig := trunc(ln(number+0.1)/ln(10))+1
      else dig := 1;

      if dig > width then begin
         writeln(output,'stringinteger: number width too small');
         writeln(output,dig:1,' digit number (',number:1,')');
         writeln(output,'does not fit in ',width:1,' characters');
         halt
      end;
      if leadingzeros
      then bigdigit := length + 1 (* no sign if leading zeros *)
      else begin
         bigdigit := length + width - dig + 1;

         if (bigdigit <= length) and (sign < 0) then begin
            writeln(output,'stringinteger: no room for sign');
            halt
         end;
      end;
      if sign < 0 then letters[bigdigit-1] := '-';

      for place := length + width downto bigdigit do begin
         case (number mod 10) of
            0: letters[place] := '0';
            1: letters[place] := '1';
            2: letters[place] := '2';
            3: letters[place] := '3';
            4: letters[place] := '4';
            5: letters[place] := '5';
            6: letters[place] := '6';
            7: letters[place] := '7';
            8: letters[place] := '8';
            9: letters[place] := '9';
         end;
         number := number div 10;
      end;
      length := length + width;
   end
end;
(* end module pic.stringinteger version = 2.63; (@ of dops.p 1991 November 2 *)

(* begin module pic.stringreal *)
procedure stringreal(number: real; var name: string;
                     width, decimal: integer);
(* make the string from the real number, start putting characters in
at the start point. use width characters and decimal characters
after the decimal place *)
   (* note that the rounding operation to get the digits below zero
      must be done first.  then the digits above zero can be lopped off.
      this makes 99.99 come out correctly to 100.0 (to 1 decimal place)
      otherwise, 99.99 -> 0.99 -> 1.0 (rounded) -> 10 (print with 1 decimal
      place), and stringinteger won't be happy about that.
    2003 Aug 28.  corrected missing minus sign for -1 < number <= 0.
*)
var
   abovezero: integer; (* the number shifted above the decimal place, to
      'decimal' positions (and rounded) *)
   shift: integer; (* power of ten used to shift a number around
      relative to the decimal point *)
   sign: integer; (* the sign of the number *)
   thedecimal: integer; (* integer version of the decimal part of the number *)
   theupper: integer; (* integer version of the upper part of the number *)
   signspot: integer; (* of the spot the sign will go. *)
begin
   (* sanity check: *)
   if name.length + width > maxstring then begin
      writeln(output,'real number =',number:1,
                     ' would exceed maxstring = ',maxstring:1);
      halt;
   end;
   if number < 0 then sign := -1
                 else sign := +1;

   number := abs(number); (* make positive *)

   (* the amount to shift the number above zero *)
   shift := round(exp(decimal*ln(10))); (* amount to move above zero *)
   abovezero := round(number*shift); (* move above zero, round off *)
   theupper := trunc(abovezero/shift);
   thedecimal := abovezero - shift*theupper;
{
writeln(output,'            stringreal:  number  = ',number:pwid:pdec);
writeln(output,'            stringreal:    sign  = ',sign:pwid);
writeln(output,'            stringreal: theupper = ',theupper:pwid);
}

   (* create the actual real number *)
   (* before decimal point *)
   (* provide a space for the sign in the resulting string: *)

   with name do begin (* put in the decimal point *)
      (* force a space for the sign by making the number negative *)
      signspot := length + 1; (* take note of the spot the sign will go. *)
      stringinteger(sign*theupper,name,width-decimal-1,false);

      (* 2003 Aug 28
         There is a very special case, known as bug1992.
         when the (number > -1) and (number < 0) the upper part
         of the number is zero (theupper = 0) BUT as an integer
         the sign cannot be passed to stringinteger, since
         -0 is of course 0 (usually, or sometimes).
         SO we have to handle that case and put a minus sign
         in 'by hand' here. *)
      if (sign < 1) and (theupper = 0)
         (* if number is tic = -0.000000000000000055511151231257827021181583405
            (a real example!!) then we would get -0.0 on rounding.
            SO round to the number of decimal places.
            The number of decimal places:
                  10^decimal = exp(ln(10^decimal)) = exp(decimal*ln(10))
         *)
         and (round(exp(decimal*ln(10))*number) <> 0)
      then begin
         (* ok, starting at signspot, move to the right until
            we are snug up against the number *)
         while letters[signspot+1]=' ' do signspot := succ(signspot);
         letters[signspot]:='-';
      end;

{
write(output,
   'stringinteger(',sign*theupper:1,',"');
writestring(output, name);
write(output, '",', width-decimal-1:1,',',false);
writeln(output,')');
}

      (* put in the decimal point *)
      length := length + 1;
      letters[length] := '.';
   end;
   stringinteger(thedecimal,name,decimal,true); (* after decimal point *)
end;
(* end module pic.stringreal clean.stringreal.p 2003 Aug 28 *)

(* begin module pic.picnumber *)
procedure picnumber(var afile: text;
                    dx, dy, number: real; width, decimal: integer;
                    justification: char);
(* Supply graphic commands for a 'number' whose center is at the relative point
(dx, dy) from the current point, 'width' characters wide and 'decimal'
characters beyond the decimal point.
If the width is zero, no number is produced.
procedure stringnumber(number: integer; start: integer; var name: string);
the location after the call is the same as before the call.
The string is optionally justified: left, centered or right: lcr. *)
var
   name: string; (* the string to pack the number into for shipping out *)
begin
   if width > 0 then begin
      mover(afile,dx,dy);

      clearstring(name);

      if decimal>0
      then stringreal(number,name,width,decimal)
      else stringinteger(round(number),name,width,false);

      graphstring(afile, name, justification);
      mover(afile,-dx,-dy);
   end
end;
(* end module pic.picnumber version = 2.63; (@ of dops.p 1991 November 2 *)

(* begin module pic.xtic *)
procedure xtic(var afile: text; length, dx, dy, number: real;
               width, decimal: integer;
               logxnormal: boolean;
               logxbase: real);
(* produce a tic mark for the x axis of "length" long.
Supply a number whose center is at the relative point (dx, dy)
from the end to the tick, 'width' characters wide and 'decimal'
characters beyond the decimal point.
If the width is zero, no number is produced.
the location after the call is the same as before the call.
If logxnormal is true, then raise the number
to logxbase. *)
begin
   liner(afile,0.0,-length);
   if logxnormal
   then picnumber(afile,dx,dy,exp(number*logxbase),width,decimal,'c')
   else picnumber(afile,dx,dy,number,width,decimal,'c');
   mover(afile,0.0,length);
end;
(* end module pic.xtic version = 2.63; (@ of dops.p 1991 November 2 *)

(* begin module pic.ytic *)
procedure ytic(var afile: text; length, dx, dy: real;
               number: real;
               width, decimal: integer;
               logynormal: boolean;
               logybase: real);
(* produce a tic mark for the y axis of "length" long.
Supply a number whose right side is started at the relative point (dx, dy)
from the end to the tick, 'width' characters wide and 'decimal'
characters beyond the decimal point.
If the width is zero, no number is produced.
the location after the call is the same as before the call.
If logynormal is true, then raise the number
to logybase. *)
begin
   liner(afile,-length,0.0);

   (* convert the number if we are doing logynormal: *)
   if logynormal then number := exp(number*logybase);

   picnumber(afile,dx,dy,number,width,decimal,'r');
   mover(afile,length,0.0);
end;
(* end module pic.ytic version = 2.63; (@ of dops.p 1991 November 2 *)

(* begin module pic.doaxis *)
procedure doaxis(var afile: text;
         theaxis: char;
         doaxisline: boolean;
         alength,fromtic,interval,totic: real;
         subintervals: real;
         length, dx, dy: real;
         width, decimal: integer;
         logscale, lognormal: boolean;
         logbase: real);
(* draw an axis starting from the current position.
Which axis it is is defined by theaxis, 'x' (horizontal) or 'y' (vertical).
Combining the code for both axes into one procedure is a little
slower, but drawing the axis does note ever take significant time,
and this allows improvements to be made on both axes simultaneously.
The length of the axis is alength.
If doaxisline is true then the axis line is drawn.
The axis is labeled with numbers starting with fromtic
at intervals given up to totic.
The remaining variables describe the form of the tic marks as in ytic.
If the width is zero, no number is produced.
the location after the call is the same as before the call.
If logscale and lognormal is true, then raise the tic numbers to logbase.
 *)
var
   half: real; (* half the jump interval.  By adding this to the while loops,
                  we assure that the very last tic gets done, and isn't lost
                  due to roundoff *)
   jump: real; (* the space to move on the graph between tic marks *)
   jumpdistance: real; (* the total jumps made.  this may not be
      a simple function of the input variables since they may
      not work out to an exact number of jumps *)
   tic: real; (* the numerical value of the tic label *)

   dosubtics: boolean; (* do sub tics *)
   subtic: real; (* the numerical value of the (unlabeled) subtic *)
   subinterval: real; (* the numerical interval between subtics *)
   subjump: real; (* the space to move on the graph between subtic marks *)
   halfsubinterval: real; (* half a subjump, see half *)

   currentspot: real; (* current graphing spot *)
   oldspot: real; (* previous graphing spot *)
   axisscale: real; (* axis scaling factor *)

begin
   writeln(afile,'gsave');

{
writeln(output,'In doaxis');
writeln(output,'interval=',interval:10:4);
writeln(output,'subintervals=',subintervals:10:4);
writeln(output,'logbase=',logbase:10:4);
}

   if theaxis = 'x' then begin
      if doaxisline then begin
         liner(afile,+alength,0.0);
         mover(afile,-alength,0.0);
      end;
   end
   else begin
      if doaxisline then begin
         liner(afile,0.0,+alength);
         mover(afile,0.0,-alength);
      end;
   end;

   if totic = fromtic then begin
      writeln(output,'doaxis: ',theaxis,' axis fromtic and totic',
                     ' cannot be equal');
      halt;
   end;

   if (alength = 0.0) or (interval = 0.0) then begin
      writeln(output,'doaxis: neither ',
                     theaxis,' axis length nor interval can be zero');

      halt;
   end;

   axisscale := alength / (totic - fromtic);
   jump := axisscale * interval;
   jumpdistance := 0;
   half := interval / 2.0;

   if subintervals > 1 then begin
      dosubtics := true;
      subinterval := interval/subintervals;
      halfsubinterval := subinterval / 2.0;
      subjump := jump/subintervals;
   end
   else begin
      dosubtics := false;
      subinterval := 0;
      halfsubinterval := 0;
      subjump := 0;
   end;
{
writeln(output,'fromtic  = ',fromtic:10:4);
writeln(output,'totic    = ',totic:10:4);
writeln(output,'interval = ',interval:10:4);
writeln(output,'half     = ',half:10:4);
}

   tic := fromtic;
   if interval > 0.0 then while tic <= totic+interval do begin
{
writeln(output,'* tic=',tic:10:4);
}
      (* adding interval/2 makes sure that the y axis does include
         the last point but that the x axis doesn't get an extra
         from (presumably) rounding effects *)
      if tic <= totic+(interval/2) then begin
         if theaxis = 'x'
         then xtic(afile,length,dx,dy,tic,width,decimal,lognormal,logbase)
         else ytic(afile,length,dx,dy,tic,width,decimal,lognormal,logbase);
      end;

      (* 2007 Aug 30 the extra interval makes subtics go to the end
         of the graph rather than ending at the last tic mark *)
      if tic <= totic + interval then begin
{
writeln(output,'totic+half + interval=',totic+half + interval:10:4);
writeln(output,'totic+interval=',totic+interval:10:4);
      if tic <= totic+half + 1 then begin
writeln(output,'TIC=',tic:10:4);
writeln(afile,'% tic=',tic:10:4);
mover(afile,0.05,0.0);
}
         if dosubtics then begin (* do subtic marks *)
            if logscale then begin (* do subtic marks on log scale *)
               (* subtic starts as a "normal" number (ie, no log taken) at tic: *)
{
writeln(output,'2^tic=',exp(tic*logbase):10:4);
writeln(output,'2^(tic+interval)=',exp((tic+interval)*logbase):10:4);
}
               subtic := exp(tic*logbase);
               (* subtic will proceed to the same but at tic+interval.
                  We divide that into the subintervals. *)
{
writeln(output,'halfsubinterval=',halfsubinterval:10:4,' original');
}
               subinterval := (exp((tic+interval)*logbase) - subtic)/subintervals;
               halfsubinterval := subinterval/2.0;

{
writeln(output,'subtic=         ',subtic:10:4);
writeln(output,'subinterval=    ',subinterval:10:4);
writeln(output,'halfsubinterval=',halfsubinterval:10:4);
}

               oldspot := axisscale * tic;
               while subtic < exp(logbase*(tic+interval)) - halfsubinterval
               do begin

                  (* although tic is on a log scale,
                     we have to have subtic on the regular scale
                     to alter the positions of the subtics *)

                  (* if subinterval is constant,
                     the following makes linearly spaced marks: *)
                  subtic := subtic + subinterval;

                  (* the actual jumps have to be in the log form: *)
                  currentspot := axisscale*ln(subtic)/logbase;
                  subjump := currentspot - oldspot;
{
writeln(output,'       SUBTIC=',subtic:10:4);
writeln(output,'       ln(SUBTIC)/logbase=',ln(subtic)/logbase:10:4);
writeln(output,'       currentspot=',currentspot:10:4);
writeln(output,'       subjump=',subjump:10:4);
writeln(output,'       oldspot=',oldspot:10:4);
writeln(afile,'% subtic=',subtic:10:4);
}
                  oldspot := currentspot;


                  if theaxis = 'x' then begin
                     xtic(afile,length/2,dx,dy,0,0,0,lognormal,logbase);
                     mover(afile,subjump,0.0);
                  end
                  else begin
                     ytic(afile,length/2,dx,dy,0,0,0,lognormal,logbase);
                     mover(afile,0.0,subjump);
                  end;

                  jumpdistance := jumpdistance + subjump;
               end
            end
            else begin (* do subtic marks on regular scale *)
               subtic := tic;
               while subtic < tic+interval-halfsubinterval do begin
                  subtic := subtic + subinterval;

                  if theaxis = 'x' then begin
                     mover(afile,subjump,0.0);
                     if subtic <= totic
                     then xtic(afile,length/2,dx,dy,0,0,0,lognormal,logbase);
                  end
                  else begin
                     mover(afile,0.0,subjump);
                     if subtic <= totic
                     then ytic(afile,length/2,dx,dy,0,0,0,lognormal,logbase);
                  end;

                  jumpdistance := jumpdistance + subjump;
               end
            end
         end
         else begin (* do regular tic marks *)
            if theaxis = 'x' then mover(afile,jump,0.0)
                             else mover(afile,0.0,jump);
            jumpdistance := jumpdistance + jump
         end
      end;
      tic := tic + interval;
   end
   else if interval < 0.0 then while tic >= totic-half do begin
      if dosubtics then writeln(output,
                        'Sorry, no subtics with negative scales');

      if theaxis = 'x'
      then xtic(afile,length,dx,dy,tic,width,decimal,lognormal,logbase)
      else ytic(afile,length,dx,dy,tic,width,decimal,lognormal,logbase);

      tic := tic + interval;
      if tic >= totic-half then begin
         if theaxis = 'x' then mover(afile,jump,0.0)
                          else mover(afile,0.0,jump);
         jumpdistance := jumpdistance + jump
      end
   end;

   if theaxis = 'x' then mover(afile,-jumpdistance,0.0)
                    else mover(afile,0.0,-jumpdistance);
   writeln(afile,'grestore');
end;
(* end module pic.doaxis version = 2.75; (@ of dops.p 2007 Aug 30 *)

(* begin module pic.xaxis *)
procedure xaxis(var afile: text;
         doaxisline: boolean;  (* line on axis is plotted *)
         axlength,fromtic,interval,totic: real;
         xsubintervals: real;
         length, dx, dy: real;
         width, decimal: integer;
         logxscale, logxnormal: boolean;
         logxbase: real);
(* draw an x axis starting from the current position. *)
begin
   doaxis(afile,
         'x',
         doaxisline,
         axlength,fromtic,interval,totic,
         xsubintervals,
         length, dx, dy,
         width, decimal,
         logxscale, logxnormal,
         logxbase)
end;
(* end module pic.xaxis version = 2.70; (@ of dops.p 1996 March 21 *)

(* begin module pic.yaxis *)
procedure yaxis(var afile: text;
         doaxisline: boolean;  (* line on axis is plotted *)
         aylength,fromtic,interval,totic: real;
         ysubintervals: real;
         length, dx, dy: real;
         width, decimal: integer;
         logyscale, logynormal: boolean;
         logybase: real);
(* draw an y axis starting from the current position. *)
begin
   doaxis(afile,
         'y',
         doaxisline,
         aylength,fromtic,interval,totic,
         ysubintervals,
         length, dx, dy,
         width, decimal,
         logyscale, logynormal,
         logybase)
end;
(* end module pic.yaxis version = 2.70; (@ of dops.p 1996 March 21 *)

(* ********************************************************************** *)

(* begin module pic.dotr *)
procedure dotr(var afile: text);
(* draw a dot at the current position *)
begin
   (* a zero length line won't work here, since some systems
      really believe it and draw nothing.  (eg. PostScript on the
      Applewriter draws nothing, but PostScript on the Sun puts one pixle *)
   drawr(afile, dotfactor,0.0,'l',0.0);
end;
(* end module pic.dotr version = 9.22; (@ of xyplo.p 2010 Nov 02 *)

(* begin module pic.boxr *)
procedure boxr(var afile: text; width, height: real);
(* make a box to file afile with width in the x direction
and height in the y direction as given.
the box goes toward the positive x and y directions.
the box is relative to the current position, so it
returns to original position afterwards *)
begin
   liner(afile,0.0,height);
   liner(afile,width,0.0);
   liner(afile,0.0,-height);
   liner(afile,-width,0.0)
end;
(* end module pic.boxr version = 4.80; (@ of piclib 1985 dec 26 *)

(* begin module pic.cboxr *)
procedure cboxr(var afile: text; width, height: real);
(* make a box to file afile with width in the x direction
and height in the y direction as given.
the box is centered at the current position.
the box is relative to the current position, so it
returns to original position afterwards *)
var h2,w2: real; (* height and width over 2 *)
begin
   h2 := height/2;
   w2 := width/2;
   mover(afile,-w2,-h2);
   liner(afile,0.0,height);
   liner(afile,width,0.0);
   liner(afile,0.0,-height);
   liner(afile,-width,0.0);
   mover(afile,w2,h2);
end;
(* end module pic.cboxr version = 9.22; (@ of xyplo.p 2010 Nov 02 *)

(* begin module pic.polrec *)
procedure polrec(r,theta: real; var x,y: real);
(* convert polar to rectangular coordinates,
theta is in radians *)
begin
   x := r*cos(theta);
   y := r*sin(theta)
end;
(* end module pic.polrec version = 9.22; (@ of xyplo.p 2010 Nov 02 *)

(* begin module pic.degtorad *)
function degtorad(angle: real):real;
(* convert angle in degrees to radians *)
begin
   degtorad := (angle / 360) * 2 * pi
end;
(* end module pic.degtorad version = 9.22; (@ of xyplo.p 2010 Nov 02 *)

(* begin module pic.spiral *)
procedure spiral(var afile: text; thickness: real; steps: integer;
                 radius: real);
(* make a spiral into file afile, at the current position,
with a certain thickness and using a certain number of steps at
whose largest radius is 'radius'.  return to same position afterward.
If steps has a negative value, then the spiral is drawn clockwise,
otherwise it is drawn counterclockwise. *)
var
   dr: real; (* change in r *)
   dtheta: real; (* change in theta *)
   r: real; (* radius of the current position *)
   theta: real; (* angle of the current position *)
   x: real; (* the x coordinate *)
   xpos: real; (* to remember the center of the spiral *)
   y: real; (* the y coordinate *)
   ypos: real; (* to remember the center of the spiral *)
begin
   if steps <> 0 then begin (* avoid explosion *)
      xpos := picxglobal;
      ypos := picyglobal;
      r := 0;
      theta := 0;
      dr := abs(thickness/steps);
      dtheta := 2 * pi / steps;

      while r < abs(radius) do begin
         r := r + dr;
         theta := theta + dtheta;
         polrec(r,theta,x,y);
         linea(afile,x+xpos,y+ypos)
      end;
      movea(afile,xpos,ypos)
   end;
end;
(* end module pic.spiral version = 4.80; (@ of piclib 1985 dec 26 *)

(* begin module pic.movepolar *)
procedure movepolar(var afile: text; angle, distance: real);
(* move relative to the current position by placing
the appropriate pic commands into afile.
the angle is in degrees, the distance is in cm.*)
var
   dx: real; (* change in x *)
   dy: real; (* change in y *)
begin
   polrec(distance, degtorad(angle) ,dx,dy);
   mover(afile,dx,dy)
end;
(* end module pic.movepolar version = 4.80; (@ of piclib 1985 dec 26 *)

(* begin module pic.boxintercept *)
procedure boxintercept(xmin,ymin,xmax,ymax,m,b: real;
                       var intercept: boolean; var x1,y1,x2,y2: real);
(* does the line y=m*x+b intercept the box defined by the points
(xmin,ymin) and (xmax,ymax)?  if so, intercept is true and
the intercept points are given by (x1,y1) and (x2,y2) *)
var
   count: integer; (* number of sides the line intersects *)
      xlo,xhi,ylo,yhi: boolean; (* whether the line intersects the
         box at the low value of x, etc *)
function fny(x: real):real;
(* calculate the y value given the x *)
begin fny := m*x+b end;

function fnx(y: real):real;
(* calculate the x value given the y *)
begin fnx := (y-b)/m end;

function between(a,b,c: real):boolean;
(* is b between a and c? *)
begin between:=(a<=b) and (b<=c) end;

procedure normalcases;
(* analyze for the usual cases when the slope m is not zero *)
begin (* normalcases *)
(* writeln(output,'m=',m:5:4,' b=',b:5:4,' in normalcases');*)
(*  ymin, ymax, xmin and xmax are the coordinates of the box;
    xlo, yhi, xhi, and ylo are which sides are intersected

                             |
     ymax +----yhi----+      |
          |           |      |
         xlo         xhi     |
          |           |      |
     ymin +----ylo----+      |
        xmin        xmax     |
                             |

*)
      xlo := between(ymin,fny(xmin),ymax);
      xhi := between(ymin,fny(xmax),ymax);
      ylo := between(xmin,fnx(ymin),xmax);
      yhi := between(xmin,fnx(ymax),xmax);
(*
writeln(output,'xlo = ',xlo);
writeln(output,'xhi = ',xhi);
writeln(output,'ylo = ',ylo);
writeln(output,'yhi = ',yhi);
*)

      intercept := true; (* optimistic *)

      (* simplify cases which intersect corners.  These are the ones
         where more than two side intersections are true. *)
      count := 0;
      if xlo then count := succ(count);
      if xhi then count := succ(count);
      if ylo then count := succ(count);
      if yhi then count := succ(count);

      if count > 2 then begin
	 (* one corner must be intersected.  Simplify by preserving the
	    opposing intersections.  When there are 4 intersections, only one
	    pair is preserved.  The calculations will be correct either way. *)
         if xlo and xhi then begin
            yhi := false;
            ylo := false;
         end
         else if  ylo and yhi then begin
            xhi := false;
            xlo := false;
         end
         else begin
            writeln(output,'error in between count!');
            halt
         end
      end;

           if xlo and xhi then begin x1 := xmin;      x2 := xmax        end
      else if xlo and ylo then begin x1 := xmin;      x2 := fnx(ymin)   end
      else if xlo and yhi then begin x1 := xmin;      x2 := fnx(ymax)   end
      else if xhi and ylo then begin x1 := xmax;      x2 := fnx(ymin)   end
      else if xhi and yhi then begin x1 := xmax;      x2 := fnx(ymax)   end
      else if ylo and yhi then begin x1 := fnx(ymin); x2 := fnx(ymax)   end
      else intercept := false;

      if intercept then begin
         y1 := fny(x1);
         y2 := fny(x2)
      end

end; (* normalcases *)
begin (* boxintercept *)
      (* note: abs(m) is required to protect against negative zero... *)
      if abs(m) = 0.0
      then begin
         intercept := between(ymin,b,ymax);
         if intercept then begin
            x1 := xmin; y1 := b;
            x2 := xmax; y2 := b;
         end
      end
      else normalcases
end; (* boxintercept *)
(* end module pic.boxintercept version = 9.22; (@ of xyplo.p 2010 Nov 02 *)

(* begin module pic.plusr *)
procedure plusr(var afile: text; width, height: real);
(* make a plus sign to file afile with width in the x direction
and height in the y direction as given.
the box is centered at the current position.
the box is relative to the current position, so it
returns to original position afterwards *)
var h2,w2: real; (* height and width over 2 *)
begin
   h2 := height/2;
   w2 := width/2;
   mover(afile,-w2,0);
   liner(afile,width,0.0);
   mover(afile,-w2,h2);
   liner(afile,0.0,-height);
   mover(afile,0,h2);
end;
(* end module pic.plusr version = 9.22; (@ of xyplo.p 2010 Nov 02 *)

(* begin module pic.xr *)
procedure xr(var afile: text; width, height: real);
(* make an x to file afile with width in the x direction
and height in the y direction as given.
the box is centered at the current position.
the box is relative to the current position, so it
returns to original position afterwards *)
var h2,w2: real; (* height and width over 2 *)
begin
   h2 := height/2;
   w2 := width/2;
   mover(afile,-w2,-h2);
   liner(afile,width,height);
        mover(afile,0,-height);
   liner(afile,-width,height);
   mover(afile,w2,-h2);
end;
(* end module pic.xr version = 9.22; (@ of xyplo.p 2010 Nov 02 *)

(* begin module pic.arc *)
procedure arc(var thefile: text; angle1, angle2, radius: real;
              steps: integer);
(* create an arc in thefile going from angle1 to angle2 (degrees) in the
positive direction of angle, with the given radius.
use the given number of steps to make it.  return to the same position
as before the arc was drawn. *)
var
   dtheta: real; (* change in theta *)
(* s: integer; (@ index to the steps *)
   theta: real; (* current angle *)
   x,y: real; (* coordinates around starting point *)
(* zerox,zeroy: real; (@ starting location, center of curve *)
begin
(*
   zerox := picxglobal;
   zeroy := picyglobal;
*)
   theta := degtorad(angle1);
   dtheta := degtorad( abs(angle2-angle1)/steps );
   polrec(radius,theta,x,y);
(* can't do this for postscript arc:
   movea(thefile,zerox+x,zeroy+y);
      ' ',scale*zerox:  picwidth:picdecim,
      ' ',scale*zeroy:  picwidth:picdecim,
*)

   (* NONSTANDARD postscript, much faster *)
   writeln(thefile,
   (* 'stroke newpath', *)(* force there to be no current point *)
   (* ' 0 0', *)
      'a',
      ' ',scale*radius: picwidth:picdecim,
      ' ',angle1:       picwidth:picdecim,
      ' ',angle2:       picwidth:picdecim);
   write(thefile, 'arc');
   if angle2 < angle1 then write(thefile,'n'); (* for negative draws *)
   (* origin move: writeln(thefile,' stroke newpath 0 0 moveto'); *)
   writeln(thefile,' n');
   (* the moveto puts us back to the origin *)

(*
   for s := 1 to steps do begin
      theta := theta + dtheta;
      polrec(radius,theta, x,y);
      linea(thefile,zerox+x,zeroy+y);
   end;

   movea(thefile,zerox,zeroy)
*)
end;
(* end module pic.arc version = 9.22; (@ of xyplo.p 2010 Nov 02 *)

(* begin module pic.circler *)
procedure circler(var afile: text; radius: real);
(* make a circle at the current position of some radius. *)
var
   steps: integer; (* number of steps to make the circle *)
begin
   (* number of segments increases with diameter, but the constant
      still should be a function of how good it looks on a particular
      graphic system, I'm afraid.  However, there should be a lower
      bound on the number of steps, so even small circles look good *)
   if radius < 1.0 then steps := 25
                   else steps := round(radius*25);

   arc(afile,0.0,360.0,radius,steps);
end;
(* end module pic.circler version = 9.22; (@ of xyplo.p 2010 Nov 02 *)

(* begin module pic.ibeam *)
procedure ibeam(var afile: text; width, height: real);
(* Make an ibeam shaped symbol to file afile with width in the x direction
and height in the y direction.  Center it at the current position.
Put a circle at the center, with radius 1/4th the width
Return to original position afterwards. *)
var h2,w2: real; (* height and width over 2 *)
    r: real; (* the radius of the circle *)
begin
   h2 := height/2;
   w2 := width/2;
   mover(afile,-w2,-h2);
   liner(afile,width,0.0);
   mover(afile,-width,height);
   liner(afile,width,0.0);
   mover(afile,-w2,0.0);
   liner(afile,0.0,-height);
   mover(afile,0.0,h2);
   r := width/8;
{ this is silly for the new printers
   if r < 0.025 then r := 0.025; (* small circles do not come out well *)
}
   circler(afile,r);
end;
(* end module pic.ibeam version = 9.22; (@ of xyplo.p 2010 Nov 02 *)

(* begin module pic.setgray *)
procedure setgray(var afile: text;
                      brightness: real);
(* set the gray scale to the requested one.
Range of the variables is 0 to 1. *)
const
   colfield = 8; (* width of numbers printed to the file *)
   colwidth = 4; (* number of decimal places for numbers *)
   (* PostScript on a Sun 4 cannot handle 5 decimal places!  Use 4 or less *)
begin
   write  (afile,'n'); (* be sure it's started cleanly *)
   write  (afile,' ',brightness:colfield:colwidth);
   writeln(afile,' setgray');
end;
(* end module pic.setgray version = 9.22; (@ of xyplo.p 2010 Nov 02 *)

(* begin module pic.setcolor *)
procedure setcolor(var afile: text;
                       hue, saturation, brightness: real);
(* Set the color to the requested one.  Range of the variables is 0 to 1.
   Colors in PostScript are defined with hue, saturation and brightness with
   the sethsbcolor function.  The standard hue runs from red at 0 to red
   at 1 with the Roy G. Biv in between.  The famous physicist Roy G. Biv
   stands for: red, orange, yellow, green, blue, indigo, violet, the colors
   of the spectrum.  On 1992 September 16 I realized that, amazingly, the
   wavelength in nm match quite nicely!

      700 red
      750 orange
      600 yellow
      650 green
      500 blue
      550 indego
      400 violet

   Since the hue runs in a circle, the spectrum is not generated from the hue
   range 0 to 1.  The routine setcolor therefore converts the input numbers
   by multiplying by 0.85.  This gives the color range from red through
   violet corresponding to values 0 to 1.  (Note:  adding 0.15 would give the
   range from yellow through red, but that is not a spectrum.)
*)
const
   colfield = 8; (* width of numbers printed to the file *)
   colwidth = 4; (* number of decimal places for numbers *)
   (* PostScript on a Sun 4 cannot handle 5 decimal places!  Use 4 or less *)
   huefactor = 0.85; (* number of decimal places for numbers *)
begin
   write  (afile,'n'); (* be sure it's started cleanly *)
   write  (afile,' ',hue*huefactor:       colfield:colwidth);
   write  (afile,' ',saturation:colfield:colwidth);
   write  (afile,' ',brightness:colfield:colwidth);
   writeln(afile,' sethsbcolor');
end;
(* end module pic.setcolor *)

(* begin module pic.rectinit *)
procedure rectinit(var outfile: text);
(* create the definition of a rectangle.  Rectsize must be called
to initialize and to change the size of the rectangle. *)
begin
   writeln(outfile,'/rct');
   writeln(outfile,'{gsave');
   writeln(outfile,'  newpath');
   writeln(outfile,'  0 0 moveto');
   writeln(outfile,'  xs 0 lineto');  (* xs is the x side length *)
   writeln(outfile,'  xs ys lineto'); (* ys is the x side length *)
   writeln(outfile,'  0 ys lineto');
   writeln(outfile,'  closepath fill');
   writeln(outfile,'grestore} def');
end;
(* end module pic.rectinit version = 9.22; (@ of xyplo.p 2010 Nov 02 *)

(* begin module pic.rectdo *)
procedure rectdo(var afile: text);
(* Make a rectangle with lower right hand corner at (x,y) and the
given side (position and color are determined by earlier calls
to setcolor and move).  Return to original position afterwards. *)
begin
   writeln(afile, ' rct');
end;
(* end module pic.rectdo version = 9.22; (@ of xyplo.p 2010 Nov 02 *)

(* begin module pic.rectsize *)
procedure rectsize(var afile: text;
                   var xsideold, ysideold, xside, yside: real);
(* determine if the values of xside and yside have changed from xsideold and
ysideold.  If either has changed, write the instructions to change the size of
the rectangle into afile. *)
begin
   if xside <> xsideold then begin
      writeln(afile,'/xs ',
                    (xside*scale):picwidth:picdecim,
                    ' def'); (* xs is the x side length used in rectinit *)
      xsideold := xside;
   end;

   if yside <> ysideold then begin
      writeln(afile,'/ys ',
                    (yside*scale):picwidth:picdecim,
                    ' def'); (* ys is the y side length used in rectinit *)
      ysideold := yside;
   end
end;
(* end module pic.rectsize version = 9.22; (@ of xyplo.p 2010 Nov 02 *)

(* ********************************************************************** *)
(* ********************************************************************** *)
(* ********************************************************************** *)

(* begin module pic.3d.determinant *)
function determinant(a: tbtarray): real;
(* compute the determinant of a *)
begin
   determinant := +a[1,1] * (a[2,2]*a[3,3] - a[3,2]*a[2,3])
                  -a[1,2] * (a[2,1]*a[3,3] - a[3,1]*a[2,3])
                  +a[1,3] * (a[2,1]*a[3,2] - a[3,1]*a[2,2])
end;
(* end module pic.3d.determinant *)

(* begin module pic.3d.d32 *)
procedure d32(o, a, b, c, v: threevector;  var xloc,yloc: real);
(* convert from 3d to 2d.  the players are:
o: the coordinate of the object point to be converted to 2d
a,b,c: define the position of the window (screen):
a: center of screen
b: screen x coordinate direction
c: screen y coordinate direction
v: the position of the viewer
xloc,yloc: the resulting image vector in screen coordinates.
   The method of graphics is to project the object (o) toward the viewer
(v) and to determine the interception of this line with the screen
as defined by a,b and c.  the result is expressed in the coordinate system
of the screen, and so can be plotted on a 2d plotting device.
   When one works through the vector math, it turns out that to find
the screen coordinates requires solving a set of linear equations.
This is done using Cramer's rule and determinants. *)
var
   ov,oa: real; (* for partial calculation *)
   j: integer; (* index to the arrays *)
   d,x,y: tbtarray;
begin
   (* define the coefficients of the equations in d,x and y *)
   for j:=1 to 3 do begin
      ov := o[j]-v[j];
      d[j,1]:=b[j];
      d[j,2]:=c[j];
      d[j,3]:=ov;

      oa:=o[j]-a[j];
      x[j,1]:=oa;
      x[j,2]:=c[j];
      x[j,3]:=ov;

      y[j,1]:=b[j];
      y[j,2]:=oa;
      y[j,3]:=ov;
   end;

   (* use cramer's rule to find the solution *)
   xloc:=determinant(x)/determinant(d);
   yloc:=determinant(y)/determinant(d);
end;
(* end module pic.3d.d32 *)

(* begin module pic.3d.view *)
procedure view(v: threevector; var gaze: threevector; smag: real;
               var a,b,c: threevector);
(* this routine converts a viewing position (v) and a viewing
direction (gaze), into the a,b,c values of a vertically oriented screen
(ie, the screen is right side up).  a is the center of the screen,
b is the x axis, c is the y axis on the screen.  This saves the user
the trouble to make sure that b, c and the direction of viewing are
orthogonal.

one may magnify the view by making smag greater than one, or one may
shrink the view by making smag less than one.

if the viewing direction vector is not large enough,
then the program halts.

note: gaze is automatically converted to a unit vector. *)
var
   db: real; (* magnitude of db *)
   dgaze: real; (* magnitude of gaze *)
   j: integer; (* index to the arrays *)
begin
   (* first check out the gaze direction *)
   dgaze := sqrt(gaze[1]*gaze[1] + gaze[2]*gaze[2] + gaze[3]*gaze[3]);
   if smag = 0.0 then begin
      writeln(output,'screen magnitude cannot be zero');
      halt
   end;
   if dgaze <= 0.001 then begin
      writeln(output,'gaze magnitude (',dgaze:5:3,') is too small');
      halt
   end;

   (* make gaze a unit vector and set up the a vector as the
      viewing point plus the gaze vector *)
   for j := 1 to 3 do begin
      gaze[j] := gaze[j]/dgaze;
      a[j] := v[j] + gaze[j]
   end;

   (* the x axis of the screen, the b vector, is horizontal and
      orthogonal to the gaze *)
   b[1] := +gaze[2];
   b[2] := -gaze[1];
   b[3] := 0;
   db := sqrt(b[1]*b[1] + b[2]*b[2] + b[3]*b[3]);

   (* check for top view case and correct if so: *)
   if db = 0.0 then begin
      db := 1;
      b[1] := 1;
      b[2] := 0;
   (* b[3] := 0; already from above *)
   end
   else for j := 1 to 3 do b[j] := b[j]/db; (* make b a unit vector *)

   (* now that the gaze is a unit vector, and we have constructed
      the x axis in the b vector also as a unit vector, the cross
      product of these two will generate the y axis as a unit
      vector, c: *)
   c[1] := +(b[2]*gaze[3] - gaze[2]*b[3]);
   c[2] := -(b[1]*gaze[3] - gaze[1]*b[3]);
   c[3] := +(b[1]*gaze[2] - gaze[1]*b[2]);

   (* now normalize both b and c vectors to be of size 1/smag *)
   for j := 1 to 3 do begin
      b[j] := b[j]/smag;
      c[j] := c[j]/smag;
   end
end;
(* end module pic.3d.view *)

(* begin module pic.3d.makescreen *)
procedure makescreen(vx,vy,vz, gx,gy,gz, smagnitude: real; var s: screen);
(* create the screen s based on the viewing location (vx,vy,vz)
and the direction of gaze (gz,gy,gz).  The screen size is scaled by
smagnitude; doubling smagnitude should double the size of the scene. *)

(* This routine makes creation of the screen very simple for the user.
One need not look at the view routine. *)
begin
   s.v[1] := vx;
   s.v[2] := vy;
   s.v[3] := vz;
   s.g[1] := gx;
   s.g[2] := gy;
   s.g[3] := gz;
   with s do view(v,g,smagnitude, a,b,c);
   s.smag := smagnitude;
   s.range := 1/smagnitude
end;
(* end module pic.3d.makescreen *)

(* begin module pic.3d.project3d *)
procedure project3d(x,y,z: real; s: screen; var xscreen,yscreen: real);
(* project the point (x,y,z) onto the screen s, to find the screen
coordinates (xscreen and yscreen). *)

(* This routine simplifies the projection function for the user. *)
var
   o: threevector; (* for passing the values to d32 *)
begin
   o[1] := x;
   o[2] := y;
   o[3] := z;
   with s do d32(o,a,b,c,v,xscreen,yscreen);
end;
(* end module pic.3d.project3d *)

(* begin module pic.3d.test.fun *)
function fun(r: real): real;
(* a function to plot *)
begin
   fun := 3/(1+r*r/2)
end;
(* end module pic.3d.test.fun *)

(* begin module pic.3d.test.test3d *)
procedure test3d(var afile: text);
(* test three dimensional graphics *)
var
   s: screen; (* the screen on which to project the 3d image *)

   xscreen, yscreen: real; (* location on the screen corresponding to the
      projection of o onto the screen defined by v,a,b,c *)
   xold,yold: real; (*  the previous valuse of xscreen and yscreen *)

   (* definition of a spiral *)
   dr: real; (* change in r *)
   dtheta: real; (* change in theta *)
   r: real; (* radius of the current position *)
   radius: real; (* the radius of the spiral *)
   theta: real; (* angle of the current position *)
   thickness: real; (* spacing between spiral arms *)
   steps: real; (* number of steps around a circle of the spiral *)

   x,y,z: real; (* the location in three space *)
begin

   makescreen(5.0,5.0,5.0, -1.0,-1.0,-1.0, 5.0, s);

   r := 0;
   theta := 0;
   steps := 15;
   thickness := 0.1;
   radius := 2.0;
   dr := thickness/steps;
   dtheta := 2 * pi / steps;

   x := 0;
   y := 0;
   z := fun(r);
   project3d(x,y,z, s, xold,yold);
   mover(afile,xold,yold); (* premove to the startpoint of the graph *)

   while r < radius do begin
      r := r + dr;
      theta := theta + dtheta;
      polrec(r,theta,x,y);

      z := fun(r);
      project3d(x,y,z, s, xscreen,yscreen);

      (* draw a line from where we where to the new place *)
      liner(afile, xscreen - xold, yscreen - yold);
      xold := xscreen;
      yold := yscreen;
   end;
end;
(* end module pic.3d.test.test3d *)

(* ********************************************************************** *)
(* ********************************************************************** *)
(* ********************************************************************** *)

(* begin module skipblanks *)

(* 2003 July 31: tab is considered a blank character *)

function isblank(c: char): boolean;
(* is the character c blank or tab? *)
const
   tab = 9; (* tab character *)
begin
   isblank := (c = ' ') or (ord(c) = tab)
end;

procedure skipblanks(var thefile: text);
(* skip over blanks until a non-blank, or end of line, is found *)
begin
   while isblank(thefile^) and not eoln(thefile) do get(thefile);
end;

procedure skipnonblanks(var thefile: text);
(* skip over nonblanks until a blank, or end of line, is found *)
begin
   while (not isblank(thefile^)) and not eoln(thefile) do get(thefile);
end;

procedure skipcolumn(var thefile: text);
(* skip over a data column *)
begin
   skipblanks(thefile); skipnonblanks(thefile)
end;
(* end module skipblanks version = 9.22; (@ of xyplo.p 2010 Nov 02 *)
(* ********************************************************************** *)

(* begin module dops.readchar *)
procedure readchar(var a: text; var c: char);
(* read from file a the character c by first
skipping preceding blanks and then skipping other non-blanks after *)
begin
      skipblanks(a);
      read(a,c);
      skipnonblanks(a)
end;
(* end module dops.readchar *)
(* ********************************************************************** *)

(* begin module dops.mkhalt *)
procedure mkhalt(var outfile: text);
(* generate the call to halt *)
begin
      write(outfile,'dops ');
      halt
end;
(* end module dops.mkhalt *)

(* begin module dops.testblank *)
procedure testblank(var infile,outfile: text);
(* Test for blank as the next character.  If it is not, terminate
the program.  If this is not done, reads may bomb on badly formed
input.  Example:  boxrz will bomb on the attempt to read the number
because it turns out to be a z *)
procedure die;
begin (* die *)
   writeln(outfile,'(badly formed graphics instruction) show');
   stoppic(outfile); (* close what we have *)
   mkhalt(outfile);
end; (* die *)
begin
      if eoln(infile)
      then die
      else if infile^<>' ' then die
end;
(* end module dops.testblank *)

(* begin module dops.mkstartpic *)
procedure mkstartpic(var infile,outfile: text);
(* generate the call to startpic *)
var scale, x, y: real; (* scale factor, and coordinate to start with *)
begin
      readln(infile,scale,x,y);
      startpic(outfile,scale,x,y,'t')
end;
(* end module dops.mkstartpic *)

(* begin module dops.mkstoppic *)
procedure mkstoppic(var infile,outfile: text);
(* generate the call to stoppic *)
begin
      readln(infile);
      stoppic(outfile)
end;
(* end module dops.mkstoppic *)

(* begin module dops.mkdrawr *)
procedure mkdrawr(var infile,outfile: text);
(* generate the call to drawr *)
var dx,dy: real; visibility: char; spacing: real;
begin
      testblank(infile,outfile);
      read(infile,dx,dy);
      readchar(infile,visibility);
      readln(infile,spacing);
      drawr(outfile,dx,dy,visibility,spacing)
end;
(* end module dops.mkdrawr *)

(* begin module dops.mkmover *)
procedure mkmover(var infile,outfile: text);
(* generate the call to mover *)
var dx,dy: real;
begin
      testblank(infile,outfile);
      readln(infile,dx,dy);
      mover(outfile,dx,dy)
end;
(* end module dops.mkmover *)

(* begin module dops.mkliner *)
procedure mkliner(var infile,outfile: text);
(* generate the call to liner *)
var dx,dy: real;
begin
      testblank(infile,outfile);
      readln(infile,dx,dy);
      liner(outfile,dx,dy)
end;
(* end module dops.mkliner *)

(* begin module dops.mkdrawa *)
procedure mkdrawa(var infile,outfile: text);
(* generate the call to drawa *)
var x,y: real; visibility: char; spacing: real;
begin
      testblank(infile,outfile);
      read(infile,x,y);
      readchar(infile,visibility);
      readln(infile,spacing);
      drawa(outfile,x,y,visibility,spacing)
end;
(* end module dops.mkdrawa *)

(* begin module dops.mkmovea *)
procedure mkmovea(var infile,outfile: text);
(* generate the call to movea *)
var x,y: real;
begin
      testblank(infile,outfile);
      readln(infile,x,y);
      movea(outfile,x,y)
end;
(* end module dops.mkmovea *)

(* begin module dops.mklinea *)
procedure mklinea(var infile,outfile: text);
(* generate the call to linea *)
var x,y: real;
begin
      testblank(infile,outfile);
      readln(infile,x,y);
      linea(outfile,x,y)
end;
(* end module dops.mklinea *)

(* begin module dops.mkdotr *)
procedure mkdotr(var infile,outfile: text);
(* generate the call to dotr *)
begin
      (* note that no testblank is needed because there are no arguments *)
      readln(infile);
      dotr(outfile)
end;
(* end module dops.mkdotr *)

(* begin module dops.mkpicnumber *)
procedure mkpicnumber(var infile,outfile: text);
(* generate the call to picnumber *)
var dx, dy, number: real; width, decimal: integer;
    centered: boolean;
begin
      testblank(infile,outfile);
      read(infile, dx, dy, number, width, decimal);
      skipblanks(infile);
      centered := (infile^='t'); (* a t means true *)
      readln(infile); (* skip past the line *)
      picnumber(outfile, dx, dy, number, width, decimal,'c');
end;
(* end module dops.mkpicnumber *)

(* begin module dops.mkxtic *)
procedure mkxtic(var infile,outfile: text);
(* generate the call to xtic *)
var length, dx, dy, number: real; width, decimal: integer;
   logxnormal: boolean;
   clogxnormal: char;
   logxbase: real;
begin
      testblank(infile,outfile);
      readln(infile, length, dx, dy, number, width, decimal,
           clogxnormal,logxbase);
      logxnormal := (clogxnormal = 't');
      xtic(outfile, length, dx, dy, number, width, decimal,
           logxnormal,logxbase);
{
procedure xtic(var afile: text; length, dx, dy, number: real;
               width, decimal: integer;
               logxnormal: boolean;
               logxbase: real);
}
end;
(* end module dops.mkxtic *)

(* begin module dops.mkytic *)
procedure mkytic(var infile,outfile: text);
(* generate the call to ytic *)
var length, dx, dy, number: real; width, decimal: integer;
   logynormal: boolean;
   clogynormal: char;
   logybase: real;
begin
      testblank(infile,outfile);
      readln(infile, length, dx, dy, number, width, decimal,
           clogynormal,logybase);
      logynormal := (clogynormal = 't');
      ytic(outfile, length, dx, dy, number, width, decimal,
           logynormal,logybase);
end;
(* end module dops.mkytic *)

(* begin module dops.mkxaxis *)
procedure mkxaxis(var infile,outfile: text);
(* generate the call to xaxis *)
var
   doaxisline: boolean;  (* line on axis is plotted *)
   axlength,fromtic,interval,totic: real;
   xsubintervals: real;
   length, dx, dy: real;
   width, decimal: integer;
   logxscale, logxnormal: boolean;
   clogxscale, clogxnormal: char;
   logxbase: real;
begin
      testblank(infile,outfile);
      readln(infile,axlength,fromtic,interval,totic,
             xsubintervals,
             length, dx, dy, width, decimal,
             clogxscale, clogxnormal,
             logxbase);
      logxnormal := (clogxnormal = 't');
      logxscale := (clogxscale = 't');
      doaxisline := true;
      xaxis(outfile,doaxisline, axlength,fromtic,interval,totic,
                    xsubintervals,
                    length, dx, dy, width, decimal,
                    logxscale, logxnormal,
                    logxbase);
end;
(* end module dops.mkxaxis *)

(* begin module dops.mkyaxis *)
procedure mkyaxis(var infile,outfile: text);
(* generate the call to yaxis *)
var
   doaxisline: boolean;  (* line on axis is plotted *)
   aylength,fromtic,interval,totic: real;
   ysubintervals: real;
   length, dx, dy: real;
   width, decimal: integer;
   logyscale, logynormal: boolean;
   clogyscale, clogynormal: char;
   logybase: real;
begin
      testblank(infile,outfile);
      readln(infile,aylength,fromtic,interval,totic,
             ysubintervals,
             length, dx, dy, width, decimal,
             clogyscale, clogynormal,
             logybase);
      logynormal := (clogynormal = 't');
      logyscale := (clogyscale = 't');
      doaxisline := true;
      yaxis(outfile,doaxisline, aylength,fromtic,interval,totic,
                    ysubintervals,
                    length, dx, dy, width, decimal,
                    logyscale, logynormal,
                    logybase);
end;
(* end module dops.mkyaxis *)

(* begin module dops.mkboxr *)
procedure mkboxr(var infile, outfile: text);
(* generate the call to the boxr routine *)
var
      width, height: real;
begin
      testblank(infile,outfile);
(*debug writeln(outfile,'in boxr!');*)
      readln(infile,width,height);
      boxr(outfile,width,height)
end;
(* end module dops.mkboxr *)

(* begin module dops.mkcboxr *)
procedure mkcboxr(var infile, outfile: text);
(* generate the call to the cboxr routine *)
var
      width, height: real;
begin
      testblank(infile,outfile);
(*writeln(outfile,'in cboxr');debug*)
      readln(infile,width,height);
(*writeln(outfile,'width height=',width:4:2,height:4:2);debug*)
      cboxr(outfile,width,height)
end;
(* end module dops.mkcboxr *)

(* begin module dops.mkibeam *)
procedure mkibeam(var infile, outfile: text);
(* generate the call to the ibeam routine *)
var
      width, height: real;
begin
      testblank(infile,outfile);
(*writeln(outfile,'in ibeam');debug*)
      readln(infile,width,height);
(*writeln(outfile,'width height=',width:4:2,height:4:2);debug*)
      ibeam(outfile,width,height)
end;
(* end module dops.mkibeam *)

(* begin module dops.mkcircler *)
procedure mkcircler(var infile, outfile: text);
(* generate the call to the circler routine *)
var
      radius: real;
begin
      testblank(infile,outfile);
      readln(infile,radius);
      circler(outfile,radius)
end;
(* end module dops.mkcircler *)

(* begin module dops.mkspiral *)
procedure mkspiral(var infile,outfile: text);
(* generate the call to spiral *)
var thickness: real; steps: integer; radius: real;
begin
      testblank(infile,outfile);
      readln(infile, thickness, steps, radius);
      spiral(outfile, thickness, steps, radius)
end;
(* end module dops.mkspiral *)

(* begin module dops.mkmovepolar *)
procedure mkmovepolar(var infile,outfile: text);
(* generate the call to movepolar *)
var angle, distance: real;
begin
      testblank(infile,outfile);
      readln(infile, angle, distance);
      movepolar(outfile, angle, distance)
end;
(* end module dops.mkmovepolar *)

(* begin module dops.mkarc *)
procedure mkarc(var infile,outfile: text);
(* generate the call to arc *)
var angle1, angle2, radius: real; steps: integer;
begin
      testblank(infile,outfile);
      readln(infile, angle1, angle2, radius, steps);
      arc(outfile, angle1, angle2, radius, steps);
end;
(* end module dops.mkarc *)

(* begin module dops.mkplusr *)
procedure mkplusr(var infile,outfile: text);
(* generate the call to plusr *)
var width, height: real;
begin
      testblank(infile,outfile);
      readln(infile, width, height);
      plusr(outfile, width, height)
end;
(* end module dops.mkplusr *)

(* begin module dops.mkxr *)
procedure mkxr(var infile,outfile: text);
(* generate the call to xr *)
var width, height: real;
begin
      testblank(infile,outfile);
      readln(infile, width, height);
      xr(outfile, width, height)
end;
(* end module dops.mkxr *)

(* begin module dops.mktest3d *)
procedure mktest3d(var outfile: text);
(* generate the call to test3d *)
begin
      (* note that no testblank is needed because there are no arguments *)
      test3d(outfile)
end;
(* end module dops.mktest3d *)

(* begin module dops.mksetgray *)
procedure mksetgray(var infile,outfile: text);
(* generate the call to setgray *)
var brightness: real;
begin
   testblank(infile,outfile);
   readln(infile, brightness);
   setgray(outfile, brightness);
end;
(* end module dops.mksetgray *)

(* begin module dops.mksetcolor *)
procedure mksetcolor(var infile,outfile: text);
(* generate the call to setcolor *)
var hue, saturation, brightness: real;
begin
   testblank(infile,outfile);
   readln(infile, hue, saturation, brightness);
   setcolor(outfile,  hue, saturation, brightness);
end;
(* end module dops.mksetcolor *)

(* begin module dops.mkrectinit *)
procedure mkrectinit(var outfile: text);
(* generate the call to rectinit *)
(* NOTE: xsideold and ysideold are GLOBALS *)
begin
      (* note that no testblank is needed because there are no arguments *)
   xsideold := 0.0;
   ysideold := 0.0;
   rectinit(outfile);
end;
(* end module dops.mkrectinit *)

(* begin module dops.mkrectdo *)
procedure mkrectdo(var outfile: text);
(* generate the call to rectdo *)
begin
   rectdo(outfile);
end;
(* end module dops.mkrectdo *)

(* begin module dops.mkrectsize *)
procedure mkrectsize(var infile,outfile: text);
(* generate the call to rectsize *)
(* NOTE: xsideold and ysideold are GLOBALS *)
var xside, yside: real;
begin
   testblank(infile,outfile);
   readln(infile, xside, yside);
   rectsize(outfile, xsideold, ysideold, xside, yside)
end;
(* end module dops.mkrectsize *)

(* begin module dops.translate *)
procedure translate(var demofile, infile, outfile: text;
                    depth: integer);
(* translate functions found in infile to pure pic input at outfile,
up to the picture end. use the file demo rather than input if
the command 'demo' is in infile.  depth keeps track of how deeply
the procedure has recursed in demonstration. *)
var
      buffer: string; (* part of a line of text from the source *)
      ch: char; (* a character read from infile *)
      go: boolean; (* continue testing characters on this line *)
      index: integer; (* a position in buffer *)
      pe: trigger; (* a trigger for the picture end *)
      ps: trigger; (* a trigger for the picture start *)

(* functions which are looked for: *)
      halt,
      demo,
      drawr, mover, liner,
      drawa, movea, linea,
      picnumber, xtic, ytic, xaxis, yaxis,
      dotr, boxr, cboxr, ibeam, circler, spiral, movepolar,
      arc, plusr, xr, test3d,
      setgray, setcolor,
      rectinit, rectsize, rectdo
         : trigger;

(* NOTE ON MAKING NEW NAMES:  A pair of names, were one name is a subset of the
other name, like "rect" and "rectsize", will cause trouble with the trigger
mechanism because the trigger for the subset ("rect") will fire first.  Then
the testblank procedure will be called when arguments are read, and the program
will halt.  To prevent this, do not allow any of the names to be subsets.  The
solution in the case of "rect" was to rename it "rectdo". *)

procedure fill;
(* fill up all the triggers *)
begin
(*                                       1         2 *)
(*                              12345678901234567890 *)
filltrigger(ps                ,'.PS                 ');
filltrigger(pe                ,'.PE                 ');

filltrigger(halt              ,'halt                ');
filltrigger(demo              ,'demo                ');
filltrigger(drawr             ,'drawr               ');
filltrigger(mover             ,'mover               ');
filltrigger(liner             ,'liner               ');
filltrigger(drawa             ,'drawa               ');
filltrigger(movea             ,'movea               ');
filltrigger(linea             ,'linea               ');
filltrigger(dotr              ,'dotr                ');
filltrigger(picnumber         ,'picnumber           ');
filltrigger(xtic              ,'xtic                ');
filltrigger(ytic              ,'ytic                ');
filltrigger(xaxis             ,'xaxis               ');
filltrigger(yaxis             ,'yaxis               ');
filltrigger(boxr              ,'boxr                ');
filltrigger(cboxr             ,'cboxr               ');
filltrigger(ibeam             ,'ibeam               ');
filltrigger(circler           ,'circler             ');
filltrigger(spiral            ,'spiral              ');
filltrigger(movepolar         ,'movepolar           ');
filltrigger(arc               ,'arc                 ');
filltrigger(plusr             ,'plusr               ');
filltrigger(xr                ,'xr                  ');
filltrigger(test3d            ,'test3d              ');
filltrigger(setgray           ,'setgray             ');
filltrigger(setcolor          ,'setcolor            ');
filltrigger(rectinit          ,'rectinit            ');
filltrigger(rectsize          ,'rectsize            ');
filltrigger(rectdo            ,'rectdo              ');
end;

procedure resetall;
(* reset all the triggers searched for *)
begin
      resettrigger(ps);
      resettrigger(pe);

      resettrigger(halt);
      resettrigger(demo);
      resettrigger(drawr);
      resettrigger(mover);
      resettrigger(liner);
      resettrigger(drawa);
      resettrigger(movea);
      resettrigger(linea);
      resettrigger(dotr);
      resettrigger(picnumber);
      resettrigger(xtic);
      resettrigger(ytic);
      resettrigger(xaxis);
      resettrigger(yaxis);
      resettrigger(boxr);
      resettrigger(cboxr);
      resettrigger(ibeam);
      resettrigger(circler);
      resettrigger(spiral);
      resettrigger(movepolar);
      resettrigger(arc);
      resettrigger(plusr);
      resettrigger(xr);
      resettrigger(test3d);
      resettrigger(setgray);
      resettrigger(setcolor);
      resettrigger(rectinit);
      resettrigger(rectsize);
      resettrigger(rectdo);
end;
procedure tests;
(* test for the functions.  if any function finds out what the
line is, it is responsible for completing the line by doing
the appropriate reading and readln'ing *)
begin
      if inpicture then begin
         if go then begin
            testfortrigger(ch,pe);
            if pe.found then begin
               if depth > 0
               then begin
                  writeln(output,'% .PE ignored in demo file');
                  while not eof(demofile)
                     do readln(demofile) (* skip remaining lines *)
               end
               else if not inpicture
                    then writeln(output,'% .PE ignored: not in picture')
                    else mkstoppic(infile,outfile);
               go := false
            end;
         end;

         if go then begin
            testfortrigger(ch,halt);
            if halt.found then begin
               mkhalt(outfile);
               go := false
            end
         end;
         if go then begin
            testfortrigger(ch,demo);
            if demo.found then begin
               if depth >= 1 then begin
                  writeln(output,'% RECURSIVE DEMONSTRATION REFUSED');
               end
               else begin
                  writeln(output,'% DEMONSTRATION BEGINS');
                  reset(demofile);
                  if inpicture then begin (* skip to ps in demo *)
                     (* this avoids use of pe.  if i called stoppic,
                        it would stops program (current incarnation
                        with no mouse response in await) *)
                     resettrigger(ps);
                     while (not ps.found) and (not eof(demofile)) do begin
                        read(demofile,ch);
                        testfortrigger(ch,ps);
                        while eoln(demofile) and (not eof(demofile))
                           do readln(demofile)
                     end;
                  end;
                  translate(demofile,demofile,outfile,depth + 1);
                  if not inpicture then mkstartpic(infile,outfile);
                  writeln(output,'% DEMONSTRATION ENDS');
               end;
               go := false
            end
         end;
         if go then begin
            testfortrigger(ch,drawr);
            if drawr.found then begin
               mkdrawr(infile,outfile);
               go := false
            end
         end;
         if go then begin
            testfortrigger(ch,mover);
            if mover.found then begin
               mkmover(infile,outfile);
               go := false
            end
         end;
         if go then begin
            testfortrigger(ch,liner);
            if liner.found then begin
               mkliner(infile,outfile);
               go := false
            end
         end;
         if go then begin
            testfortrigger(ch,drawa);
            if drawa.found then begin
               mkdrawa(infile,outfile);
               go := false
            end
         end;
         if go then begin
            testfortrigger(ch,movea);
            if movea.found then begin
               mkmovea(infile,outfile);
               go := false
            end
         end;
         if go then begin
            testfortrigger(ch,linea);
            if linea.found then begin
               mklinea(infile,outfile);
               go := false
            end
         end;
         if go then begin
            testfortrigger(ch,dotr);
            if dotr.found then begin
               mkdotr(infile,outfile);
               go := false
            end
         end;
         if go then begin
            testfortrigger(ch,picnumber);
            if picnumber.found then begin
               mkpicnumber(infile,outfile);
               go := false
            end
         end;
         if go then begin
            testfortrigger(ch,xtic);
            if xtic.found then begin
               mkxtic(infile,outfile);
               go := false
            end
         end;
         if go then begin
            testfortrigger(ch,ytic);
            if ytic.found then begin
               mkytic(infile,outfile);
               go := false
            end
         end;
         if go then begin
            testfortrigger(ch,xaxis);
            if xaxis.found then begin
               mkxaxis(infile,outfile);
               go := false
            end
         end;
         if go then begin
            testfortrigger(ch,yaxis);
            if yaxis.found then begin
               mkyaxis(infile,outfile);
               go := false
            end
         end;
         if go then begin
            testfortrigger(ch,cboxr);
            if cboxr.found then begin
               mkcboxr(infile,outfile);
               go := false
            end
         end;
         if go then begin
            testfortrigger(ch,ibeam);
            if ibeam.found then begin
               mkibeam(infile,outfile);
               go := false
            end
         end;
         if go then begin (* note that boxr will compete with cboxr for
               triggering since they will both trigger at the same time.
               we want cboxr to win and suppress boxr, so
               the test for boxr must follow that of cboxr *)
            testfortrigger(ch,boxr);
            if boxr.found then begin
               mkboxr(infile,outfile);
               go := false
            end
         end;
         if go then begin
            testfortrigger(ch,circler);
            if circler.found then begin
               mkcircler(infile,outfile);
               go := false
            end
         end;
         if go then begin
            testfortrigger(ch,spiral);
            if spiral.found then begin
               mkspiral(infile,outfile);
               go := false
            end
         end;
         if go then begin
            testfortrigger(ch,movepolar);
            if movepolar.found then begin
               mkmovepolar(infile,outfile);
               go := false
            end
         end;
         if go then begin
            testfortrigger(ch,arc);
            if arc.found then begin
               mkarc(infile,outfile);
               go := false
            end
         end;
         if go then begin
            testfortrigger(ch,plusr);
            if plusr.found then begin
               mkplusr(infile,outfile);
               go := false
            end
         end;
         if go then begin
            testfortrigger(ch,xr);
            if xr.found then begin
               mkxr(infile,outfile);
               go := false
            end
         end;
         if go then begin
            testfortrigger(ch,test3d);
            if test3d.found then begin
               mktest3d(outfile);
               go := false
            end
         end;

         if go then begin
            testfortrigger(ch,setgray);
            if setgray.found then begin
               mksetgray(infile,outfile);
               go := false
            end
         end;

         if go then begin
            testfortrigger(ch,setcolor);
            if setcolor.found then begin
               mksetcolor(infile,outfile);
               go := false
            end
         end;

         if go then begin
            testfortrigger(ch,rectinit);
            if rectinit.found then begin
               mkrectinit(outfile);
               go := false
            end
         end;

         if go then begin
            testfortrigger(ch,rectdo);
            if rectdo.found then begin
               mkrectdo(outfile);
               go := false
            end
         end;
         if go then begin
            testfortrigger(ch,rectsize);
            if rectsize.found then begin
               mkrectsize(infile,outfile);
               go := false
            end
         end;

         if go then begin
            testfortrigger(ch,ps);
            if ps.found then begin
               writeln(output,'% .PS ignored: already in picture');
               go := false
            end
         end;

         if go then begin
            if ch = ' ' then begin
               if not eoln(infile) then with buffer do begin
                  while (not eoln(infile)) and (length < maxstring) do begin
                     length := length + 1;
                     read(infile,letters[length]);
                  end;
                  graphstring(outfile,buffer,'l');
               end;
               go := false
            end;
         end
      end (* corresponds to: if inpicture then begin *)
      else begin (* not in a picture yet *)
         testfortrigger(ch,ps);
         if ps.found then begin
            mkstartpic(infile,outfile);
            go := false
         end
      end;
end;
begin
      fill;

      (* look at each line at a time: *)
      while not eof(infile) do begin
         resetall;
         clearstring(buffer);
         index := 0;
         go := true;
         while go do begin
            if eoln(infile) then begin
               (* nothing was recognized in the tests, so just dump: *)
               readln(infile);
               if inpicture
               then graphstring(outfile,buffer,'l')
               else begin
                  (* NONSTANDARD postscript comment *)
                  write(outfile,'%');
                  writestring(outfile,buffer);
                  writeln(outfile);
               end;
               go := false
            end
            else begin
               if index < maxstring then begin
                  read(infile, ch);

                  index := succ(index);
                  buffer.letters[index] := ch;
                  buffer.length := index;

                  tests
               end
               else begin
                  writeln(outfile);
                  writeln(outfile,'translate: line too long');
                  mkhalt(outfile)
               end
            end
         end
      end
end;
(* end module dops.translate *)

(* begin module dops.themain *)
procedure themain(var demo, fromfile, tofile: text);
(* the main procedure of the program *)
begin
      writeln(output,'%!'); (* this allows use of psview *)
      writeln(output,'% dops ',version:4:2);
      (* prestart as a convenience *)
(*
      startpic(tofile,defscale,5,5);
      writeln(tofile,'% picture has been started');
*)
      translate(demo,fromfile,tofile,0)
end;
(* end module dops.themain *)

begin
      themain(demo, input, output);
1: end.
