program sorth(hlist, shlist, list, sorthp, output);
(* sorth: sort helix list
   thomas schneider

module libraries required: delman, prgmods, delmods, auxmods
*)

label 1; (* end of program *)

const
(* begin module version *)
version = 2.42; (* of sorth.p 1995 July 26
origin: 1984 feb 4
(* end module version *)

(* begin module describe.sorth *)
(*
name
      sorth: sort helix list

synopsis
      sorth(hlist: in, shlist: out, list: out, sorthp: in, output: out)

files
      hlist: a list of helixes generated from program helix.
      shlist: a list of helixes, where the longest or strongest helix
         has been chosen from each piece to piece comparison ('set').
      list: progress of the program.
      sorthp: parameters to control the program.
         1. characters on the first line of the file determine the priority
         order for sorting the helixes.  all commands must end with 'a' to
         indicate 'ambiguous'.  the commands are:
            ea - sort on energies (see technical notes)
            la - sort on lengths (see technical notes)
            ela - sort first on energies then on lengths.
            lea - sort first on lengths then on energies.
         2. the second line of the file must contain one integer,
         'top'.  up to 'top' of the strongest helixes will be
         written to shlist. if 'top' = 1, then any set of helixes that are
         ambiguous are not copied to the shlist.  this allows
         one to find the strongest unambiguous helix in each set.
         3. the third line is the minimum length or maximum energy
         of helixes to be sorted.
      output: messages to the user.

description
      the strongest helixes in hlist are sorted and copied to shlist.
      the user can sort on energy, length, energy then length, or
      length then energy.  the user may chose more than one helix
      to be output (eg, the top 10).

see also
      helix.p

author
      thomas dana schneider

bugs
      none known

technical notes
      when only one variable is sorted on, the order of the other
      variable will not be meaningful because it is determined by
      the way the sort algorithm works.
      the constant 'maxhelix' determines the maximum number of
      helixes that can be sorted.
*)
(* end module describe.sorth *)

(* begin module sorth.const *)
      prime = '"'; (* single quote mark *)
      lowpriority = 3; (* the number of priority levels: 1 to this.
         see the operations type *)
      maxhelix = 1000; (* one more than the largest number of helixes the
         program can handle. *)
(* end module sorth.const *)

type
(* begin module sorth.type *)
      operations = record
         (* the order that the sorting operations should take place.
         the priority array gives the order.  e = energy, l = length,
         a = ambiguous *)
         priority: array[1..lowpriority] of char;
         low: integer; (* the lowest priority used in the priority array *)
         getenergy: boolean; (* whether or not energy is to be gotten *)
         top: integer; (* the top 'top' helixes are printed to shlist *)
         minlenmaxene: real; (* the minimum length or maximum energy of
            helixes to be accepted from hlist *)
         doinglength: boolean; (* true when minlenmaxene is positive *)
      end;

      position = 0..maxhelix; (* the positions of helixes *)
      helix = record (* information about a helix *)
         x,y: integer; (* the coordinates of a helix *)
         length: integer; (* the length of a helix *)
         energy: real; (* the energy of a helix *)
      end;
(* end module sorth.type *)

var
(* begin module sorth.var *)
      hlist, (* helix list from helix *)
      shlist, (* sorted list of helixes *)
      list, (* progress of the program *)
      sorthp: (* parameters of the program *)
         text;

      (* globals used by the lessthan, swap and doset *)
      ordering: array [1..maxhelix] of position; (* the order of the helixes *)
      helixes: array [1..maxhelix] of helix; (* the helixes in the order that
         they were read in *)
      sorton: char; (* what to sort the helixes on (by reordering the
         ordering array) *)
      both: boolean; (* true means that sorting should be done
         on both length and energy, using the value of
         sorton to tell which to do first. *)
(* end module sorth.var *)

(* 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 = 'delmod 6.51 85 apr 17 tds/gds' *)

(* 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 = 'delmod 6.51 85 apr 17 tds/gds' *)

(* 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 = 'delmod 6.51 85 apr 17 tds/gds' *)

(* ************************************************************************ *)
(* end module package.primitive version = 'delmod 6.51 85 apr 17 tds/gds' *)

(* begin module findcolon *)
procedure findcolon(var thefile: text);
(* move the file to the characters just past ': ' (colon, space) *)
var
      found: boolean;
procedure die;
begin writeln(output,' no helix list data'); halt end;
begin
      found:=false;
      while not found do begin
         if eof(thefile)
            then die
            else if thefile^ = ':'
               then begin
                  get(thefile);
                  if eof(thefile)
                     then die
                     else if thefile^ = ' '
                        then begin
                           get(thefile);
                           found:=true
                        end
               end
               else get(thefile)
      end
end;
(* end module findcolon version = 'auxmod 1.37 85 apr 4 gds/tds'; *)

(* begin module gethelix *)
procedure gethelix(var hlist: text;
                   var x,y: integer;
                   var length: integer;
                       getenergy: boolean;
                   var energy: real;
                   var done: boolean);
(* get the (x,y,length) info on some helix listed in hlist.
if the x,y pair changes, then done is true and no values are
returned.  checknames may be used at this point to start the
next pair.  if getenergy is true, then the routine expects that there
will be an energy on the helix line *)
begin (* gethelix *)
      if eof(hlist)
         then done:=true
         else begin
            get(hlist); (* skip the blank *)
            if hlist^ <> ' '
               then done:=true
               else begin
                  done:=false;
                  findcolon(hlist); read(hlist,x);
                  findcolon(hlist); read(hlist,y);
                  findcolon(hlist); read(hlist,length);
                  if getenergy then read(hlist, energy)
                               else energy := 0;
                  readln(hlist)
               end
         end
end; (* gethelix *)
(* end module gethelix version = 'auxmod 1.37 85 apr 4 gds/tds'; *)

(* begin module sorth.sortfunctions *)
function lessthan(a, b: position): boolean;
(* see quicksort *)
var
      ea, eb: real; (* the absolute value of the energy values
         at positions a and b *)
      la, lb: integer; (* the length values at a and b *)
function byenergy: boolean;
(* give the function less than using energies *)
begin (* byenergy *)
      ea :=   abs(helixes[ordering[a]].energy);
      eb :=   abs(helixes[ordering[b]].energy);
      byenergy := ea < eb
end; (* byenergy *)
function bylength: boolean;
(* give the function less than using lengths *)
begin (* bylength *)
      la := helixes[ordering[a]].length;
      lb := helixes[ordering[b]].length;
      bylength := la < lb
end; (* bylength *)
begin (* lessthan *)
      case sorton of
         'e': if byenergy
              then lessthan := true
              else if both
                 then if ea = eb
                    then lessthan := bylength
                    else lessthan := false
                 else lessthan := false;
         'l': if bylength
              then lessthan := true
              else if both
                 then if la = lb
                    then lessthan := byenergy
                    else lessthan := false
                 else lessthan := false
      end
end; (* lessthan *)

procedure swap(a, b: position);
(* see quicksort *)
var
      hold: position; (* for swapping *)
begin (* swap *)
      hold := ordering[a];
      ordering[a] := ordering[b];
      ordering[b] := hold
end; (* swap *)
(* end module sorth.sortfunctions *)

(* begin module quicksort *)

procedure quicksort(left, right: position);
(* quick sort a list between positions left and right, into ascending order.
   a position is simply a scalar of the form 0..max.
   the array to be sorted is dimensioned 1..max.
(the difference in the ranges is important to the correct operation
of the sort...)
   two external routines are used:
      function lessthan(a, b: position): boolean is a generalized test for
          value-at-a < value-at-b.
      procedure swap(a, b: position) switches the items at positions a and b.
   since these routines are external, the procedure is general.

   this procedure taken from the book
   'algorithms + data structures = programs' by niklaus wirth,
   prentice-hall, inc., englewood cliffs, n.j.(1976), pp. 76-82 *)
var
      lower, upper: position; (* the positions looked at currently *)
      center: position; (* the rough center of the region being sorted *)
begin
      lower  := left;
      center := (left + right) div 2;
      upper  := right;

      repeat
         while lessthan(lower, center) do lower := succ(lower);
         while lessthan(center, upper) do upper := pred(upper);

         if lower <= upper then begin
            (* keep track of the center through the map: *)
                 if lower = center then center:=upper
            else if upper = center then center:=lower;
            swap(lower, upper);
            lower := succ(lower);
            upper := pred(upper)
         end
      until lower > upper;

      if  left < upper then quicksort(left,  upper);
      if lower < right then quicksort(lower, right)
end;

(* end module quicksort version = 'prgmod 3.97  85 may 5 tds'; *)

(* begin module sorth.hcopy *)
procedure hcopy(var hlist, shlist, list: text);
(* hlist copy.  copy the header of hlist to shlist and list. *)
var
      dummy: integer; (* to capture output of copylines *)
begin (* hcopy *)
      reset(hlist);
      rewrite(shlist);
      write(shlist,' sorth ',version:4:2,' sorted helixes from ');
      (* do a minimal check that the hlist is reasonable: *)
      if copylines(hlist,shlist,8) <> 8 then begin
         writeln(output, 'hlist does not have a header');
         halt
      end;

      reset(hlist);
      rewrite(list);
      write(list,' sorth ',version:4:2,' sorted helixes from ');
      dummy := copylines(hlist,list,7);
      readln(hlist) (* skip first header line for gethelix *)
end; (* hcopy *)
(* end module sorth.hcopy *)

(* begin module sorth.getparams *)
procedure getparams(var hlist: text;
                    var sorthp: text; var order: operations;
                    var list: text);
(* get whether energies are in the hlist, and learn from sorthp
what priority order should be used.  learn the top number of helixes
to print.  report the results to list. *)
var
      i: integer; (* general index *)
      doe, dol: integer; (* counts of the priorities 'e' and 'l' *)
begin (* getparams *)
      (* first pick up data from hlist *)
      reset(hlist);
      for i:=1 to 5 do readln(hlist); (* move to energy data *)
      get(hlist); (* skip the blank *)
      (* determine if one has the potential to sorton energies: *)
      if hlist^ = 'e' then order.getenergy := true
                      else order.getenergy := false;
      (* move to first set of helixes *)
      for i:=1 to 3 do readln(hlist);

      (* read the order and see if it is consistant with getenergy *)
      reset(sorthp);
      if eof(sorthp) then begin
         writeln(output,'sorthp parameter file is empty');
         halt
      end;
      with order do begin
         low := 0; doe := 0; dol := 0;
         while (not eoln(sorthp)) and (low < lowpriority) do begin
            if sorthp^=' '
            then get(sorthp) (* skip spaces *)
            else begin
               low := low + 1;
               read(sorthp,priority[low]);
               if priority[low] = 'e' then doe := doe + 1
               else if priority[low] = 'l' then dol := dol + 1
               else if priority[low] <> 'a' then begin
                  writeln(output,'priority "',priority[low],
                                 '" is not allowed');
                  halt
               end
            end
         end;
         if (doe > 1) or (dol > 1) then begin
            writeln(output,'eea and lla are not allowed parameters');
            halt
         end;
         if low < 2 then begin
            writeln(output,'priority list is too short in file sorthp');
            halt
         end;
         if priority[low] <> 'a' then begin
            writeln(output,'priority list in file sorthp must end with an "a"');

            halt
         end;
            (* check that priority does not depend on energy when
               energy is not in the hlist *)
         if not getenergy then begin
            i := 1;
            while (i < low) and (priority[i] <> 'e') do i := i + 1;
            if priority[i] = 'e' then begin
               writeln(output,'sorthp requests sorting by energy');
               writeln(output,'but there are no energies in hlist');
               halt
            end
         end;
         readln(sorthp);

         (* determine top *)
         if eof(sorthp)
         then begin
            top := 1;
            minlenmaxene := 0.0
         end
         else begin
            readln(sorthp, top);
            if top < 1 then begin
               writeln(output,' "top" parameter must be > 1');
               halt
            end;
            if eof(sorthp)
            then minlenmaxene := 0.0
            else readln(sorthp, minlenmaxene)
         end;
         doinglength := (minlenmaxene > 0.0) or (not order.getenergy);

         (* display priorities to list file *)
         write(list,' ');
         for i := 1 to low do write(list,priority[i]);
         write(list,' priority order of sorting:');
         for i := 1 to low-1 do begin
            case priority[i] of
               'e': write(list,' energy');
               'l': write(list,' length')
            end;
            write(list,' then')
         end;
         writeln(list,' ambiguous.');
         write(list,' ',top:4,' helixes or fewer are written to shlist.  ');
         if top > maxhelix
         then writeln(list,'all helixes guaranteed to be written.')
         else writeln(list,'some helixes may be removed.');

         if doinglength
         then writeln(list,' ',trunc(minlenmaxene):1,
            ' is the minimum length helix recorded')
         else writeln(list,' ',minlenmaxene:7:2,
            ' is the maximum energy helix recorded');

         (* set the global "both" *)
         if low = 3
         then both := true
         else both := false
      end
end; (* getparams *)
(* end module sorth.getparams *)

(* begin module sorth.hrecord *)
procedure hrecord(var hlist: text; x,y,length: integer;
                  doenergy: boolean; energy: real);
(* record the helix information in hlist, just like hrecord in program
helix does.  do energy if doenergy is true. *)
begin (* hrecord *)
      write(hlist,'    x5',prime,': ',x:6,
                  '    y5',prime,': ',y:6,
                  '    length: ',length:3);
      if doenergy then write(hlist,'  ', energy:7:2, ' kcal');
      writeln(hlist)
end; (* hrecord *)
(* end module sorth.hrecord *)

(* begin module sorth.doset *)
procedure doset(var hlist,shlist: text;
                order: operations;
                var n, u, e, l, a: integer);
(* do a set of helixes from hlist (one set is the helixes from the
comparision of two pieces).  use the order given to select those
helixes to put into shlist.  record the results in the totals: n, u, e, l, a
(none, unique, energy, length or ambiguous cases).
note the references to globals:
      ordering, helixes, sorton and both
these are necessary because quicksort must use lessthan and swap
which in turn must do the dirty work.  see sorth.var. *)
var
      done: boolean; (* when there are no more helixes in a set *)
      h: integer; (* which helix is being read,
         then the number of helixes *)
      s: integer; (* output sorting index *)
begin (* doset *)
      h := 0;
      with order do repeat
         (* advance to next open position of helixes and ordering array *)
         h := h + 1;
         if h > maxhelix then begin
            writeln(output,' there are more helixes than can be handled.');
            writeln(output,' (there are at least ',maxhelix:1,' of them.)',
               '  increase constant maxhelix.');
            halt
         end;
         with helixes[h] do gethelix(hlist,x,y,length,getenergy,energy,done);

         if done
         then h := h - 1 (* the last h is not used *)
         else begin
            ordering[h] := h; (* initial ordering is 1 to 1 *)

            (* reject unwanted helixes by decreasing h *)
            case doinglength of (* using two cases is faster *)
               true:  if helixes[h].length < minlenmaxene then h := h - 1;
               false: if helixes[h].energy > minlenmaxene then h := h - 1
            end
         end
      until done;

      if h = 0
      then begin (* no helixes in this set *)
         sorton := 'n';
      end
      else if h = 1 then begin (* got a unique helix *)
         with helixes[ordering[h]],order
            do hrecord(shlist,x,y,length,getenergy,energy);
         sorton := 'u'
      end
      else begin (* more than one helix to do *)
         sorton := order.priority[1]; (* first priority sort follows *)
         quicksort(1,h);
         if (not lessthan(h-1,h)) and (order.top = 1)
         then sorton := 'a' (* it is ambiguous... *)
         else begin (* choose the top helix *)
            s := h;
            with order do while (s > 0) and (s > h - top)
            do begin
               with helixes[ordering[s]]
               do hrecord(shlist,x,y,length,getenergy,energy);
               s := s - 1
            end
         end
      end;

      (* display the findings to list.  make sure that no
      top is specified if none was found *)
      write(list,' helixes: ',h:3);
      if (sorton = 'a') or (sorton = 'n')
      then write(list,'      ',' ':3)
      else write(list,' top: ',ordering[h]:3);
      write(list,' by: ',sorton);

      case sorton of
         'n': n := n + 1; (* no helixes *)
         'u': u := u + 1; (* unique helix *)
         'e': e := e + 1; (* energy *)
         'l': l := l + 1; (* length *)
         'a': a := a + 1; (* ambigous *)
      end;

      (* copy the next id line if not end of file *)
      if not eof(hlist) then begin
         write(shlist,' '); (* replace space lost to gethelix procedure *)
         copyaline(hlist,shlist)
      end
end; (* doset *)
(* end module sorth.doset *)

(* begin module sorth.themain *)
procedure themain(var hlist, shlist, list, sorthp: text);
(* the main procedure of sorth.  only the strongest helixes of hlist are
copied to shlist, according to the parameter rules sorthp.  report in list. *)
var
      order: operations; (* the order of sorting operations *)
      (* totals for the various possible cases for each set:
      none, unique, energy, length and ambiguous. *)
      n, u, e, l, a: integer;
      sets: integer; (* set number *)
      endofline: boolean; (* triggers end of lines in the list file *)
begin (* themain *)
      writeln(output,' sorth ',version:4:2);

      hcopy(hlist,shlist,list);
      getparams(hlist, sorthp, order, list);

      (* pick up sets of helixes and chose one. *)
      sets := 0; n := 0; u := 0; e := 0; l := 0; a := 0;
      while not eof(hlist) do begin
         sets := sets + 1;
         endofline := ((sets mod 2) = 1);
         if endofline then writeln(list);
         write(list,' set: ',sets:3);

         doset(hlist,shlist,order,n,u,e,l,a);

         if endofline then write(list,'.');
      end;

      writeln(list);
      writeln(list);
      writeln(list,' classification of sets in hlist',
         ' by number of helixes:');
      writeln(list,' none:      ',n:4,' no helixes');
      writeln(list,' unique:    ',u:4,' one helix');
      writeln(list,' energy:    ',e:4,' sorted by energy first');
      writeln(list,' length:    ',l:4,' sorted by length first');
      writeln(list,' ambiguous: ',a:4,' sorted, with no single',
         ' unambiguous strongest helix');
      writeln(list,' total:     ',(n+u+e+l+a):4);

      writeln(output,' ',sets:1,' sets of helixes sorted');

      if (n+u+e+l+a) <> sets then begin
         writeln(output,'sorth.themain, program error');
         halt
      end
end; (* themain *)
(* end module sorth.themain *)

begin (* sorth *)
      themain(hlist, shlist, list, sorthp);
1: end. (* sorth *)