program denav(data, denavp, bfile, output);
(* denav: density average

  Dr. Thomas D. Schneider
  National Institutes of Health
  National Cancer Institute
  Center for Cancer Research Nanobiology Program
  Molecular Information Theory Group
  Frederick, Maryland  21702-1201
  toms@ncifcrf.gov
  permanent email: toms@alum.mit.edu (use only if first address fails)
  http://www.ccrnp.ncifcrf.gov/~toms/

 *)

label 1; (* end of program *)

const
(* begin module version *)
   version = 1.00; (* of denav.p 2009 Apr 14
2009 Apr 14, 1.00: origin
*)
updateversion = 1.00; (* defines lowest acceptable current parameter file *)
(* end module version *)

(* begin module describe.denav *)
(*
name
   denav: density average

synopsis
   denav(data: in, denavp: in, bfile: out, output: out)

files

   data:  Scan data file prepared for density plotting with denplo

   bfile:  (not used)

   denavp:  parameters to control the program.  The file must contain the
      following parameters, one per line:

      parameterversion: The version number of the program.  This allows the
         user to be warned if an old parameter file is used.

      piece coordinate column in data

      Ri evaluation (bits per site)

   output: messages to the user

description

   When an individual information scan has been set up with zero
   coordinates, it can be plotted using denplo.  What is the average
   information for the negative, zero and positive coordinates?

examples

1.00  version of denav that this parameter file is designed for.
4     piece coordinate column in data
6     Ri evaluation (bits per site)

documentation

see also

   {scan program to generate data:} scan.p
   {density plotting program:} denplo.p
   {example parameter file for current scan data file:} denavp

author

   Thomas Dana Schneider

bugs

technical notes

*)
(* end module describe.denav *)

var
   data, (* file used by this program *)
   denavp, (* file used by this program *)
   bfile: text; (* file used by this program *)

(* 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 = 5.40; (@ of prgmod.p 2008 Aug 06 *)

(* begin module skipblanks *)

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

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

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

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

procedure skipcolumn(var thefile: text);
(* skip over a data column *)
begin
   skipblanks(thefile); skipnonblanks(thefile)
end;
(* end module skipblanks version = 5.40; (@ of prgmod.p 2008 Aug 06 *)

(* begin module denav.themain *)
procedure themain(var data, denavp, bfile: text);
(* the main procedure of the program *)
const
   debugging = true; (* create output for debugging? *)
var
   parameterversion: real; (* parameter version number *)
   column: integer; (* current column being inspected in the data file *)
   coordinatecolumn: integer; (* piece coordinate column in data *)
   ricolumn: integer; (* Ri evaluation (bits per site) column *)
   coordinate: integer; (*piece coordinate column in data *)
   ri: real; (* Ri evaluation (bits per site) *)
   risumneg: real; (* Ri sum from negative coordinate values *)
   risumzer: real; (* Ri sum from zero     coordinate values *)
   risumpos: real; (* Ri sum from positive coordinate values *)
   negatives: integer; (* count of negative coordinate cases *)
   zeros:     integer; (* count of zero     coordinate cases *)
   positives: integer; (* count of positive coordinate cases *)
begin
   writeln(output,'denav ',version:4:2);
   reset(denavp);
   readln(denavp, parameterversion);
   if round(100*parameterversion) < round(100*updateversion) then begin
      writeln(output, 'You have an old parameter file!');
      halt
   end;

   readln(denavp, coordinatecolumn);
   readln(denavp, ricolumn);

   risumneg := 0.0;
   risumpos := 0.0;
   risumzer := 0.0;
   negatives := 0;
   zeros     := 0;
   positives := 0;

   reset(data);
   while not eof(data) do begin
      if (data^ = '*') or (eoln(data)) then begin
         readln(data); (* skip blank lines and comment lines *)
      end
      else begin
         column := 1;
         while column < coordinatecolumn do begin
            column := succ(column);
            skipcolumn(data);
         end;

         read(data, coordinate);
         column := succ(column);

         while column < ricolumn do begin
            column := succ(column);
            skipcolumn(data);
         end;
         read(data, ri);

         readln(data);

{
         if debugging then writeln(output, coordinate:1, ' ',ri:10:8);
}

         (* at this point we have the coordinate and the ri value *)

         if coordinate < 0
         then begin
            risumneg := ri + risumneg;
            negatives := succ(negatives);
{
if debugging then writeln(output, '= ', negatives:1, ' ', risumneg:1);
}
         end
         else if coordinate > 0
         then begin
            risumpos := ri + risumpos;
            positives := succ(positives);
{
if debugging then writeln(output, '= ', positives:1, ' ', risumpos:1);
}
         end
         else begin (* coordinate = 0 *)
            risumzer := ri + risumzer;
            zeros := succ(zeros);
         end;

      end;

   end;

   (* now compute the averages! *)
   if negatives > 0 then begin
{
      writeln(output, 'negatives: ', negatives:1, ' ', risumneg:10);
}
      writeln(output,
            'Average of ',negatives:7,' negative coordinate Ri values: ',
            (risumneg / negatives): 10:6, ' bits per site');
   end
   else begin
      writeln(output, '                   No negative coordinates');
   end;

   if zeros > 0 then begin
{
      writeln(output, 'zeros: ', zeros:1, ' ', risumzer:10);
}
      writeln(output,
            'Average of ',zeros:7,    ' zero     coordinate Ri values: ',
            (risumzer / zeros): 10:6, ' bits per site');
   end
   else begin
      writeln(output, '                   No zero coordinates');
   end;

   if positives > 0 then begin
{
      writeln(output, 'positives: ', positives:1, ' ', risumpos:10);
}
      writeln(output,
            'Average of ',positives:7,' positive coordinate Ri values: ',
            (risumpos / positives): 10:6, ' bits per site');
   end
   else begin
      writeln(output, '                   No positive coordinates');
   end;

   if negatives+zeros+positives > 0 then begin
      writeln(output,
            'Average of ',positives+zeros+negatives:7,
                       ' total    coordinate Ri values: ',
            (risumneg + risumzer + risumpos)
             /
            (negatives+zeros+positives): 10:6, ' bits per site');
   end
   else begin
      writeln(output, '                   No coordinates');
   end;

end;
(* end module denav.themain *)

begin
   themain(data, denavp, bfile);
1: end.