program search(book, searchinst, result, input, searchfeatures, output);
(* search: search a book for strings
      by Thomas Schneider, modified by Gary Stormo
      some concepts from Jon McCabe and Stephen O'Haire

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

module libraries: delman, delmods, prgmods *)

label 1; (* end of program *)

const
(* begin module version *)
   version = 7.03; (* of search.p 2017 Oct 30
2017 Oct 30, 7.03: clean up
2017 Sep 06, 7.02: insttype=2 fail - fix in delmod
2017 Jul 14, 7.01: namelength/linelength to 200/210
2017 Jul 12, 7.00: namelength/linelength to 200
2017 Feb 09, 6.99: checking marked mutations = passed!
2017 Feb 09, 6.98: fix lost mutation numbers
2017 Feb 09, 6.97: fix lines for negative mutations
2016 Oct 08, 6.96: remove extra " for quotestring
2011 May 08, 6.95: bug: org/chr missing
2011 Apr 05, 6.94: create proper link for demo.colors.search.jpg
2010 Nov 25, 6.93: remove debug 'dofeatures ====='; update address
2008 Apr 21, 6.92: first base of mutation counting is zero of feature
2008 Apr 21, 6.91: duplication removal by lister when mutations identical
2008 Apr 21, 6.90: cleanup dofeatures
2008 Apr 21, 6.88: hoining dofeatures
2008 Apr 21, 6.87: negparts need adjusting ... done?
2008 Apr 21, 6.86: negparts need adjusting
2008 Apr 21, 6.85: posparts and negparts need adjusting
2008 Apr 21, 6.84: located the bug, fixing ... in dofeatures
2008 Apr 21, 6.83: search string inversion bug
2008 Apr 19, 6.82: avoid infinite loop when bad command - quit!
2008 Apr 16, 6.81: clean up
2008 Apr 16, 6.80: R and I commands
2006 May 09, 6.78: fix bug: "name" was ""name"" from
                   procedure display at
                   "name the piece into the delila instructions"
2005 Sep 15, 6.77: Infinite loop if no final 'q' command!
                   This is hard to solve because the program needs to work both
                   interactively and with 'search < searchp'.
                   Recorded as a known bug.
2004 Sep  8, 6.76: fix bug in getcolorbackground (wrong call for newTparam)
2004 Aug 25, 6.75: searchp example built, delila subsystem takes comments
2004 Aug  5, 6.74: demo.colors.search.jpg
2004 Aug  5, 6.73: demo.colors.search directory (as tarball)
2004 Jul 21, 6.72: implement petal height control; document C command
2004 Jul 21, 6.71: implement petal height control; document C command
2004 Jul 18, 6.70: B command features written properly
2004 Jul 18, 6.69: B command reads ok.
2004 Jul 18, 6.68: module upgrades completed; new B and C commands
2004 Jul 18, 6.67: spelling; begin color named objects in searchfeatures
6.66 2002 Aug 28: add to wish list - Delila format mutations
6.65 2001 Sep  9: bug in ascii feature placement when complement:
                  featureinserts
6.64 2001 Sep  7: work on spiechange
6.63 2001 Sep  7: fix bug: complement should accept (/) symbols
6.62 2001 Aug 23: variables for both book and searchinst (not debugged).
6.61 2000 Nov  6: set dnamax to 10,000,000 for speed (it had been reset)
6.60 2000 Jun 16: spell improve, doc upgrade
6.59 2000 Jun 26: absorb upgrades to delmod
6.57 2000 Jun 3: allow [ and ] in search patterns
6.50 1999 Dec 15: named searches go into searchinst as name "whatever";
6.49 1999 Dec 13: delmod/y2k upgrade
6.44 1999 June 3: mark mismatches on features
6.42 1999 April 29: type 3 instructions
6.38 1999 April 28: delila instructions use 'same', fixed and in delmod
6.33 1999 April 8: delila instructions use 'same'
6.30 1999 April 7: "A" command resets objects to be arrows
6.29 1998 June 24: inst file renamed searchinst so as not to kill previous inst!!
6.25 1998 Jan 4: upgrade dnamax to 10 million base sequences for search of genomes
6.24 1997 Oct 10: strings made 2000 characters long
6.16 1997 May 8: relationship uses $ now, not ^
6.12 1997 May 7: upgrade to searchfeatures
origin before 1983 *)
(* end module version *)

(* begin module describe.search *)
(*
name
      search: search a book for strings

synopsis
      search(book: in, searchinst: out, result: out, input: intty,
             searchfeatures: out, output: out)

files
      book: any book from the Delila system

      searchinst: Delila instructions of the form 'get from 56 -5 to same
	 +5;' that define the location of found strings.  One must turn on
	 printing to the searchinst file to obtain these (see below).  If
	 there are instructions with names inside double quotes, then these
	 will be put out as Delila name instructions.  See the searchfeatures
	 file; this is turned on at the same time.  See examples.

      result: a transcript of the results seen on the output file.
         Lines not containing numerical data begin with an '*' so that
         they can be ignored by other programs such as genhis and xyplo.

      input: typed input from the user, or a file of rules.

      searchfeatures: features for the lister program.
         To start the file, simply provide a name inside double quotes
         (eg "EcoRI").
         Subsequent searches (eg gaattc) will be labeled with that name.
         To turn off the features, use an empty quote string, as "".
         The searchfeatures file can be concatenated with other features
         to create the features file for lister.

      output: messages, results and prompts to the user.

description
      (note: in the following examples, do not type the quote marks.)
      the search program allows one to look for simple patterns in a book.
      the patterns can be like 'ggag', that is, with particular bases
      (always written 5' to 3') or it can include unknown 'spacing' bases,
      as in 'ggagnnnnnnnnnatg'.  any base will be allowed in the n positions.
      one can shorten the instruction: 'ggag9natg', and one can make some of
      the spacing 'extentable' as in 'ggag5e4natg' which allows a 5 to 9
      spacing between the two elements.  one can obtain Delila instructions
      for the strings found by turning on printing, setting 'from' and 'to'
      values and searching.  for example: 'd p f -5 t +10 q gga6e3n#atg'
      sets up printing, with from=-5, to=+10.  the search will result in
      instructions for strings centered on the a of the atg (by the # symbol).
      the form '(a/g)ct' means to search for both 'act' and 'gct'.
      you may specify numbers of mismatches, and control how much is printed.
      you can type many commands on one line, separated by spaces.
      you can also search for relations between bases.  currently the
      allowed relations are: identity, non-identity, complementarity and
      non-complementarity.  see delman.use.search or type 'help' while inside
      the program to get more information.

      NOTE: Many commands now are capital letters to avoid
      confusion between commands.  See the help function (H, ?)
      for details.

      2004 July 21: There is a new command, B, which
      defines colors on seach features as displayed by the
      lister program.

examples

      If one is working with an odd binding site (one with an odd number
      of bases) one should use the # symbol to obtain Delila instructions.
      The complement sequence will continue to number the central base.
          gaa#nttc complemented becomes gaa#nttc

      If one is working with an even binding site (one with an even number
      of bases) one should use the % symbol to obtain Delila instructions.
      The complement sequence will continue to number the following base.
          ga%attc complemented becomes ga#attc

      The program is pretty smart about writing the Delila instructions.  If
      one searches for the complementary sequence, the instructions are
      automatically written to extract the complementary pattern found.  Thus
      if one searches for #gtt in ex0bk (the example book file), there is one
      found in the positive direction of the fragment.  Then if one takes the
      complement with the "~" command, one is searching for aa#c.  Two of these
      are located in the piece.  The instructions are written so that the gtt's
      all line up, as is easily checked by extracting the fragments with delila
      and looking with alist.

      To create searchfeatures, define the name of a search string by typing
      the name inside quotes (as: "EcoRI") and then search.  Vertical bars or
      carets (| or ^) in the search string (as: g|aattc) will carry over to
      the feature.

If you have search instructions:

D
from -200
to   +200
q
"aceB"
agttatcaagtatttttaattaaaatggaaattgtttttgattttgcattttaaatgagtagtcttagtt#n
q

then the resulting Delila instructions in searchinst will look like:

organism E.coli; chromosome E.coli;
piece U00096;
name "aceB";
get from 4212981 -200 to same +200 direction +;

documentation
   delman.use.search

see also
   {An example search parameter, searchp file:} searchp
   {A search parameter file for most restriction enzymes:} enzyme

   {Program to create input DNA sequences:} delila.p
   {Example book for input:} ex0bk
   {Program that takes the bookand the scan features as input
   to create a map of the DNA sequence:} lister.p

   {Example colors:}
   demo.colors.search.jpg

   {Example of colored boxes in search results:}
   demo-colors-search.zip

   {A second colored box example:} search-bar-color.zip

author
   Thomas D. Schneider, modified by Gary Stormo

bugs
   There is overlap between the letters used as commands to the program and
   letters used as ambiguous bases.  For instance, h can mean (a/c/t) or it
   can mean 'help'.  The best way to avoid confusion is to always start
   search strings with either a,c,g,t,n or (.  Warning: if you use a file for
   input, be sure that the rules include a quit command and have no errors in
   them.  It is possible that errors will lead to an infinite loop, though
   this has never been observed.   (This may be a general problem with
   interactive i/o in pascal on your computer.)

   The searchfeatures only work with rigid strings since only the first
   definition will be accepted by lister.  It is not clear how to handle
   variable sizes (like gat3etag).

   A search for "incorrect.symmetry" cc#gg will be located correctly, but
   following that by ~ = (to get the complement and search again) will not
   give a correct display.  Using cc%gg will because the % matches the
   symmetry of the site.  Likewise, "ok" tc#nga ~ = "not.ok" tc%nga ~ =
   It is not yet clear whether the user can be protected against this.

   The variables orgchange, chrchange and piechange do not work properly.
   Probably they should be done for both the book being read AND for the
   delila instructions being written.
   HOWEVER it seems that since the book has a certain structure,
   the output instructions can follow it to some extent.

   1999 June 2:  It would be nice to record mismatches on the feature, but
   this would require changing the definition for every found piece!

   2002 Aug 28: add to wish list - Delila format mutations
---
* piece: U24170, #1, configuration: linear, direction: +, begin: 3158, end: 3238
* "atcctgggaatttctgggaa"
  "agactgggcatgtctgggca"
  " xx     x  x      x " 5 mismatche(s)
   ^ 3194 *  
---
Have search spit out the "with"-type Delila instructions for these.

   2005 Sep 15:  The program will go into an infinite loop if one calls it
   non-interactively like this:

     search < searchp

   if there is no final 'q' command.  It is not clear how to solve
   this because interactive input relys on detecting the eof.  So just
   make sure that your searchp files end with 'q'!

*)
(* end module describe.search *)

(* no longer pointed to:
   demo.colors.search.tar.Z
   {(See untarball in the toolkit,}
   https://alum.mit.edu/www/toms/toolkit.html
   {for unpackaging the tar file.)}
*)

(****************************************************
*     programming note
* search uses a period to end search strings rather
* than the string length.  this is a historical
* artifact that must be lived with until a more
* advanced search is written (preferably based on
* a bnf).
*****************************************************)

(* more constants: *)
(* BEGIN module interact.const *)
(* lock in place for larger strings ... *)
      maxstring = 2000; (* the maximum string *)
(* END module interact.const version = 'prgmod 3.96  85 mar 18 tds'; *)

   maxquotestring = 2000; (* maximum length of a string inside "" *)
   answersmax = 10; (* maximum number of times to try to get an answer
     from input.  If the input is a file, we NEVER will get an answer
     and go into an infinite loop.  To avoid that, count the number
     of answers.  If it exceeds this constant, the program stops. *)

(* LOCK begin module book.const *)
(* constants needed for book manipulations *)
      dnamax = 10000000;  (* length of dna arrays *)
      namelength = 200; (* maximum key name length *)
      linelength = 210; (* maximum line readable in book *)
(* LOCK end module book.const version = 7.47; {of delmod.p 2000 June 26} *)

(* begin module datetime.const *)
   datetimearraylength = 19; (* length of dataarray for dates,
       It is just long enough to include the 4 digit year - solving the
       year 2000 problem:
1980/06/09 18:49:11
123456789 123456789 
         1         2
*)
(* end module datetime.const version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module info.const *)
   infofield = 10; (* size of field for printing information in bits *)
   infodecim = 6; (* number of decimal places for printing information *)
       (* these are used for conlist only *)
   nfield = 8; (* size of field for printing n, the number of sites *)
(* end module info.const version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

(* begin module ribl.const *)
   minribl = -2000; (* lowest  ribl matrix from allowed *)
   maxribl = +2000; (* highest ribl matrix to   allowed *)
   defnegativeinfinity = -1000; (* default for negative infinity
                                   for the Ri(b,l) table *)
(* end module ribl.const version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

type

(* begin module datetime.type *)
   (* array for dates *)
   datetimearray = packed array[1..datetimearraylength] of char;
(* end module datetime.type version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.type *)
(* types needed for book manipulations *)

      chset = set of 'a'..'z';

      (* types defined in book definition *)

      alpha = packed array[1..namelength] of char; (* this is not alfa *)

      (* name is a left justified string with blanks following the
         characters *)
      name = record
         letters: alpha;
         length: 0..namelength (* zero means an unspecified structure *)
      end;

      lineptr = ^line;
      line = record (* a line of characters *)
         letters: packed array [1..linelength] of char;
         length: 0..linelength;
         next: lineptr
      end;

      direction = (plus, minus, dircomplement, dirhomologous);
      configuration = (linear, circular);
      state = (on, off);

      header = record (* header of key *)
         keynam: name; (* key name of structure *)
         fulnam: lineptr; (* full name of structure *)
         note:   lineptr (* note key *)
      end;

(* begin module base.type *)
   (* define the four nucleotide bases *)
   base = (a,c,g,t);
(* end module base.type version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

      (* sequence types *)
      dnaptr = ^dnastring;
      dnarange = 0..dnamax;
      seq = packed array[1..dnamax] of base;
      dnastring = record
         part: seq;
         length: dnarange;
         next: dnaptr
      end;


      orgkey = record (* organism key *)
         hea: header;
         mapunit: lineptr (* genetic map units *)
      end;

      chrkey = record (* chromosome key *)
         hea: header;
         mapbeg: real; (* number of genetic map beginning *)
         mapend: real (* number of genetic map ending *)
      end;

      pieceptr = ^piece;
      piekey = record (* piece key *)
         hea: header;
         mapbeg: real; (* genetic map beginning *)
         coocon: configuration; (* configruation (circular/linear) *)
         coodir: direction; (* direction (+/-) relative to genetic map *)
         coobeg: integer; (* beginning nucleotide *)
         cooend: integer; (* ending nucleotide *)
         piecon: configuration; (* configruation (circular/linear) *)
         piedir: direction; (* direction (+/-) relative to coordinates *)
         piebeg: integer; (* beginning nucleotide *)
         pieend: integer; (* ending nucleotide *)
      end;


      piece = record
         key: piekey;
         dna: dnaptr
      end;

      reference = record
         pienam : name;  (* name of piece referred to *)
         mapbeg : real;  (* genetic map beginning *)
         refdir : direction;  (* direction relative to coordinates *)
         refbeg : integer;  (* beginning nucleotide *)
         refend : integer;  (* ending nucleotide *)
      end;

      genkey = record  (* gene key *)
         hea : header;
         ref : reference;
      end;

      trakey = record  (* transcript key *)
         hea : header;
         ref : reference;
      end;

      markerptr = ^marker;
      markey = record (* marker key *)
         hea : header;
         ref : reference;
         sta : state;
         phenotype : lineptr;
         next : markerptr;
      end;

      marker = record
         key : markey;
         dna : dnaptr;
      end;

(* end module book.type version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* 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 = 4.86; (@ of prgmod.p 2004 Sep 8 *)
(* end module interact.type version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

{
   ribltypeptr = integer; (* dummy to satisfy compiler NO GOOD!!!! *)
}

(* begin module wave.type *)
      waveptr = ^waveparam; (* to link several wave definitions together *)
      waveparam = record (* parameters to define a cosine wave *)
         (* define a cosine wave: *)
         extreme: char; (* h or l, the high or low extreme to be defined *)
         wavelocation: real; (* the location in bases of the extreme *)
         wavebit: real; (* the location in bits of the extreme *)
         waveamplitude: real; (* the amplitude of the wave in bits *)
         wavelength: real; (* the wave length of the wave in bases *)
{ 2000 Nov 1.  This is difficult because the length of a curve is not
simply multiplied by a scale factor.
         BASEdashon: real; (* dashon interval in bases *)
         BASEdashoff: real; (* dashon interval in bases *)
         BASEdashoffset: real; (* dashon interval in bases *)
}
         dashon: real; (* the size of on dashes in cm.
                          dashon <= 0 means no dashes *)
         dashoff: real; (* the size of on dashes in cm *)
         dashoffset: real; (* the size of off dashes in cm *)
         thickness: real; (* thickness of the wave in cm.  <= 0 means default *)
         next: waveptr; (* the next wave *)
      end;
(* end module wave.type version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

(* begin module scan.type *)
   (* define the ribl matrix.  Source library is prgmod.p
      REQUIRES ribl.const module
      REQUIRES base.type or book.type module
      REQUIRES wave.type module for wave definition! *)
   rblarray = array[a..t, minribl..maxribl] of real; (* real(B,L) *)
   ribltype = record
      riblname: string; (* name of the weight matrix *)
      riblheader: stringptr; (* misc header information, starting with '*' *)
      data: rblarray; (* real(B,L) *)
      numbers: array[a..t, minribl..maxribl] of integer; (* n(b,l) *)
      frombase, tobase: integer; (* range of the matrix *)
      mean, stdev: real; (* mean and standard deviation of distribution *)
      consensus: real; (* Ri value of consensus sequence *)
      anticonsensus: real; (* Ri value of anticonsensus sequence *)
      averageRi: real; (* average Ri for random sequence *)
      n: integer; (* number of sequences used to create matrix *)
      symmetry: char; (* the symmetry of the Riblmatrix *)
      waves: waveptr; (* cosine wave definitions for the matrix *)
      cmperbase, cmperbit: real; (* required for wave dash definition *)
      Ribound: real; (* the Ri boundary for this definition *) 
      Zbound: real; (* the Z boundary for this definition *) 
      Pbound: real; (* the probability boundary for this definition *) 
   end;
   ribltypeptr = ^ribltype; (* pointer to ribltype so it does
      not need to be created unless actually used *)
(* end module scan.type version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

(* begin module definition.feature.type *)

   (* Features and their definitions for the lister program *)
   petaltype = record
      (* define parts of a colored rectangle around a walker *)
      el: integer; (* edgelinewidth: edge linewidth (integer) *)
      c : char;    (* color kind: how color is defined: r for RGB, h for HSB *)
      (* ---------------------------------------------- *)
      eh: real;    (* edgeh: edge hue        OR red     *)
      es: real;    (* edges: edge saturation OR blue    *)
      eb: real;    (* edgeb: edge brightness OR green   *)
      (* ---------------------------------------------- *)
      fh: real;    (* fillh: fill hue        OR red     *)
      fs: real;    (* fills: fill saturation OR blue    *)
      fb: real;    (* fillb: fill brightness OR green   *)
      (* ---------------------------------------------- *)
      spare: real;    (* spare variable *)
   end;

   definitiontypeptr = ^definetype; (* pointer to definetype *)
   definetype = record (* store a definition *)
(*
         define "Site" "-" "<v0*]" "[*0^>" -7 -3.5 0 +3.0 +5
*)
      nametag: string; (* name of the feature *)
      background: string; (* background padding characters *)
      negparts: string; (* parts of the feature display *)
      posparts: string; (* parts of the feature display *)
      locations: array [1..maxstring] of real; (* locations of the marks
         of the feature *)
      (* ---------------------------- *)
      (* computed values: *)
      marks: integer; (* the number of marks used in this feature *)
      min: real; (* one end of the feature (minimum location) *)
      max: real; (* the other end of the feature (maximum location) *)
      number: integer; (* the number of the definition read in *)
(*[[*)
      matrix: ribltypeptr; (* an ribl matrix if available *)
      (* The following boundaries are DIFFERENT from the originally defined
         boundaries of the matrix.  This allows the matrix to keep its bounds
         (eg for default values for scanning AND it allows the lister program
         to define separate boundaries for displaying green versus red on the
         zero bar.  *)
      RiboundFeature: real; (* the Ri boundary for this feature *)
      ZboundFeature: real; (* the Z boundary for this feature *)
      PboundFeature: real; (* the probability boundary for this feature *)
      petal: petaltype; (* color rectangle background *)
(*]]*)
      next: definitiontypeptr; (* the next definition recorded *)
   end;

(* ---------------------------- *)
   featuretypeptr = ^featuretype; (* pointer to featuretype *)
   featuretype = record (* store a feature *)
(*
         @ K01789 270 +1 "Site" "Positive 270"
         @ K01789 270 -1 "Site" "Negative 270"
         @ K01789 280 +1 "Site" "Positive 280"
[[
         @ K01789 229.0 -1 "dnaA" "+12.2 bits " 12.200338 -0.473212 0.318031
]]
*)
      id: string; (* name of the sequence piece the feature is on *)
      coordinate: real; (* location of zero base *)
      orientation: integer; (* -1 or +1 is in direction of piece *)
      nametag: string; (* the name tag - name of the feature *)
      othertag: string; (* the other tag - another tag for misc use *)
      (* ---------------------------- *)
      definition: definitiontypeptr; (* the definition of the feature *)
(*[[*)
      (* ---------------------------- *)
      (* the following three values are only defined when there
         is a definition for the feature matrix *)
      Ri: real; (* individual information of this site *)
      Z: real; (* Z score for this site *)
      probability: real; (* probability for this site *)

      (* Three parameters replace Ri, Z and p when there is no matrix,
         Aparam, Bparam and Cparam.  If there is a matrix, Tparam
         will be blank.  Otherwise: *)
      Tparam: char; (* a parameter that defines the kind of parameter
                       that the remainder parameters are. *)
      Aparam: real; (* a parameter *)
      Bparam: real; (* a parameter *)
      Cparam: real; (* a parameter *)
      Dparam: real; (* a parameter *)
         (* Since these parameters
         are now (2004 July 17) able to substitute for Ri, Z and
         probability, the Tparam cannot be placed in front of the
         three parameters in feature files.  Tparamvalues:
            h = hsb feature colors: hue, saturation, brightness
            r = rgb feature colors: red, green, blue
            b = bit range for the feature
            ' ' = defaults to sequence walker
         Summary:
         Tparam  Aparam     Bparam     Cparam     Dparam     For
          ' '    Ri         Z          P          (none)     walker
          'h'    hue        saturation brightness (none)     ascii lines
          'r'    red        green      blue       (none)     ascii lines
          'H'    hue        saturation brightness thickness  color bar
          'R'    red        green      blue       thickness  color bar
         In this scheme one can control the color of the bar as an
         ascii string OR one can control the colored rectangle behind
         the bar but NOT BOTH.
         *)
{zzz}
(*]]*)
      (* ---------------------------- *)
      evencoordinate: real; (* location of zero base for asymmetric and odd
         symmetry sites, with 0.0 as the decimal.  For even symmetry sites
         it's the point of symmetry, with 0.5 as the decimal. *)
      (* ---------------------------- *)
      (* computed values: *)
      unsatisfied: boolean; (* true when the feature has not yet been used
                               during printing a line *)
      fromrange: real; (* the lower part of the feature *)
        torange: real; (* the higher part of the feature *)
      number: integer; (* the number of the feature read in *)
      desiredline: integer; (* The line the feature wants to be printed on.
                            Zero means any value.  Line increases down the
                            page. *)
      (* ---------------------------- *)
      next: featuretypeptr; (* the next feature recorded *)
   end;
(* end module definition.feature.type version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

         (* for keeping track of related bases *)
      intarray = array[1..maxstring] of integer;

var
      book, (* the book being searched *)
      searchinst, (* Delila instructions *)
      result, (* for hardcopy output *)
      searchfeatures: (* search features for lister program *)
      (* input, search instructions *)
      (* output, display to screen *)
         text;

      datetime: datetimearray; (* the date and time of this search *)
      buffer: string; (* a rule or typed line from input *)
      gotten: boolean; (* whether a token desired was obtained *)
      command: char; (* a user's command *)

      (* user settable varibles *)
      (* variables defined in viewhelp: *)
      viewt,viewe,viewi,viewf,views,viewp,viewd,viewm,viewb: boolean;
      viewnothing: boolean; (* a flag that the user is not viewing
         anything.  if true, the user is prevented from doing a search
         for which no results will be shown *)

      mismatches: integer; (* maximum number of mismatches allowed *)
      printinginstructions: boolean; (* true means that when a pattern
         is found, instructions for getting to it in the library are
         printed into searchinst *)
      fromplace,toplace: integer; (* what to get around a found pattern
         #'d base (Delila print command) *)
      insttype: integer; (* instruction type.
          insttype=1 means the form   get from p  -/+f to same +/-t dir +/-;
          insttype=2 means the form   get from p1 -/+f to p2   +/-t dir +/-;
          where p or p2 are the current location, p1 is the previous
          location.
          insttype=3 means the form   get from p1 -/+f to p2   +/-t dir +/-;
          where p1 is the current location 5' end, p2 is the 3' end
      *)
      oktosearch: boolean; (* if true, searches are performed when
          the = command is given or when the pattern is typed.  If false,
          searches are performed only when the = command is given *)

      firstinput: boolean; (* this variable helps the user the first time,
          when, under unix operating systems, the first line is ignored *)

      (* user defined search start point and shift between search steps: *)
      startpoint, shift: integer;

      ok: boolean; (* the pattern typed in was ok *)

      pattern, (* the pattern to search for *)
      exppattern, (* the pattern expanded *)
      lastpattern: (* the pattern previously searched for *)
         string;

      bookmatches: integer; (* total matches in the book by a pattern *)

         (* relatedbase is an array of integer that tells what position
            in the pattern this base is related to *)
      relatedbase: intarray;
      relational: boolean;  (* true if there is a relation in the pattern *)

      flip: boolean; (* if true, reverse the direction of written Delila
         instructions.  This is used when the fromposition is greater
         than the to position. *)

      reusethepiece: boolean; (* if true, don't read the book anymore
         because there is only one sequence in the book!  This should
         speed up searches. *)

      (* variables for searchfeatures *)
      dosearchfeatures: boolean; (* whether to create the searchfeatures *)
      featuredefwritten: boolean; (* whether the current feature definition
          has been written to searchfeatures *)
      thedefinition: definitiontypeptr; (* the current definition *)
      thefeature: featuretypeptr; (* the current feature *)
      doarrow: boolean; (* the next feature will be an arrow *)
      featureinserts: boolean; (* inserts to the feature defined or not.
         For historical reasons symbols like # are inserted into the
         feature at the final match stage.  But if one searches with a
         string and takes the complement, then they must NOT be inserted
         again during the second search.  I suppose the same thing
         could happen with mismatches.  So this boolean makes sure
         that the symbols are inserted only one time. *)
      poundsign: boolean; (* true if a # was used in the search string.
         Normally the search feature ends are just specified from
         the left edge.  However, when a user requests #, it means
         that the feature is specified from the point given by the #. *)

(* ooo *)
      (* 2001 Aug 23: I will DEFINE the following to be tracking
         the book state.  See below for the new variables. *)
      orgchange, chrchange, piechange: boolean; (* that any changed *)
      orgopen, chropen, pieopen: boolean; (* used by getocp *)
      org: orgkey; (* the current organism *)
      chr: chrkey; (* the current chromosome *)
      pie: pieceptr; (* the current piece *)
      theline: integer; (* tracking line number in book *)

      (* 2001 Aug 23: these variables track the state of the searchinst *)
      sorgchange, schrchange, spiechange: boolean; (* that any changed *)
      sorgopen, schropen, spieopen: boolean; (* which are open *)

      patterncomplement: boolean; (* if false, then the pattern is as typed
         in.  If true, then the pattern is the complement of that typed in.
         By keeping track of this, the orientation of the Delila instructions
         can be made correct. *)

   answers: integer; (* current number of answers to a question *)

(* begin module book.var *)
(* ************************************************************************ *)
(* global variables needed for book manipulations *)

      (* free storage: *)
      freeline: lineptr; (* unused lines *)
      freedna: dnaptr; (* unused dnas *)

      readnumber: boolean; (* whether to read a number from the notes, or
                              to read in the notes *)
      number: integer; (* the number of the item just read *)
      numbered: boolean; (* true when the item just read is numbered *)
      skipunnum: boolean; (* a control variable to allow skipping of
                             un-numbered items in the book *)

(* ************************************************************************ *)
(* end module book.var version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* procedure modules *****************************************************)

(* begin module package.primitive *)
(* ************************************************************************ *)
(* 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 = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* 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 = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module copylines *)
function copylines(var fin, fout: text; n: integer): integer;
(* copy n lines of file fin to file fout.
the actual number of lines copied is returned. *)
var
      index: integer; (* the current line number *)
begin (* copylines *)
      index := 0;
      while (not eof(fin)) and (index < n) do begin
         copyaline(fin, fout);
         index := succ(index)
      end;

      copylines := index
end; (* copylines *)
(* end module copylines version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* ************************************************************************ *)
(* end module package.primitive version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module package.getocp *)
(* ************************************************************************ *)
(* begin module package.brpiece *)
(* ************************************************************************ *)
(* begin module book.basis *)
(* procedures needed for book manipulations *)

(* get procedures should be used for all linked lists of records *)

procedure getline(var l: lineptr);
(* obtain a line from the free line list or by making a new one *)
begin
      if freeline<>nil
         then begin
            l:=freeline;
            freeline:=freeline^.next
         end
         else new(l);
      l^.length:=0;
      l^.next:=nil
end;

procedure getdna(var l: dnaptr);
begin
      if freedna<>nil
         then begin
            l:=freedna;
            freedna:=freedna^.next
         end
         else new(l);
      l^.length:=0;
      l^.next:=nil
end;

(* clear procedures should be called each time the records are no longer needed
         failure to do this may result in a stack overflow. *)

procedure clearline(var l: lineptr);
(* return a line to the free line list *)
var   lptr: lineptr;
begin
      if l<>nil then begin
         lptr:=l;
         l:=l^.next;
         lptr^.next:=freeline;
         freeline:=lptr
      end
end;

procedure writeline(var afile: text; l: lineptr;
                        carriagereturn: boolean);
(* write a line to a file, with carriage return if
carriagereturn is true. *)
var
   index: integer; (* index to characters in l *)
begin
   with l^ do begin
      for index := 1 to length
      do write(afile, letters[index]);
   end;
   if carriagereturn then writeln(afile);
end;

procedure showfreedna;
(* show the freedna list *)
var
   counter: integer; (* count of freedna list *)
   l: dnaptr; (* pointer into freedna list *)
begin
   l := freedna;
   counter := 0;
   while l <> nil do begin
      counter := succ(counter);
      write(output,counter:1);
      write(output, ', length = ',l^.length:1);
{
This is illegal according to gpc because one cannot
write a pointer to a text file.  It can be unearthed
for debugging.
      write(output, ', pointer id: ',l:1);
}
      writeln(output);
      l := l^.next
   end;
end;

procedure cleardna(var l: dnaptr);
(* clear the dna strutures to the free list *)
var   lptr: dnaptr;
begin
      if l<>nil then begin
         lptr:=l;
         l:=l^.next;
         lptr^.next:=freedna;
         freedna:=lptr
      end
end;

procedure clearheader(var h: header);
(* clear the header h (remove lines to free storage) *)
begin
      with h do begin
         clearline(fulnam);
         while note<>nil do clearline(note)
      end
end;

procedure clearpiece(var p: pieceptr);
(* clear the dna of the piece *)
begin
      while p^.dna<>nil do cleardna(p^.dna);
      clearheader(p^.key.hea)
end;

function chartobase(ch:char):base;
(* convert a character into a base *)
begin
      case ch of
         'a': chartobase:=a;
         'c': chartobase:=c;
         'g': chartobase:=g;
         't': chartobase:=t
      end
end;

function basetochar(ba:base):char;
(* convert a base into a character *)
begin
      case ba of
         a: basetochar:='a';
         c: basetochar:='c';
         g: basetochar:='g';
         t: basetochar:='t';
      end
end;

function complement(ba:base):base;
(* take the complement of ba *)
begin
      case ba of
         a: complement:=t;
         c: complement:=g;
         g: complement:=c;
         t: complement:=a;
      end
end;

function chomplement(b: char): char;
(* create the character complement of base b.  I must be getting hungry! *)
begin
  chomplement := basetochar(complement(chartobase(b)));
end;

function pietoint(p: integer; pie: pieceptr): integer;
(* p is a coordinate on the piece.
   we want to transform p into a number
   from 1 to n: an internal coordinate system for
   easy manipulation of piece coordinates *)
(* Note: the dirhomologous and dircomplement are treated as
plus and minus directions, which MIGHT NOT BE RIGHT! *)
var   i: integer; (* an intermediate value *)
begin
      with pie^.key do begin
         case piedir of
            dirhomologous,
            plus:  if p>=piebeg
                   then i:=p-piebeg+1
                   else i:=(p-coobeg)+(cooend-piebeg)+2;
            dircomplement,
            minus: if p<=piebeg
                   then i:=piebeg-p+1
                   else i:=(cooend-p)+(piebeg-coobeg)+2
         end;
         pietoint:=i
      end
end;

function inttopie(i: integer; pie: pieceptr):integer;
(* i is in the range 1 to some maximum.  it is an internal coordinate
   system for the program.  we want to do a
   coordinate transformation to obtain
   a value in the range of the piece called pie:
   i=1 corresponds to piebeg and
   i=its maximum corresponds to pieend *)
(* Note: the dirhomologous and dircomplement are treated as
plus and minus directions, which MIGHT NOT BE RIGHT! *)
var   p: integer; (* an intermediate value *)
begin
      with pie^.key do begin
         case piedir of
            dirhomologous,
            plus: begin
               p:=piebeg+(i-1);
               if p>cooend
                  then if coocon=circular
                     then p:=p-(cooend-coobeg+1)
            end;
            dircomplement,
            minus: begin
               p:=piebeg-(i-1);
               if p<coobeg
                  then if coocon=circular
                     then p:=p+(cooend-coobeg+1)
            end
         end;
         inttopie:=p
      end
end;

function piecelength(pie: pieceptr): integer;
(* return the length of the dna in pie *)
begin
      piecelength:=pietoint(pie^.key.pieend,pie)
end;

(* end module book.basis version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.getto *)
function getto(var thefile: text; var theline: integer; ch: chset): char;
(* search the file for a character in the first line which is a
   member of the set ch.
Note:  on 1999 March 10 the definition of this function was cleaned
up.  Instead of putting thefile on the line AFTER the charcter ch
has been found, it puts thefile ON the line.  Other routines like
brdna and brpiece have to move to the next line themselves.  This makes
getto give the OBJECT.  *)
var
   achar: char; (* a character in thefile *)
   done: boolean; (* done finding achar *)
begin
   achar:=' ';
   done := false;
   while not done do begin
      if eof(thefile)
      then done := true
      else begin
         achar := thefile^;
         if achar in ch
         then done := true
         else begin
            readln(thefile);
            theline := succ(theline)
         end
      end
   end;
   if (achar in ch) then getto:=achar
                    else getto:=' '
{ The old method -
      while (not(achar in ch)) and (not eof(thefile))
         do begin
            readln(thefile,achar);
            theline := succ(theline)
         end;
   if (achar in ch) then getto:=achar
                    else getto:=' '
}
end;
(* end module book.getto version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.skipstar *)
procedure skipstar(var thefile: text);
(* skip start of line (or star = '*'). *)
begin (* skipstar *)
   if eof(thefile) then begin
      writeln(output,' procedure skipstar: end of book found');
      halt
   end
   else begin
      if thefile^ <> '*' then begin
         writeln(output,' procedure skipstar: bad book');
         writeln(output,' "*" expected as first character on the line, but "',
                        thefile^,'" was found');
         halt
      end;
      get(thefile); (* skip the star *)
      if thefile^ <> ' ' then begin
         writeln(output,' procedure skipstar: bad book');
         writeln(output,' "* " expected on a line but "*',
                        thefile^,'" was found');
         halt
      end;
      get(thefile) (* skip the blank *)
   end
end; (* skipstar *)
(* end module book.skipstar version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.brreanum *)
procedure brreanum(var thefile: text; var theline: integer; var reanum: real);
(* read a real number from the file *)
begin
      skipstar(thefile);
      readln(thefile,reanum);
      theline := succ(theline)
end;
(* end module book.brreanum version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.brnumber *)
procedure brnumber(var thefile: text; var theline: integer; var num: integer);
(* read a number from the file *)
begin
      skipstar(thefile);
      readln(thefile,num);
      theline := succ(theline)
end;
(* end module book.brnumber version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.brname *)
procedure brname(var thefile: text; var theline: integer; var nam: name);
(* read a name from the file *)
var   i: integer; (* an index to the name *)
      c: char; (* a character read *)
begin (* brname *)
      skipstar(thefile);
      with nam do begin
         length:=0;
         repeat
            length:=succ(length);
            read(thefile,c);
            letters[length] := c
         until (eoln(thefile)) or
               (length>=namelength) or
               (letters[length]=' ');
         if letters[length]=' ' then length:=length-1;
         if length<namelength
            then for i:=length+1 to namelength do letters[i]:=' '
      end;
      readln(thefile);
      theline := succ(theline)
end; (* brname *)
(* end module book.brname version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.brline *)
procedure brline(var thefile: text; var theline: integer; var l: lineptr);
(* read a line from the file *)
var
      i: integer;
      acharacter: char;
begin
      skipstar(thefile);
      i:=0;
      (* protection added 2005 Sep 15 *)
      while (not eoln(thefile)) and (i < linelength)  do begin
         i:=succ(i);
         read(thefile,acharacter);
         l^.letters[i]:=acharacter
      end;
      (* protection added 2005 Sep 15 *)
      if not eoln(thefile) then begin
         writeln(output,'***********************************************');
         writeln(output,'* WARNING: brline: book line length exceeded');
         writeln(output,'* linelength > ',linelength:1,' characters');
         writeln(output,'* Only ',linelength:1,' characters read from book');
         writeln(output,'***********************************************');
      end;
      l^.length:=i;
      l^.next:=nil;
      readln(thefile);
      theline := succ(theline)
end;
(* end module book.brline version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.brdirect *)
procedure brdirect(var thefile: text; var theline: integer; var direct: direction);
(* read a direction *)
var   ch: char;
begin
      skipstar(thefile);
      readln(thefile,ch);
      theline := succ(theline);
      if ch='+' then direct:=plus
                else direct:=minus
end;
(* end module book.brdirect version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.brconfig *)
procedure brconfig(var thefile: text; var theline: integer; var config: configuration);
(* read a configuration *)
var   ch: char;
begin
      skipstar(thefile);
      readln(thefile,ch);
      theline := succ(theline);
      if ch='l' then config:=linear
                else config:=circular
end;
(* end module book.brconfig version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.brnotenumber *)
procedure brnotenumber(var thefile: text; var theline: integer; var note: lineptr);
(* book note reading to obtain the number of the object.
the procedure returns the value of the number as a global.
(this is not such a good practice, but we are stuck with it for now.) *)
begin (* brnotenumber *)
      note:=nil;
      numbered := false;
      number := 0; (* force number to zero if there
         is no number at all *)
      (* the next character is n or * depending on whether there are notes *)
      if thefile^ = 'n' then begin
         readln(thefile);
         theline := succ(theline);
         if thefile^ <> 'n' then begin
            skipstar(thefile);
            if not eoln(thefile) then begin
               if thefile^ = '#' then begin
                  numbered := true;
                  get(thefile); (* move past the number symbol *)
                  read(thefile,number);
               end
            end;
            repeat
               readln(thefile);
               theline := succ(theline)
            until thefile^ = 'n';
            readln(thefile);
            theline := succ(theline)
         end
         else begin
            readln(thefile);
            theline := succ(theline)
         end
      end
end; (* brnotenumber *)
(* end module book.brnotenumber version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.brnote *)
procedure brnote(var thefile: text; var theline: integer; var note: lineptr);
(* read note key *)
var
      newnote: lineptr; (* the new note *)
      previousnote: lineptr; (* the last line of the notes *)
begin (* brnote *)
      note:=nil;
      if thefile^ = 'n' then begin (* enter note *)
         readln(thefile);
         theline := succ(theline);
         if thefile^ <> 'n' then begin (* abort null note (n/n) *)
            getline(note);
            newnote:=note;
            while thefile^ <> 'n' do begin (* wait until end of note *)
               brline(thefile,theline,newnote);
               previousnote:=newnote;
               (* get next note *)
               getline(newnote^.next);
               newnote:=newnote^.next;
            end;
            (* last note was not used, so: *)
            clearline(newnote);
            previousnote^.next:=nil;
            readln(thefile);
            theline := succ(theline);
         end
         else begin
            readln(thefile);
            theline := succ(theline);
         end;
      end
end; (* brnote *)
(* end module book.brnote version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.brheader *)
procedure brheader(var thefile: text; var theline: integer; var hea: header);
(* read the header of a key. *)
begin
      with hea do begin
         readln(thefile); (* move past the object name - new definition 1999 Mar 13 *)
         theline := succ(theline);
{bbb}
         (* read key name *)
         brname(thefile,theline,keynam);

         (* read full name *)
         getline(fulnam);
         brline(thefile,theline,fulnam);

         (* read note key *)
         if readnumber then brnotenumber(thefile,theline,note)
                       else brnote(thefile,theline,note)
      end
end;
(* end module book.brheader version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.copyheader *)
procedure copyheader(fromhea: header; var tohea: header);
(* copy the header fromhea into tohea.  Note that the
linked objects are NOT copied, but merely pointed to. *)
begin
   tohea.keynam.letters := fromhea.keynam.letters;
   tohea.keynam.length  := fromhea.keynam.length;
   tohea.note           := fromhea.note;
   tohea.fulnam         := fromhea.fulnam;
end;
(* end module book.copyheader version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.brpiekey *)
procedure brpiekey(var thefile: text; var theline: integer; var pie: piekey);
(* read piece key, track the line number *)
begin
      with pie do begin
         brheader(thefile,theline,hea);
         brreanum(thefile,theline,mapbeg);
         brconfig(thefile,theline,coocon);
         brdirect(thefile,theline,coodir);
         brnumber(thefile,theline,coobeg);
         brnumber(thefile,theline,cooend);
         brconfig(thefile,theline,piecon);
         brdirect(thefile,theline,piedir);
         brnumber(thefile,theline,piebeg);
         brnumber(thefile,theline,pieend);
      end
end;
(* end module book.brpiekey version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.brdna *)
procedure brdna(var thefile: text; var theline: integer; var dna: dnaptr);
(* read in dna from thefile, track the line *)
(* note: if the dna were circularized, by linking the last dnastring
to the first, then the cleardna routine could not clear properly,
and would loop forever... there is no reason to do that, since a simple
mod function will allow one to access the circle. *)
var
      ch: char;
      workdna: dnaptr;
begin
      getdna(dna);
      workdna:=dna;
      ch:=getto(thefile,theline,['d']);
      readln(thefile);
      theline := succ(theline);
      read(thefile,ch); (* skipstar *)
      while (ch = '*') do
      begin
         read(thefile,ch); (* skip blank *)
         repeat
            read(thefile,ch);
            if ch in ['a','c','g','t'] then begin
               if workdna^.length=dnamax then begin
                  getdna(workdna^.next);
                  workdna:=workdna^.next
               end;
               workdna^.length:=succ(workdna^.length);
               workdna^.part[workdna^.length]:=chartobase(ch)
            end
         until eoln(thefile);
         readln(thefile); (* go to next line *)
         theline := succ(theline);
         read(thefile,ch); (* ch is either '*' or 'd' *)
      end;
      readln(thefile); (* read past the d *)
      theline := succ(theline);
end;
(* end module book.brdna version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.brpiece *)
procedure brpiece(var thefile: text; var theline: integer; var pie: pieceptr);
(* read in a piece, change theline to reflect the lines traversed *)
begin
{
   readln(thefile); (* move past the word 'piece' - new definition 1999 Mar 13 *)
   theline := succ(theline); (* BUG: was below! *)
bbb}
   brpiekey(thefile,theline,pie^.key);
   if numbered or (not skipunnum)
   then brdna(thefile,theline,pie^.dna);
   readln(thefile); (* move past the word 'piece' - new definition 1999 Mar 13 *)
   theline := succ(theline);
end;
(* end module book.brpiece version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.brinit *)
procedure brinit(var book: text; var theline: integer);
(* check that the book is ok to read, and
   set up the global variables for br routines *)
begin (* brinit *)
      (* halt if the book is bad (first word is 'halt') or the first
         character is not * *)
      reset(book);
      if not eof(book) then begin
         (* check for the date line *)
         if book^ <> '*' then begin
            if book^ <> 'h'
            then writeln(output, ' this is not the first line of a book:')
            else writeln(output, ' bad book:');
            write(output, ' ');

            while not (eoln(book) or eof(book)) do begin
               write(output, book^);
               get(book)
            end;
            writeln(output);
            halt
         end
      end
      else begin
         writeln(output, ' book is empty');
         halt
      end;

      (* initialize free storage *)
      freeline:=nil;
      freedna:=nil;

      readnumber:=true; (* usually we read in numbers for items *)
      number:=0; (* arbitrary value *)
      numbered:=false; (* the piece has no number (none yet read in) *)
      skipunnum:=false;
   theline := 1;
end; (* brinit *)
(* end module book.brinit version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* ************************************************************************ *)
(* end module package.brpiece version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.brorgkey *)
procedure brorgkey(var thefile: text; var theline: integer; var org: orgkey);
(* read organism key *)
begin
      with org do begin
{bbb}
         brheader(thefile,theline,hea);
         getline(mapunit);
         brline(thefile,theline,mapunit);
      end
end;
(* end module book.brorgkey version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.brchrkey *)
procedure brchrkey(var thefile: text; var theline: integer; var chr: chrkey);
(* read chromosome key *)
begin
      with chr do begin
{bbb}
         brheader(thefile,theline,hea);
         brreanum(thefile,theline,mapbeg);
         brreanum(thefile,theline,mapend);
      end
end;
(* end module book.brchrkey version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.getocp *)
procedure getocp(var thefile: text; var theline: integer;
                 var org: orgkey;   var orgchange, orgopen: boolean;
                 var chr: chrkey;   var chrchange, chropen: boolean;
                 var pie: pieceptr; var piechange, pieopen: boolean);
(* Get the next piece and its organism and chromosome keys.  The three change
variables indicate whether or not a new organism, chromosome or piece name
was found.  If a piece is not found, then pieopen will be false.  orgopen,
chropen and pieopen are used by getocp to tell when it has entered an
organism, chromosome or piece.  All booleans should be set to false
initially.  There should be one triplet for each book read.

It is important to initialize ALL variables, including pie:

   orgchange := false; orgopen := false;
   chrchange := false; chropen := false;
   piechange := false; pieopen := false;
   pie := nil;
   theline := 0;

1999 June 2  The book reading routines now treat data objects more
precisely.  Rather than test for eof, the endo of book occurs when pieopen is
returned as false.  A book reading loop now looks like this:

   repeat
      getocp(book, theline,
                   org, orgchange, orgopen,
                   chr, chrchange, chropen,
                   pie, piechange, pieopen);
      writeln(output,'pieopen: ',pieopen);
      if pieopen then begin
         writeln(output,'piece at line: ',theline:1);
      end;
   until not pieopen;
*)

var
      ch: char;
      newchr: chrkey;
      neworg: orgkey;
      newpie: pieceptr;
begin
      ch:='a';
      while not (ch in [' ','p']) do begin
         ch:=getto(thefile,theline,['o','c','p']);
         if ch <> ' ' then begin
            case ch of
               'o': if orgopen
                  then begin
                     readln(thefile);
                     (* move past the word 'organism'
                        - new definition 1999 Mar 13 *)
                     orgopen:=false (* close organism *)
                  end
                  else begin
                     brorgkey(thefile,theline,neworg);
                     if (neworg.hea.keynam.letters <> org.hea.keynam.letters) and
                        (neworg.hea.keynam.length <> org.hea.keynam.length)
                     then begin
{
writeln(output,'--------orgchanged!');
write  (output,'--------old org:"',   org.hea.keynam.letters);
writeln(output,              '" ',    org.hea.keynam.length:1);
write  (output,'--------new org:"',neworg.hea.keynam.letters);
writeln(output,              '" ',neworg.hea.keynam.length:1);
}
(*ccc*)
                        orgchange:=true;
                        copyheader(neworg.hea,org.hea);
                        (* move the mapunit over to the org! *)
                        org.mapunit := neworg.mapunit;
                        clearline(neworg.mapunit);

                     end
                     else orgchange:=false;
                     orgopen:=true;
                  end;
               'c': if chropen
                  then begin
                     readln(thefile);
                     (* move past the word 'chromosome'
                        - new definition 1999 Mar 13 *)
                     chropen:=false (* close chromosome *)
                  end
                  else begin
                     brchrkey(thefile,theline,newchr);
                     if (newchr.hea.keynam.letters <> chr.hea.keynam.letters) and
                        (newchr.hea.keynam.length <> chr.hea.keynam.length)
                     then begin
{
writeln(output,'--------chrchanged!');
write  (output,'--------old chr:"',   chr.hea.keynam.letters);
writeln(output,              '" ',    chr.hea.keynam.length:1);
write  (output,'--------new chr:"',newchr.hea.keynam.letters);
writeln(output,              '" ',newchr.hea.keynam.length:1);
}
                        chrchange:=true;
                        copyheader(newchr.hea,chr.hea);
                        (* move the map range over to the chr! *)
                        chr.mapbeg := newchr.mapbeg;
                        chr.mapend := newchr.mapend;
                     end
                     else chrchange:=false;
                     chropen:=true;
                  end;
               'p': if pieopen
                  then begin
                     pieopen:=false; (* close last piece *)
                     ch:='a' (* prevent falling out of the loop *)
                  end
                  else begin
                     new(newpie);
                     brpiece(thefile,theline,newpie);
                     if pie = nil
                     then piechange := true
                     else begin
                        if (newpie^.key.hea.keynam.letters
                         <> pie^.key.hea.keynam.letters) and
                           (newpie^.key.hea.keynam.length
                         <> pie^.key.hea.keynam.length)
                        then begin
                           piechange:=true;
                        end
                        else piechange:=false;
                     end;
                     pieopen:=true;
                     (* we always have to switch over to the new piece,
                        because although the name may be the same,
                        the DNA sequence could be different.  That is,
                        the book may contain two pieces with the same
                        name, and we want to be sure to search the
                        new one, not the old one. *)
                     if pie <> nil then begin
                        clearpiece(pie); (* save the links *)
                        dispose(pie); (* close up shop *)
                     end;
                     pie := newpie;
                  end
            end
         end
         else begin
            pieopen := false
         end
      end
end;
(* origin: search version = 6.39 *)
(* end module book.getocp version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* ************************************************************************ *)
(* end module package.getocp version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.getbase *)
function getbase(position: integer; pie: pieceptr):base;
(* Get a base from the position (internal coordinates) of the
   piece.  Protection is made against positions outside the piece.
In the case of circles it would be convenient to wrap around when
requests are off the end.  So the routine will do a modular
wrap for positions outside the range 1 to the length.
This is a new feature as of 2000 March 22.
*)
var
   workdna: dnaptr; (* pointer to the dna part of pie *)
   p: integer; (* current count of bases into the workdna *)
   spot: integer; (* the last base of the dna part *)
   thelength: integer; (* the length of the piece *)
begin
{
writeln(output,'NEW getbase: position=',position:1,'^^^^^^^^^^^^^^^^^^^^');
}
   (* handle cases of position out of range by circular wrapping *)
   thelength := piecelength(pie);
   while position < 1         do position := position + thelength;
   while position > thelength do position := position - thelength;

   workdna:=pie^.dna;
   p:=workdna^.length;
   while position > p do begin
{
writeln(output,' workdna^.length=',workdna^.length:1);
}
      workdna := workdna^.next;
      if workdna = nil then begin
        writeln(output,'error in function getbase!');
        halt
      end;
      p := p + workdna^.length;
   end;
{
writeln(output,'p=',p:1);
}
   if workdna = nil then begin
      writeln(output,'error in getbase: request off end of piece');
      halt
   end
   else begin
      spot := workdna^.length - (p-position);
{
writeln(output,'spot=',spot:1);
showdnasegment(output,workdna, spot);
}
      if (spot <= 0) then begin
         writeln(output,'error in getbase, spot (= ',spot:1, ') must be positive');
         halt
      end;
      if (spot > workdna^.length) then begin
         writeln(output,'error in getbase, spot (=',spot:1,
                        ') must be less than length (=',workdna^.length:1,')');
         halt
      end;
{
writeln(output,'base = ', workdna^.part[spot]);
}
      getbase:=workdna^.part[spot]
   end
end;
(* end module book.getbase version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.stepbase *)
function stepbase(startdna: dnaptr; var dna: dnaptr; var d: dnarange): base;
(* advance d by one base in dna and then return the base at the new d.
   (this means that one should initialize d to zero)
   if we go past the last base, we restart at startdna.
   note: d is not the number of the base... it is used as a
   record for stepbase.  do not mess with it, and do not use it to find
   out what base you are on.  use a separate counter.  *)
begin
      if (d=dnamax) or (d=dna^.length) then begin
         d:=1;
         dna:=dna^.next;
         if dna=nil then dna:=startdna
      end
      else d:=succ(d);
      stepbase:=dna^.part[d]
end;
(* end module book.stepbase version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.iwname *)
procedure iwname(var thefile: text; thename: name);
(* write the name to the file *)
var   c: integer;
begin
      for c:=1 to thename.length
         do write(thefile, thename.letters[c])
end;
(* end module book.iwname version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.iworg *)
procedure iworg(var afile: text; org: orgkey);
(* write an organism specification.  no writeln is done to allow
   write orgchr to do this. *)
begin
      write(afile,'organism ');
      iwname(afile,org.hea.keynam);
      write(afile,';');
end;
(* end module book.iworg version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.iwchr *)
procedure iwchr(var afile: text; chr: chrkey);
(* write an chromosome specification.  no writeln is done to allow
   write orgchr to do this. *)
begin
      write(afile,'chromosome ');
      iwname(afile,chr.hea.keynam);
      write(afile,';');
end;
(* end module book.iwchr version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.iwpie *)
procedure iwpie(var afile: text; pie: piekey);
(* write a piece specification *)
begin
      write(afile,'piece ');
      iwname(afile,pie.hea.keynam);
      writeln(afile,';');
end;
(* end module book.iwpie version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.iworgchr *)
procedure iworgchr(var afile: text;
                       org: orgkey; orgchange, orgopen: boolean;
                       chr: chrkey; chrchange, chropen: boolean);
(* write both organism and chromosome specifications, based
on whether the organism or chromosome changed (orgchange and chrchange)
and whether they are currently open (orgopen, chropen).
See getocp in the br routines. *)
begin
   if orgchange and orgopen then iworg(afile,org);
   if orgchange and chrchange and orgopen and chropen then write(afile,' ');
   if chrchange and chropen then iwchr(afile,chr);
   if (orgchange and orgopen) or (chrchange and chropen) then writeln(afile)
end;
(* end module book.iworgchr version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.iwget *)
procedure iwget(var afile: text; pie: pieceptr;
                    fromplace, pieceplace, toplace: integer;
                    flip: boolean; insttype: integer;
                    carriagereturn: boolean);
(* print a get delila instruction in the orientation of pie,
   from fromplace to toplace pieceplace. +/- direction;
If flip is false, the piece direction is as on the piece,
if it is true, the it is the opposite direction.
insttype: instruction type.
          insttype=1 means the form   get from p  -/+f to same +/-t dir +/-;
          insttype=2 means the form   get from p1 -/+f to p2   +/-t dir +/-;
          where p, p1 and p2 are locations
carriagereturn: if true, add a carriage return to the end of the line.
 *)
procedure iwposition(relative: integer; sameallowed: boolean);
(* write the <POSITION> *)
procedure iwrelative(relative: integer);
begin
           if relative>=0 then write(afile,' +', relative:1)
      else if relative< 0 then write(afile,' ',relative:1)
end;
begin (* iwposition *)
      if (insttype = 1) and sameallowed
      then write(afile,' same')
      else write(afile,' ',pieceplace:1);
      case pie^.key.piedir of
         plus:  iwrelative(+relative);
         minus: iwrelative(-relative);
      end
end;
begin (* iwget *)
      write(afile,'get from');
      iwposition(fromplace, false);
      write(afile,' to');
      iwposition(toplace, true);
      write(afile,' direction');
      case pie^.key.piedir of
                  dirhomologous, (* handle case, may not be right *)
         plus:  case flip of
                  false: write(afile,' +');
                   true: write(afile,' -');
                end;
                  dircomplement, (* handle case, may not be right *)
         minus: case flip of
                  false: write(afile,' -');
                   true: write(afile,' +');
                end
      end;         
      write(afile,';');
      if carriagereturn then writeln(afile);
end;
(* end module book.iwget version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.iwget2 *)
procedure iwget2(var afile: text; pie: pieceptr;
                    fromplace, place1, toplace, place2: integer;
                    flip: boolean;
                    carriagereturn: boolean);
(* print a get Delila instruction in the orientation of pie,
The form of the instructions is:
  get from place1 +/-fromplace
        to place2 +/-toplace direction +/-;
If flip is false, the piece direction is as on the piece,
if it is true, the it is the opposite direction.
carriagereturn: if true, add a carriage return to the end of the line.
*)
procedure iwposition(place, relative: integer);
   procedure iwrelative(relative: integer);
      begin
              if relative>=0 then write(afile,' +', relative:1)
         else if relative< 0 then write(afile,' ',relative:1)
      end;
   begin (* iwposition *)
         write(afile,' ',place:1);
         case pie^.key.piedir of
            plus:  iwrelative(+relative);
            minus: iwrelative(-relative);
         end
   end;
begin (* iwget2 *)
      write(afile,'get from');
      iwposition(place1,fromplace);
      write(afile,' to');
      iwposition(place2,toplace);
      write(afile,' direction');
      case pie^.key.piedir of
                  dirhomologous, (* handle case, may not be right *)
         plus:  case flip of
                  false: write(afile,' +');
                   true: write(afile,' -');
                end;
                  dircomplement, (* handle case, may not be right *)
         minus: case flip of
                  false: write(afile,' -');
                   true: write(afile,' +');
                end
      end;
      write(afile,';');
      if carriagereturn then writeln(afile);
end;
(* end module book.iwget2 version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module book.iwcombk *)
procedure iwcombk(var book, afile: text);
(* make a comment in the file that says the name of the book *)
begin
      write(afile,'(* ');
      reset(book);
      if copylines(book, afile, 1) = 0 then begin
         writeln(output,
            ' book is empty, can not write comment for instructions');
         halt
      end;
      writeln(afile,'*)');
end;
(* end module book.iwcombk version = 7.81; {of delmod.p 2017 Sep 14} (@ *)

(* begin module package.datetime *)
(* ************************************************************************ *)
(* begin module getdatetime *)
procedure getdatetime(var adatetime: datetimearray);
(* get the date and time into a single array from the system clock.
   adatetime contains the date:
       1980/06/09 18:49:11
         ye mo da ho mi se
(year, month, day, hour, minute, second).
As of 2000 February 18, the Sun Pascal compiler requires a formatting
statement.  This statement allows the date to be generated in this
standard Delila format in a single call.  Information about the
formatting statement is available on the manual page for date in Unix.
If a computer does not have this method, see the 'oldgetdatetime' routine
in delmod.p (https://alum.mit.edu/www/toms/delila/delmod.html)
for some conversion code.

GPC Functions:
function  GetUnixTime (var MicroSecond : Integer) : UnixTimeType;

http://agnes.dida.physik.uni-essen.de/~gnu-pascal/gpc_109.html#SEC109

7.10.8 Date And Time Routines 

procedure GetTimeStamp (var t : TimeStamp); 
function Date (t : TimeStamp) : packed array [1 .. DateLength] of Char; 
function Time (t : TimeStamp) : packed array [1 .. TimeLength] of Char; 

DateLength and TimeLength are implementation dependent constants. 

GetTimeStamp (t) fills the record `t' with values. If they are valid, the Boolean
flags are set to True. 

TimeStamp is a predefined type in the Extended Pascal standard. It may be
extended in an implementation, and is indeed extended in GPC. For the full
definition of `TimeStamp', see section 8.255 TimeStamp. 

*)
var
    t: TimeStamp; 

(* begin module pluckdigit *)
function pluckdigit(number, logplace:integer): char;
(* return the digit at the place value ('logplace') position of number.
example:
      pluckdigit(13625, 3) = 3
      pluckdigit(13625, 4) = 1
This routine was taken from module numberdigit in prgmod.p, but is
modified so as not to give the sign.  Instead it gives zeros
above the digits.  'myabsolute' replaced 'absolute', which
is apparently a keyword for GPC.  The name is kept for
to keep the code looking similar to its origin.
*)
var
      place: integer; (* the exponent of logplace *)
      count: integer; (* used to make place *)
      myabsolute: integer; (* the absolute value of number *)
      acharacter: char; (* the character to be returned *)
procedure digit;
(* extract a digit at the place position *)
var
      tenplace: integer; (* ten times place *)
      z: integer; (* an intermediate value *)
      d: integer; (* the digit extracted *)
begin (* digit *)
      tenplace:=10*place;
      z:=myabsolute-((myabsolute div tenplace)*tenplace);
      if place = 1
         then d:=z
         else d:= z div place;
      case d of
         0: acharacter:='0';
         1: acharacter:='1';
         2: acharacter:='2';
         3: acharacter:='3';
         4: acharacter:='4';
         5: acharacter:='5';
         6: acharacter:='6';
         7: acharacter:='7';
         8: acharacter:='8';
         9: acharacter:='9';
      end
end; (* digit *)
begin (* pluckdigit *)
      place:=1;
      for count:=1 to logplace do place:=10*place;
      if number=0 then begin
         acharacter:='0'
      end
      else begin
         myabsolute:=number;
         if myabsolute >= place
         then digit
         else acharacter := '0'
      end;
      pluckdigit:=acharacter
end; (* pluckdigit *)
(* end module pluckdigit *)

begin

 (* according to:
http://agnes.dida.physik.uni-essen.de/~gnu-pascal/gpc_109.html#SEC109
*)
   GetTimeStamp(t); 

(* Predefined time stamp:
http://agnes.dida.physik.uni-essen.de/~gnu-pascal/gpc_389.html#SEC389
  TimeStamp = {@@packed} record
      DateValid,
      TimeValid   : Boolean;
      Year        : Integer;
      Month       : 1 .. 12;
      Day         : 1 .. 31;
      DayOfWeek   : 0 .. 6;  { 0 means Sunday }
      Hour        : 0 .. 23;
      Minute      : 0 .. 59;
      Second      : 0 .. 61; { to allow for leap seconds }
      MicroSecond : 0 .. 999999
    end;
*)

  with t do begin
    if TimeValid then begin
{
       writeln(output,'valid time');
       writeln(output,'year   =',year:4);
       writeln(output,'month  =',month:2);
       writeln(output,'day    =',day:2);
       writeln(output,'hour   =',hour:2);
       writeln(output,'minute =',minute:2);
       writeln(output,'second =',second:2);
}
       adatetime := 'year/mm/dd hh:mm:ss';
       adatetime[ 1] := pluckdigit(Year,3);
       adatetime[ 2] := pluckdigit(Year,2);
       adatetime[ 3] := pluckdigit(Year,1);
       adatetime[ 4] := pluckdigit(Year,0);
       adatetime[ 6] := pluckdigit(Month,1);
       adatetime[ 7] := pluckdigit(Month,0);
       adatetime[ 9] := pluckdigit(Day,1);
       adatetime[10] := pluckdigit(Day,0);
       adatetime[12] := pluckdigit(Hour,1);
       adatetime[13] := pluckdigit(Hour,0);
       adatetime[15] := pluckdigit(Minute,1);
       adatetime[16] := pluckdigit(Minute,0);
       adatetime[18] := pluckdigit(Second,1);
       adatetime[19] := pluckdigit(Second,0);
    end
    else begin
       writeln(output,'getdatetime: invalid time!');
       halt;
    end
  end;
{
   Sun compiler method:
   date(adatetime,'%Y/%m/%d %H:%M:%S');
}
end;
(* end module getdatetime version = 1.17; (@ of timegpc.p 2002 Oct 9 *)

(* begin module readdatetime *)
procedure readdatetime (var thefile: text; var adatetime: datetimearray);
(* read the date and time from the file.  It must have this format:
123456789 123456789
         1         
1980/06/09 18:49:11
 *)
(*
2001 Jan 2: thefile is not written to!  Force 4 digit years.
2000 Oct 11: upgraded so that the p2c compiler does not object
to writing out the adatetime; added checks for the date.
*)
var
   index: integer; (* to the udatetime *)
   (* the following is an unpacked date time array, to avoid reading
   into a packed array.  reading into a packed array is not transportable *)
   udatetime: array[1..datetimearraylength] of char;
begin
   for index:=1 to datetimearraylength do read(thefile,udatetime[index]);
   pack(udatetime, 1, adatetime);
   if (adatetime[3]='/') and (adatetime[12]=':') then begin
      writeln(output,'You have an old datetime (only 2 year digits): ');
      for index:=1 to datetimearraylength do write(output,adatetime[index]);
      writeln(output);
      writeln(output,'Convert your database to 4 digit years.');
      halt;
   end;
   (* check the adatetime format.  Note that further checks for the
      other positions in the array could be done to be sure that
      they are numbers.  But this should be pretty good. *)
   if (adatetime[ 5]<>'/') or
      (adatetime[ 8]<>'/') or
      (adatetime[14]<>':') or
      (adatetime[17]<>':')
   then begin
      writeln(output,'readdatetime: bad date time read:');
      for index:=1 to datetimearraylength do write(output,adatetime[index]);
      writeln(output);
      halt
   end;
end;
(* end module readdatetime version = 1.17; (@ of timegpc.p 2002 Oct 9 *)

(* begin module writedatetime *)
procedure writedatetime(var thefile: text; adatetime: datetimearray);
(* expand the date and time out and print in the file *)
var
      index: integer; (* index of datetime *)
begin
      for index:=1 to datetimearraylength do write(thefile,adatetime[index])
end;
(* end module writedatetime version = 1.17; (@ of timegpc.p 2002 Oct 9 *)

(* begin module timeseed *)    
(* Read the computer date and time.  Reverse the order of the digits and put
a decimal point in front.  This gives a fraction between zero and one that
varies quite quickly, and is always unique (if the computer has sufficient
accuracy).  It is to be used as a seed to a random number generator.  This
has the nice property that the seed changes every second and does not repeat
for thousands of years!  *)
procedure addtoseed(var seed, power: real; c: char);
(* add the digit represented by c to the seed at the power position *)
var
      n: integer; (* the character represented by c *)
begin (* addtoseed *)
      power := power/10;
{
writeln(output,'addtoseed, c = ',c);
writeln(output,'addtoseed, ord(c) = ',ord(c));
}
   n := ord(c) - ord('0');
   if (n < 0) or (n > 9) then begin
      writeln(output,'timeseed: error in datetime');
      writeln(output,'it contains "',c,'" which is not a number.');
      writeln(output,'The getdatetime routine must be fixed.');
      halt;
   end;
   seed := seed + power*n
end; (* addtoseed *)

procedure makeseed(adatetime: datetimearray; var seed: real);
(* convert adatetime to a real number in seed, reversed order
Here is the standard adatetime format:
123456789 123456789
         1         
1980/06/09 18:49:11
*)
var
   power: real; (* a digit of the seed such as 0.01 *)
begin
   seed := 0.0;
   power := 1.0;
   addtoseed(seed, power, adatetime[19]);
   addtoseed(seed, power, adatetime[18]);
   (* : *)
   addtoseed(seed, power, adatetime[16]);
   addtoseed(seed, power, adatetime[15]);
   (* : *)
   addtoseed(seed, power, adatetime[13]);
   addtoseed(seed, power, adatetime[12]);
   (*   *)
   addtoseed(seed, power, adatetime[10]);
   addtoseed(seed, power, adatetime[ 9]);
   (* / *)
   addtoseed(seed, power, adatetime[ 7]);
   addtoseed(seed, power, adatetime[ 6]);
   (* / *)
   addtoseed(seed, power, adatetime[ 4]);
   addtoseed(seed, power, adatetime[ 3]);
   addtoseed(seed, power, adatetime[ 2]);
   addtoseed(seed, power, adatetime[ 1]);
end;

procedure orderseed(adatetime: datetimearray; var seed: real);
(* convert adatetime to a real number in seed, normal order *)
var
   power: real; (* a digit of the seed such as 0.01 *)
begin
   seed := 0.0;
   power := 1.0;
   addtoseed(seed, power, adatetime[ 3]);
   addtoseed(seed, power, adatetime[ 4]);
   addtoseed(seed, power, adatetime[ 6]);
   addtoseed(seed, power, adatetime[ 7]);
   (* / *)
   addtoseed(seed, power, adatetime[ 9]);
   addtoseed(seed, power, adatetime[10]);
   (* / *)
   addtoseed(seed, power, adatetime[12]);
   addtoseed(seed, power, adatetime[13]);
   (*   *)
   addtoseed(seed, power, adatetime[15]);
   addtoseed(seed, power, adatetime[16]);
   (* : *)
   addtoseed(seed, power, adatetime[18]);
   addtoseed(seed, power, adatetime[19]);
end;

procedure timeseed(var seed: real);
(* read the computer date and time.  reverse the order of the digits
and put a decimal point in front.  this gives a fraction between
zero and one that varies quite quickly, and is always unique (if the
computer has sufficient accuracy).  it is to be used as a seed to
a random number generator. *)
var
   adatetime: datetimearray; (* a date and time *)
begin (* timeseed *)
   getdatetime(adatetime);
{
writeln(output,'timeseed: adatetime: ',adatetime);
}
   makeseed(adatetime, seed);
end; (* timeseed *)
(* end module timeseed version = 1.17; (@ of timegpc.p 2002 Oct 9 *)
 
(*[[*)
(* begin module limitdate *)
procedure limitdate(a,b,c,d: char; limitdatetime: datetimearray);
(* test whether the current time is before the limit.
If it is later, halt the program *)
var
   adatetime: datetimearray; (* a date and time *)
   Dday: real; (* the critical day *)
   now: real;  (* this very moment *)
begin
   getdatetime(adatetime);
{
writeln(output,'adatetime:',adatetime);
writeln(output,'adatetime[1]: ', adatetime[1]);
writeln(output,'adatetime[2]: ', adatetime[2]);
writeln(output,'adatetime[3]: ', adatetime[3]);
writeln(output,'adatetime[4]: ', adatetime[4]);
writeln(output,'adatetime[5]: ', adatetime[5]);
}

   orderseed(adatetime, now);

   if
      (limitdatetime[1] <> ' ') or
      (limitdatetime[2] <> ' ') or
      (limitdatetime[3] <> ' ') or
      (limitdatetime[4] <> ' ') then halt;

   limitdatetime[1] := a;
   limitdatetime[2] := b;
   limitdatetime[3] := c;
   limitdatetime[4] := d;

   orderseed(limitdatetime, Dday);

{
   writeln(output,'now:  ',now:20:10);
   writeln(output,'Dday: ',Dday:20:10);
}
   if now > Dday then begin
       writeln(output,'This program expired on ',limitdatetime);
       writeln(output,'See: https://alum.mit.edu/www/toms/walker/contacts.html');
       halt
   end
end;
(* end module limitdate version = 1.17; (@ of timegpc.p 2002 Oct 9 *)
(*]]*)

(* ************************************************************************ *)
(* end module package.datetime version = 10.57; (@ of lister.p 2004 Jul 22 *)

(* begin module package.interact *)
(* ************************************************************************ *)
(* begin module interact.prompt *)
procedure prompt(var afile: text);
(* prompt a file.
   the prompt is sent to the output file, and
   a line is read into the pascal line buffer.
   (for the cyber system this means to readln afile.)
   guarantee no bomb *)
begin (* prompt *)
      if eof(afile) then reset(afile);
      readln(afile);
end; (* prompt *)
(* end module interact.prompt version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

(* begin module interact.readchar *)
procedure readchar(var afile:text; var ch: char);
(* read a character from afile, guarantee no bomb *)
begin (* readchar *)
      if eof(afile) then prompt(afile);
      read(afile,ch);
(*writeln(output,'"',ch,'"') *)
end; (* readchar *)
(* end module interact.readchar version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

(* begin module clearstring *)
procedure clearstring(var ribbon: string);
(* empty the string *)
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; (* 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. *)
begin (* initializestring *)
   clearstring(ribbon);
   ribbon.next := nil;
end; (* initializestring *)
(* end module clearstring version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

(* begin module interact.readstring *)
procedure readstring(var afile: text; var line: string);
(* read in a string from afile, protect against bombing *)
var
      index: integer; (* for line *)
      cha: char; (* a character read in *)
      done: boolean; (* used for removing trailing blanks from the line *)
      acceptable: boolean; (* was the line typed short enough? *)
begin (* readstring *)
      with line do begin
         repeat
            clearstring(line);
            prompt(afile);
            index := 0; (* we now count characters *)
            while (not eoln(afile)) and (index < maxstring) do begin
               index:=succ(index);
               readchar(afile,cha);
               line.letters[index]:=cha
            end;
            if not eoln(afile)
            then begin
               writeln(output, 'type lines shorter than ', (maxstring+1): 1,
                               ' characters.  please retype the line...');
               acceptable := false
            end
            else acceptable := true
         until acceptable;

         length := index;
         if length > 0 then begin
            done := false;
            repeat (* remove blanks from the line.  note that a while loop
                  can not be used because one must avoid letters[0], since that
                  position does not exist... *)
               if letters[length] = ' '
               then length := pred(length)
               else done := true;

               if length = 0 then done := true
            until done
         end;

         if length > 0 then current := 1
                       else current := 0
      end
end; (* readstring *)
(* end module interact.readstring version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

(* begin module interact.figurestring *)
procedure figurestring(
      var line: string; (* a string of characters to figure out *)
      var first: integer; (* first found non-blank character in the line *)
      var last:  integer; (* last character before a blank after first *)
      var whzat: char; (* what the token is *)
      var c:     char;  (* the first character of the token *)
      var i:     integer; (* integer value of token if it is integer; or 0 *)
      var r:     real); (* the real value if it is real; or 0.0 *)
(* figurestring figures out the tokens in a string.
 it recognizes words, integers, reals and poorly formed numbers.
 you can easily use it to parse lines.

our goal is to figure out what thing is on a string.
start looking at the current place on the line.
first and last are the first 'token' in line after start.
the current place is updated to the letter after last.
the thing found is described by the value of whzat:
      'c': character (when the token does not begin with a digit, '+', or '-')
      'i': integer
      'r': real
      ' ': blank line
      'g': garbage, cannot figure it out
and the value of the thing found is the appropriate variable *)
var
      numbers: set of '0'..'9';
      sign: integer; (* sign of a number *)
      numberstart: integer; (* the point a number starts,
                               beyond its sign, if any *)
      point: integer; (* location of decimal point *)
      power: integer; (* of 10 representing a place value in the number *)

      l: integer; (* an index for dissecting numbers *)

function figureinteger(first,last:integer):integer;
(* figure the integer in the token *)
var
      i: integer; (* index *)
      sum, increment: integer;
begin (* figureinteger *)
      power:=1; (* start at ones place *)
      sum:=0; (* start sum at zero *)
      for i:=last downto first do begin
         case line.letters[i] of
            '0': increment:=0;
            '1': increment:=1;
            '2': increment:=2;
            '3': increment:=3;
            '4': increment:=4;
            '5': increment:=5;
            '6': increment:=6;
            '7': increment:=7;
            '8': increment:=8;
            '9': increment:=9
         end;
         sum:=sum+power*increment;
         power:=power*10
      end;
      figureinteger:=sum
end; (* figureinteger *)

begin (* figurestring *)
      numbers:=['0','1','2','3','4','5','6','7','8','9'];
(* c:=' '; i:=0; r:=0.0; do not affect these variables unless necessary *)
      point:=0;
      whzat := '.'; (* assume that we have someting to work on *)
      (* now to see if that is true: *)

      with line do if (length = 0) or (current < 1) or (current > length)
      then whzat := ' '
      else begin
         (* figure out where the first token is in the line *)
         first:=line.current;
         while (line.letters[first]=' ') and (first < line.length)
            do first:=succ(first);
         if (first = line.length) and (line.letters[first] = ' ')
            then whzat := ' ';
      end;

      if whzat <> ' ' then begin
         last:=first;
         while (line.letters[last]<>' ') and (last < line.length)
            do last:=succ(last);
         if line.letters[last] = ' ' then last := pred(last);

         (* the token is between inclusive first and last *)

         c:=line.letters[first];
         if (c in numbers) or (c in ['+','-']) then begin
            if c in ['+','-'] then begin
               case c of
                  '+': sign:=+1;
                  '-': sign:=-1;
               end;
               numberstart:=succ(first)
            end
            else begin
               sign:=+1;
               numberstart:=first
            end;

            whzat:='i';
            for l:=numberstart to last do begin
               if not(line.letters[l] in numbers)
               then if line.letters[l]='.' (* we found a period *)
                  then if whzat='i' (* if so far it is numbers *)
                     then begin
                        whzat:='r'; (* it is actually real *)
                        point:=l
                     end
                     else whzat:='g' (* it is a second '.', ie garbage *)
                  else whzat:='g' (* it is garbage *)
            end;
            (* if it is only numbers, it is integer *)

            (* build number *)
            (* if it ends in a period, it is integer *)
            if (whzat = 'r') and (point = last) then whzat:='i';

            if whzat = 'i' then begin
               if point = last (* had an ending decimal point *)
               then i:=sign * figureinteger(numberstart,pred(last))
               else i:=sign * figureinteger(numberstart,last);
               r:=i
            end
            else if whzat = 'r' then begin
               i:=figureinteger(numberstart,point-1);
               r:=sign * (i+figureinteger(point+1,last)/power);
               i:=sign * i
            end
         end
         else begin
            whzat:='c';
         end;

         (* move the start to just beyond the last character of the token *)
         line.current:=succ(last)
      end
end; (* figurestring *)

(* end module interact.figurestring version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

(* ************************************************************************ *)
(* end module package.interact version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

(* begin module copystring *)
procedure copystring(a: string; var b: string);
(* copy string a to b *)
var l: integer; (* index to the string *)
begin
   b.length := a.length;
   for l := 1 to a.length do b.letters[l] := a.letters[l]
end;
(* end module copystring version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

(* begin module nametostring *)
procedure nametostring(a: name; var b: string);
(* copy delila type name a to string b *)
var l: integer; (* index to the strings *)
begin
   b.length := a.length;
   for l := 1 to a.length do b.letters[l] := a.letters[l]
end;
(* end module nametostring *)

(* begin module package.interact.gets *)
(* ************************************************************************ *)
(* begin module interact.nostring *)
function nostring(var buffer: string): boolean;
(* true if there are no characters in the rest of the buffer;
   false if there are characters.
   also, if there is no buffer, then buffer.length is set to 0 *)
var   answer: boolean; (* the answer returned *)
procedure kill;
(* destroy the line *)
begin (* kill *)
      answer := true; (* blood and gore *)
      clearstring(buffer) (* total death *)
end; (* kill *)
begin (* nostring *)
      with buffer do begin
         if length > 0 then begin
            if length < maxstring
            then while (letters[current] = ' ') and (current < length)
               do current := succ(current);

            if current <= maxstring
            then if letters[current] = ' '
                 then kill
                 else answer := false
            else kill
         end
         else kill
      end;
      nostring := answer
end; (* nostring *)
(* end module interact.nostring version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

(* 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 = 4.73; (@ of prgmod.p 2004 Jul 22 *)

(* begin module interact.flagstring *)
procedure flagstring(var afile: text; var buffer: string);
(* flag an error in the buffer at the current place, and clear the buffer *)
begin (* flagstring *)
      with buffer do length := current; (* chop off the rest of the buffer *)
      writestring(afile, buffer); (* show the buffer *)
      write(afile,'? ');
      clearstring(buffer)
end; (* flagstring *)
(* end module interact.flagstring version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

(* begin module interact.getchar *)
procedure getchar(var afile: text; var buffer: string;
                  var cha: char; var gotten: boolean);
(* get a character from the buffer, or refill the buffer and
   let the calling program figure out whether the buffer has non
   blank characters in it. *)
var
      (* variables for calling figurestring: *)
      first: integer;
      last: integer;
      what: char;
      int: integer;
      rea: real;
begin (* getchar *)
      if buffer.length = 0 then begin
         gotten:=false;
         readstring(afile,buffer)
      end
      else begin
         figurestring(buffer,first,last,what,cha,int,rea);
         gotten := (what <> ' ')
      end
end; (* getchar *)
(* end module interact.getchar version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

(* begin module interact.getinteger *)
procedure getinteger(var afile: text; var buffer: string;
                     var int: integer; var gotten: boolean);
(* get the integer int from the buffer or interactive file afile *)
var
      (* variables for calling figurestring: *)
      first: integer;
      last: integer;
      what: char;
      cha: char;
      rea: real;
begin (* getinteger *)
      if buffer.length = 0
      then begin
         gotten:=false;
         readstring(afile,buffer)
      end
      else begin
         figurestring(buffer,first,last,what,cha,int,rea);

         if what <> 'i'
         then begin
            flagstring(output,buffer);
            writeln(output,' please type an integer');
            gotten:=false
         end
         else gotten:=true
      end
end; (* getinteger *)
(* end module interact.getinteger version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

(* begin module interact.getreal *)
procedure getreal(var afile: text; var buffer: string;
                  var rea: real; var gotten: boolean);
(* get the real rea from the buffer or interactive file afile
integer values are also accepted. *)
var
      (* variables for calling figurestring: *)
      first: integer;
      last: integer;
      what: char;
      cha: char;
      int: integer;
begin (* getreal *)
      if buffer.length = 0
      then begin
         gotten:=false;
         readstring(afile,buffer)
      end
      else begin
         figurestring(buffer,first,last,what,cha,int,rea);

         if not (what in ['r', 'i'])
         then begin
            flagstring(output,buffer);
            writeln(output,' please type a real number');
            gotten:=false
         end
         else gotten:=true
      end;

      (* handle integers *)
      if what = 'i' then rea := int
end; (* getreal *)
(* end module interact.getreal version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

(* begin module interact.token *)
procedure token(var buffer, atoken: string; var gotten: boolean);
(* get a token from the buffer *)
var
      (* variables for calling figurestring: *)
      first: integer;
      last: integer;
      what: char;
      cha: char;
      int: integer;
      rea: real;
      index: integer; (* to the buffer *)
begin
      figurestring(buffer,first,last,what,cha,int,rea);

      if what = ' '
      then gotten := false
      else begin
         clearstring(atoken);
         for index := first to last
            do atoken.letters[index-first+1] := buffer.letters[index];
         atoken.length := last - first + 1;
         atoken.current := 1;

         gotten:=true
      end
end;
(* end module interact.token version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

(* begin module interact.gettoken *)
procedure gettoken(var afile: text; var buffer: string;
                   var atoken: string; var gotten: boolean);
(* get a token from the buffer or interactive file afile *)
begin (* gettoken *)
      if buffer.length = 0
      then begin
         gotten:=false;
         readstring(afile,buffer)
      end
      else token(buffer,atoken,gotten)
end; (* gettoken *)
(* end module interact.gettoken version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

(* ************************************************************************ *)
(* end module package.interact.gets version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

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

(* begin module writename *)
procedure writename(var afile: text; s: name);
(* write the string s to afile *)
var i: integer; (* index to the string s *)
begin
   with s do for i := 1 to length do write(afile,letters[i]);
end;
(* end module writename version = 8.47; (@ of lister.p 1997 May 7 *)

{
(* begin module clearquotestring *)
procedure clearquotestring(var s: quotestring);
(* clear the quotestring s.  If there is a program
error, an exclamation mark will appear on the output. *)
var i: integer; (* index to s *)
begin
   for i := 1 to maxquotestring do s.string[i] := '!';
   s.length := 0;
end;
(* end module clearquotestring version = 8.47; (@ of lister.p 1997 May 7 *)
}

(* begin module detabstring *)
procedure detabstring(var s: string);
(* Tabs are converted to blanks.  This allows the tab to be a special "solid"
kind of blank, but the string written won't vary with position on the line
depending on the tabstops. *)
const
   tabcharacter = 9; (* the ordinal of the tab character *)
var
   c: char; (* a character write out *)
   i: integer; (* index to the string s *)
begin
   with s do for i := 1 to length do begin
      c := letters[i];
      if ord(c) = tabcharacter then c := ' '; (* tab converted to blank *)
      letters[i] := c
   end;
end;
(* end module detabstring version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

(* begin module writequotestring *)
procedure writequotestring(var afile: text; s: string);
(* Write the string s to afile, surrounded by double quote marks.  Tabs are
converted to blanks.  This allows the tab to be a special "solid" kind of
blank, but the string written won't vary with position on the line depending
on the tabstops. *)
begin
   detabstring(s);
   write(afile, '"');
   writestring(afile, s);
   write(afile, '"');
end;
(* end module writequotestring version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

(* begin module writeadefinition *)
procedure writeadefinition(var afile: text; d: definitiontypeptr);
(* write the definition d to file afile *)
var i: integer; (* index to the string s *)
begin
   with d^ do begin
      write(afile,'define');
      write(afile,' '); writequotestring(afile,nametag);
      write(afile,' '); writequotestring(afile,background);
      write(afile,' '); writequotestring(afile,negparts);
      write(afile,' '); writequotestring(afile,posparts);
      for i := 1 to marks do write(afile,' ',locations[i]:1:1);
      writeln(afile);
{
      (* don't do this here - it's awkward if one wants
         to just write the parts above.  Just rember to write the
      matrix after the definition! *)
      writeln(afile);
      putriblmatrix(afile, matrix);
}
   end;
end;
(* end module writeadefinition version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

(* begin module clearfeature *)
procedure clearfeature(var f: featuretypeptr);
(* clear out the feature f *)
begin
   with f^ do begin
      initializestring(id);
      coordinate := 0;
      orientation := +1;
      initializestring(nametag);
      initializestring(othertag);
      definition := nil;
(*[[*)
      Ri := 0.0;
      Z := 0.0;
      probability := 0.0;
(*]]*)
      Tparam := ' ';
      Aparam := 0.0;
      Bparam := 0.0;
      Cparam := 0.0;
      Dparam := 0.0;
      evencoordinate := 0.0;
      unsatisfied := true;
      fromrange := 0;
      torange := 0;
      number := 0;
      desiredline := 0;
      next := nil;
   end
end;
(* end module clearfeature version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

(* begin module writeafeature *)
procedure writeafeature(var afile: text; f: featuretypeptr);
(* write the feature f to file afile in the format:
@ K00000 36 +1 "nametag" "othertag" 12.200338 -0.473212 0.318031
where the last three numbers are the Ri, Z and p.
Do not carriage return.
*)
begin
   with f^ do begin
      write(afile,'@ ');
      writestring(afile,id);
      write(afile,' ',coordinate:1:1);
      write(afile,' ');
      if orientation > 0 then write(afile,'+',orientation:1)
                         else write(afile,    orientation:1);
      write(afile,' ');
      writequotestring(afile,nametag);
      write(afile,' ');
      writequotestring(afile,othertag);
(*[[*)
{zzz}
      if definition^.matrix <> nil then begin
         write(afile,' ',Ri:infofield:infodecim);
         write(afile,' ',Z:infofield:infodecim);
         write(afile,' ',probability:infofield:infodecim);
      end
      else begin (* parameters, new as of 2004 July 17 *)
         write(afile,' ',Tparam);
         write(afile,' ',Aparam:infofield:infodecim);
         write(afile,' ',Bparam:infofield:infodecim);
         write(afile,' ',Cparam:infofield:infodecim);
         write(afile,' ',Dparam:infofield:infodecim);
      end

(* 2004 July 17: I went ahead and removed this.  As far as I know
it is NOT used in feature files:
      if evencoordinate <> coordinate
      then write(afile,' ',evencoordinate:1:1);
*)
(*]]*)
   end
end;
(* originally from 9.29; (@ of lister.p 1999 July 5 *)
(* end module writeafeature version = 4.73; (@ of prgmod.p 2004 Jul 22 *)

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

(* end of procedure modules **********************************************)

procedure addtoquotestring(var q: string; c: char);
(* add the character c to the end of q *)
(* 2004 Jul 18: no longer using quotestring. *)
{zzz}
begin
   with q do begin
      if length < maxquotestring then begin
         length := succ(length);
         letters[length] := c;
      end
      else begin
         writeln(output,'addtoquotestring: string cannot be extended:');
         write(output,' '); {bug 2006 May 09}
         writequotestring(output,q);
      end
   end;
end;
{zzzFFF} {zzzfff}

procedure numtoquotestring(var q: string; i: integer);
(* add the number i to the end of q.  Handle negative numbers too!! *)
var
{
   digit: integer; (* one digit of i as a number *)
}
   digit: 0..9; (* one digit of i as a number *)
   newdigit: 0..127; (* ord(digit) *)
   c: integer; (* index to hold *)
   count: integer; (* how many digits were inserted *)
   chardigit: char; (* one digit of i as a character *)
   hold: string; (* for switching the digits around *)
   sign: integer; (* sign of i to keep track of the symbol *)
begin
{
writeln(searchfeatures,'numtoquotestring: i = ',i:1);
   clearquotestring(hold);
}
   clearstring(hold);
   if i < 0 then begin
      i := -i;
      sign := -1;
   end
   else sign := +1;
   with hold do begin
      while i > 0 do begin
         digit := i mod 10;

         (* gpc: warning: type casts are a Borland Pascal extension *)
{
         chardigit := char(digit + ord('0'));
}

         newdigit := digit + ord('0');

         (* gpc: warning: type casts are a Borland Pascal extension *)
         chardigit := char(newdigit);
{
2001 Sep 10: it is rather puzzling why this is being objected to ...
}

         addtoquotestring(hold,chardigit);
         count := succ(count);
         i := i div 10;
{
writeln(searchfeatures,'i = ',i:1);
}
      end;
   end;
   if sign < 0 then addtoquotestring(hold,'-');
   (* We got the digits from smallest to largest, so reverse them: *)
   for c := hold.length downto 1
   do addtoquotestring(q,hold.letters[c]);
end;
{zzzFFF} {zzzfff}

procedure dofeatures(var searchfeatures: text;
                     mismatches: integer;
                        pattern: string;
                     mispattern: string;
                     thefeature: featuretypeptr;
                     thedefinition: definitiontypeptr;
                     patterncomplement: boolean;
                     theplace: integer);
(* write out the features.  If there are mismatches, define a new
feature with them marked.
2008 Apr 22: rebuild the logic.  The original code did not
allow for additional marks besides the mutation.  That is,
aacta|tactgtatataaaaacagta|tcaat ~ =
has '|' whch was not being placed correctly.
theplace is the place of # or % in the string *)
const
   debug = false; (* turn on to test code *)
var
   holdNPend: char; (* hold negative part end *)
   holdPPend: char; (* hold negative part end *)
   holdLOCATIONend: real; (* hold location end *)
   m: integer; (* index to mismatches *)
   stuffatend: boolean; (* there is a change at the end of the feature.  In this
      case we do not want to overwrite the end with the original symbol. *)
   tmpdefinition: definitiontypeptr; (* the current definition with x changes *)
   tmpfeature: featuretypeptr; (* the current feature with x changes *)
   reversespot: integer; (* index to mismatches, the complement of m *)

{yyy}
   done: boolean; (* done moving stuff *)
   mm: integer; (* the location of the place to put the base, moved
      across the sequence until it is in the right spot. *)
   thelocation: integer; (* location of a mutation to record in the feature.
      The first base is now defined to be zero! *)
begin
   if mismatches = 0
   then begin
      writeafeature(searchfeatures,thefeature);
      writeln(searchfeatures);
{
writeln(output,'WRITEAFEATURE ---------- BEGIN');
writeafeature(output,thefeature);
writeln(output);
writeln(output,'WRITEAFEATURE ---------- END');
}
   end
   else begin
      write(searchfeatures,'* ',mismatches:1,' mismatch');
      if mismatches <> 1
      then write(searchfeatures,'s');
      writeln(searchfeatures);
      write(searchfeatures,'* ');
      writestring(searchfeatures,pattern);
      writeln(searchfeatures);
      write(searchfeatures,'* ');
      writestring(searchfeatures,mispattern);
      writeln(searchfeatures);
      new(tmpfeature);
      clearfeature(tmpfeature);
      new(tmpdefinition);
      (* copy thefeature to tmpfeature so that
         the changes don't affect the next result *)
      tmpfeature^ := thefeature^;
      tmpdefinition^ := thedefinition^;
(*
write(output,'y-thedefinition '); {yyy}
writeadefinition(output,thedefinition); {yyy}
*)

      (* Handle the ends of the feature *)
{ yyy neutralize old method 2008 apr 22
      with tmpdefinition^ do begin
         (* Do the end first: remove to make space for change marks *)
         if mispattern.letters[mispattern.length] = 'x' then begin
            (* throw away end character *)
            stuffatend := true;
         end
         else begin
            (* take off the end and hold it *)
            stuffatend := false;
            holdNPend := negparts.letters[2];
            holdPPend := posparts.letters[2];
            holdLOCATIONend := locations[2];
         end;
         marks := pred(marks);
         posparts.length := marks;
         negparts.length := marks;

         (* Now do the begin *)
         if mispattern.letters[1] = 'x' then begin
            (* remove now to never add later since there is X *)
            marks := pred(marks);
            posparts.length := marks;
            negparts.length := marks;
         end
      end;
}

      with tmpfeature^ do begin
         definition := tmpdefinition;
      end;

{
if debug then writeln(searchfeatures,'* holdLOCATIONend ',holdLOCATIONend:1:1);
if debug then writeln(searchfeatures,'* stuffatend ',stuffatend);
}

{QQQ}
if patterncomplement then writeln(searchfeatures, '* patterncomplement!!');

      for m := 1 to mispattern.length do begin
if debug then writeln(output,'*yyy ',m:2,' ',pattern.letters[m],mispattern.letters[m]);
         if mispattern.letters[m] = 'x' then begin
{zzzqqq}
{writeln(output,'*yyy ',m:2,' ',pattern.letters[m],mispattern.letters[m]); zq}
	    (* the first feature position is the zero coordinate, not 1: *)
            with tmpfeature^ do begin
{QQQ}
               if patterncomplement (* first base is zero! *)
               then thelocation := mispattern.length - m
               else thelocation := m - 1;
{zq}
               (* 2017 Feb 09 bug fix! correct thelocation by theplace! *)
               thelocation := thelocation - (theplace-1);
               addtoquotestring(   tmpfeature^.nametag,'.');
               numtoquotestring(   tmpfeature^.nametag,thelocation);
            end;
            with tmpdefinition^ do begin
               addtoquotestring(nametag,'.');
               numtoquotestring(nametag,thelocation);
               marks := succ(marks);
               if marks > maxstring then begin (* 2008 apr 22 *)
                  writeln(output, 'dofeatures: too many parts to a feature');
                  halt
               end;
	       (* first feature position is the zero coordinate, not 1: *)
               negparts.length := marks;
               posparts.length := marks;
               posparts.letters[posparts.length]:='!'; {yyy}

               (* move the negative parts to the right to create space *)
               for mm := marks downto 2
               do negparts.letters[mm]:= negparts.letters[mm-1];
               negparts.letters[1]:='!';

               (* new 2008 Apr 22: move the locations up *)
               (* mm is the hole location *)
               mm := marks-1;
if debug then writeln(output,'initial mm=',mm:3);
               done := false;
if debug then write(output,'A-tmpdefinition '); {yyy}
if debug then writeadefinition(output,tmpdefinition); {yyy}
               (* initialize upper location: *)
               locations[marks] := locations[marks-1];
if debug then write(output,'a-tmpdefinition '); {yyy}
if debug then writeadefinition(output,tmpdefinition); {yyy}
               while not done do begin 
if debug then write(output,'m=',m:3);
if debug then write(output,' mm=',mm:3);
if debug then writeln(output,' locations[mm]=',locations[mm]:3:1);
                  if (locations[mm] <= m) or (mm < 1) then begin

if debug then writeln(searchfeatures,'*    mm = ', mm:1);
if debug then writeln(searchfeatures,'*     m = ', m :1);
if debug then writeln(searchfeatures,'* marks = ', marks :1);
if debug then writeln(searchfeatures,'* mispattern.length = ', mispattern.length :1);
if debug then writeln(output,'(marks-mm)+1=', (marks-mm)+1:1); {yyy}

                     if not patterncomplement then begin
                        posparts.letters[mm+1] := pattern.letters[m];
         {REVERSE!}     negparts.letters[marks-mm] := chomplement(pattern.letters[m]);
                        locations[mm+1] := thelocation; (* first base is zero! *)
{QQQ}
                     end
                     else begin
                        posparts.letters[mm+1] := chomplement(pattern.letters[m]);
         {REVERSE!}     negparts.letters[marks-mm] := pattern.letters[m];
                        locations[mm+1] := thelocation; (* first base is zero! *)
                     end;
                     done := true;
                  end
                  else begin
if debug then if mm+1>marks then writeln(output,'mm+1>marks',(mm+1):1,'>',marks:1);
                     posparts.letters[mm+1] := posparts.letters[mm];
if debug then writeln(output,'zzz(marks-mm)+1=', (marks-mm)+1:1); {yyy}
if debug then writeln(output,'zzz(marks-mm)=', (marks-mm):1); {yyy}
   {REVERSE!!}       negparts.letters[(marks-mm)]
                                      :=negparts.letters[(marks-mm)+1];
                     locations[mm+1] := locations[mm];
                     locations[mm] := -1; (* track the open spot *)
                     posparts.letters[mm] := '_';
   {REVERSE!!}       negparts.letters[(marks-mm)+1] := '_';
                     mm := mm - 1;
                  end;


if debug then write(output,'B-tmpdefinition '); {yyy}
if debug then writeadefinition(output,tmpdefinition); {yyy}
               end;
if debug then write(output,'C-tmpdefinition '); {yyy}
if debug then writeadefinition(output,tmpdefinition); {yyy}
if debug then writeln(output,'negparts.length=', negparts.length:1); {yyy}
if debug then writeln(output,'posparts.length=', posparts.length:1); {yyy}
{halt;}

{QQQ}
{ original old code:
               if not patterncomplement then begin
                  negparts.letters[marks] := pattern.letters[m];
               end
               else begin
                  reversespot := mispattern.length - m + 1;
                  posparts.letters[marks] := pattern.letters[reversespot];
                  negparts.letters[marks] := pattern.letters[reversespot];
               end
}
            end
         end
;if debug then writeln(searchfeatures);
      end;

{zzzfff} {zzzFFF}
if debug then write(output,'x-tmpdefinition '); {yyy}
if debug then writeadefinition(output,tmpdefinition); {yyy}
{yyy 2008 apr 22 old code removed:
      (* If needed, put the end back on *)
      if not stuffatend then with tmpdefinition^ do begin
         addtoquotestring(negparts,holdNPend);
         addtoquotestring(posparts,holdPPend);
         marks := succ(marks);
         locations[marks] := holdLOCATIONend;
      end;
}

      (* write out the new definition and feature *)
{yyy}
if debug then writeln(output,'tmpdefinition marks = ',tmpdefinition^.marks:1); {yyy}
if debug then write(output,'y-tmpdefinition '); {yyy}
if debug then writeadefinition(output,tmpdefinition); {yyy}
      writeadefinition(searchfeatures,tmpdefinition);
      writeln(searchfeatures);
      writeafeature(searchfeatures,tmpfeature);
      writeln(searchfeatures);

{
(* note: pattern.length INCLUDES the '#' =59; sequence is only 58 zq *)
writeln(output); (*zq*)
writeln(output,'locrecord= ',locrecord:1); (*zq only in proc match*)
writeln(output,'theplace= ',theplace:1); (*zq only in proc match*)
writeln(output,'theplace= ',theplace:1); (*zq in proc match; passed now*)
writeln(output,'pattern.length= ',pattern.length:1); (*zq*)
writeln(output,'tmpdefinition:'); (* zzzqqq *)
writeadefinition(output,tmpdefinition); (* zzzqqq *)
writeln(output); (* zzzqqq *)
writeln(output,'tmpfeature:'); (* zzzqqq *)
writeafeature(output,tmpfeature); (* zzzqqq *)
writeln(output); (* zzzqqq *)
}

      (* clean up memory *)
      dispose(tmpfeature);
      dispose(tmpdefinition);
(*
halt; {yyy}
*)
   end;
end;

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

procedure answercheck(var answers: integer);
(* if the number of answers exceeds answersmax, we are probably
in an infinite loop from a file input so we halt. *)
begin
   answers := succ(answers);
   if answers > answersmax then begin
      writeln(output,'too many tries, quitting');
      halt;
   end;
end; (* answercheck *)

procedure searchhelp(var f: text);
(* write out all the commands to file f *)
(* for extra cirricular services contact minerva baldwin *)
var   answer: char;
begin (* searchhelp *)
writeln(f,'commands:');
writeln(f,'A: set the next named feature to be an arrows like [--->');
{
writeln(f,'   This must be done following the " command.');
}
writeln(f,'C: Define colors for the features in lister.');
writeln(f,'   The form is:');
writeln(f,'      C Tparam Aparam Bparam Cparam Dparam');
writeln(f,'   where');
writeln(f,'     Tparam is one of hHrR');
writeln(f,'       h color the letters using HSB colors');
writeln(f,'       r color the letters using RGB colors');
writeln(f,'       H put a colored rectangle behind the letter, HSB');
writeln(f,'       R put a colored rectangle behind the letter, RGB');
writeln(f,'     Aparam is either Hue        or Red');
writeln(f,'     Bparam is either Saturation or Green');
writeln(f,'     Cparam is either Brightness or Blue');
writeln(f,'     Dparam is the scale to multiply the rectangle height by');
writeln(f,'   Tparam determines the state used for Aparam, Bparam and Cparam');
{ppp}

write  (f,'D: enter Delila subsystem for writing instructions');
writeln(f,' (help is available there)');
writeln(f,'H: help (this list)');
writeln(f,'?: help (this list)');
writeln(f,'L: the single letters that correspond to the (x/y) forms');
writeln(f,'M: set maximum number of mismatches allowed (currently ',
              mismatches:1,')');

writeln(f,'N: a toggle switch:');
writeln(f,'   If true, searches are performed when');
writeln(f,'   the = command is given or when the pattern is typed.');
writeln(f,'   If false, searches are performed only when the =',
            ' command is given.');
write  (f,'   (currently ');
if oktosearch = true
then writeln(f,'true).')
else writeln(f,'false).');

writeln(f,'P: set start point (phase) and shift between search steps');
writeln(f,'q: quit');
writeln(f,'Q: quit');
writeln(f,'V: view: set viewing parameters (help is available there)');
write  (f,'*: the rest of the line is a comment sent to');
writeln(f,  ' the searchinst and result files');
writeln(f,'=: search for the last typed search pattern');
writeln(f,'~: complement the search pattern');
writeln(f,'   This reverses the order AND complements each character.');

writeln(f,'": create searchfeatures file');

writeln(f,'   * To create searchfeatures, define the name of a search string');
writeln(f,'   by typing the name inside double quotes (as: "EcoRI") and then');
writeln(f,'   do the search (eg type: gaattc).');

writeln(f,'   * Vertical bars or carets (| or ^) in the search string (as:');
writeln(f,'   g|aattc) will carry over to the feature.');

writeln(f,'   * All subsequent searches will be labeled with the same name');
writeln(f,'   until you give a new one.');

writeln(f,'   * To turn off the features, use an empty quote string, as "".');
writeln(f,'   * Because space is used to parse, no spaces are allowed in');
writeln(f,'   * the names.  However, the lister program will treat a tab.');
writeln(f,'   * as a space, so this can be used instead.');

writeln(f,'   * The searchfeatures file can be concatenated with other');
writeln(f,'     features to create the features file for the lister program.');

writeln(f,'   * Before using the name command, use the "A" command to ');
writeln(f,'     have the feature be an arrow.');
writeln(f,'I: invert: complement the characters in the search string.');
writeln(f,'   This DOES NOT change the order of the characters.');
writeln(f,'R: reverse the order of the characters in the search string.');
writeln(f,'   This reverses the order BUT DOES NOT complement characters.');
{qqq}

{zzz}

writeln(f,'one can type several commands on one line separated by spaces.');
writeln(f,'example (without quotes):  "M 1 V n V s V d V p gaattc * ecori"');
writeln(f);
      answers := 0;
      repeat
         if nostring(buffer) then
            writeln(output,'do you want to see the patterns allowed? (y/n)');
         getchar(input, buffer, answer, gotten);
         gotten := (answer = 'y') or (answer = 'n');
         answercheck(answers);
      until gotten;

      if answer='y' then begin
write  (f,'anything else recognizable is a pattern.');
writeln(f,' (the ''s are not needed):');
writeln(f,'the simplest pattern contains only "a", "c", "g", or "t".');
writeln(f,'one can also ask for purines ("r") or pyrimidines ("y").');
writeln(f,'the form "(x/y)" or "(x/y/z)" allows search for "x or y" or');
writeln(f,'"x or y or z".  x, y, and z must be "a", "c", "g", or "t".');
writeln(f);
writeln(f,'"n" means any base.');
writeln(f,'"e" means an extension: the scan is made with and without');
writeln(f,'    the e.  during the scan, e is treated like an "n".');
writeln(f,'one may ask for several letters, by stating the number. ');
writeln(f,'    so "5a" is the same as "aaaaa".');
writeln(f,'"#" the number of the next base is returned on a match.');
writeln(f,'    default: the first base beyond n''s and e''s');
writeln(f,'    Position of # is changed when complement (~) is taken.');
writeln(f,'"%" the number of the base before the current one is given.');
writeln(f,'    % is like #, but does NOT change position when ~ is taken.');
writeln(f,'"<" allow no mismatches until ">"');
writeln(f,'">" allow mismatches until "<" (the default)');
writeln(f,'    "a<cag>a" makes "cag" searched for exactly.');
writeln(f,'"string" Double quoted string is a name.');
writeln(f,'    It starts searchfeatures.  Empty string turns it off.');
writeln(f,'"|" symbol that will be put into searchfeatures BETWEEN bases');
writeln(f,'"^" symbol that will be put into searchfeatures BETWEEN bases');
writeln(f,'"[" symbol that will be put into searchfeatures BETWEEN bases');
writeln(f,'"]" symbol that will be put into searchfeatures BETWEEN bases');
writeln(f,'example: "acag(t/c/a)5e#2a<3r3y>".');
writeln(f);
      answers := 0;
      repeat
         answer := ' ';
         if nostring(buffer) then
            writeln(output,'do you want to see the relational ',
                     'patterns allowed? (y/n)');
         getchar(input, buffer, answer, gotten);
         gotten := (answer = 'y') or (answer = 'n');
         answercheck(answers);
      until gotten;

      if answer='y' then begin
writeln(f,'"$" specifies a relation between two bases. the format is:');
writeln(f,'    "$#r", where # is a position in the pattern (preceding');
writeln(f,'    this position) that is related to this position by r.');
writeln(f,'    r must be one of the following: i (identity);');
writeln(f,'    ni (non-identity); c (complement); nc (non-complement).');
writeln(f,'    w (complement including g-t pair); nw (non-w).');
writeln(f);
writeln(f,'examples: "n$1i" would search for aa, cc, gg and tt.');
writeln(f,'          "5n$1c" would search for complementary bases ');
writeln(f,'          separated by 4 unspecified bases.');
writeln(f,'    NOTE: this was formerly [before 1997] the ^ symbol.');
         end
      end;
      clearstring(buffer); (* make sure the buffer does not propagate *)
end; (* searchhelp *)

procedure showletters(var f: text);
(* show to file f the letter conversions of forms like "(x/y)" *)
begin (* showletters *)
      writeln(f,'a is a   ');
      writeln(f,'c is  c  ');
      writeln(f,'m is ac  ');
      writeln(f,'g is   g ');
      writeln(f,'r is a g ');
      writeln(f,'s is  cg ');
      writeln(f,'v is acg ');
      writeln(f,'t is    t');
      writeln(f,'w is a  t');
      writeln(f,'y is  c t');
      writeln(f,'h is ac t');
      writeln(f,'k is   gt');
      writeln(f,'d is a gt');
      writeln(f,'b is  cgt');
      writeln(f,'n is acgt');
      writeln(f,'e is acgt or nothing');
writeln(f,'warning: some of these ambiguous nucleotide symbols overlap ');
writeln(f,'         some of the search command symbols.  for instance, h');
writeln(f,'         can mean "a or c or t" or it can mean "help".');
writeln(f,'         to avoid problems, begin search patterns with a, c, g,');
writeln(f,'         t or n, or use the (x/y/z) form for ambiguous bases.');
end; (* showletters *)

procedure Delilahelp(var f: text);
(* write out commands in Delila subsystem to file f *)
begin (* Delilahelp *)
writeln(f,'this is a search subsystem for creating Delila');
writeln(f,'instructions based on search results.');
writeln(f);
writeln(f,'q: quit to search.');
writeln(f,'p: print Delila instructions on/off switch.');
writeln(f,'   when on, instructions will be printed in searchinst for each');
writeln(f,'   pattern found.  if n is the value of the numbered (#''d)');
writeln(f,'   base, f and t are the values from and to, then the');
writeln(f,'   instruction is a get of the form:');
writeln(f,'      get from n +f to n +t;');
writeln(f,'   organism, chromosome and piece specifications are');
writeln(f,'   automatically taken care of.');
writeln(f);
writeln(f,'f: set from: number of bases before #''d base to get from');
writeln(f,'t: set to: number of bases after #''d base to get to');
writeln(f,'   both from and to are with respect to the piece,');
writeln(f,'   so you do not need to worry about orientations.');
writeln(f,'i  set instruction type: one or two places in get.');
writeln(f,'   insttype=1 means the form   get from p  -/+f to same +/-t dir +/-;');
writeln(f,'   insttype=2 means the form   get from p1 -/+f to p2   +/-t dir +/-;');
writeln(f,'   where p or p2 are the current location, p1 is the previous,');
writeln(f,'   f is the from range and t is the to range.');
writeln(f,'   insttype=3 means the form   get from p1 -/+f to p2   +/-t dir +/-;');
writeln(f,'   where p1 is the current location 5'' end, p2 is the 3'' end,');
writeln(f,'w: to write your own Delila instructions.');
writeln(f,'   warning: they are not checked at all.');
writeln(f,'   q to quit.  common uses: title and defaults.');
writeln(f,'s: show all current Delila instructions');
writeln(f,'k: kill all current Delila instructions');
writeln(f,'h: help (this list).');
write  (f,'*: the rest of the line is a comment sent to');
writeln(f,  ' the searchinst and result files');
end; (* Delilahelp *)

procedure viewhelp(var f: text);
(* write to file f the sub-commands for the view command *)
procedure state(a: boolean);
(* write the state of a out *)
begin (* state *)
      if a then write(f,' on ')
           else write(f,'    ');
end; (* state *)
begin (* viewhelp *)
writeln(f,'the view command sets the amount of search results made visible.');
writeln(f,'the parts correspond to the order of printed results.');
writeln(f,'the current state of the toggle switches available');
writeln(f,'is shown below, along with the commands:');

state(viewt); writeln(f,'t: typed pattern');
state(viewe); writeln(f,'e: expanded pattern');
state(viewi); writeln(f,'i: identification information (organism, etc)');
state(viewf); writeln(f,'f: forms of patterns (with variable e)');
state(views); writeln(f,'s: sequence found');
state(viewp); writeln(f,'p: position of matches on piece');
state(viewd); writeln(f,'d: distance between matches');
state(viewm); writeln(f,'m: matches per piece');
state(viewb); writeln(f,'b: matches per book');
writeln(f);
state(false); writeln(f,'a: all results viewed');
state(false); writeln(f,'n: nothing viewed (even total matches)');
state(false); writeln(f,'q: quit view process (do nothing)');
state(false); writeln(f,'h: help (this list)')
end; (* viewhelp *)

function patternlength(pattern: string): integer;
(* give the length of the pattern. *)
var   p: integer; (* index to pattern *)
begin (* showpattern *)
   p:=1;
   while pattern.letters[p]<>'.' do begin
{
      write(afile,pattern.letters[p]);
}
      p:=succ(p)
   end;
   patternlength := p-1;
end; (* patternlength *)

procedure showpattern(var afile:text; pattern: string; ln: boolean);
(* write the pattern to the file.  do line feed depending on ln. *)
var   p: integer; (* index to pattern *)
begin (* showpattern *)
      write(afile,'"');
      p:=1;
      while pattern.letters[p]<>'.' do begin
         write(afile,pattern.letters[p]);
         p:=succ(p)
      end;
      write(afile,'"');
      if ln then writeln(afile)
end; (* showpattern *)

procedure expandpattern(var p:string; var ok: boolean; var pp: integer);
(* make p into an expanded pattern:
   all or's and numbers become characters.
   if ok is false then pp is the position of an error in the pattern
   p, otherwise it is an aribitrary position in p. *)
var
      t: string; (* temporary string for building absolute format *)
      tp: integer; (* the position in t *)
      rbase: integer;  (* position of related base *)
      number: integer; (* the number of times to repeat the unit *)
      numbers: set of '0'..'9';
      letter: char; (* the letter to put into t *)
      binary: packed array[1..4] of char; (* represents "(x/y)" forms *)
      expecting: boolean; (* a letter is expected after "(" or "/" *)

procedure getnumber(p: string; var pp, number: integer);
(* pull out a number from p starting at pp.
   move pp past the number *)
var
      n, (* a second position in p *)
      power, (* power of 10 *)
      increment (* used to calculate the number *)
         : integer;
      done: boolean; (* are we done yet? *)
begin (* getnumber *)
      number:=0;
      while (p.letters[pp] in numbers) do pp:=succ(pp);
      (* pp is now past the number *)
      n:=pred(pp); (* n is at the lowest digit *)

      power:=1;
      done:=false;

      repeat
         case p.letters[n] of
            '0': increment:=0;
            '1': increment:=1;
            '2': increment:=2;
            '3': increment:=3;
            '4': increment:=4;
            '5': increment:=5;
            '6': increment:=6;
            '7': increment:=7;
            '8': increment:=8;
            '9': increment:=9;
         end;
         number:=number+power*increment;
         power:=power*10;
         n:=pred(n);
         if n=0 then done:=true
                else done:=not(p.letters[n] in numbers)
      until done;
      if p.letters[pp]='.' then ok:=false;
      if number=0 then ok:=false
end; (* getnumber *)

procedure getrelation;
(* we now have a pattern symbol ('$') that specifies a relation between
   the base at this position and another base.  we will now determine
   what that relation is and to what position it refers.  *)
var   non: boolean;  (* means the negation of a specific relation *)

begin   (* getrelation *)
      relational := true;
      non := false;
      pp := succ(pp);  (* get the next letter *)
      if (p.letters[pp] in numbers) then getnumber(p,pp,rbase)
      else begin
         ok := false;
         writeln(output,' number required after relation symbol');
      end;
      if (p.letters[pp] = 'n') then begin
         non := true;  (* the negation of the relation *)
         pp := succ(pp);
      end;
      if (p.letters[pp] = 'i') then
            (* identity *)
         if non then letter := 'p' else letter := 'l'
      else if (p.letters[pp] = 'c') then
            (* complementarity *)
         if non then letter := 'q' else letter := 'j'
      else if (p.letters[pp] = 'w') then
            (* complementarity including g-t pairs *)
         if non then letter := 'z' else letter := 'x'
      else begin
         ok := false;
         if (p.letters[pp] = '.') then
            writeln(output,'early end to pattern, relation expected')
            else writeln(output,'illegal relation specified: ',p.letters[pp]);
      end;
      end;

begin (* expandpattern *)
      numbers:=['0','1','2','3','4','5','6','7','8','9'];
      pp:=1;
      tp:=0;
      ok:=true;
      relational := false;  (* until we see a relation in the pattern *)
      while (p.letters[pp]<>'.') and ok do begin
         if p.letters[pp] in numbers then getnumber(p,pp,number)
                             else number:=1;
         if ok then begin
            if p.letters[pp]='(' then begin (* begin "(x/y..." *)
               expecting:=true;
               pp:=succ(pp);
               if p.letters[pp]<>'.' then begin
                  binary:='    ';
                  while (p.letters[pp]<>')') and ok do begin
                     expecting:=false;
                     if p.letters[pp] in ['a','c','g','t']
                     then begin
                        case p.letters[pp] of
                           'a': binary[1]:=p.letters[pp];
                           'c': binary[2]:=p.letters[pp];
                           'g': binary[3]:=p.letters[pp];
                           't': binary[4]:=p.letters[pp];
                        end;
                     end
                     else if p.letters[pp]='.' then begin
                        expecting:=true;
                        ok:=false
                     end
                     else begin
                        ok:=false;
                        writeln(output,'character not in acgt');
                     end;
                     if ok then begin
                        pp:=succ(pp);
                        if p.letters[pp] <> ')'
                           then if p.letters[pp] = '/'
                              then begin
                                 pp:=succ(pp);
                                 expecting:=true; (* need char after / *)
                              end
                              else begin
                                 ok:=false;
                                 writeln(output,'i expected / or )');
                              end;
                     end
                  end
               end
               else begin
                  ok:=false;
                  writeln(output,'unclosed (')
               end;
               if expecting then begin
                 ok:=false;
                 writeln(output,'i expect more after ( or /');
                 pp:=pred(pp)
               end;

               if ok then
                       if binary='    ' then ok:=false
                  else if binary='a   ' then letter:='a'
                  else if binary=' c  ' then letter:='c'
                  else if binary='ac  ' then letter:='m'
                  else if binary='  g ' then letter:='g'
                  else if binary='a g ' then letter:='r'
                  else if binary=' cg ' then letter:='s'
                  else if binary='acg ' then letter:='v'
                  else if binary='   t' then letter:='t'
                  else if binary='a  t' then letter:='w'
                  else if binary=' c t' then letter:='y'
                  else if binary='ac t' then letter:='h'
                  else if binary='  gt' then letter:='k'
                  else if binary='a gt' then letter:='d'
                  else if binary=' cgt' then letter:='b'
                  else if binary='acgt' then letter:='n'
            end (* close "(x/y...)" *)
            else begin
               if (p.letters[pp] in
                  ['a','c','g','t','n','e','#','%','|','^','[',']','$',
                   'm','r','w','s','y','k','v','h','d','b'])
                   or (p.letters[pp]='<') or (p.letters[pp]='>')
                   or (p.letters[pp]='^')
               then letter:=p.letters[pp]
               else begin
                  ok:=false;
                  writeln(output,'funny character');
               end;
            end;
               (* if letter = '$' we need to find what relation is specified
                  and to what base it applies *)
            if (letter = '$') then getrelation
               else rbase := 0;
            if ok then begin (* we have the letter now *)
               while (number<>0) and ok do begin
                  tp:=succ(tp);
                  if tp>=maxstring then begin
                     ok:=false;
                     writeln(output,'expanded pattern too long')
                  end
                  else begin (* store in t *)
                     t.letters[tp]:=letter;
                     if (tp > rbase) then relatedbase[tp] := rbase
                        else begin
                           writeln(output,' related base must',
                                          ' preceed relation');
                           ok := false;
                        end;
                     number:=pred(number)
                  end
               end;
               if ok then pp:=succ(pp)
            end
         end
      end;

      if ok then begin
         p:=t; (* replace pattern from temporary pattern *)
         p.letters[tp+1]:='.' (* must end on period *)
      end
end; (* expandpattern *)

function complementbase(c: char): char;
(* complement the base *)
begin
   if not (c in [
         'a','c','m','g','r','s','v','t','w','y','h','k','d','b','n','%',
         '<','>','e','#','[',']',
         '|','^',
         '(','/',')',
         '0','1','2','3','4','5','6','7','8','9'])
   then begin
      writeln(output,'bad character "',c,'" found during complement');
   end
   else begin
      case c of
         'a': complementbase := 't'; (* is a    *)
         'c': complementbase := 'g'; (* is  c   *)
         'm': complementbase := 'k'; (* is ac   *)
         'g': complementbase := 'c'; (* is   g  *)
         'r': complementbase := 'y'; (* is a g  *)
         's': complementbase := 'w'; (* is  cg  *)
         'v': complementbase := 'b'; (* is acg  *)
         't': complementbase := 'a'; (* is    t *)
         'w': complementbase := 's'; (* is a  t *)
         'y': complementbase := 'r'; (* is  c t *)
         'h': complementbase := 'd'; (* is ac t *)
         'k': complementbase := 'm'; (* is   gt *)
         'd': complementbase := 'h'; (* is a gt *)
         'b': complementbase := 'v'; (* is  cgt *)
         'n': complementbase := 'n'; (* is acgt *)
         '<': complementbase := '>'; (* start of match *)
         '>': complementbase := '<'; (* end of match *)
         'e': complementbase := 'e'; (* expand *)
         '(': complementbase := ')'; (* parens reverse *)
         ')': complementbase := '('; (* parens reverse *)
         '%','#','|','^','0','1','2','3','4','5','6','7','8','9',
         '[',']', (* bracket does not change *)
         '/'      (* slash does not change *)
{zzz}
            : complementbase := c
      end
   end
end; (* complementbase *)

procedure stringcomplement(var p: string);
(* create the complement of the string p *)
(* when numbered strings of the form 32n or #5 (not %) are found,
the substring (32n or #5) must not be reversed in order.  The
machinery below does this by keeping track of the start and stop
points of every string found. *)
var
   h: string; (* for holding the p string *)
   i: integer; (* index to the strings *)
   j: integer; (* index to the strings *)
   endofflip: boolean; (* true if the end of a numbered region is found *)
   flipnumber: boolean; (* true if there is a number to be flipped *)
   upper: integer; (* a constant to make things a bit faster *)
   start: integer; (* start of a flipped region *)
   stop: integer; (* stop of a flipped region *)
begin
   clearstring(h);
   h.length := p.length;
   flipnumber := false;
   upper := p.length+1;
   for i := 1 to p.length do begin
      endofflip := true;
      stop := i;
      h.letters[upper-i] := complementbase(p.letters[i]);
      if p.letters[i] in ['#','0','1','2','3','4','5','6','7','8','9']
      then begin
         endofflip := false;
         if not flipnumber then begin
            flipnumber := true;
            start := i
         end
      end;
{
write(output,'before flip: "');writestring(output,h); writeln(output,'"');
}
      if endofflip and flipnumber
      then begin (* move the entire number and the next base
                    without turning it around *)
{
writeln(output,'endofflip');
}
         for j := start to stop
         do h.letters[(upper-stop)+(j-start)] := complementbase(p.letters[j]);
         flipnumber := false
      end
{
;write(output,' after flip: "');writestring(output,h); writeln(output,'"');
}
   end;
   copystring(h,p);
   write(output,'* the complementary search string is now: "');
   writestring(output,p);
   writeln(output,'"');
   write(result,'* the complementary search string is now: "');
   writestring(result,p);
   writeln(result,'"');
   if printinginstructions then begin
      write(searchinst,'(* the complementary search string is now: "');
      writestring(searchinst,p);
      writeln(searchinst,'" *)');
   end;
{yyy}
{
halt;
}
end;

{qqq}
procedure stringreverse(var p: string);
(* Reverse the order of the string p WITHOUT taking
the complement. No symbols are altered. *)
var
   upper: integer; (* a constant to make things a bit faster *)
   i: integer; (* index to the string *)
   midpoint: integer; (* half of the string length *)
   hold: char; (* a character in the string *)
begin
   upper := p.length+1;
   midpoint := p.length div 2;
   for i := 1 to midpoint do begin
      hold := p.letters[i];
      p.letters[i] := p.letters[upper-i];
      p.letters[upper-i] := hold;
   end;
   write(output,'* the reversed search string is now: "');
   writestring(output,p);
   writeln(output,'"');
   write(result,'* the reversed search string is now: "');
   writestring(result,p);
   writeln(result,'"');
   if printinginstructions then begin
      write(searchinst,'(* the reversed search string is now: "');
      writestring(searchinst,p);
      writeln(searchinst,'" *)');
   end;
end; (* stringreverse *)

{qqq}
procedure stringinvert(var p: string);
(* complement each letter in the string.
The order is not changed. *)
var
   i: integer; (* index to the string *)
begin
   for i := 1 to p.length do begin
      p.letters[i] := complementbase(p.letters[i]);
   end;
   write(output,'* the inverted (complement only) search string is now: "');
   writestring(output,p);
   writeln(output,'"');
   write(result,'* the inverted (complement only) search string is now: "');
   writestring(result,p);
   writeln(result,'"');
   if printinginstructions then begin
      write(searchinst,'(* the inverted (complement only)',
                       ' search string is now: "');
      writestring(searchinst,p);
      writeln(searchinst,'" *)');
   end;
end; (* stringinvert *)

procedure analysepattern(var p,e: string; var ok: boolean);
(* the pattern is forced to end in a blank, and it
   is expanded into e.  ok is true when no errors were found. *)

var
      n: integer; (* position in p *)
      elate: boolean; (* e on end of pattern *)
      numofnums: integer; (* number of # symbols *)
      error: string; (* used to show an error *)
begin(* analysepattern *)
      clearstring(error);
      ok:=true;

      (* locate periods in pattern *)
      for n:=1 to p.length
         do if p.letters[n]='.' then begin
            ok := false;
            (* The original convention in search.p was to define the end of the 
            string with '.'.  Make it possible to change/upgrade this by
             supplying the length *)
            p.length := n - 1;
         end;
{zzzsss}

      if not ok
         then writeln(output,'no periods (.) allowed in the search pattern');

      if p.length+1 <= maxstring
      then with p do letters[length+1] := '.'
      else begin
         ok:= false;
         writeln(output,'search string is too long')
      end;

      (* begin expansion *)
      e:=p;

      (* find last possible e *)
      (* do it now since we know end: length *)
      n:=p.length;
      while (n<>1) and ((e.letters[n] in
         ['0','1','2','3','4','5','6','7','8','9']) or
         (e.letters[n]='#') or (e.letters[n]='%') or 
         (e.letters[n]='<') or (e.letters[n]='>') or
         (e.letters[n]='|') or (e.letters[n]='^') or
         (e.letters[n]='[') or (e.letters[n]=']') )
         do n:=pred(n);
      if e.letters[n]='e' then begin
         ok:=false;
         error.letters[n]:='e';
         elate:=true
      end
      else elate:=false;

      if ok then begin
         expandpattern(e,ok,n);
         if not ok then error.letters[n]:='^'
      end;

      (* find first possible e.  this is done here since expandpattern
         removes numbers, making it easy to spot a prefixed e *)
      n:=1;
      while ((e.letters[n]='#') or (e.letters[n]='%') or
             (e.letters[n]='<') or (e.letters[n]='>') or
             (e.letters[n]='|') or (e.letters[n]='^') or
             (e.letters[n]='[') or (e.letters[n]=']') or
             (e.letters[n]='(') or (e.letters[n]=')') or
             (e.letters[n]='$') or (e.letters[n]='/')
             ) {zzzddd}
                and (n<maxstring)
         do n:=succ(n);
      if e.letters[n]='e' then error.letters[n]:='e';
      if (e.letters[n]='e') or elate then begin
         ok:=false;
         writeln(output,'no e''s allowed on ends')
      end;

      (* check that # is correct *)
      numofnums:=0;
      n:=1;
      while e.letters[n] <> '.' do begin
         if (e.letters[n]='#') or (e.letters[n]='%')
         then numofnums:=succ(numofnums);
         n:=succ(n)
      end;
      if numofnums>1 then begin
         writeln(output,'extra # or % symbols');
         ok:=false;
         error.letters[n-1]:='#';
      end;
      if (e.letters[pred(n)]='#') or (e.letters[pred(n)]='%')
      then begin
{
         writeln(output,'# symbol not allowed at end');
         ok:=false;
         error.letters[n-1]:='#';
}
         writeln(output,'# or % symbol found at end: adding n to end');
         if succ(p.length) > maxstring then begin
            writeln(output,'... that makes the string too long');
            ok:=false;
            error.letters[n-1]:='#';
         end
         else begin
            p.length := succ(p.length);
            p.letters[p.length] := 'n';
            p.letters[p.length+1] := '.';
            e.length := p.length;
            e.letters[e.length] := 'n';
            e.letters[e.length+1] := '.';
         end
      end;

      (* make sure it is not all n's, e's, <'s, >'s, #'s and %'s *)
      n:=1;
      while (e.letters[n] in ['n','e']) or
            (e.letters[n]='$') or {zzzddd}
            (e.letters[n]='<') or
            (e.letters[n]='>') or
            (e.letters[n]='#') or
            (e.letters[n]='%') or
            (e.letters[n]='|') or
            (e.letters[n]='^') or
            (e.letters[n]='[') or
            (e.letters[n]=']') or
            (e.letters[n]='(') or
            (e.letters[n]=')') or
            (e.letters[n]='/')
            do n:=succ(n);
      if e.letters[n]='.' then begin
         ok:=false;
         error.letters[n-1]:='&';
         writeln(output,'You must have other characters in addition to:',
                        ' ne$<>#%|^[]');
      end;

      if ok then begin
         if viewt then begin
            writeln(output);
            write(output,'   typed pattern: ');
            showpattern(output,p,true);

            write(result,'*    typed pattern: ');
            showpattern(result,p,true);
         end;
         if viewe then begin
            write(output,'expanded pattern: ');
            showpattern(output,e,true);

            write(result,'* expanded pattern: ');
            showpattern(result,e,true);
         end;

         if printinginstructions then begin
            writeln(searchinst);
            write(searchinst,'(*   typed pattern: ');
            showpattern(searchinst,p,false);
            writeln(searchinst,' *)')
         end
      end
      else begin (* write a message out for errors *)
         (* get . for error from pattern *)
         n:=1;
         while e.letters[n]<>'.' do n:=succ(n);
         error.letters[n]:='.';

         showpattern(output,e,true);
         showpattern(output,error,false);
         writeln(output,' unusable pattern');

         (* ignore further commands on the line: *)
         clearstring(buffer)
      end
end; (* analysepattern *)

procedure compress(var p: string; s: integer; var rb: intarray);
(* remove from p the character pointed to by s
   and shift the rest of the pattern to the left *)
var
      ss: integer; (* successor of s *)
begin (* compress *)
   while p.letters[s] <> '.' do begin
      ss:=succ(s);
      p.letters[s]:=p.letters[ss];
      rb[s] := rb[ss];
      s:=ss
   end;
   p.length := ss;
end;(* compress *)

function related(b1, b2: base; r : char): boolean;
(* true if b1 and b2 are related as specified by r *)
begin
      case r of
            (* identity *)
         'l' : related := (b1 = b2);
            (* non-identity *)
         'p' : related := (b1 <> b2);
            (* complementarity *)
         'j' : related := (b1 = complement(b2));
            (* non-complementary *)
         'q' : related := (b1 <> complement(b2));
            (* complementarity including g-t pairs *)
         'x' : related := (  (b1 = complement(b2)) 
                       or   ((b1 = g) and (b2 = t))
                       or   ((b1 = t) and (b2 = g))  );
            (* non-complementarity, including g-t pairs *)
         'z' : related := (  (b1 <> complement(b2)) 
                       and not ((b1 = g) and (b2 = t))
                       and not ((b1 = t) and (b2 = g))  );
      end;
end;

procedure match(pattern: string; mismax: integer;
                relatedbase: intarray;
                var org: orgkey; var orgchange: boolean;
                var chr: chrkey; var chrchange: boolean;
                var pie: pieceptr; var piechange: boolean;
                patterncomplement: boolean); (* if true, reverse inst's *)
(* match this exact pattern in pie, allow mismax mismatches *)
var
      (* for scanning piece *)
      i: integer; (* internal position on piece *)
      length: integer; (* length of this piece *)
      searchlength: integer; (* the actual length of dna to search *)
      ba: base; (* a base in the piece *)
      equal: boolean; (* the pattern is equal to the piece at this place *)
      startdna, (* the first dnastring in the piece *)
      currentdna: dnaptr; (* where we are *)
      cd: dnarange; (* index to currentdna *)
      subcurrentdna: dnaptr; (* the place we are checking *)
      subcd: dnarange; (* index to subcurrentdna *)

      (* pointers in pattern *)
      patstart,
      patstop, (* the pattern with n's and e's removed from ends *)
      patspot: (* the location in between (inclusive) patstart and patstop
                  or a spot in general on the pattern *)
         integer;

      (* for matches *)
      patmatches: integer; (* number of pattern matches on this piece *)
      j: integer; (* internal coordinate position of first letter in pattern *)
      blank: boolean; (* put a blank at that spot (when outside the piece) *)
      realpattern: string; (* the dna in the piece corresponding to pattern *)
      circle: boolean; (* configuration of the piece (for speed) *)
{zzzqqq}
      locrecord, (* the location in pattern to record: # and % command *)
      theplace, (* the actual place on the pattern recorded # and % *)
      pieceplace, (* like theplace, but on the piece *)
      internalplace, (* internal coordinate location for pieceplace *)
      firstinternalplace, (* the internal coordinate of the first match *)
      lastinternalplace, (* the internalplace matched before this one *)
      distance: (* the distance to last match or -1 *)
         integer;

      (* for mismatches *)
      fail: boolean; (* fail to match the piece *)
      mismatches: integer; (* the actual number of mismatches *)
      mispattern: string; (* the places mismatches occured *)
      blankpattern: string; (* for clearing mispattern rapidly *)

      (* where mismatches are allowed: *)
      allowmismatch: packed array[1..maxstring] of boolean;
      allowed: boolean; (* used to build allowmismatch *)

      firsttime: boolean; (* the first time a match is found, the
         piece name instruction is written; otherwise be silent *)

      r: integer; (* index for reversing feature string *)

procedure display(var afile: text);
(* display the match.
   global variables used:
   realpattern, mismatches, theposition, pieceplace *)
var   p: integer;
begin (* display *)
      if views then begin
         write(afile,'  ');
         showpattern(afile,realpattern,true);

         (* show mismatches *)
         if mismax<>0 then begin
            write(afile,'  ');
            showpattern(afile,mispattern,false);
{
            writeln(afile,' ',mismatches:1,' mismatche(s)')
}
            write(afile,' ',mismatches:1,' mismatch');
            if mismatches <> 1 
            then write(afile,'s');
            writeln(afile);
         end;

         if viewp then begin
            (* mark location *)
            write(afile,'  ');
            for p:=1 to theplace do write(afile,' '); (* 0->1 1990 Sep 4 *)
            write(afile,'^ ',pieceplace:1,' ')
{yyy}
         end
      end
      else if viewp then write(afile,' ',pieceplace:6);

      if viewd
      then if distance > 0
         then write(afile,' ',distance:6,' bases from last match ');

      if views or viewd or viewp then writeln(afile,'*')
end; (* display *)

begin (* match *)
      firsttime := true;

      (* In this part we check the special case in which a search
      string begins with [ or ends with ].  In these cases,
      remove the given so that it can be replaced BETWEEN letters *)
{
      showpattern(output,pattern,false);
zzzsss this seems unnecessary!
}
      if (pattern.letters[1]='[')
      then begin
         with thedefinition^ do begin
            posparts.letters[1] := '-';
         end
      end;
{zzzNNN}
      (* Do the same thing for the right side *)
      (* find the end of the string in patspot *)
      patspot := 1;
      while pattern.letters[patspot]<>'.'
      do patspot := succ(patspot);
      if (pattern.letters[patspot-1] = ']')
      then begin
         with thedefinition^ do begin
            posparts.letters[2] := '-';
         end
      end;

{zzzsss}
{
writeln(output);
write(output,'BEFORE ---> "');
writeadefinition(output,thedefinition);
writeln(output,'"');
}

      (* extract # and % command and mismatch areas *)
      locrecord := 1;
      allowed := true;
      patspot := 1;
      poundsign := false;
      while pattern.letters[patspot]<>'.' do begin

         if (pattern.letters[patspot]='<') or (pattern.letters[patspot]='>') or
            (pattern.letters[patspot]='#') or (pattern.letters[patspot]='%') or
            (pattern.letters[patspot]='|') or (pattern.letters[patspot]='^') or
            (pattern.letters[patspot]='[') or (pattern.letters[patspot]=']')
         then begin

                 if (pattern.letters[patspot]='#') then poundsign := true;
{if poundsign then writeln(output,'POUNDSIGN FOUND') 
             else writeln(output,'no POUNDSIGN found'); zzzsss}
                 if pattern.letters[patspot]='#' then locrecord:=patspot
            else if pattern.letters[patspot]='%' then locrecord:=patspot-1
            else if pattern.letters[patspot]='<' then allowed:=false
            else if pattern.letters[patspot]='>' then allowed:=true;

            if dosearchfeatures then with thedefinition^ do begin
               if (pattern.letters[patspot]='#') or
                  (pattern.letters[patspot]='%') or
                  (pattern.letters[patspot]='|') or
                  (pattern.letters[patspot]='^') or
                  (pattern.letters[patspot]='[') or
                  (pattern.letters[patspot]=']')
               then if not featureinserts then with thedefinition^ do begin
                  marks := succ(marks);

                  with posparts do begin (* move end mark out *)
                     length := succ(length);
                     letters[length] := letters[length-1];
                     letters[length-1] := pattern.letters[patspot];
                  end;

                  (* the location is relative to the current one
                     For '%', the location to put the % is between bases
                     hence there is a subtraction of 0.5. *)
                  if (pattern.letters[patspot]='%') or
                     (pattern.letters[patspot]='[') or
                     (pattern.letters[patspot]=']')
                  then locations[posparts.length-1] := patspot - locrecord - 0.5
                  else locations[posparts.length-1] := patspot - locrecord;

                  (* Note:
                     "EcoRI0"   #gaattc
                     "EcoRI1"   g#aattc
                     "EcoRI2"   ga#attc
                     "EcoRI3"   gaa#ttc
                     "EcoRI4"   gaat#tc
                     "EcoRI5"   gaatt#c
                     "EcoRI6"   gaattc#
                     The last two of this series are handled a bit oddly.
                     Ecori5 gives no # and 6 extends the length.
                     This may never matter... *)

                  (* vertical bar is between bases *)
                  if (pattern.letters[patspot]='|') or
                     (pattern.letters[patspot]='^')
                  then locations[posparts.length-1] :=
                       locations[posparts.length-1] - 0.5;

                  (* do negparts here *)
                  with negparts do begin (* move end mark out *)
                     length := succ(length);
                     letters[length] := letters[length-1];
                     letters[length-1] := pattern.letters[patspot];
                  end;

                  (* the middle part of the string has to be reversed *)
                  for r := 2 to posparts.length-1 do begin
                     negparts.letters[(posparts.length-1)-r + 2]
                        := posparts.letters[r];
                  end

               end;
            end;
            compress(pattern,patspot,relatedbase)
         end
         else begin
            allowmismatch[patspot]:=allowed;
            patspot:=succ(patspot)
         end
      end;
      pattern.length := patspot; (* that's the current length *)

{zzzsss}
{
write(output,'AFTER ---> "');
writeadefinition(output,thedefinition);
writeln(output,'"');
}

      if dosearchfeatures then with thedefinition^ do begin
         (* pattern.length includes ".", the definition runs from zero up,
            so subtract 1: *)
{ zzzsss}
{
writeln(output,'****> locrecord: ',locrecord:1);
writeln(output,'****> locations[1]: ',locations[1]:1:1);
writeln(output,'****> marks: ',marks:1);
writeln(output,'****> pattern.length: ',pattern.length:1);
}
         (* only do this once per feature! *)
         if not featureinserts then begin
{
writeln(output,'***** DONE');
}
            locations[1] := 1 - locrecord;
            locations[marks] := (pattern.length-1) - locrecord;
         end;
         featureinserts := true;
{
}
{zzzfff}
{
writeln(searchfeatures,'in match');
writeln(searchfeatures,'pattern.length = ',pattern.length:1);
writeln(searchfeatures,'length = ',length:1);
writeln(searchfeatures);
writeln(searchfeatures,'* THE CURRENT FEATURE DEFINITION at match IS:');
write(output,'DOUBLE AFTER ---> "');
writeadefinition(output,thedefinition);
writeln(output,'"');
}
      end;

      (* display pattern form *)
      if viewf then begin
         write(output,'  '); (* spaces to match the comment in results *)
         showpattern(output,pattern,true);
         write(result,'* ');
         showpattern(result,pattern,true)
      end;

      (* find pattern w/o n's or e's on ends *)
      (* find start of pattern *)
      patstart:=1;
         (* if this search contains a relation, we must use all of it,
            even beginning n's, otherwise we can skip them *)
      if not relational then
      while pattern.letters[patstart] in ['n','e']
         do patstart:=succ(patstart);
      (* find end of pattern *)
      patstop:=patstart;
      while pattern.letters[patstop]<>'.' do patstop:=succ(patstop);
      (* find end of pattern without n's or e's *)
      while pattern.letters[patstop] in ['.','n','e']
         do patstop:=pred(patstop);


      (* calculate length of dna *)
if pie = nil then begin
writeln(output,'PROGRAM ERROR: empty DNA POINTER!');
halt;
end;
      length:=pietoint(pie^.key.pieend,pie);
      (* calculate length of dna to search *)
      case pie^.key.piecon of
         linear: searchlength:=length-(patstop-patstart);
         circular: searchlength:=length
      end;

      (* set up dna for stepbase *)
      startdna:=pie^.dna;
      currentdna:=startdna;
      cd:=0;

      (* decide configuration *)
      circle:=(pie^.key.piecon=circular);

      (* clear blank pattern *)
      patspot:=1;
      while pattern.letters[patspot]<>'.' do begin
         blankpattern.letters[patspot]:=' ';
         patspot:=succ(patspot);
      end;
      blankpattern.letters[patspot]:='.'; (* last one *)
      blankpattern.length:=patspot-1;
{FFF}

      (* search the piece for the pattern *)
      patmatches:=0;
      lastinternalplace:=0; (* start nowhere *)
      distance := -1; (* no previous match to calculate distance with *)

      (* now move to the correct start *)
      for i := 2 to startpoint
         do ba := stepbase(startdna,currentdna,cd);
      i := startpoint; (* the user defined start *)
      while i <= searchlength do begin
         patspot:=patstart;
         subcurrentdna:=currentdna;
         subcd:=cd;
         fail:=false;
         mismatches:=0;
         mispattern:=blankpattern;
         realpattern := blankpattern;
         repeat (* until the pattern matches or doesn't *)
            ba:=stepbase(startdna,subcurrentdna,subcd);
            if relational then realpattern.letters[patspot] := basetochar(ba);
            with pattern do
            if (letters[patspot] in ['l','p','j','q','x','z']) then
               equal := related(ba,
                  chartobase(realpattern.letters[relatedbase[patspot]]),
                  letters[patspot])
            else case ba of
               a: equal:=
                  letters[patspot] in ['a','m','r','v','w','h','d','n','e'];
               c: equal:=
                  letters[patspot] in ['c','m','s','v','y','h','b','n','e'];
               g: equal:=
                  letters[patspot] in ['g','r','s','v','k','d','b','n','e'];
               t: equal:=
                  letters[patspot] in ['t','w','y','h','k','d','b','n','e'];
            end;
            if not equal then begin
               mismatches:=succ(mismatches);
               mispattern.letters[patspot]:='x';
{
mispattern.length := patspot;
zzzFFF
}
               if allowmismatch[patspot]
                  then fail:= (mismatches>mismax)
                  else fail:= true
            end;
            patspot:=succ(patspot)
         until fail or (patspot=succ(patstop));

         if not fail then begin (* yahoo got un *)
            patmatches:=succ(patmatches); (* tally up *)
            theplace:=locrecord;
{writeln(output,'****** theplace: ',theplace:1); (*zq*)}
            patspot:=1;
            (* pull out the real pattern *)
            repeat
               (* internal coordinate position at a letter in pattern *)
               j:=i-(patstart-patspot);
               blank:=false;

               (* take care of pattern off end of piece *)
               if j<1
               then if circle
                  then while j<1 do j:=j+length
                  else blank:=true
               else if j>length
               then if circle
                  then while j>length do j:=j-length
                  else blank:=true;

                  (* get the real pattern *)
               if blank
               then realpattern.letters[patspot]:=' '
               else realpattern.letters[patspot]:=basetochar(getbase(j,pie));

               if patspot=locrecord then begin
                  theplace:=patspot;
                  if blank
                  then if patspot<patstart
                          then theplace:=patstart
                          else theplace:=patstop;
                  if viewd or (insttype=2)
                  then lastinternalplace := internalplace;
                  internalplace:=i-(patstart-theplace);
                  pieceplace:=inttopie(internalplace,pie)
               end;

               patspot:=succ(patspot)
            until pattern.letters[patspot]='.';
            realpattern.letters[patspot]:='.'; (* close pattern *)

            if viewd or (insttype=2)
            then if distance < 0
               then begin
                  distance := 0; (* no distance on first match *)
                  firstinternalplace := internalplace (* for circular case *)
               end
               else distance := internalplace - lastinternalplace;
            display(output);
            display(result);
{
writeln(output,'printinginstructions: ', printinginstructions);
writeln(output,'distance: ', distance:1);
}
            if printinginstructions then begin

{zzz I can't believe that the iworgchr did NOT do the write
here - the reason must be that the important thing is that
when delila instructions are turned on, orgchange and chrchange should
be set true!}
                     (*force these for all instructions: *)
{zzziii 1999 apr 29}

{1999 dec 24 zzzccc
christmas changes:
                     orgchange := true;
                     chrchange := true;

                     iworgchr(searchinst,org,orgchange,chr,chrchange);
                     (* at this point, so far as we are concerned,
                     the organism and chromosome changes are not
                     important because we have written instructions
                     out about them.  So we set them false. *)
                     orgchange := false;
                     chrchange := false;
}

{zzzttt}
                     (* 2000 June 26 I added orgopen, chropen here.
                        This may well solve the problems! *)
                     if orgchange or chrchange
                     then iworgchr(searchinst,org,orgchange,orgopen,
                                              chr,chrchange,chropen);

                     if firsttime or piechange
                     then iwpie(searchinst,pie^.key);
                     if firsttime then firsttime := false;

                     (* name the piece into the delila instructions *)
                     if dosearchfeatures then begin
                        {write(searchinst,'name "'); bug 2006 May 09}
                        write(searchinst,'name ');
{zzznnn}
                        (* NOTE: the quotes are give by the
                           writequotestring procedure! *)
                        writequotestring(searchinst,thedefinition^.nametag);
                        {writeln(searchinst,'";'); bug 2006 May 09}
                        writeln(searchinst,';');
                     end;

{zzz} {
writeln(searchinst,'                (* theeeespot *)');
if piechange then writeln(searchinst,'                (* piechange = true *)');
if firsttime then writeln(searchinst,'                (* firsttime = true *)');
writeln(output,'Bubba 0 in case insttype = 2, distance = ', distance:1);
}
               case insttype of
               1: begin

                     if patterncomplement
                     then begin
                     iwget(searchinst,pie,-fromplace,pieceplace,-toplace,
                           not flip, insttype, true);

{ writeln(output,'first iwget2 not flip'); }
                     end
                     else begin
                     iwget(searchinst,pie,fromplace,pieceplace,toplace,
                           flip, insttype, true);

{ writeln(output,'first iwget2 flip'); }
                     end

                  end;
               2: if distance > 0 then begin
{
writeln(output,'Bubba 1 in case insttype = 2');
}
                  (* this instruction type gives distances between sites *)

                     if patterncomplement
                     then begin
                         iwget2(searchinst,pie,
                         -fromplace,inttopie(lastinternalplace,pie),
                         -toplace,  pieceplace, not flip, true);
{ writeln(output,'second iwget2 not flip'); }
                     end
                     else begin
                         iwget2(searchinst,pie,
                         fromplace,inttopie(lastinternalplace,pie),
                         toplace,  pieceplace, flip, true);
{ writeln(output,'second iwget2 flip'); }
                     end

{
                     if patterncomplement
                     then iwget2(searchinst,pie,
                         fromplace,inttopie(lastinternalplace,pie),
                         toplace,  pieceplace, not flip, true)
                     else iwget2(searchinst,pie,
                         fromplace,inttopie(lastinternalplace,pie),
                         toplace,  pieceplace, flip, true);
}

                  end;
               3: begin
                  (* this instruction type gives the current search string *)

                     if patterncomplement
                     then begin
{zzzbbb to be done}
{
                         iwget2(searchinst,pie,
                         -fromplace,inttopie(lastinternalplace,pie),
                         -toplace,  pieceplace, not flip);
}

                         iwget2(searchinst,pie,
-fromplace,
pieceplace,
-toplace,
inttopie (
pietoint(pieceplace,pie) + patternlength(pattern) - 1,
          pie),
not flip, true);

{writeln(output,'second iwget2 not flip');}
                     end
                     else begin
{
writeln(output,'patternlength(e) = ',patternlength(e):1);
writeln(output,'currentpatternlength = ',currentpatternlength:1);
writeln(output,'patternlength(pattern) = ',patternlength(pattern):1);
}

                         iwget2(searchinst,pie,
                         fromplace, pieceplace,
                         toplace,
inttopie (
pietoint(pieceplace,pie) + patternlength(pattern) - 1,
          pie),

flip, true);
{ writeln(output,'second iwget2 flip'); }
                     end
                  end;
               end;
            end;

            if dosearchfeatures then begin
               (* only put out named features *)
               if thedefinition^.nametag.length > 0 then begin

                  (* put out the definition if this is the first time *)
                  if not featuredefwritten then begin
                     writeln(searchfeatures);
                     (* show the string *)
                     write(searchfeatures,'* ');
                     writestring(searchfeatures, pattern);
                     writeln(searchfeatures);
                     writeadefinition(searchfeatures, thedefinition);
                     writeln(searchfeatures);
                     featuredefwritten := true;
                  end;

                  (* put out the searchfeatures file *)
                  nametostring(pie^.key.hea.keynam, thefeature^.id);
{
2004 July 18:
This compiled before gpcc, but the types are different!!
                  thefeature^.id := pie^.key.hea.keynam;
}

                  (* The orientation of the feature as searched is the same
                     as the orientation of the piece
                     UNLESS WE FLIPPED TO COMPLEMENT!! *)
                  if patterncomplement
                  then begin
                     if pie^.key.piedir = plus
                     then thefeature^.orientation := -1
                     else thefeature^.orientation := +1
                  end
                  else begin
                     if pie^.key.piedir = plus
                     then thefeature^.orientation := +1
                     else thefeature^.orientation := -1
                  end;

                  (* we have to keep track of whether the feature
                     was in the same or a different orientation from
                     the piece so that its end coordinate can be
                     set correctly *)

{
original code:
                  if ((thefeature^.orientation = +1) and
                      (pie^.key.piedir = minus))
                      or
                     ((thefeature^.orientation = -1) and
                      (pie^.key.piedir = plus))
if poundsign then writeln(output,'POUNDSIGN FOUND')
             else writeln(output,'no POUNDSIGN found');
}

                  if (((thefeature^.orientation = +1) and
                       (pie^.key.piedir = minus))
                       or
                      ((thefeature^.orientation = -1) and
                       (pie^.key.piedir = plus))
                     ) and (not poundsign)
                  then begin (* shift it so it covers the right region *)
             (* THERE IS NO NEED TO SHIFT IF THE # MECHANISM IS IN PLAY!
                In that case the exact location has already been determined. *)
{zzzaaa}

{
writeln(output,'~~~~~~thefeature^.definition^.marks:',
thefeature^.definition^.marks:1);
writeln(output,'~~~~~~thefeature^.definition^.locations[1]:',
                      thefeature^.definition^.locations[1]:1:1);
writeln(output,'~~~~~~thefeature^.definition^.locations[thefeature^.definition^.marks]:',
                      thefeature^.definition^.locations[thefeature^.definition^.marks]:1:1);
}
{zzzsss}

                     thefeature^.coordinate := inttopie(internalplace
                          + round(- thefeature^.definition^.locations[1]
                                  + thefeature^.definition^.locations
                                     [thefeature^.definition^.marks])
                            ,pie);
{
writeln(output,'~~~~~~thefeature^.coordinate:',
                      thefeature^.coordinate:1:1);
}


{
writeln(output,'=================>flip');
}
                  end
                  else thefeature^.coordinate := inttopie(internalplace,pie);
{
write(output,'thedefinition:::::::::');
writeadefinition(output, thedefinition);
writeln(output);
---
writeln(output,'thefeature^.coordinate = ', thefeature^.coordinate:1:1);
writeln(output,'thefeature^.definition^.min = ', thefeature^.definition^.min:1:1);
writeln(output,'thefeature^.definition^.max = ', thefeature^.definition^.max:1:1);
writeln(output,'thefeature^.definition^.locations[1] = ',
                thefeature^.definition^.locations[1]:1:1);
writeln(output,'thefeature^.definition^.locations[2] = ',
                thefeature^.definition^.locations[2]:1:1);
writeln(output,'thedefinition^.min = ', thedefinition^.min:1:1);
writeln(output,'thedefinition^.max = ', thedefinition^.max:1:1);
writeln(output,'thedefinition:');
writeadefinition(output, thedefinition);
writeln(output);
writeln(output);
}

{
if patterncomplement then
   thefeature^.orientation := -thefeature^.orientation;
the above is not the right code to use
}

{zq}
		  dofeatures(searchfeatures, mismatches, pattern, mispattern,
                             thefeature, thedefinition, patterncomplement,
                             theplace);
(* writeln(output,'dofeatures ==========================='); {yyy} *)

               end;
            end;

         end;
         (* move to next position in piece *)
         i := i + shift;
         for j := 1 to shift (* use j as an index to advance the position *)
            do ba:=stepbase(startdna,currentdna,cd)
      end;

      bookmatches := bookmatches + patmatches;

      (* handle circular case for distances *)
      if viewd then if (pie^.key.piecon = circular) and (patmatches > 0)
      then begin
         distance := (piecelength(pie) - internalplace) + firstinternalplace;
         writeln(output,'circle closing distance: ',distance:1);
         writeln(result,'* circle closing distance: ',distance:1);
      end;

      if viewm then begin
         write(output,' ',patmatches:1);
         if (patmatches = 1) then writeln(output,' match in piece')
                             else writeln(output,' matches in piece');
         writeln(output);
         write(result,' ',patmatches:1);
         if (patmatches = 1) then writeln(result,' match in piece')
                             else writeln(result,' matches in piece');
         writeln(result,'*************')
      end
end; (* match *)

procedure multi(pattern: string; mismax: integer;
                relatedbase: intarray;
                var org: orgkey; var orgchange: boolean;
                var chr: chrkey; var chrchange: boolean;
                var apiece: pieceptr; var piechange: boolean;
                patterncomplement: boolean); (* if true, reverse delila inst's *)
(* make multiple calls to match, allow mismax mismatches *)

procedure reduce(pattern: string; place: integer; relatedbase: intarray);
(* reduce the pattern by some e's, and recurse *)
var
      subplace: integer;
begin (* reduce *)
      while not (pattern.letters[place] in ['e','.']) do place:=succ(place);
      if pattern.letters[place]='.'
         then match(pattern,mismax,relatedbase, (* no more e's *)
                    org,orgchange,chr,chrchange,apiece,piechange,
                    patterncomplement)
         else begin (* at an e *)
            repeat (* removing e's and recursing *)
               subplace:=place;
               (* find end of e's *)
               while pattern.letters[subplace]='e'
                  do subplace:=succ(subplace);

               (* do everything to the right *)
               reduce(pattern,subplace,relatedbase);

               (* remove one e *)
               compress(pattern,place,relatedbase)
            until pattern.letters[place]<>'e';
            reduce(pattern,place,relatedbase); (* do pattern without e's *)
         end
end; (* reduce *)
begin (* multi *)
      reduce(pattern,1,relatedbase)
end; (* multi *)

procedure multimatch(pattern: string; var book: text; mismax: integer;
      var org: orgkey; (* the current organism *)
      var chr: chrkey; (* the current chromosome *)
      var pie: pieceptr; (* the current piece *)
      var orgchange, chrchange, piechange: boolean; (* that either changed *)
      var orgopen, chropen, pieopen: boolean; (* used by getocp *)
      var sorgchange, schrchange, spiechange: boolean; (* that either changed *)
      var sorgopen, schropen, spieopen: boolean; (* who is open in searchinst *)
          patterncomplement: boolean); (* if true, reverse delila inst's *)
(* try to match the pattern to each piece in book allow mismax mismatches *)
var
   piecesinbook: integer; (* count the number of pieces in the book.
      If there is only one, then set reusethepiece true so that
      we no longer read the book *)
(* ooo *)
{
      org: orgkey; (* the current organism *)
      chr: chrkey; (* the current chromosome *)
      orgchange, chrchange: boolean; (* that either changed *)
      orgopen, chropen, pieopen: boolean; (* used by getocp *)
      pie: pieceptr;
}


procedure showpiece(var afile: text; pie: pieceptr);
(* show the information for the piece *)
begin (* showpiece *)
      with pie^.key do begin
         with hea.keynam do write(afile,' piece: ',letters:length,',');
         if numbered then write(afile,' #',number:1,',');
         write(afile,' configuration: ');
         if piecon=linear then write(afile,'linear,')
                          else write(afile,'circular,');
         write(afile,' direction: ');
         if piedir=plus then write(afile,'+,')
                        else write(afile,'-,');
         write(afile,' begin: ',piebeg: 1,',',
                     ' end: ', pieend:1);
         writeln(afile)
      end
end; (* showpiece *)

begin (* multimatch *)
      reset(book);
{ooo
      new(pie);
      orgopen:=false;
      chropen:=false;
      pieopen:=false;
}
      sorgopen:=false;
      schropen:=false;
      spieopen:=false;

      bookmatches := 0;
      theline := 1;

      piecesinbook := 0;
      while ((not reusethepiece) and (not eof(book))) or
            (reusethepiece and (piecesinbook < 1))
      do begin
{writeln(output,'multimatch 1: piechange=',piechange);}
         if not reusethepiece
         then getocp(book,theline,org,orgchange,orgopen,
                     chr,chrchange,chropen,
                     pie,piechange,pieopen)
         else begin
{writeln(output,'all change variables MADE FALSE'); zzz}
{christmas zzzccc
the effect of deleting these is that the org/chr/pie
insts are now given when not needed!
}
            orgchange := false;
            chrchange := false;
            piechange := false
         end;
{zzz writeln(output,'multimatch 2: piechange=',piechange);}
         if not eof(book) then begin
            piecesinbook := succ(piecesinbook);

            if viewi then begin
               (* show organism and chromosome *)

{OOO}
               iworgchr(output,org,sorgchange,sorgopen,
                               chr,schrchange,schropen);
               if sorgchange or schrchange then write(result,'* ');

{
               (* 2001 Sep 7 *)
               if sorgchange or schrchange
               then spiechange := true
               else spiechange := false;
}

               iworgchr(result,org,sorgchange,sorgopen,
                               chr,schrchange,schropen);

               (* show piece *)
               showpiece(output,pie);
               writeln(result,'*'); (* put space from the previous one *)
               write(result,'*'); (* removed space 1990 Sep 4 *)
               showpiece(result,pie);
            end;

(* write the org, chr and pie instructions ONLY once the first
one has been found!!!  So the following is not used: ...
2012 May 08 add it back and update it.
*)
            if printinginstructions then begin
{
writeln(searchinst,'(','* iworgchr being called from multimatch zzz *',')');
}

(*
writeln(searchinst,'{sorgchange: ', sorgchange ,' BUBBA}');
writeln(searchinst,'{schrchange: ', schrchange ,' BUBBA}');
writeln(searchinst,'{sorgopen  : ', sorgopen   ,' BUBBA}');
writeln(searchinst,'{schropen  : ', schropen   ,' BUBBA}');
*)
               sorgopen := not sorgopen;
               schropen := not schropen;
               iworgchr(searchinst,org,sorgchange,sorgopen,
                                   chr,schrchange,schropen);
               if piechange then iwpie(searchinst,pie^.key);
            end;

            if orgchange then clearheader(org.hea);
            if chrchange then clearheader(chr.hea);

            multi(pattern,mismax,relatedbase,
                  org,orgchange,chr,chrchange,
                  pie,piechange, (* pass to allow silence when not found *)
                  patterncomplement);

(* the old location: clearpiece(pie) *)
         end
      end;

      if not reusethepiece then if piecesinbook = 1 then begin
         reusethepiece := true;
         writeln(output,'(reusing the single piece in the book for efficiency)');
         reset(book); (* fool it into doing the loop above ... *)
      end;

(* the old location: dispose(pie); *)
      if viewb then begin
         write(output,' ',bookmatches:1);
         if (bookmatches = 1) then writeln(output,' match in book')
                              else writeln(output,' matches in book');
         writeln(output);
         write(result,'* ',bookmatches:1);
         if (bookmatches = 1) then writeln(result,' match in book')
                              else writeln(result,' matches in book');
         writeln(result,'**********************') (* space result file a bit *)
      end
end; (* multimatch *)

procedure instheader;
(* make header for instruction file *)
begin (* instheader *)
   rewrite(searchinst);
   write(searchinst,'title "');
   writedatetime(searchinst,datetime);
   writeln(searchinst,' search ',version:4:2,'";');
   iwcombk(book,searchinst);

   (* no longer needed - these are defaults for delila
   writeln(searchinst,'default numbering piece;');
   writeln(searchinst,'default numbering 1;');
   *)
   writeln(searchinst,'set out-of-range reduce-range;');
end; (* instheader *)

procedure copyfile(var fromfile, tofile: text; space: boolean);
(* copy one file into the other.
   no resets or rewrites are done.
   if space is true, then put an '* ' before each line. *)
begin (* copyfile *)
      while not eof(fromfile) do begin
         if space then write(tofile,'* ');
         copyaline(fromfile, tofile)
      end
end; (* copyfile *)

procedure getmismatches;
(* obtain a new value of the global variable, mismatches *)
var   newmismatches: integer; (* the value the person typed *)
begin (* getmismatches *)
      answers := 0; 
      repeat
         if nostring(buffer) then writeln(output,'mismatches were ',
                                  mismatches:1,', type new:');
         getinteger(input,buffer,newmismatches,gotten);
         if gotten then if newmismatches < 0 then begin
            writeln(output,'type zero or a positive integer');
            clearstring(buffer);
            gotten := false
         end;
         answercheck(answers);
      until gotten;
      mismatches := newmismatches;
      writeln(result,'* maximum number of mismatches now allowed: ',
                       mismatches:1)
end; (* getmismatches *)

{ppp}
procedure getcolorbackground;
(* obtain a new value of the color background for letters.
This is called by the C command *)
var newTparam: char; (* the value the person typed *)
    newparam: real; (* the value the person typed *)
begin (* getcolorbackground *)
   with thefeature^ do begin

      (* Tparam *)
      answers := 0; 
      repeat
         if nostring(buffer) then writeln(output,'Tparam was "',
                                  Tparam:1,'", type new (allowed: "hrHR"):');
         getchar(input,buffer,newTparam,gotten);
         if gotten then if not (newTparam in ['h','r','H','R'])
          then begin
            writeln(output,'Tparam must be one of " hrHR"');
            clearstring(buffer);
            gotten := false
         end;
         answercheck(answers);
      until gotten;
      Tparam := newTparam;
      writeln(result,'* Tparam colorbackground: ',
                       Tparam:1);

      (* Aparam *)
      answers := 0;
      repeat
         if nostring(buffer) then writeln(output,'Aparam was ',
                                  Aparam:infofield:infodecim,', type new:');
         getreal(input,buffer,newparam,gotten);
         if gotten then if newparam < 0 then begin
            writeln(output,'type zero or a positive integer');
            clearstring(buffer);
            gotten := false
         end;
         answercheck(answers);
      until gotten;
      Aparam := newparam;
      writeln(result,'* Aparam colorbackground: ',
                       Aparam:infofield:infodecim);

      (* Bparam *)
      answers := 0;
      repeat
         if nostring(buffer) then writeln(output,'Bparam was ',
                                  Bparam:infofield:infodecim,', type new:');
         getreal(input,buffer,newparam,gotten);
         if gotten then if newparam < 0 then begin
            writeln(output,'type zero or a positive integer');
            clearstring(buffer);
            gotten := false
         end;
         answercheck(answers);
      until gotten;
      Bparam := newparam;
      writeln(result,'* Bparam colorbackground: ',
                       Bparam:infofield:infodecim);

      (* Cparam *)
      answers := 0;
      repeat
         if nostring(buffer) then writeln(output,'Cparam was ',
                                  Cparam:infofield:infodecim,', type new:');
         getreal(input,buffer,newparam,gotten);
         if gotten then if newparam < 0 then begin
            writeln(output,'type zero or a positive integer');
            clearstring(buffer);
            gotten := false
         end;
         answercheck(answers);
      until gotten;
      Cparam := newparam;
      writeln(result,'* Cparam colorbackground: ',
                       Cparam:infofield:infodecim);

      (* Dparam *)
      answers := 0;
      repeat
         if nostring(buffer) then writeln(output,'Dparam was ',
                                  Dparam:infofield:infodecim,', type new:');
         getreal(input,buffer,newparam,gotten);
         if gotten then if newparam < 0 then begin
            writeln(output,'type zero or a positive integer');
            clearstring(buffer);
            gotten := false
         end;
         answercheck(answers);
      until gotten;
      Dparam := newparam;
      writeln(result,'* Dparam colorbackground: ',
                       Dparam:infofield:infodecim);

   end;
end; (* getcolorbackground *)

{
un built procedure - could dump 2004 jul 18
procedure getcolorletters;
(* obtain a new value of the color letters for letters.
This is called by the C command *)
var   newcolorletters: integer; (* the value the person typed *)
begin (* getcolorletters *)
      answers := 0;
      repeat
         if nostring(buffer) then writeln(output,'colorletters were ',
                                  colorletters:1,', type new:');
         getinteger(input,buffer,newmismatches,gotten);
         if gotten then if newmismatches < 0 then begin
            writeln(output,'type zero or a positive integer');
            clearstring(buffer);
            gotten := false
         end;
         answercheck(answers);
      until gotten;
      colorletters := newcolorletters;
      writeln(result,'* maximum number of colorletters now allowed: ',
                       colorletters:1)
end; (* getcolorletters *)
}

procedure toggleoktosearch;
(* toggle the oktosearch variable *)
begin
   oktosearch := not oktosearch;
   if oktosearch
   then begin
      writeln(result,'* searches will proceed when pattern is input',
                     ' or with the = command');
      writeln(output,'* searches will proceed when pattern is input',
                     ' or with the = command');
   end
   else begin
      writeln(result,'* searches will proceed only when = is input');
      writeln(output,'* searches will proceed only when = is input');
   end
end;

procedure acomment;
(* the rest of the input buffer is a commment.  the entire line
must be written to the result, and perhaps the searchinst files *)
begin (* acomment *)
      (* send it to result file *)
      write(result,'* ');
      writestring(result,buffer);
      writeln(result);

      if printinginstructions then begin
         (* send it to the searchinst file *)
         write(searchinst,'(* ');
         writestring(searchinst,buffer);
         writeln(searchinst,' *)')
      end;

      (* remove the rest of the line *)
      clearstring(buffer)
end; (* acomment *)

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

(* Routines for handline features *)

procedure setfeaturename;
(* set the feature name: define the feature and prepare for
writing to file searchfeatures. *)
var
   i: integer; (* location in a string *)
begin
   if not dosearchfeatures then begin
      rewrite(searchfeatures);
      writeln(searchfeatures,'* search ',version:4:2);
      dosearchfeatures := true;
   end;
   featuredefwritten := false;

   with thedefinition^ do begin
      (* transfer the string from pattern *)
      (* remove final quote (if there is one): *)
      if pattern.letters[pattern.length] = '"'
      then pattern.length := pattern.length-1;

      for i := 2 to pattern.length
      do nametag.letters[i-1] := pattern.letters[i];
      nametag.length := pattern.length -1;

      (* tell the user about this *)
      if nametag.length > 0 then begin
         {write(output,'feature: "'); bug 2006 May 09}
         write(output,'feature: ');
         writequotestring(output,nametag);
         write(output,' will be written to searchfeatures');
         if doarrow then writeln(output,' as an arrow.');
         writeln(output);
         write(result,'* feature: '); {bug 2006 May 09}
         writequotestring(result,nametag);
         write(result,' will be written to searchfeatures'); {bug 2006 May 09}
         if doarrow then writeln(result,' as an arrow.');
         writeln(result);
      end
      else begin
         writeln(output,'No more features will be written to searchfeatures');
         writeln(result,'* No more features will be written to searchfeatures');
      end;

      if doarrow then begin
         background.letters[1] := '-';
         background.length := 1;
         negparts.length := 2;
         negparts.letters[1] := '<';
         negparts.letters[2] := ']';
         posparts.length := 2;
         posparts.letters[1] := '[';
         posparts.letters[2] := '>';
         doarrow := false (* reset *)
{zzzaaa}
{zzzfff}
      end
      else begin
         background.letters[1] := '-';
         background.length := 1;
         negparts.length := 2;
         negparts.letters[1] := '[';
         negparts.letters[2] := ']';
         posparts.length := 2;
         posparts.letters[1] := '[';
         posparts.letters[2] := ']';
{zzzNNN}
      end;

      locations[1] := 0;
      locations[2] := 1;
      marks := 2;
      min := 0.0; {zzz?}
      max := 0.0; {zzz?}
      number := 1;
(*[[*)
      matrix := nil;
{
      Ribound := 0;
      Zbound := 0;
      Pbound := 0;
}
(*]]*)
      next := nil;

   end;
{
}
{zzzaaa}
   writeln(searchfeatures);
   writeln(searchfeatures,'* THE CURRENT FEATURE DEFINITION IS:');
   write(searchfeatures,'* ');
   writeadefinition(searchfeatures, thedefinition);
   writeln(searchfeatures);

(* 2004 Jul 18: I THINK this is a good place to clear the feature! *)
{
   clearfeature(thefeature);
}
   with thefeature^ do begin
      with id do begin
         letters[1] := '!';
         length := 1;
      end;
{ This is handled in clearfeature:
      coordinate := 0.0;
      orientation := +1;
}
      thefeature^.nametag := thedefinition^.nametag;
{ This is handled in clearfeature:
      with othertag do begin
         letters[1] := '?';
         length := 0;
      end;
}
      definition := thedefinition;
{ This is handled in clearfeature:
(*[[*)
      Ri := 0.0;
      Z := 0.0;
      probability := 0.0;
(*]]*)
      unsatisfied := true;
      fromrange := 0.0;
        torange := 0.0;
      number := 1;
      desiredline := 1;
      next := nil;
}
      number := 1;
      desiredline := 1;
   end;

{
writeln(searchfeatures);
write(searchfeatures,'* INITIAL FEATURE: ');
write(searchfeatures,'* ');
writeafeature(searchfeatures,thefeature);
writeln(searchfeatures);
}

    featureinserts := false; (* inserts to the feature not defined yet *)
{zzzfff}
end;

procedure setarrowfeature;
(* reset the feature to be an arrow *)
begin
   doarrow := not doarrow;
{zzzaaa}
{
   if doarrow
   then writeln(output,'Doing arrows for features')
   else writeln(output,'Doing rectangles for features')
}
end;

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

procedure writeDelilainst;
(* a search subsystem for creating Delila instructions based
   on search results.  see Delilahelp. *)
var   copy: text;  (* internal file for showing searchinst *)
begin (* writeDelilainst *)
      writeln(output,'Entering Delila sub system');
      if not printinginstructions then begin
         (* if we were not printing instructions, then we must make
            sure that instructions printed from this point have
            the org, chr or pie information.  We force it by
            resetting the variables. *)
{hhh}

         orgchange := true;
         chrchange := true;
         piechange := true;

         instheader;
         writeln(output,'Any previous instructions were destroyed.');
         writeln(output,'A title has been provided.');
         writeln(output,'Defaults have been set.');
      end
      else begin
         writeln(output,'The searchinst file is not empty.');
         writeln(output,'use the show command to see them and');
         writeln(output,'use the kill command to restart them');
      end;

      printinginstructions:=true;
      writeln(output,'Delila instruction printing',
                     ' has been automatically turned ON.');

      writeln(output,'from place is ', fromplace:1);
      writeln(output,'  to place is ', toplace:1);

      repeat
         answers := 0;
         repeat
            if nostring(buffer) then write('Delila ');
            getchar(input,buffer,command,gotten);
            answercheck(answers);
         until gotten;

         if not (command in ['q','p','f','t','w','s','k','h','i','*'])
         then begin
            flagstring(output, buffer);
            answers := 0;
            repeat
               if nostring(buffer)
                  then writeln(output,'do you need help? (y/n)');
               getchar(input,buffer,command,gotten);
               answercheck(answers);
            until gotten;

            if command<>'n' then Delilahelp(output)
         end
         else if command='h' then Delilahelp(output)

         else if command='i' then begin
            answers := 0;
            repeat
               if nostring(buffer)
               then writeln(output,'instruction type was ',
                                   insttype:1,', type new:');
               getinteger(input,buffer,insttype,gotten);
               if gotten then if
                  ((insttype<1) and (insttype>3))
               then writeln(output,' instruction type must be',
                                   ' 1, 2 or 3 ');
               answercheck(answers);
            until gotten and ((insttype=1) or (insttype=2)
                                           or (insttype=3));
            writeln(result,'* instruction type: ',insttype:1);
            writeln(output,'* instruction type: ',insttype:1);
         end

         else if command='p' then begin
            printinginstructions:=not printinginstructions;
            write(output,'now');
            if not printinginstructions then write(output,' NOT');
            writeln(output,' printing Delila instructions');
            write(result,'* now');
            if not printinginstructions then write(result,' NOT');
            writeln(result,' printing Delila instructions')
         end
         else if command='f' then begin
            answers := 0;
            repeat
               if nostring(buffer)
               then writeln(output,'from place was ',
                                   fromplace:1,', type new:');
               getinteger(input,buffer,fromplace,gotten);
               answercheck(answers);
            until gotten;
            writeln(result,'* bases from: ',fromplace:1);
         end
         else if command='t' then begin
            answers := 0;
            repeat
               if nostring(buffer)
               then writeln(output,'to place was ',
                                   toplace:1,', type new:');
               getinteger(input,buffer,toplace,gotten);
               answercheck(answers);
            until gotten;
            writeln(result,'* bases to: ',toplace:1);
         end
         else if command='w' then begin
            writeln(output,'write Delila instructions (q to quit):');
            repeat
               readstring(input, buffer); (* pick up a fresh line *)
               command := buffer.letters[1];

               if command<>'q' then begin
                  writestring(searchinst,buffer);
                  writeln(searchinst);

                  write(result,'* ');
                  writestring(result,buffer);
                  writeln(result)
               end
            until command='q';
            gettoken(input,buffer,pattern,gotten);
               (* that gets us past the 'quit', pattern is not used so we
                  are safe, and the user may type commands after the quit *)
            command:='w' (* prevent falling out of Delila *)
         end
         else if command='s' then begin
            writeln(output,'current Delila instructions:');
            writeln(result,'* current Delila commands:');
            reset(searchinst);
            rewrite(copy);
            copyfile(searchinst,copy,false);

            (* get back searchinst: *)
            rewrite(searchinst);
            reset(copy);
            copyfile(copy,searchinst,false);

            (* show copy: *)
            reset(copy);
            copyfile(copy,output,false);

            reset(copy);
            copyfile(copy,result,true);

            command:='s'; (* prevent falling out *)
         end
         else if command='k' then begin
            answers := 0; 
            repeat
               if nostring(buffer) then
                  writeln(output,'shall I kill all Delila instructions?');
               getchar(input,buffer,command,gotten);
               answercheck(answers);
            until gotten;
            if command='y' then begin
               instheader;
               writeln(result,'* all previous Delila instructions killed');
               writeln(output,'aaaARRRRGGgg!!')
            end
            else writeln(output,'Scardy Cat!')
         end
         else if command='*' then acomment;

         if command = 'q' then begin
            if fromplace > toplace then flip := true
                                   else flip := false
         end
      until command='q';
      command:='d'; (* prevent falling out of search *)
      write(output,'Delila instruction printing is ');
      if printinginstructions
      then writeln(output,'ON')
      else writeln(output,'OFF')
end; (* writeDelilainst *)

procedure setstartshift;
(* set the search startpoint and shift between search positions *)
begin (* setstartshift *)
      answers := 0;
      repeat
         if nostring(buffer)
         then writeln(output,'start point was ',startpoint:1,
                  ' type a new one:');

         getinteger(input,buffer,startpoint,gotten);
         if gotten then if startpoint <= 0 then begin
            writeln(output,'give a positive integer');
            gotten := false
         end;
         answercheck(answers);
      until gotten;

      answers := 0;
      repeat
         if nostring(buffer)
         then writeln(output,'shift was ',shift:1,
                  ' type a new one:');

         getinteger(input,buffer,shift,gotten);
         if gotten then if shift <= 0 then begin
            writeln(output,'give a positive integer');
            gotten := false
         end;
         answercheck(answers);
      until gotten;

      writeln(result,'* start (phase): ',startpoint:1,
                     ' shift: ',shift:1)
end; (* setstartshift *)

procedure viewvariables;
(* set the characteristics one wants to see *)
var
      ch: char; (* user typed command character *)
begin (* view *)
      answers := 0;
      repeat
         if nostring(buffer)
         then writeln(output,'view (type h for help)');
         getchar(input,buffer,ch,gotten);
         if gotten then gotten :=
            (ch in ['t','e','i','f','s','p','d','m','b','h','n','a','q']);
         answercheck(answers);
      until gotten;

           if ch = 'h' then viewhelp(output)
      else if ch = 'a' then begin
         viewt := true; viewe := true; viewi := true;
         viewf := true; views := true; viewp := true;
         viewd := true; viewm := true; viewb := true;
      end
      else if ch = 'n' then begin
         viewt := false; viewe := false; viewi := false;
         viewf := false; views := false; viewp := false;
         viewd := false; viewm := false; viewb := false;
      end
      else case ch of
         't': viewt := not viewt;
         'e': viewe := not viewe;
         'i': viewi := not viewi;
         'f': viewf := not viewf;
         's': views := not views;
         'p': viewp := not viewp;
         'd': viewd := not viewd;
         'm': viewm := not viewm;
         'b': viewb := not viewb;
         'q': ;
      end;
      viewnothing := not (viewt or viewe or viewi
                       or viewf or views or viewp
                       or viewd or viewm or viewb)
end; (* view *)

begin (* search *)

      rewrite(result);

      (* get the date *)
      getdatetime(datetime);

      write(output,'* ');
      writedatetime(output,datetime);
      writeln(output,' search ',version:4:2);
      write(result,'* ');
      writedatetime(result,datetime);
      writeln(result,' search ',version:4:2);

      brinit(book, theline);
      reset(book);
      copyaline(book,output);
      reset(book);
      copyaline(book,result);

(* don't touch the instructions until Delila subsystem starts!!!
      instheader;
*)

      clearstring(pattern);
      clearstring(lastpattern);

      (* defaults for user settable variables *)
      viewt := true; viewe := true; viewi := true;
      viewf := true; views := true; viewp := true;
      viewd := true; viewm := true; viewb := true;
      viewnothing:=false;
      mismatches:=0;
      oktosearch := true;
      startpoint := 1;
      shift := 1;
      printinginstructions:=false;
      fromplace:=-100;
      toplace:=+100;
      insttype:=1;
      firstinput := true;
      dosearchfeatures := false;
      featuredefwritten := false;
      flip := false;

      patterncomplement := false;

      (* ooo *)
      org.hea.keynam.length := 0;
      chr.hea.keynam.length := 0;
      pie := nil;
{OOO}
      orgopen:=false;
      chropen:=false;
      pieopen:=false;
      orgchange:=true;
      chrchange:=true;
      piechange:=true;

      sorgopen:=false;
      schropen:=false;
      spieopen:=false;
      sorgchange:=true;
      schrchange:=true;
      spiechange:=true;

      theline:=0;

      (* initialize feature names *)
      new(thefeature);
      new(thedefinition);
      clearfeature(thefeature);

      reusethepiece := false;

      clearstring(buffer);
      writeln(output,' NOTE: Your first typed line will be ignored.');
      repeat
         answers := 0;
         repeat
            if firstinput then begin
               write(output,'(type a carriage return now.',
                            '  Your next line will be listened to.)');
               firstinput := false;
               reset(input);
            end
            else if nostring(buffer) then write(output,'search ');
            gettoken(input, buffer, pattern, gotten);
            answercheck(answers);
         until gotten;
         command := pattern.letters[1];

              if command='D' then writeDelilainst
{ppp}
         else if command='C' then getcolorbackground
         else if command='H' then searchhelp(output)
         else if command='?' then searchhelp(output)
         else if command='L' then showletters(output)
         else if command='M' then getmismatches
         else if command='N' then toggleoktosearch
         else if command='P' then setstartshift
         else if command='Q' then writeln(output,'quit')
         else if command='q' then writeln(output,'quit')
         else if command='V' then viewvariables
         else if command='"' then setfeaturename
         else if command='A' then setarrowfeature
{zzzaaa}
        else if command='*' then acomment
         else if viewnothing then begin
            flagstring(output,buffer);
            writeln(output);
            writeln(output,' your search command was ignored because');
            writeln(output,' no results would be displayed.');
            writeln(output,' use the view command to turn displays on.');
         end

         else if command='~' then begin
            if lastpattern.length <> 0
            then begin
               stringcomplement(lastpattern);
               patterncomplement := not patterncomplement;
               if patterncomplement
               then writeln(output,'now complement to typed pattern')
               else writeln(output,'now homologous to typed pattern')
            end
            else writeln(output,'no previous pattern')
         end

         else if command='R' then begin
            if lastpattern.length <> 0
            then begin
               stringreverse(lastpattern);
{qqq}
            end
            else writeln(output,'no previous pattern')
         end

         else if command='I' then begin
            if lastpattern.length <> 0
            then begin
               stringinvert(lastpattern);
{qqq}
            end
            else writeln(output,'no previous pattern')
         end

         else if command='=' then begin
            if lastpattern.length <> 0
            then begin
               analysepattern(lastpattern, exppattern, ok);
               if ok then multimatch(exppattern,book,mismatches,
                                     org, chr, pie,
                                     orgchange, chrchange, piechange,
                                     orgopen, chropen, pieopen,
                                     sorgchange, schrchange, spiechange,
                                     sorgopen, schropen, spieopen,
                                     patterncomplement)
            end
            else writeln(output,'no previous pattern')

            (* these are the other allowed first letters.
               note that some commands above exclude their use *)
         end
         else if not (
              (command in
                [ '(',
                  '0','1','2','3','4','5','6','7','8','9',
                  'a','c','g','t',
                  'm','r','s','v','w','y','h','k','d','b','n','e' ]) or
                (command='#') or (command='%') or
                (command='<') or (command='>') or
                (command='|') or (command='^') or
                (command='[') or (command=']')
              )
         then begin
            flagstring(output, buffer);
            answers := 0;
            repeat
               if nostring(buffer) then
                  writeln(output,'do you need help? (y/n)');
               getchar(input,buffer,command,gotten);
               answercheck(answers);
            until gotten;
            if command<>'n' then searchhelp(output);
            command := ' '; (* this prevents the program
               from stopping if the answer was 'q' *)
         end

         else begin
            analysepattern(pattern,exppattern,ok);

            if ok then begin
               patterncomplement := false; (* reset at type in *)
               lastpattern := pattern;
               if oktosearch then multimatch(exppattern,book,mismatches,
                          org, chr, pie,
                          orgchange, chrchange, piechange,
                          orgopen, chropen, pieopen,
                          sorgchange, schrchange, spiechange,
                          sorgopen, schropen, spieopen,
                          patterncomplement)
            end;
         end
      until command='q';
1: end.  (* search *)
