BSMA3 control program
BSMA3A.PAS
version 20061115
Hannes.Tammet@ut.ee

Selected bookmarks:
Parameters of simulator
Internal calibration constants
External calibration constants
Small global procedures
Function Mobility
Function ADCoriginal
Setcontrols and Startposition
Read BSMA3A.INI
Read BSMA3.DLF
TEST
   local procedures
   main routine
MEASUREMENT
   local routines: message, timereserve, midnight, showkeys
   init flowrate, meteo and balance
   estimate inlet penetration
   wait for starttime
   window settings
   write diagnostics
   calibration
   balance bridge
   mark time, correct clock, and measure flow zero
   scan
   calculate mobility distribution
   calculate size distribution
   show data
   store data
   save diagram tables
   save data
   save details
   condition data
   run cycle
   failure pause
   MAIN ROUTINE
Extratest
Welcome
MAIN PROGRAM



Program BSMA3A;
{BSMA3 control program: test BSMA3, record mobility and size distributions.
 Requires MS DOS and should be compilated by Turbo Pascal 7}
{$M 65000, 0, 60000}
{$R+}
Uses DOS, CRT;

CONST version = 'BSMA3A version HT20061115';
    developer = false; {true for code developer and false for user}
    simulator = false; {true when testing the program without BSMA,
                         false when compiling the program for measurement}

{Output is saved in tab-separated monthly files 3Ayymm00.* and may be saved
 in diurnal files 3Ayymmdd.* and diagram table files. * marks the extension
 (usually xl), which can be selected in BSMA3A.INI. Monthly files are saved
 onto internal disk, the destination of diurnal files and diagram tables is
 to be selected by the user. See the manual for the detailed explanations.

 A line of a main output file begins with a 7-number heading:

    yymmdd hhmm day_of_year temperature humidity pressure noise

 Time hhmm marks the center of the period of time averaging.

 The heading is followed by 10+10+16+16 values of fraction concentrations (cm-3)
 of + and - size distribution and + and - mobility distribution.

 At the end of the file is an diagnostic appendix consisting of 5 values:

    powervoltage filtervoltage flowrate number_of_scans balanceindex

 There are 64 values in a line.

 Limits of distribution fractions follow the "decade to eight" scheme

 11 limits of 10 size distribution fractions are:
     0.422, 0.562, 0.75, 1.00, 1.33, 1.78, 2.37, 3.16, 4.22, 5.62, 7.50 nm

 17 limits of 16 mobility distribution fractions are:
     0.0316, 0.0422, 0.0562, 0.0750, 0.100, 0.133. 0.178, 0.237,
     0.316, 0.422, 0.562, 0.75, 1.00, 1.33, 1.78, 2.37, 3.16 cm2V-1s-1.

 Concentration of particles less than 0.42 nm is negligible, thus the first
 size fraction 0.422-0.562 nm can be considered as a fraction of 0-0.562 nm.

 The last mobility fraction is always empty in the nature. It controls the
 broadening of the transfer function.

 Times of measuring operations and delays are expressed in millisconds
 in the program.

 NB! The active folder must consist of a calibration file BSMA3A.INI
     and 4 subfolders MONTHS, DAYS, TABLES and DETAILS.}

{PARAMETERS OF SIMULATOR}
CONST
  simulatordelay = 40; {ms per measurement, natural ADC delay is 100 ms}
  acceleration = 5; {compared with a natural measurement, up to 5}
  simupos : array [1..16] of integer =
      (40,35,30,25,20,15,10,05,00,00,70,120,200,120,30,00);
  simuneg : array [1..16] of integer =
      (70,60,50,45,40,30,20,10,05,05,10,20,100,200,150,50);
  simunoisefactor = 1; {allows to increase the noise, normal value = 1}
VAR simulatormark : integer; {time mark for simulator}

{INTERNAL CALIBRATION CONSTANTS SPECIFIC FOR BSMA3}
CONST
  asymmetry = 1; {divider for negative and multiplier for positive ion n}
  voltagefactor = 222; {at standardflowrate and balanceresistor = 0}
  concentrationfactor = 222; {at standardflowrate, electrometergain = 1,
                              and neglecting inlet losses}
  grid_adsorptioncoefficient      = 0.06;  {for Z = 1 cm2V-1cm-1, 0 C, 1013 mb}
  longplate_adsorptioncoefficient = 0.055; {for long filterplates only at
                                                Z = 1 cm2V-1cm-1, 0 C, 1013 mb}
  adsorptionratio = 0.867; {ratio of short and long filterplate adsorption}
  edgeloss = 0.043; {loss of ions on the filter due to the edge effect}
  noisefactor = 0.75; {experimental value}
  c_powervoltage = 0.005062;  {powervoltage / ADC_counts}
  c_filtervoltage = 0.1982; {filter voltage / ADC_counts}
  standardpowervoltage = 24; {V, corresponds to the standardflowrate}
  standardflowrate = 40; {dm3/s}
  electrometerdelay = 222; {ms, time constant of the input circuit}
  dlftolerance = 0.01; {permitted relative uncertainty of dlf verification}
  flowzerocycle = 5; {primary adjustment of flow sensor zero}
  flowzerohours : set of byte = [10, 22]; {hours of flow sensor zero adjustment}
 {Critical values}
  minpower = 21;
  maxpower = 26;
  minflow = 32;
  maxflow = 42;
  minfilter = 450;
  maxfilter = 520;
  minbias = -5;
  maxbias = 5;
  minbalance = 1;
  maxbalance = 14;
  minfractiontime = 6;
  maxoverload = 2; {allowed number of overload events during one scan}

VAR
{EXTERNAL CALIBRATION CONSTANTS, explanations in BSMA3A.INI}
  clockcorrection : integer;
  flowsensor : boolean;
  cycletime,
  c_voltage,
  c_concentration,
  cpressurea, cpressureb,
  ctemperaturea, ctemperatureb,
  chumiditya, chumidityb,
  cflowrate,
  balanceconstant,
  electrometergain : single;
  extrapath : string;     {for external or network disk}
  extension : string [3]; {for output file}
  initialcontrols : string [4]; {initial setting of output controls}

{DLF CONSTANTS}
  ADCballasttime_a, {100 - zero_voltage_measurement_time, ms}
  ADCballasttime_b, {extra measurement time per ADC count, ms}
  delayfactor : single;

{LPT CONTROLS, 0 = off or open, 1 = on or closed, balance = 0..15}
  fan_on,
  HV_on,
  plusiongate_open,
  minusiongate_open,
  balance : integer;

{ONE-DAY DATA STORAGE}
         storage : array [1..144, 1..64] of single;
  storagecounter,
      oldcounter : integer;

Procedure Clearkeyboard;
  Var c : char;
  Begin
  while keypressed do c := readkey;
  end;

Procedure Blankscreen;
  Begin
     textmode (3);
     textbackground (blue);
     textcolor (yellow);
     clrscr;
  End;

Function Sec100 : integer; {centisecond from sharp minute,
                            clock step is about 5 cs}
  Var h, m, s, t : word;
  Begin gettime (h, m, s, t); sec100 := 100 * s + t End {of sec100};

Procedure Pause (ms : single); {up to 60000 ms}
  Begin delay (round (delayfactor * ms)) End;

Procedure Checktime;
   {Check and adjust the computer clock}
   function bad : boolean;
   var c : char;
       s : string [3];
   begin
     repeat
        readln (s);
        if s = '' then c := 'Y' else c := upcase (s [1]);
     until c in ['Y', 'N'];
     bad := c = 'N';
     writeln;
   end;
Var yy, mm, dd, tt, hh, ss, s0 : word;
Begin
   clrscr;
   writeln;
   writeln ('A hint: when asked a Y/N question ',
            'then immediate ENTER means Yes.');
   writeln;
   {Setting of date and time}
   getdate (yy, mm, dd, tt);
   write ('Computer date is: year = ', yy,
                          ' month = ', mm,
                            ' day = ', dd);
   write ('. Is this correct (Y/N)? ');
   if bad then begin
      writeln ('Please tell the correct numbers:');
      write   ('    year = '); readln (yy);
      write   ('   month = '); readln (mm);
      write   ('     day = '); readln (dd);
      setdate (yy, mm, dd);
   end;
   s0 := 99;
   repeat
      gettime (hh, mm, ss, tt);
      if ss <> s0 then begin
         gotoxy (1, 6); clreol;
         write ('Computer time is: hour = ', hh,
                               ' minute = ', mm,
                               ' second = ', ss);
         write ('. Is this correct (Y/N)? ');
      end;
      s0 := ss;
   until keypressed;
   if bad then begin
      writeln ('Please tell the correct numbers (minute at sharp minute):');
      write   ('    hour = '); readln (hh);
      write   ('  minute = '); readln (mm);
      settime (hh, mm, 0, 0);
   end;
End {of checktime};

function Mobility_a {cm2 V-1 s-1}
        (pressure {mb},
         Celsius,
         massdiameter {nm} : double) : double;
  const {air parameters}
    GasMass = 28.96;
    Polarizability = 0.00171;
    VisCon1 = 0.3036;
    VisCon2 = 44;
    VisCon3 = 0.8;
    {nanometer particles expected:}
    particledensity = 2.08;
    ParticleCharge = 1;

  function Omega11 (x : double) : double;  {*(1,1)*(T*) for (*-4) potential}
    var p, q : double;                   {and elastic-specular collisions}
    begin
      if x > 1 then Omega11 := 1 + 0.106 / x + 0.263 / exp ((4/3) * ln (x))
      else begin p := sqrt (x); q := sqrt (p);
        Omega11 := 1.4691 / p - 0.341 / q + 0.181 * x * q + 0.059 end;
    end;

  const a = 1.2; b = 0.5; c = 1; {the slip factor coefficients}
        ExtraDistance = 0.115 {nm}; TransitionDiameter = 2.48 {nm};
  var   GasDiameter, MeanVelocity, Viscosity, FreePath, DipolEffect,
        DeltaTemperature, CheckMark, ParticleMass, CollisionDistance,
        Kn, Omega, s, x, y, Temperature {K} : double;

  begin
    Temperature := Celsius + 273.15;
    Viscosity {microPa s} := 0.02713 * sqrt (GasMass * Temperature) /
      sqr (VisCon1 * (1 + exp (VisCon3 * ln (VisCon2 / Temperature))));
    MeanVelocity {m/s} := 145.5 * sqrt (Temperature / GasMass);
    FreePath  {nm} := (166251 * Viscosity * Temperature) /
                      (GasMass * Pressure * MeanVelocity);
    ParticleMass {amu} := 315.3 * ParticleDensity *
                                  exp (3 * ln (MassDiameter));
    DeltaTemperature := Temperature;
    repeat
      CheckMark := DeltaTemperature;
      GasDiameter {nm} := VisCon1 *
        (1 + exp (VisCon3 * ln (VisCon2 / DeltaTemperature)));
      CollisionDistance {nm} := MassDiameter / 2 + ExtraDistance +
                                 GasDiameter / 2;
      DipolEffect := 8355 * ParticleCharge * Polarizability /
                            sqr (sqr (CollisionDistance));
      DeltaTemperature := Temperature + DipolEffect;
    until abs (CheckMark - DeltaTemperature) < 0.01;
    if ParticleCharge = 0 then Omega := 1
                          else Omega := Omega11 (Temperature / DipolEffect);
    Kn := FreePath / CollisionDistance;
    if Kn < 0.03 {underflow safe} then y := 0 else y := exp (- c / Kn);
                  {NB! erratum in JAS (y := 1) corrected!}
    x := (273.15 / DeltaTemperature) *
         exp (3 * ln (TransitionDiameter / MassDiameter));
    if x > 30 {overflow safe} then s := 1
      else if x > 0.001
      then s := 1 + exp (x) * sqr (x / (exp (x) - 1)) * (2.25 / (a + b) - 1)
      else {underflow safe}  s := 1 + (2.25 / (a + b) - 1);
    Mobility_a := 1.602 * ((2.25 / (a + b)) / (Omega + s - 1))  *
                sqrt (1 + GasMass / ParticleMass) *
                (1 + Kn * (a + b * y)) /
                (6 * PI * Viscosity * CollisionDistance);
  end;

Function ADCsimulator (i : integer) : integer;
Var h, m, s, f : word;
             k : integer;
          t, x : single;
Begin
   pause (round (simulatordelay / acceleration));
         {equivalent of ADC measurement time}
   gettime (h, m, s, f);
   k := 100 * s + f; k := k - simulatormark;
   if k < 0 then k := k + 6000;
   t := acceleration * k / 300; {unit = RC time constant in cs}
   case i of
        1 : x := 0; {zero check channel}
        2 : begin {simulated error}
            x := 6;
            for k := 1 to 12 do x := x - random;
            if random < 0.005 then x := x + 20 * random;
            x := simunoisefactor * x * electrometergain / 10;
            end;
        3 : x := 8100 * exp (-t);
        4 : x := 2810; {flowrate}
        5 : x := 6700; {pressure}
        6 : x := 4000; {humidity}
        7 : x := 700;  {temperature}
        8 : x := 4700; {fan voltage}
   else x := 0;
   end;
   if abs (x) > 16383 then x := 16383 * x / abs (x);
   adcsimulator := round (x);
End {of ADCsimulator};

Function ADCoriginal (channel : integer) : longint;
  {Specific code for control or measurement with PICO ADC16 in BSMA
   Channel -1 = close ADC16,
   Channel  0 = open ADC16 connected to COM1 port,
   Channel  # = measure channel # with resolution 13 bits (time 55-78 ms)
               (actual time is 53.5 + 0.003*|counts| ms
                with a variation less than 1 ms)
   If an overflow or an failure occurred then result will be 16383
   NB - the codes are specific for BSMA}
  CONST
  I8250_DLAB            = $80;
  I8250_DATA_READY      = $01;
  I8250_RTS             = $02;
  RS232_base  : word    = 0;
  port_addresses  : array [1..4] of word = ($3F8,$2F8,$3E8,$2E8);
  VAR
  resolution : integer;
  bios_timer : longint absolute$0040:$006C;
  controlstr : string;
  divisor    : longint;
  sign    : char;
  upper   : word;
  lower   : word;
  value   : longint;
  timeout : longint;
  control : byte;
  BEGIN
  if channel < 0 then begin channel := 0; resolution := 0 end
  else if channel = 0 then resolution := 1
  else resolution := 13;
  value := 0;
  if channel = 0 then begin
    if resolution = 0 then port [RS232_base+4] := $00
    else if resolution in [1..4] then begin
      RS232_base := port_addresses [resolution];
      divisor := 115200 div 9600;
      port [RS232_base+3] := I8250_DLAB;
      port [RS232_base+1] := (divisor shr 8) and $FF;
      port [RS232_base]  :=  divisor and $FF;
      port [RS232_base+3] := 3;
      port [RS232_base+4] := I8250_RTS;
    end
  end
  else if (channel in [1..8]) and (resolution in [8..16]) then begin
    while port [RS232_base+5] and I8250_DATA_READY <> 0 do
      lower := port [RS232_base];
      control := ((channel-1) shl 5) + ((resolution-1) shl 1) + 1;
    port [RS232_base] := control;
    timeout := bios_timer + 20;
    while (bios_timer < timeout)
      and (port [RS232_base+5] and I8250_DATA_READY = 0) do;
    sign := chr (port [RS232_base]);
    while (bios_timer < timeout)
      and (port [RS232_base+5] and I8250_DATA_READY = 0) do;
    upper := port [RS232_base];
    while (bios_timer < timeout)
      and (port [RS232_base+5] and I8250_DATA_READY = 0) do;
    lower := port [RS232_base];
    if bios_timer < timeout then begin
       value := (upper shl 8) + lower;
       if sign = '-' then value := -value;
    end else value := 16383;
  end;
  adcoriginal := value;
  End {of ADCoriginal};

Function ADC (channel : integer) : integer;
   {if not simulator then 100 ms independent of value}
   Var n : integer;
       ballasttime : single;
   Begin
   if simulator then n := adcsimulator (channel)
   else begin
      n := adcoriginal (channel);
      ballasttime := ADCballasttime_a - ADCballasttime_b * abs (n);
      if ballasttime > 1 then delay (round (delayfactor * ballasttime));
   end;
   adc := n;
   End {of ADC};

Procedure Setcontrols (t {milliseconds after settings} : integer);
  Var x : byte;
  Begin
  x := fan_on + 2 * HV_on + 4 * (1 - plusiongate_open)
     + 8 * (1 - minusiongate_open) + 16 * (balance mod 2)
     + 32 * ((balance mod 4) div 2) + 64 * ((balance mod 8) div 4)
     + 128 * (balance div 8);
  if simulator then simulatormark := sec100 else port [888] {LPT1} := x;
  pause (t);
  End {of setcontrols};

Procedure Startposition; begin
  fan_on := 0;
 {switch off voltages}
  HV_on := 0;
  plusiongate_open := 1;
  minusiongate_open := 1;
  setcontrols (0);
end {of startposition};

Procedure Read_ini;
  {Read BSMA3A.INI and evaluate calibration variables}
  Var       p : integer;
           OK : boolean;
    rg, re, x : single;
            f : text;
      s, a, b : string;
  procedure failure (prompt : string); {in read_ini}
    begin
      writeln;
      writeln ('FAILURE: ', prompt, ',');
      writeln ('press ENTER to escape!');
      readln; halt;
    end; {of failure}

  Begin
  {control values}
  x := -999; p := 999;
  extension := '***';
  extrapath := '***';
  initialcontrols := '***';
  clockcorrection := p;
  flowsensor := false;
  cycletime := x;
  rg := x;
  re := x;
  cpressurea := x;
  cpressureb := x;
  ctemperaturea := x;
  ctemperatureb := x;
  chumiditya := x;
  chumidityb := x;
  cflowrate := x;
  balanceconstant := x;
  write (' Reading BSMA3A.INI ... ');
  assign (f, 'BSMA3A.INI');
  {$I-} reset (f); {$I+}
  if IOresult <> 0 then
    failure ('cannot find BSMA3A.INI in the active folder');
  while not eof (f) do begin
    readln (f, s);
    if s = '' then s := ' '; {to allow empty lines}
    if s [1] <> ' ' then begin
       {a line that begins with space is considered as a comment}
       p := pos ('=', s);
       if p = 0 then
          failure ('an assignment line does not consist of = in BSMA3A.INI');
       a := copy (s, 1, p-1);
       b := copy (s, p+1, length (s) - p);
       p := pos (' ', b);
       if p > 0 then b := copy (b, 1, p - 1);
           {allow space-delimited end of comments}
       if a = 'extension' then extension := b else
       if a = 'extrapath' then extrapath := b else
       if a = 'controls'  then initialcontrols := b else begin
          val (b, x, p);
          if p <> 0 then failure ('wrong number presentation in BSMA3A.INI');
          if a = 'clockcorrection' then clockcorrection := round (x) else
          if a = 'flowsensor' then flowsensor := (x = 1) else
          if a = 'cycletime' then cycletime := round (x) else
          if a = 'gainresistor' then rg := x else
          if a = 'balanceresistor' then re := x else
          if a = 'cpressurea' then cpressurea := x else
          if a = 'cpressureb' then cpressureb := x else
          if a = 'ctemperaturea' then ctemperaturea := x else
          if a = 'ctemperatureb' then ctemperatureb := x else
          if a = 'chumiditya' then chumiditya := x else
          if a = 'chumidityb' then chumidityb := x else
          if a = 'cflowrate' then cflowrate := x else
          if a = 'balanceconstant' then balanceconstant := x else
          failure ('wrong keyword in BSMA3A.INI');
       end {of a = 'controls'};
    end {of s [1] <> ' '};
  end {of while not eof (f)};
  close (f);
  {check control values}
  if abs (clockcorrection) > 25 then
     failure ('clockcorrection exceeds +-25');
  if (cycletime <> 10) and (cycletime <> 15) then
     failure ('cycletime illegal or not initiated');
  x := -999;
  if rg = x then failure ('gainresistor not initiated');
  if re = x then failure ('balanceresistor not initiated');
  if cpressurea = x then failure ('cpressurea not initiated');
  if cpressureb = x then failure ('cpressureb not initiated');
  if ctemperaturea = x then failure ('ctemperaturea not initiated');
  if ctemperatureb = x then failure ('ctemperatureb not initiated');
  if chumiditya = x then failure ('chumiditya not initiated');
  if chumidityb = x then failure ('chumidityb not initiated');
  if cflowrate = x then failure ('cflowrate not initiated');
  if balanceconstant = x then failure ('balanceconstant not initiated');
  if extrapath = '***' then failure ('extrapath not initiated');
  if (length (extrapath) < 2) and (extrapath <> '-') then
     failure ('illegal extrapath in ini-file');
  if extension = '***' then failure ('extension not initiated');
  if length (initialcontrols) <> 4 then
     failure ('controls not initiated or too short');
  OK := true;
  for p := 1 to 4 do OK := OK and (initialcontrols [p] in ['0', '1']);
  if not OK then failure ('wrong symbol in controls');
  if rg = 0 then electrometergain := 1
            else electrometergain := 1 + 50000 / rg;
  c_voltage := voltagefactor * 25 / (25 + re / (2000 + re));
  c_concentration := concentrationfactor / electrometergain;
  p := adc (0); {open ADC}
  writeln ('done.');
  End {of initialization};

Procedure Read_dlf;
  {Read and verify or create timing calibration file BSMA3.DLF}
  procedure verifydlf (n : integer; var c1, c2, c3 : single);
                      {cycles, d_factor, ballast_a, ballast_b}
    var     i, j, k : integer;
     u1 ,v1, u2, v2,
            x, y, z : single;
    begin
    writeln (' Initial value of the delayfactor = ', delayfactor:5:3);
    writeln (' Measured delayfactor:');
    y := 0;
    for i := 1 to n do begin
       k := sec100;
       for j := 1 to 100 do delay (100);
       k := sec100 - k;
       if k < 0 then k := k + 6000;
       x := 1000 / k;
       y := y + x / n;
       write (x:8:3);
    end;
    c1 := y;
    if n > 5 then delayfactor := y;
    if n mod 10 <> 0 then writeln;
    writeln (' Average measured delayfactor = ', y:5:3);
    writeln (' ADC time (should be 100 ms) at ', adc (3), ' counts:');
    v1 := 0; u1 := 0;
    for i := 1 to n do begin {ADC3 about 0 counts}
       z := 0;
       k := sec100;
       for j := 1 to 100 do z := z + abs (adc (3));
       k := sec100 - k;
       if k < 0 then k := k + 6000;
       x := k / 10; {as converted from cs to ms}
       z := z / 100;
       write (x:8:1);
       v1 := v1 + z / n; {counts}
       u1 := u1 + x / n; {time : ms}
    end;
    if n mod 10 <> 0 then writeln;
    writeln (' Average ADC counts ', round (v1),
             ', average time ', u1:4:2, ' ms');
    writeln (' ADC time (should be 100 ms) at ', adc (5), ' counts:');
    v2 := 0; u2 := 0;
    for i := 1 to n do begin {ADC5 about 6700 counts}
       z := 0;
       k := sec100;
       for j := 1 to 100 do z := z + abs (adc (5));
       k := sec100 - k;
       if k < 0 then k := k + 6000;
       x := k / 10; {as converted from cs to ms}
       z := z / 100;
       write (x:8:1);
       v2 := v2 + z / n; {counts}
       u2 := u2 + x / n; {time : ms}
    end;
    if n mod 10 <> 0 then writeln;
    writeln (' Average ADC counts ', round (v2),
             ', average time ', u2:4:2, ' ms');
    if v1 = v2 then y := 0
               else y := (u2 - u1) / (v2 - v1); {extra time per ADC count}
    x := u1 - y * v1; {time at 0 ADC count (should be 100 ms)}
    c2 := ADCballasttime_a - (x - 100); {if 0-time too long then decrease a}
    c3 := ADCballasttime_b + y; {if time increases with counts, increase b}
    end; {of verify in Read_dlf}

  Var dlf0, dlf1, dlf2, dlf3 : single;
              yy, mm, dd, ff : word;
                        k, t : longint;
                  OK, bypass : boolean;
                           f : text;
  Begin {Read_dlf}
  {set default values}
  delayfactor := 1.5;
  ADCballasttime_a := 46.5;
  ADCballasttime_b := 0.003;
  if simulator then exit;
  write (' Reading BSMA3.DLF ... ');
  assign (f, 'BSMA3.DLF');
  {$I-} reset (f); {$I+}
  OK := (IOresult = 0);
  if not OK then writeln (' File BSMA3.DLF not available.');
  if OK then begin
     readln (f, delayfactor, ADCballasttime_a, ADCballasttime_b);
     close (f);
     writeln ('Press any key if you wish to bypass the DLF-calibration!');
     clearkeyboard;
     k := sec100;
     repeat
        t := sec100 - k; if t < 0 then t := t + 6000;
     until (t > 1000) or keypressed;
     if keypressed then begin clearkeyboard; exit end;
     dlf0 := delayfactor;
     writeln; writeln (' Verifying DLF constants (about 2 minutes) ...');
     verifydlf (4, dlf1, dlf2, dlf3);
     OK := (abs (dlf1 - dlf0) < (dlftolerance * dlf1)) and
           (abs (dlf2 - ADCballasttime_a) < (dlftolerance * dlf2)) and
           (abs (dlf3 - ADCballasttime_b) < (4 * dlftolerance * dlf3));
     if not OK then writeln ('DLF constants not exact enough');
  end;
  if not OK then begin
     dlf0 := delayfactor;
     dlf2 := ADCballasttime_a;
     dlf3 := ADCballasttime_b;
     writeln (' Old values of DLF constants:');
     writeln ('         Delayfactor =', dlf0:6:3);
     writeln ('    ADCballasttime_a =', dlf2:6:2);
     writeln ('    ADCballasttime_b =', dlf3:8:5);
     writeln (' Measuring of justified values (about 10 minutes) ...');
     verifydlf (20, delayfactor, ADCballasttime_a, ADCballasttime_b);
     rewrite (f);
     writeln (f, delayfactor:6:4, #09,
                 ADCballasttime_a:5:3, #09,
                 ADCballasttime_b:8:6);
     close (f);
     writeln (' New values of DLF constants:');
     writeln ('         Delayfactor =', delayfactor:6:3);
     writeln ('    ADCballasttime_a =', ADCballasttime_a:6:2);
     writeln ('    ADCballasttime_b =', ADCballasttime_b:8:5);
     getdate (yy, mm, dd, ff); k := yy;
     k := 10000 * k + 100 * mm + dd; {yyyymmdd}
     assign (f, 'TIMING.LOG');
     {$I-} append (f); {$I+}
     if IOresult <> 0 then rewrite (f);
     writeln (f, 'Timing calibration changed ', k);
     writeln (f, '   Old constants:', dlf0:8:4,
                                      dlf2:8:3,
                                      dlf3:10:6);
     writeln (f, '   New constants:', delayfactor:8:4,
                                      ADCballasttime_a:8:3,
                                      ADCballasttime_b:10:6);
     close (f);
  end;
  end; {of Read_dlf}

Procedure Test;
{Displays immediately ADC counts}
Const             nn = 30;
            testtime = 1000;
Var      c, previous : char;
                   x : array [0..9] of byte;
          z1, z2, z3 : array [1..nn] of single;
    tc, rh, mb, v, w : single;
    a, i, j, n, zero : integer;
                   p : longint;
                   y : array [0..20, 1..9] of integer;

procedure instruction; {in test}
   begin
   writeln;
   writeln (' BSMA3 test.');
   writeln (' If triple-connector cable is disconnected from the computer',
            '  then voltages on');
   writeln (' the LPT connector can be checked using a voltmeter connected',
            ' between a signal');
   writeln (' pin (one of # 2..9) and a ground pin (any of # 18..25). ',
            'Pins are numbered:');
   writeln;
   writeln ('         +----------------------------------------------------+');
   writeln ('         | 13  12  11  10  09  08  07  06  05  04  03  02  01 |');
   writeln ('         \                                                    / ');
   writeln ('          \  25  24  23  22  21  20  19  18  17  16  15  14  /');
   writeln ('           \________________________________________________/');
   writeln;
   writeln (' Pin voltages: low about 0 and high about 5 V. ',
            'Initially all pins are set low.');
   writeln (' If triple-connector cable is disconnected from the ADC ',
            'then ADC inputs can be');
   writeln (' checked using a controlled source ',
            'of DC voltage -2.5 V ... +2.5 V.  8195 ADC ');
   writeln (' counts should correspond to the voltage of 2.5 V.');
   writeln;
   writeln (' Remember test commands:');
   writeln ('    2...9  : set one LPT pin 2..9 high and display LPT settings');
   writeln ('    Q...I  : set LPT pin (see # above the letter)',
            ' low and display LPT settings');
   writeln ('  space bar: display centiseconds, readings of ADC,',
            ' and LPT settings');
   writeln ('      R    : repeat the space bar action 25 times');
   writeln ('      S    : show temperature, pressure and RH');
   writeln ('      K    : test scan with saving of the record');
   writeln ('      L    : test scan relative to saved record');
   writeln ('      A    : check ADC measurement time');
   writeln ('      B    : bridge balance test ',
            '(not stable during the first half hour)');
   writeln ('      N    : bridge noise test');
   writeln ('      P    : display explanation of LPT pins and ADC channels');
   writeln ('      X    : exit the test');
   writeln ('   another : display again the instruction above');
   end;

procedure pins_and_channels; {in test}
   begin
   writeln ('                      0 / 1');
   writeln (' Pin2 = fan          off/on        ADC1 = control zero');
   writeln (' Pin3 = HV           off/on        ADC2 = electrometer output');
   writeln (' Pin4 = +ion gate   open/closed    ADC3 = analyzer voltage');
   writeln (' Pin5 = -ion gate   open/closed    ADC4 = flow rate sensor');
   writeln (' Pin6 = balance1     off/on        ADC5 = barometric sensor');
   writeln (' Pin7 = balance2     off/on        ADC6 = humidity sensor');
   writeln (' Pin8 = balance4     off/on        ADC7 = temperature sensor');
   writeln (' Pin9 = balance8     off/on        ADC8 = power&filter voltage');
   end;

procedure setport; {in test}
   begin
   if not simulator then port [888] := x [1] + 2 * x [2] + 4 * x [3] +
      8 * x [4] + 16 * x [5] + 32 * x [6] + 64 * x [7] + 128 * x [8];
   balance := x [5] + 2 * x [6] + 4 * x [7] + 8 * x [8];
   pause (5);
   end;

procedure showport; {in test}
   var i : integer;
   begin
   write ('   LPT = ');
   for i := 1 to 8 do write (x [i]);
   writeln;
   end;

procedure electrometertransition; {can be used in the hidden Z-operation}
   var i : integer;
   begin {switch voltage on (key 3) before the measurement}
   writeln (' Zero measurement, wait for 10 seconds');
   if not simulator then pause (10000);
   zero := round ((adc (2) + adc (2) + adc (2) + adc (2) + adc (2)) / 5);
   writeln (' Electrometer transition process');
   for i:= 1 to 10 do write (adc (2) - zero:8); writeln;
   x [6] := 1; setport;
   for i:= 1 to 30 do write (adc (2) - zero:8); writeln;
   x [6] := 0; setport;
   for i:= 1 to 30 do write (adc (2) - zero:8); writeln;
   end;

procedure flowratetransition; {can be used in the hidden Z-operation}
   var i, a : integer;
   begin
   x [1] := 0; setport; pause (5000);
   for a := 1 to 9 do begin
      x [1] := 1; setport; pause (5000);
      for i := 1 to 10 do begin
         pause (900); y [i, a] := adc (4) - 3170;
      end;
      x [1] := 0; setport;
      for i := 11 to 20 do begin
         pause (900); y [i, a] := adc (4) - 3170;
      end;
   end;
   for i := 1 to 10 do begin
      write (i + 5:2, ')');
      v := 0;
      for a := 1 to 9 do v := v + y [i, a] / 9;
      write (v:6:0);
      write (i:8, ')');
      v := 0;
      for a := 1 to 9 do v := v + y [i + 10, a] / 9;
      writeln (v:6:0);
   end;
   end;

procedure sensorcalibration; {used in the hidden Z-operation}
   var i, n : integer;
          s : array [4..8] of single;
   begin
   clearkeyboard;
   writeln (' Press any key to stop the measurement and write the results!');
   n := 0;
   for i := 4 to 8 do s [i] := 0;
   repeat
      n := n + 1;
      for i := 4 to 8 do s [i] := s [i] + adc (i);
   until keypressed;
   clearkeyboard;
   writeln (' N = ', n, '   Q = ', round (s [4] / n),
                        '   P = ', round (s [5] / n),
                        '   H = ', round (s [6] / n),
                        '   T = ', round (s [7] / n),
                        '   V = ', round (s [8] / n));
   end;

Begin {test}
   for j := 0 to 20 do for i := 1 to 9 do y [j, i] := 0;
   for i := 1 to 8 do x [i] := 0;
   setport;
   previous := 'Z';
   textmode (259);
   textbackground (blue); textcolor (yellow);
   clrscr;
   instruction;
   repeat
      if previous <> ' ' then write (' Please enter command: ');
      c := upcase (readkey);
      if previous <> ' ' then writeln (c);
      case c of
       '2'..'9' : begin x [ord (c) - 49] := 1; setport; showport end;
            ' ' : begin
                  if previous <> ' ' then
                  writeln (' Time:cs   ADC1   ADC2   ADC3   ADC4   ',
                           'ADC5   ADC6   ADC7   ADC8   LPT');
                  write (sec100:8);
                  for i := 1 to 8 do write (ADC (i):7);
                  write ('   ');
                  for i := 1 to 8 do write (x [i]);
                  writeln;
                  end;
            'R' : begin
                  writeln (' Time:cs   ADC1   ADC2   ADC3   ADC4   ',
                           'ADC5   ADC6   ADC7   ADC8   LPT');
                  for j := 1 to 25 do begin
                     write (sec100:8);
                     for i := 1 to 8 do write (ADC (i):7);
                     write ('   ');
                     if (j = 1) or (j = 25)
                     then for i := 1 to 8 do write (x [i]);
                  writeln;
                  end;
                  end;
            'Q' : begin x [1] := 0; setport; showport end;
            'W' : begin x [2] := 0; setport; showport end;
            'E' : begin x [3] := 0; setport; showport end;
            'R' : begin x [4] := 0; setport; showport end;
            'T' : begin x [5] := 0; setport; showport end;
            'Y' : begin x [6] := 0; setport; showport end;
            'U' : begin x [7] := 0; setport; showport end;
            'I' : begin x [8] := 0; setport; showport end;
            'S' : begin
                  tc := ctemperaturea * adc (7) + ctemperatureb;
                  rh := (chumiditya * adc (6) + chumidityb) /
                        (1.0546 - 0.00216 * tc);
                  mb := cpressurea * adc (5) + cpressureb;
                  writeln ('          Temperature = ', tc:3:1, ' C');
                  writeln ('    Relative humidity = ', rh:3:1, ' %');
                  writeln (' Atmospheric pressure = ', mb:3:1, ' mb');
                  end;
        'K','L' : begin {test scan}
                  writeln (' Balance ', balance,
                           ' ... Charging the capacitor ...');
                  x [2] := 1; setport; pause (1000);
                  writeln ('Sec  s100  Volt  (Electrometer output)...');
                  x [2] := 0; setport;
                  simulatormark := sec100;
                  for j := 0 to 20 do begin
                     write (j:3, sec100:6, adc (3):6);
                     for i := 1 to 9 do begin
                        a := adc (2);
                        if c = 'K' then y [j, i] := a;
                        if c = 'L' then a := a - y [j, i];
                        write (a:6);
                     end;
                     writeln;
                  end;
                  end;
            'A' : if not simulator then begin {ADC time test}
                  writeln (' The voltage of the ADC channel 8 ',
                           'will be measured 100 times.');
                  writeln (' The voltage can be manipulated ',
                           'by LPT pins 4 and 5.');
                  writeln (' The average time must be 100.0+-0.5 ms ',
                           'independent of the voltage.');
                  write (' measuring...');
                  p := 0;
                  a := sec100;
                  for i := 1 to 100 do p := p + adc (8);
                  a := sec100 - a;
                  if a < 0 then a := a + 6000;
                  writeln ('   average ADC counts ', round (p / 100),
                           '   time of one measurement ', a / 10:4:1, ' ms');
                  end;
            'B' : begin {balance test}
                  x [2] := 0; setport;
                  writeln (' Zero measurement, wait for 20 seconds');
                  if not simulator then pause (20000);
                  zero := round ((adc (2) + adc (2) + adc (2)) / 3);
                  writeln (' Write the value of the balance 0..15',
                           ' and get 10 first readings ');
                  writeln (' of electrometer in a test scan.',
                           ' A value above 15 will interrupt the test');
                  repeat write (' balance = '); readln (a);
                     if a in [0..15] then begin
                        balance := a;
                        x [5] := a mod 2;
                        x [6] := (a mod 4) div 2;
                        x [7] := (a mod 8) div 4;
                        x [8] := a div 8;
                        setport; pause (100);
                        write (' charging');
                        x [2] := 1; setport; pause (1000);
                        x [2] := 0; setport;
                        simulatormark := sec100;
                        for i := 1 to 9 do write (adc (2):6);
                        writeln (adc (2):6);
                     end;
                  until not (a in [0..15]);
                  end;
            'N' : begin
                  writeln (' Bridge noise test, ', nn,
                           ' measurements, please be patient!');
                  x [3] := 1; x [4] := 1; setport;
                  for i := 0 to nn do begin
                     x [2] := 1; setport; pause (testtime);
                     x [2] := 0; setport; pause (1000);
                     if i > 0 then begin
                        z1 [i] := adc (2);
                        z2 [i] := adc (2);
                        z3 [i] := adc (2);
                        write (z2 [i]:8:0);
                     end;
                     pause (1000);
                  end;
                  v := 0;
                  for i := 1 to nn do v := v
                     + sqr (z2 [i] - (z1 [i] + z3 [i]) / 2);
                  v := sqrt (v / nn);
                  for i := 1 to nn do z2 [i] := (z1 [i] + z2 [i] + z3 [i]) / 3;
                  w := 0;
                  for i := 2 to nn - 1 do w := w
                     + sqr (z2 [i] - (z2 [i - 1] + z2 [i + 1]) / 2);
                  w := sqrt (w / (nn - 2));
                  writeln;
                  writeln (' Bridge external noise ', w:3:1,
                                  ' internal noise ', v:3:1);
                  end;
            'P' : pins_and_channels;
            'X' : begin
                  if not simulator then port [888] := 0;
                  exit;
                  end;
            'Z' : if developer then sensorcalibration; {can be modified!}
      else instruction;
      end; {of case c}
      previous := c;
   until false;
End {of test};

Procedure Measurement (standardregime : boolean);
{cycle pattern 0-0+0-0+0-0+0-0....}
Const nfraction = 16; {mobility fractions from 0.032 to 3.16}
      nscanmax  = 60; {NB: 30 scans would take 10 minutes without any reserve}
      nvoltmax = 220; {number of measurements in the calibration procedure}
      nflowmeasurements = 4; {charging time in a scan will be
                              (2 * nflowmeasurements + 3) * 100 ms}
Var datetext,
    failuremessage : string [32];
    electrometeroverload,
    nscan, {number of the scan in a cycle}
    nzero, {number of zero scans in a cycle}
    cycle, {if simulator then 1 else cycletime in minutes}
    electrometerzero {ADC counts} : integer;
    electrometerbias, {mV}
    powervoltage, {V}
    plusionfiltervoltage, {V}
    minusionfiltervoltage, {V}
    yymmdd, {date}
    hhmm, {time}
    dayofyear, {date and time}
    temperature, {C}
    humidity, {%}
    pressure, {mb}
    flowrate {l/s} : single;
    flowpressure_zero, {adc (4) at zero flowrate}
    cyclenumber : integer;
    run,
    datarecording,
    diagramtables,
    scandetails,
    extrastorage : boolean;
    oldvoltage : array [1..nvoltmax] of integer;
    sensorsum : array [4..7] of single;
    swingsum : single; {parameter of bridge balance}
    fractiontime : array [0..nfraction] of integer;
                  {how many 100 ms measurements in a fraction}
    fractionhighmobility : array [0..nfraction] of single;
             {lower mobility limit of i-fraction = fractionhighmobility [i-1]}
    inlet_standardpenetration : array [1..nfraction] of single;
    fractiontable : array [1..nscanmax, 1..nfraction] of single;
                   {to collect raw measurements during a cycle}
    displaytable : array [1..5, 1..65] of single;
                   {to collect results for display and save}
    column: integer; {display index of table 1..5}
   {Structure of table according to the second index:
    1) yymmdd, 2) hhmm, 3) dayofyear,
    4) temperature, 5) humidity, 6) pressure, 7) noise sigma,
    8-17) + size fractions, 18-27) - size fractions,
    28-43) + mobility fractions, 44-59) - mobility fractions,
    60) N+, 61) N-, 62) n+, 63) n-, 64) Z+, 65) Z-
   (Values 1..59 will be saved following the order above,
    last six values are for display only)}
    decimals : array [1..64] of integer; {in output tables}

  procedure message (s : string);
   var k : integer;
   begin
      gotoxy (31, 2); clreol;
      k := length (s);
      gotoxy (80 - k, 2); write (s);
   end;

  function timereserve : integer; {in measurement}
   var h, m, s, t : word;
   begin
      gettime (h, m, s, t);
      t := m mod cycle;
      timereserve := 60 * (cycle - t) - s - 1;
   end;

  function midnight : boolean;
   var h, m, s, t : word;
   begin
      gettime (h, m, s, t);
      t := 60 * h + m;
      midnight := (t < 5) or (t > 1435);
   end;

   procedure showkeys;
    begin
       if datarecording then begin
          gotoxy (60, 43); write ('Data storage ');
          textcolor (green); write ('ON '); textcolor (black);
          gotoxy (62, 44); write ('(Alt+D turns off)')
       end
       else begin
          gotoxy (60, 43); write ('Data storage ');
          textcolor (red); write ('0FF'); textcolor (black);
          gotoxy (62, 44); write ('(Ctrl+D turns on)')
       end;
       if diagramtables then begin
          gotoxy (60, 45); write ('Diagram tables ');
          textcolor (green); write ('ON '); textcolor (black);
          gotoxy (62, 46); write ('(Alt+T turns off)')
       end
       else begin
          gotoxy (60, 45); write ('Diagram tables ');
          textcolor (red); write ('OFF'); textcolor (black);
          gotoxy (62, 46); write ('(Ctrl+T turns on)')
       end;
       if scandetails then begin
          gotoxy (60, 47); write ('Scan details ');
          textcolor (green); write ('ON '); textcolor (black);
          gotoxy (62, 48); write ('(Alt+S turns off)')
       end
       else begin
          gotoxy (60, 47); write ('Scan details ');
          textcolor (red); write ('OFF'); textcolor (black);
          gotoxy (62, 48); write ('(Ctrl+S turns on)')
       end;
       if extrapath <> '-' then begin
          if extrastorage then begin
             gotoxy (60, 49); write ('External path ');
             textcolor (green); write ('ON '); textcolor (black);
             gotoxy (62, 50); write ('(Alt+E turns off)')
          end
          else begin
             gotoxy (60, 49); write ('External path ');
             textcolor (red); write ('OFF'); textcolor (black);
             gotoxy (62, 50); write ('(Ctrl+E turns on)')
          end;
       end
       else extrastorage := false;
    end;

   procedure checkkey;
    var a, x : integer;
    begin if keypressed then begin
       textcolor (black);
       a := 1; x := 1;
       while keypressed do begin
          a := x;
          x := (ord (upcase (readkey)));
       end;
       if a = 0 then x := -x;
       case x of
          24 : run := false;
           4 : datarecording := true;
         -32 : datarecording := false;
          20 : diagramtables := true;
         -20 : diagramtables := false;
          19 : scandetails := true;
         -31 : scandetails := false;
           5 : if extrapath <> '-' then extrastorage := true;
         -18 : extrastorage := false;
           6 : if developer then failuremessage := 'forced failure';
       end;
       showkeys;
    end; end;

  procedure init_flowrate_meteo_balance; {only before the first cycle}
   const n = 1580;
   var   i : integer;
         f : single;
   begin {requires 6 minutes}
      message ('Air flow warmup (145 s) ...   ');
      startposition;
      fan_on := 1; setcontrols (0);
      for i := 1 to 120 do pause (1000);
      message ('Measuring flowrate zero (40 s) ...   ');
      fan_on := 0; setcontrols (0);
      for i := 1 to 20 do pause (1000);
      f := 0;
      for i := 1 to 120 do f := f + adc (4);
      flowpressure_zero := round (f / 120);
      fan_on := 1; setcontrols (3000);
      message ('Bridge balance (5 s)');
      balance := 7;
      plusiongate_open := 0; minusiongate_open := 0;
      HV_on := 1; setcontrols (1100); {charge the capacitor}
      HV_on := 0; setcontrols (600);
      pause (electrometerdelay);
      f := (adc (2) + adc (2) + adc (2) + adc (2) +
            adc (2) + adc (2) + adc (2) + adc (2)) /
            (8 * electrometergain * balanceconstant);
      balance := round (balance + f);
      if balance < 0 then balance := 0;
      if balance > 15 then balance := 15;
      setcontrols (3000);
      message ('Measuring flowrate (170 s) ...   ');
      pause (12000);
      f := 0;
      for i := 1 to n do f := f + adc (4);
      f := f / n - flowpressure_zero;
      temperature := ctemperaturea * ((adc (7) + adc (7) + adc (7)) / 3)
                   + ctemperatureb;
      pressure := cpressurea * ((adc (5) + adc (5) + adc (5)) / 3)
                + cpressureb;
      flowrate := cflowrate * sqrt (f * (temperature + 273) / pressure);
      message ('Bridge balance (5 s)');
      plusiongate_open := 0; minusiongate_open := 0;
      HV_on := 1; setcontrols (1100); {charge the capacitor}
      HV_on := 0; setcontrols (600);
      pause (electrometerdelay);
      f := (adc (2) + adc (2) + adc (2) + adc (2) +
            adc (2) + adc (2) + adc (2) + adc (2)) /
            (8 * electrometergain * balanceconstant);
      balance := round (balance + f);
      if balance < 0 then balance := 0;
      if balance > 15 then balance := 15;
      setcontrols (3000);
   end;

  procedure estimate_inlet_standardpenetration;
   var                  i : integer;
           gridadsorption,
     shortplateadsorption,
      longplateadsorption,
            longplateloss,
                        c : single;
   begin
   for i := 1 to nfraction - 1 do begin
      c := exp (0.167 * ln (pressure / 1013)) *
           exp (0.333 * ln ((temperature + 273) / 273)) *
           exp (0.667 * ln (0.89 * fractionhighmobility [i]));
      gridadsorption := c * grid_adsorptioncoefficient;
      longplateadsorption := c * longplate_adsorptioncoefficient;
      shortplateadsorption := adsorptionratio * longplateadsorption;
      longplateloss := sqrt (sqr (longplateadsorption) + sqr (edgeloss));
      inlet_standardpenetration [i] := (1 - gridadsorption) *
                              (1 - shortplateadsorption - longplateloss);
   end;
   inlet_standardpenetration [nfraction] :=
                               inlet_standardpenetration [nfraction - 1];
   end;

  procedure wait_for_starttime; {in measurement}
   var a, x : integer;
   begin
      a := 0;
      message ('Waiting for sharp start time ...   ');
      repeat
         x := timereserve;
         if x <> a then begin
            gotoxy (77, 2);
            write (x:3);
         end;
         a := x;
         checkkey;
         if not run then begin startposition; exit end;
      until x = 0;
   end;

  procedure fullwindow;
   begin window (1, 1, 80, 50); end;

  procedure scanwindow;
   begin window (25, 43, 56, 50); end;

  procedure write_diagnostics;
   begin
      fullwindow;
      if (powervoltage < minpower) or (powervoltage > maxpower)
      then textcolor (red + 128) else textcolor (black);
      if developer then begin gotoxy (15, 46); write (powervoltage:5:2) end
                   else begin gotoxy (16, 46); write (powervoltage:4:1) end;
      if (minusionfiltervoltage < minfilter)
      or (minusionfiltervoltage > maxfilter)
      or (-plusionfiltervoltage < minfilter)
      or (-plusionfiltervoltage > maxfilter)
      then textcolor (red + 128) else textcolor (black);
      gotoxy ( 9, 47); write ('+', minusionfiltervoltage:3:0,
                             ' / ', plusionfiltervoltage:4:0);
      if (flowrate < minflow) or (flowrate > maxflow)
      then textcolor (red + 128) else textcolor (black);
      if not flowsensor then textcolor (yellow);
      if developer then begin gotoxy (12, 48); write (flowrate:5:2) end
                   else begin gotoxy (13, 48); write (flowrate:4:1) end;
      if (electrometerbias < minbias) or (electrometerbias > maxbias)
      then textcolor (red + 128) else textcolor (black);
      gotoxy (15, 49); write (electrometerbias:4:1);
      if (balance < minbalance) or (balance > maxbalance)
      then textcolor (red + 128) else textcolor (black);
      gotoxy (20, 50); write (balance:2);
      textcolor (black);
   end;

  procedure calibration; {in measurement}
   var   i, j, k, x : integer;
               r, y : single;
                  z : longint;
              wrong : boolean;
            voltage : array [1..nvoltmax] of integer; {ADC counts}
   begin
      if not run then exit;
      message ('Calibration procedures');
     {measure power and filter voltages}
      HV_on := 0; plusiongate_open := 1; minusiongate_open := 1;
      setcontrols (200); k := adc (8);
      powervoltage := c_powervoltage * k;
      plusiongate_open := 0;
      setcontrols (200);
      if simulator then plusionfiltervoltage := -500
      else plusionfiltervoltage := c_filtervoltage * (adc (8) - k);
      plusiongate_open := 1; minusiongate_open := 0;
      setcontrols (200);
      if simulator then minusionfiltervoltage := 500
      else minusionfiltervoltage := c_filtervoltage * (adc (8) - k);
     {measure relaxation}
      HV_on := 1; plusiongate_open := 0; minusiongate_open := 0;
      setcontrols ((2 * nflowmeasurements + 3) * 100);
      scanwindow; clrscr;
      writeln ('  HV relaxation:');
      textcolor (yellow);
      HV_on := 0;
      setcontrols (0);
      y := 0;
      if simulator then for i := 1 to nvoltmax do begin
         voltage [i] := round (7500 * exp (-i / 40));
         pause (20);
         if (i - 1) mod 10 = 0 then write (voltage [i]:8);
      end
      else for i := 1 to nvoltmax do begin
         x := abs (adc (3));
         voltage [i] := x;
         {ballast operations simulating the calculations during a scan}
         y := y + x;
         if (i - 1) mod 10 = 0 then write (x:8);
      end;
      if oldvoltage [1] <> 0 then for i := 1 to nvoltmax do begin
         voltage [i] := round ((oldvoltage [i] + voltage [i]) / 2);
         oldvoltage [i] := voltage [i];
      end;
      writeln;
     {measure electrometer zero and bias}
      z := 0;
      for i := 1 to 10 do z := z + adc (2);
      electrometerzero := z div 10;
      electrometerbias := 0.305 * (electrometerzero / electrometergain);
         {mV, allowed up to 5}
      {calculate fraction times for scan}
      j := 1;
      if flowsensor then r := flowrate / standardflowrate
                    else r := powervoltage / standardpowervoltage;
      for i := 0 to nfraction do begin
         y := c_voltage * r / fractionhighmobility [i];
         k := 0;
         if voltage [nvoltmax] < y then
            while voltage [j] > y
            do begin j := j + 1; k := k + 1 end;
         fractiontime [i] := k;
      end;
      textcolor (green);
      write ('  Fraction times [0..16]:', fractiontime [0]:7);
      for i := 1 to 16 do write (fractiontime [i]:4);
      textcolor (black);
      wrong := false;
      for i := 1 to nfraction do
         wrong := wrong or (fractiontime [i] < minfractiontime);
      if wrong then failuremessage := 'corrupted HV';
      if (minusionfiltervoltage < minfilter)
      or (minusionfiltervoltage > maxfilter)
      or (-plusionfiltervoltage < minfilter)
      or (-plusionfiltervoltage > maxfilter)
      then failuremessage := 'filter voltage';
      if failuremessage = '' then writeln ('  Electrometer ADC counts:');
      write_diagnostics;
      checkkey;
   end; {of calibration}

   procedure balancebridge (swing : single); {in measurement}
   var inc : single;
   begin
      inc := swing / (electrometergain * balanceconstant);
      if abs (inc) > 0.75 then balance := round (balance + inc);
      if balance < 0 then balance := 0;
      if balance > 15 then balance := 15;
      setcontrols (10);
   end; {of balancebridge}

  procedure mark_time_correct_clock_and_flowsensor_zero; {in measurement}
   {get date and time, prepare beginning of table,
    adjust the clock at 01:05 or 01:07,
    adjust flowpressure_zero at 22:05 or 22:07}
   const daybasis : array [1..12] of integer =
                    (0, 31, 59, 90, 120, 151, 181, 212, 243, 273, 304, 334);
   var yy, dd, hh, mm, ss, ff : word; {for date and time}
                            b : byte;
                       x, day : single;
                            i : integer;
                            s : string;
   begin
     {make beginning of the table}
      getdate (yy, mm, dd, ff); yy := yy mod 100;
      yymmdd := 10000.0 * yy + 100 * mm + dd; {yymmdd}
      day := daybasis [mm] + dd;
      if (yy div 4 = 0) and (mm > 2) then day := day + 1;
      gettime (hh, mm, ss, ff);
      hhmm := 100 * hh + mm; {hhmm}
      dayofyear := day + (hh + mm / 60) / 24; {day}
      if (hh = 1) and (mm < 10) and (not simulator) then begin {adjust clock}
         repeat gettime (hh, mm, ss, ff) until ss = 30;
         settime (hh, mm, ss + clockcorrection, ff);
      end;
      if simulator then exit;
      b := hh;
      if flowsensor and
         (((b in flowzerohours) and (mm < 10))
           or (cyclenumber = flowzerocycle)) then begin
         message ('Measuring flow pressure zero ...   ');
         {requires 10+10+15 = 35 seconds, which skips 2 scans}
         fan_on := 0; setcontrols (10000);
         x := 0; for i := 1 to 100 do x := x + adc (4) / 100;
         fan_on := 1; setcontrols (15000);
         flowpressure_zero := round ((flowpressure_zero + x) / 2);
      end;
   end;

  procedure scan; {single scan in measurement}
   {Measures, stores in fractiontable, and
    displays in scanwindow fraction readings,
    collects sensor readings and parameter of balance swing }
   var h, m, s, t : word;
          i, j, n : integer;
        a, x, sum : single;
             sign : char;
   begin
      gettime (h, m, s, t);
      gotoxy (14, 44);
      if h < 10 then write ('0'); write (h, ':');
      if m < 10 then write ('0'); write (m, ':');
      if s < 10 then write ('0'); write (s);
      if not run then exit;
      if nscan mod 4 = 0  then  plusiongate_open := 1
                          else plusiongate_open := 0;
      if (nscan + 2) mod 4 = 0  then  minusiongate_open := 1
                                else minusiongate_open := 0;
      if not standardregime then begin
         plusiongate_open := 0; minusiongate_open := 0 end;
         {permanently closed inlets in the noise test}
      gotoxy (7, 45); write (nscan:2);
      gotoxy (15, 45);
      if plusiongate_open = 1 then write (' + ions') else
      if minusiongate_open = 1 then write (' - ions') else write ('no ions');
      message ('charging the high voltage capacitor');
      HV_on := 1; setcontrols (0); {charge the capacitor}
     {collect the sensor signals}
      for i := 1 to nflowmeasurements do
         sensorsum [4] := sensorsum [4] + adc (4);
      for i := 5 to 7 do sensorsum [i] := sensorsum [i] + adc (i);
      for i := 1 to nflowmeasurements do
         sensorsum [4] := sensorsum [4] + adc (4);
      HV_on := 0; setcontrols (0); {begin the capacitor relaxation}
      message ('     scanning mobility distribution');
      scanwindow; gotoxy (1, 8);
      sign := 'Z';
      if plusiongate_open = 1 then begin textcolor (red); sign := '+'; end;
      if minusiongate_open = 1 then begin textcolor (blue); sign := '-'; end;
      electrometeroverload := 0;
      if not simulator then {wait for the lower limit of the first fraction}
         pause (electrometerdelay + 100 * fractiontime [0]);
      for i := 1 to nfraction do begin
         n := fractiontime [i]; sum := 0;
         for j := 1 to n do sum := sum + adc (2);
         x := sum / n;
         if simulator then begin
            x := x + x / i
                   + plusiongate_open  * simupos [i] / c_concentration
                   - minusiongate_open * simuneg [i] / c_concentration;
         end;
         fractiontable [nscan, i] := x;
         if i mod 4 = 1 then write (sign, x:7:0)
                        else write (x:8:0);
         if abs (x) > 8190 then electrometeroverload :=
                                electrometeroverload + 1;
      end;
      {collect balance swing parameter, zero scans only}
      if plusiongate_open + minusiongate_open = 0 then begin
         nzero := nzero + 1;
         swingsum := swingsum + fractiontable [nscan, 1]
                              - fractiontable [nscan, 16];
      end;
      {check keys and display overload message}
      fullwindow;
      textcolor (black);
      checkkey;
      gotoxy (50, 4);
      if electrometeroverload > 0 then begin
         textcolor (red + 128);
         write ('NB: overload!');
         textcolor (black);
      end
      else write ('             ');
   end;

  procedure mobilitydistribution; {in measurement}
   {Uses fractiontable [scan, fraction], calculates and writes
    fractionconcentrations into the displaytable}
   var i, j, n: integer;
       a, b, c, d, x, s,
       sumz, min, max : single;
       y : array [1..nscanmax] of single;
       pos, neg : array [1..nfraction] of single;
   begin
      message ('processing the measurements');
      sumz := 0;
      for j := 1 to nfraction do begin
         {c = factor for the ADC to cm-3 conversion}
         if flowsensor then x := standardflowrate / flowrate
                       else x := standardpowervoltage / powervoltage;
         c := c_concentration * x / inlet_standardpenetration [j];
         for i := 1 to nscan do y [i] := fractiontable [i, j];
         s := 0; n := 0; min := 1e9; max := -min;
         for i := 2 to nscan - 1 do if (i mod 4) = 2 then begin {negative}
            x := (y [i - 1] + y [i + 1]) / 2 - y [i];
            if x < min then min := x;
            if x > max then max := x;
            s := s + x;
            n := n + 1;
         end;
         neg [j] := c * (s - min - max) / (n - 2);
         s := 0; n := 0; min := 1e9; max := -min;
         for i := 4 to nscan - 1 do if (i mod 4) = 0 then begin {positive}
            x := y [i] - (y [i - 1] + y [i + 1]) / 2;
            if x < min then min := x;
            if x > max then max := x;
            s := s + x;
            n := n + 1;
         end;
         pos [j] := c * (s - min - max) / (n - 2);
         if (j > 3) and (j < 14) then begin {zero noise of 10 central fractions}
            s := 0; n := 0; min := 1e9; max := -min;
            for i := 3 to nscan - 2 do if (i mod 4) = 1 then begin
               x := y [i] - (y [i - 2] + y [i + 2]) / 2;
               if x < min then min := x;
               if x > max then max := x;
               s := s + x;
               n := n + 1;
            end;
            sumz := sumz + sqr (c * (s - min - max) / (n - 2));
         end;
      end;
      displaytable [column, 7] := noisefactor * sqrt (sumz / 10);
     {store fractions}
      for i := 1 to nfraction do begin
         displaytable [column, 27 + i] := pos [i] * asymmetry;
         displaytable [column, 43 + i] := neg [i] / asymmetry;
      end;
     {calculate concentrations of particles}
      a :=  0; c := 0;
      for i := 1 to nfraction - 7 do begin {for mobility 0.032...0.42}
         a := a + displaytable [column, 27 + i];
         c := c + displaytable [column, 43 + i];
      end;
      displaytable [column, 60] := a;
      displaytable [column, 61] := c;
     {calculate concentrations and average ln (Z) of cluster ions}
      a :=  0; b := 0; c := 0; d := 0;
      for i := nfraction - 5 to nfraction do begin {for mobility 0.56...3.2}
         x := ln (0.866 * fractionhighmobility [i]); {ln (centermobility)}
         a := a + displaytable [column, 27 + i];
         b := b + displaytable [column, 27 + i] * x;
         c := c + displaytable [column, 43 + i];
         d := d + displaytable [column, 43 + i] * x;
      end;
      displaytable [column, 62] := a;
      displaytable [column, 63] := c;
      if (b > a) or (b < -a / 2)
      then displaytable [column, 64] := 0
      else displaytable [column, 64] := exp (b / a);
      if (d > c) or (d < -c / 2)
      then displaytable [column, 65] := 0
      else displaytable [column, 65] := exp (d / c);
   end;

  procedure sizedistribution; {in measurement}
   const
   {Mobilities are numbered from low = 0.032 to high = 3.16 cm2 V-1 s-1}
      nmob = 16;
      mob0 = 0.03162278; {cm2V-1s-1, low limit of first fraction}
      moblogstep = 0.2878231; {ln (10) / 8, 8 fractions per decade}
   {Sizes are numbered from high = 7.5 to low = 0.4 nm during calculations
    but written from low to high in the output file}
      nsize = 10;
      size0 = 7.498942 {nm, high limit of first fraction};
      sizelogstep = 0.2878231; {ln (10) / 8, 8 fractions per decade}
   var
    mobilityfraction : array [0..nmob + 1] of single; {fraction concentrations
                       according to logarithmic mobility scale,
                       edge fractions are for extrapolation}
        sizefraction : array [1..nsize] of single; {fraction concentrations
                   according to logarithmic size scale}
     i, j, left, right : integer;
                moblog : array [-1..nmob + 1] of single;
            sizemoblog : array [0..nsize] of single;
   begin
     {positive particles}
      for i := 1 to nfraction do
         mobilityfraction [i] := displaytable [column, 27 + i];
      {Convert mobilityfraction [1..nmob] to sizefraction [1..nsize]}
      {Extrapolation}
      mobilityfraction [0] := 2 * mobilityfraction [1] - mobilityfraction [2];
      mobilityfraction [nmob + 1] := 0;
      {Boundaries of mobility fractions on scale of ln (mobility)}
      for i := -1 to nmob + 1 do moblog [i] := ln (mob0) + i * moblogstep;
      {Boundaries of size fractions on scale of ln (mobility)}
      for i := 0 to nsize - 1 do sizemoblog [i] :=
         ln (mobility_a (pressure, temperature,
                         exp (ln (size0) - i * sizelogstep)));
      sizemoblog [nsize] := moblog [nmob];

      {moblog               moblog               moblog               moblog
      [left-1]              [left]             [right-1]             [right]
         |  mobilityfr [left] |                    | mobilityfr [right] |
       -----------------===================================-----------------
                       |         sizefraction [i]         |
                  sizemoblog                         sizemoblog
                    [i-1]                               [i]                 }

      {debug
         for i := 0 to nsize do writeln (i:2, exp (sizemoblog [i]):9:3);
         readln; halt;}
      for i := 1 to nsize do begin
         left := 0;
         while moblog [left] < sizemoblog [i-1] do left := left + 1;
         right := 0;
         while moblog [right] < sizemoblog [i] do right := right + 1;
         if right = left then sizefraction [i] :=
            mobilityfraction [left] * (sizemoblog [i] - sizemoblog [i-1]) /
                                 moblogstep
         else sizefraction [i] :=
            mobilityfraction [left] * (moblog [left] - sizemoblog [i-1]) /
                                 moblogstep
          + mobilityfraction [right] * (sizemoblog [i] - moblog [right - 1]) /
                                 moblogstep;
         if (right - left) > 1 then for j := left + 1 to right - 1 do
            sizefraction [i] := sizefraction [i] + mobilityfraction [j];
      end;
      for i := 1 to nsize do
         displaytable [column, 7 + i] := sizefraction [1 + nsize - i];
     {negative particles, see comments above}
      for i := 1 to nfraction do
         mobilityfraction [i] := displaytable [column, 43 + i];
      mobilityfraction [0] := 2 * mobilityfraction [1] - mobilityfraction [2];
      mobilityfraction [nmob + 1] := 0;
      for i := -1 to nmob + 1 do moblog [i] := ln (mob0) + i * moblogstep;
      for i := 0 to nsize - 1 do sizemoblog [i] :=
         ln (mobility_a (pressure, temperature,
                         exp (ln (size0) - i * sizelogstep)));
      sizemoblog [nsize] := moblog [nmob];
      for i := 1 to nsize do begin
         left := 0;
         while moblog [left] < sizemoblog [i-1] do left := left + 1;
         right := 0;
         while moblog [right] < sizemoblog [i] do right := right + 1;
         if right = left then sizefraction [i] :=
            mobilityfraction [left] * (sizemoblog [i] - sizemoblog [i-1]) /
                                 moblogstep
         else sizefraction [i] :=
            mobilityfraction [left] * (moblog [left] - sizemoblog [i-1]) /
                                 moblogstep
          + mobilityfraction [right] * (sizemoblog [i] - moblog [right - 1]) /
                                 moblogstep;
         if (right - left) > 1 then for j := left + 1 to right - 1 do
            sizefraction [i] := sizefraction [i] + mobilityfraction [j];
      end;
      for i := 1 to nsize do
         displaytable [column, 17 + i] := sizefraction [1 + nsize - i];
   end; {of sizedistribution}

  procedure showdata; {in measurement}
   var         i, k : integer;
     yy, mm, dd, ff : word; {required for getdate}
           longdate : longint;
                  x : single;
                  s : string;
   begin
      fullwindow;
      clrscr;
      textcolor (black);
      gotoxy (2, 3); for i := 1 to 78 do write (#205);
      gotoxy (2, 5); for i := 1 to 78 do write (#205);
      gotoxy (2, 11); for i := 1 to 78 do write (#205);
      gotoxy (2, 13); for i := 1 to 78 do write (#205);
      gotoxy (2, 30); for i := 1 to 78 do write (#205);
      gotoxy (2, 32); for i := 1 to 78 do write (#205);
      gotoxy (2, 38); for i := 1 to 78 do write (#205);
      gotoxy (2, 42); for i := 1 to 78 do write (#205);
      for i := 4 to 41 do begin gotoxy (14, i); write (#186); end;
      gotoxy (14, 3); write (#203);
      gotoxy (14, 5); write (#206);
      gotoxy (14, 11); write (#206);
      gotoxy (14, 13); write (#206);
      gotoxy (14, 30); write (#206);
      gotoxy (14, 32); write (#206);
      gotoxy (14, 38); write (#206);
      gotoxy (14, 42); write (#202);
      gotoxy (23, 42); write (#203);
      gotoxy (58, 42); write (#203);
      for i := 43 to 50 do begin
         gotoxy (23, i); write (#186);
         gotoxy (58, i); write (#186);
      end;
      gotoxy (2, 2);
      if standardregime then write (version)
                        else begin
                           textcolor (red);
                           write ('NB: recording of noise!');
                           textcolor (black);
                        end;
      gotoxy (3, 4); write ('Parameter');
      gotoxy (16, 4); write ('Values of parameters ...');
      gotoxy (67, 4); write ('Exit = ');
      textcolor (red);
      write ('Ctrl+X');
      textcolor (black);
      gotoxy (3, 6);  write ('Time HH:MM');
      gotoxy (3, 7);  write ('     T : C');
      gotoxy (3, 8);  write ('    RH : %');
      gotoxy (3, 9);  write ('    p : mb');
      gotoxy (2, 10); write ('+- noise -+');
      gotoxy (3, 12); write ('Mobility', #25);
      gotoxy (16, 12); write ('Mobility fraction concentrations cm-3 ...');
      for i := 1 to nfraction  do begin
         gotoxy (2, 13 + i);
         write (fractionhighmobility [i - 1]:5:3, '-',
                fractionhighmobility [i]:5:3);
      end;
      gotoxy (3, 31); write ('Diameter', #25);
      gotoxy (16, 31); write ('Size fraction concentrations cm-3 ...');
      gotoxy (3, 33); write ('0.42-0.75');
      gotoxy (3, 34); write ('0.75-1.33');
      gotoxy (3, 35); write ('1.33-2.37');
      gotoxy (3, 36); write ('2.37-4.22');
      gotoxy (3, 37); write ('4.22-7.50');
      gotoxy (3, 39); write ('N-particle');
      gotoxy (3, 40); write ('n-cluster');
      gotoxy (3, 41); write ('Z-cluster');
      getdate (yy, mm, dd, ff); longdate := yy;
      gotoxy (2, 43); write ('Date        ', 10000 * longdate + 100 * mm + dd);
      gotoxy (2, 44); write ('Time');
      gotoxy (2, 45); write ('Scan');
      gotoxy (2, 46); write ('Power              V');
      gotoxy (2, 47); write ('Filter             V');
      gotoxy (2, 48); write ('Flow rate        L/s');
      gotoxy (2, 49); write ('E-meter bias      mV');
      gotoxy (2, 50); write ('Balance (0...15)    ');
      if column > 0 then begin
         write_diagnostics;
         for k := 1 to column do begin
            str (round (10000 + displaytable [k, 2]), s);
            gotoxy (13 * k + 9, 6); write (copy (s, 2, 2), ':',
                                           copy (s, 4, 2));
            gotoxy (13 * k + 9, 7); write (displaytable [k, 4]:5:1);
            gotoxy (13 * k + 9, 8); write (displaytable [k, 5]:5:1);
            gotoxy (13 * k + 8, 9); write (displaytable [k, 6]:6:1);
            gotoxy (13 * k + 9, 10);
            x := displaytable [k, 7];
            if x < 999 then  write (x:3:0) else write ('!!!');
            textcolor (red);
            gotoxy (13 * k + 8, 10); write ('+');
            for i := 1 to nfraction do begin
               gotoxy (13 * k + 3, i + 13);
               write (displaytable [k, 27 + i]:6:0);
            end;
            for i := 1 to 5 do begin
               gotoxy (13 * k + 3, i + 32);
               write (displaytable [k, 6 + 2 * i]
                    + displaytable [k, 7 + 2 * i]:6:0);
            end;
            gotoxy (13 * k + 3, 39);
            write (displaytable [k, 60]:6:0);
            gotoxy (13 * k + 3, 40);
            write (displaytable [k, 62]:6:0);
            gotoxy (13 * k + 3, 41);
            if displaytable [k, 64] = 0 then write ('    ? ')
            else write (displaytable [k, 64]:6:2);
            textcolor (blue);
            gotoxy (13 * k + 14, 10); write ('-');
            for i := 1 to nfraction do begin
               gotoxy (13 * k + 9, i + 13);
               write (displaytable [k, 43 + i]:6:0);
            end;
            for i := 1 to 5 do begin
               gotoxy (13 * k + 9, i + 32);
               write (displaytable [k, 16 + 2 * i]
                    + displaytable [k, 17 + 2 * i]:6:0);
            end; {of i}
            gotoxy (13 * k + 9, 39);
            write (displaytable [k, 61]:6:0);
            gotoxy (13 * k + 9, 40);
            write (displaytable [k, 63]:6:0);
            gotoxy (13 * k + 9, 41);
            if displaytable [k, 65] = 0 then write ('    ? ')
            else write (displaytable [k, 65]:6:2);
            textcolor (black);
         end; {of k}
      end; {of column > 0}
      checkkey;
      showkeys;
   end; {of showdata}

  procedure storedata; {in measurement}
   var i : integer;
       x : single;
   begin
      storagecounter := storagecounter + 1;
      for i := 1 to 59 do storage [storagecounter, i] :=
                          displaytable [column, i];
      storage [storagecounter, 60] := powervoltage;
      if abs (plusionfiltervoltage) < minusionfiltervoltage
      then storage [storagecounter, 61] := plusionfiltervoltage
      else storage [storagecounter, 61] := minusionfiltervoltage;
      storage [storagecounter, 62] := flowrate;
      if standardregime then storage [storagecounter, 63] := nscan
                        else storage [storagecounter, 63] := -nscan;
      storage [storagecounter, 64] := balance;
   end; {of storedata}

  procedure savediagramtables;
   var     n : integer;
       table : array [1..20, 0..144] of single;
    measured : array [0..144] of integer;
         row : array [1..27] of single;

     procedure upgradetable;
        var hhmm, m, i, j : integer;
                     a, b : single;
        begin
       {correct zero error if negative value}
        for i := 1 to 10 do begin
           a := row [7 + i];
           b := row [17 + i];
           if (a < 0) and (b > 0) then begin
              if (a + b) < 0 then begin a := a + b; b := 0 end
                             else begin b := a + b; a := 0 end
           end
           else if (b < 0) and (a > 0) then begin
              if (a + b) > 0 then begin a := a + b; b := 0 end
                             else begin b := a + b; a := 0 end
           end;
           row [7 + i] := a;
           row [17 + i] := b;
        end;
       {write row to the table}
        hhmm := round (row [2]);
        m := 60 * (hhmm div 100) + hhmm mod 100;
        j := 1 + m div cycle;
        for i := 1 to 20 do table [i, j] := 8 * row [7 + i];
                                           {8 = 1 / d(log d)}
        end;

     procedure filter; {in savediagramtables}
        var i, j : integer;
            a, b : single;
               x : array [0..144] of single;
        begin
           for j := 0 to n do begin
              measured [j] := 0;
              for i := 1 to 20 do if table [i, j] > 0
                                  then measured [j] := 1;
           end;
           for i := 1 to 20 do begin
             {interpolation}
              x [0] := table [i, 1];
              for j := 1 to n - 1 do begin
                 a := measured [j] * table [i, j] +
                      measured [j + 1] * table [i, j + 1];
                 b := measured [j] + measured [j + 1];
                 if b = 0 then x [j] := 0
                 else begin x [j] := a / b; measured [j] := 1 end;
              end;
              x [n] := table [i, n];
             {grade 1 smoothing}
              for j := 1 to n - 1 do begin
                  a := measured [j - 1] * x [j - 1] +
                       2 * measured [j] * x [j] +
                       measured [j + 1] * x [j + 1];
                  b := measured [j - 1] +
                       2 * measured [j] +
                       measured [j + 1];
                  if b > 1.5 then table [i, j] := a / b;
              end; {of j}
              for j := 0 to n do
              if table [i, j] < 0 then table [i, j] := 0;
           end; {of i}
        end; {of filter}

     procedure fillgaps; {in savediagramtables}
        {uses variables table and cycle}
        const maxjump = 6; {= largest gap + 1}
        var i, j, p, q, d : integer;
        begin
           for p := 0 to n - maxjump do {p = left border}
           if (measured [p] = 1) and (measured [p + 1] = 0) then begin
              {gap found}
              q := p + 2; {q = right border}
              while (measured [q] = 0) and (q < (p + maxjump)) do q := q + 1;
              if measured [q] = 1 then begin  {fill the gap}
                 d := q - p;
                 for j := p + 1 to q - 1 do for i := 1 to 20 do begin
                    table [i, j] :=  {linear interpolation}
                                    ((j - p) / d) * table [i, p]
                                  + ((q - j) / d) * table [i, q];
                 end;
              end;
           end;
        end;

     procedure savetable (filename : string; b : integer); {in savediagramtables}
        var i, j : integer;
               f : text;
        begin
        {$I-}
        assign (f, filename); rewrite (f);
        if IOresult <> 0 then begin
           message ('Error 1 saving diagram table');
           pause (3000);
           exit;
        end;
        for i := 1 to 10 do begin
           for j := 0 to n do begin
              write (f, round (table [b + i, j]));
              if IOresult <> 0 then begin {SI+}
                 message ('Error 2 saving diagram table');
                 pause (3000);
                 close (f);
                 exit;
              end;
              if j = n then writeln (f) else write (f, #09);
           end;
           if IOresult <> 0 then begin {SI+}
              message ('Error 3 saving diagram table');
              pause (3000);
              close (f);
              exit;
           end;
        end;
        close (f);
        {$I+}
        end; {of savetable in savediagramtables}

   var  i, j : integer;

   begin {savediagramtables}
     message ('saving of diagram tables');
     n := 1440 div cycle;
     for i := 1 to 20 do for j := 0 to 144 do table [i, j] := 0;
     for j := 1 to storagecounter do begin
        for i := 1 to 27 do row [i] := storage [j, i];
        upgradetable;
     end;
    {complete the zero column}
     for i := 1 to 20 do table [i, 0] := table [i, 1];
     filter;
     fillgaps;
     if extrastorage then savetable (extrapath + 'p' + datetext + '.xl', 0)
                     else savetable ('tables\p' + datetext + '.xl', 0);
     if extrastorage then savetable (extrapath + 'n' + datetext + '.xl', 10)
                     else savetable ('tables\n' + datetext + '.xl', 10);
   end; {of savediagramtables}

  procedure savedata (allfiles : boolean); {in measurement}
   const heading1 =
      'YYMMDD' + #09 + 'HHMM' + #09 + 'DAY' + #09 +
      'T:C' + #09 + 'RH:%' + #09 + 'p:mb' + #09 +
      'noise' + #09 +
      'd+.42-.56' + #09 + 'd+.56-.75' + #09 +
      'd+.75-1.0' + #09 + 'd+1.0-1.3' + #09 +
      'd+1.3-1.8' + #09 + 'd+1.8-2.4' + #09 +
      'd+2.4-3.2' + #09 + 'd+3.2-4.3' + #09 +
      'd+4.3-5.6' + #09 + 'd+5.6-7.5' + #09 +
      'd-.42-.56' + #09 + 'd-.56-.75' + #09 +
      'd-.75-1.0' + #09 + 'd-1.0-1.3' + #09 +
      'd-1.3-1.8' + #09 + 'd-1.8-2.4' + #09 +
      'd-2.4-3.2' + #09 + 'd-3.2-4.3' + #09 +
      'd-4.3-5.6' + #09 + 'd-5.6-7.5' + #09;
      heading2 =
      'z+.03-.04' + #09 + 'z+.04-.06' + #09 +
      'z+.06-.08' + #09 + 'z+.08-.10' + #09 +
      'z+.10-.13' + #09 + 'z+.13-.18' + #09 +
      'z+.18-.24' + #09 + 'z+.24-.32' + #09 +
      'z+.32-.42' + #09 + 'z+.42-.56' + #09 +
      'z+.56-.75' + #09 + 'z+.75-1.0' + #09 +
      'z+1.0-1.3' + #09 + 'z+1.3-1.8' + #09 +
      'z+1.8-2.4' + #09 + 'z+2.4-3.2' + #09;
      heading3 =
      'z-.03-.04' + #09 + 'z-.04-.06' + #09 +
      'z-.06-.08' + #09 + 'z-.08-.10' + #09 +
      'z-.10-.13' + #09 + 'z-.13-.18' + #09 +
      'z-.18-.24' + #09 + 'z-.24-.32' + #09 +
      'z-.32-.42' + #09 + 'z-.42-.56' + #09 +
      'z-.56-.75' + #09 + 'z-.75-1.0' + #09 +
      'z-1.0-1.3' + #09 + 'z-1.3-1.8' + #09 +
      'z-1.8-2.4' + #09 + 'z-2.4-3.2' + #09 +
      'power'  + #09 + 'filter' + #09 +
      'flowrate' + #09 + 'scans'  + #09 + 'balance';
   var     i, k : integer;
       filename,
              s : string;
              f : text;
   begin
     if storagecounter = 0 then exit;
     message ('saving of the data');
     str (round (1000000 + storage [1, 1]), s);
     datetext := copy (s, 2, 6);
     if storagecounter <> oldcounter then begin {Monthly file}
        assign (f, 'months\' + '3A' + copy (datetext, 1, 4) + '00.' + extension);
        {$I-} append (f); {$I+}
        if IOresult <> 0 then begin
           {$I-} rewrite (f);
           writeln (f, heading1, heading2, heading3);
           {$I+}
            if IOresult <> 0 then begin
               clrscr;
               writeln ('Fatal error initiating monthly file, press ENTER!');
               readln;
               halt;
            end;
        end;
        {$I-}
        for i := 1 to 63 do
          write (f, storage [storagecounter, i]:1:decimals [i], #09);
        writeln (f, storage [storagecounter, 64]:1:decimals [64]);
        oldcounter := storagecounter;
        close (f);
        {$I+}
         if IOresult <> 0 then begin
            message ('Error writing monthly file');
            pause (10000);
         end;
     end; {of monthly file}
     if not allfiles then exit;
    {Diurnal file}
     if extrastorage then filename := extrapath + datetext + '.' + extension
                     else filename := 'days\3A' + datetext + '.' + extension;
     assign (f, filename);
     {$I-} append (f); {$I+}
     if IOresult <> 0 then begin
        {$I-} rewrite (f);
        writeln (f, heading1, heading2, heading3);
        {$I+}
         if IOresult <> 0 then begin
            message ('Error writing diurnal file');
            pause (10000);
         end;
     end;
     {$I-}
     for k := 1 to storagecounter do begin
        for i := 1 to 63 do write (f, storage [k, i]:1:decimals [i], #09);
        writeln (f, storage [k, 64]:1:decimals [64]);
     end;
     close (f);
     {$I+}
     if IOresult <> 0 then begin
        message ('Error writing diurnal file');
        pause (10000);
     end;
    {Diagram tables}
     if diagramtables then savediagramtables;
     storagecounter := 0;
   end; {of savedata}

  procedure savedetails; {in measurement}
   var     i, j : integer;
       filename,
              s : string [12];
              f : text;
   begin
     message ('saving of scan details');
     str (round (1000000 + displaytable [column, 1]), s);
     filename := '3S' + copy (s, 2, 6) + '.' + extension;
     assign (f, 'details\' + filename);
     {$I-} append (f); {$I+}
     if IOresult <> 0 then rewrite (f);
     for i := 1 to 2 do write (f, displaytable [column, i]:1:0, #09);
     write (f, displaytable [column, 3]:5:3, #09);
     for i := 4 to 6 do write (f, displaytable [column, i]:3:1, #09);
     write (f, displaytable [column, 7]:1:0, #09);
     writeln (f, displaytable [column, 64]:1:0);
     for i := 1 to 35 do begin
        for j := 1 to nfraction - 1 do write (f, fractiontable [i, j]:1:0, #09);
        writeln (f, fractiontable [i, nfraction]:1:0);
     end;
     close (f);
   end; {of savedetails}

  procedure conditiondisplaydata; {in measurement}
   var   i, d : integer;
      a, b, c : single;
   begin
     for i := 8 to 60 do if i in [8..17, 28..43, 60] then begin
        if i < 20 then d := 10
        else if i < 50 then d := 16
        else d := 1;
        a := displaytable [column, i];
        b := displaytable [column, i + d];
        c := a + b;
        if a * b < 0 then begin
           if a < 0 then begin
              if c < 0 then begin a := c; b := 0 end
                       else begin a := 0; b := c end;
           end
           else begin
              if c < 0 then begin a := 0; b := c end
                       else begin a := c; b := 0 end;
           end;
           displaytable [column, i] := a;
           displaytable [column, i + d] := b;
        end;
     end;
   end; {of conditiondisplaydata}

  procedure run_cycle; {in measurement}
   var i, j, t : integer;
            pf : single;
     firsthalf : boolean;
   begin
     cyclenumber := cyclenumber + 1;
     if timereserve < (45 * cycle) then wait_for_starttime;
     if not run then exit;
     calibration;
     if failuremessage <> '' then exit;
     if simulator then t := 2 else t := 30;
     nscan := 0;
     nzero := 0;
     for i := 4 to 7 do sensorsum [i] := 0;
     swingsum := 0;
     firsthalf := true;
     while run and (timereserve > t) do begin
        nscan := nscan + 1;
        scan;
        if electrometeroverload > maxoverload then begin
           failuremessage := 'electrometer overload';
           exit;
        end;
        if firsthalf and (timereserve < 30 * cycle) then begin
           firsthalf := false;
           mark_time_correct_clock_and_flowsensor_zero;
           {corrections are made only in special hours}
        end;
     end;
     if run and (nzero > 0) then begin
        column := column + 1;
        if column > 5 then begin
           for i := 1 to 4 do for j := 1 to 65 do
           displaytable [i, j] := displaytable [i + 1, j];
           column := 5;
        end;
        if not simulator then balancebridge (swingsum / nzero);
        temperature := ctemperaturea * sensorsum [7] / nscan + ctemperatureb;
        humidity := (chumiditya * sensorsum [6] / nscan + chumidityb)
                    / (1.0546 - 0.00216 * temperature);
        pressure := cpressurea * sensorsum [5] / nscan + cpressureb;
        pf := sensorsum [4] / (2 * nflowmeasurements * nscan)
              - flowpressure_zero; {pressure in ADC counts, about 2500}
        flowrate := cflowrate *
              sqrt (abs ((pf * (temperature + 273) / pressure)));
        displaytable [column, 1] := yymmdd;
        displaytable [column, 2] := hhmm;
        displaytable [column, 3] := dayofyear;
        displaytable [column, 4] := temperature;
        displaytable [column, 5] := humidity;
        displaytable [column, 6] := pressure;
        estimate_inlet_standardpenetration;
        mobilitydistribution;
        sizedistribution;
        if datarecording then storedata;
        if scandetails then savedetails;
        conditiondisplaydata; {the displaytable data is already stored}
        showdata;
        write (#07); {beep}
     end;
     savedata (midnight);
   end; {of cycle}

  procedure failurepause; {in measurement}
   var h, m, s, t : word;
                k : integer;
         filterOK,
          restart : boolean;
   begin
   repeat
      gotoxy (50, 4);
      write ('             ');
      write_diagnostics;
      textcolor (red + 128);
      message ('Stopped (' + failuremessage + ')');
      textcolor (black);
      HV_on := 0;
      plusiongate_open := 1;
      minusiongate_open := 1; {voltages off}
      fan_on := 0; setcontrols (100);
      repeat {wait for approaching full hour}
         checkkey;
         if not run then exit;
         gettime (h, m, s, t);
         if simulator then restart := (s = 50)
         else restart := (m = 58) and (s = 0)
      until restart;
     {check filter voltages}
      k := adc (8);
      plusiongate_open := 0;
      setcontrols (200);
      if simulator then plusionfiltervoltage := -500
      else plusionfiltervoltage := c_filtervoltage * (adc (8) - k);
      plusiongate_open := 1; minusiongate_open := 0;
      setcontrols (200);
      if simulator then minusionfiltervoltage := 500
      else minusionfiltervoltage := c_filtervoltage * (adc (8) - k);
      filterOK := (minusionfiltervoltage > minfilter) and
                  (minusionfiltervoltage < maxfilter) and
                  (-plusionfiltervoltage > minfilter) and
                  (-plusionfiltervoltage < maxfilter);
      fullwindow;
      if filterOK then textcolor (black)
                  else begin textcolor (red + 128);
                             failuremessage := 'filter voltage';
                       end;
      gotoxy ( 9, 47); write ('+', minusionfiltervoltage:3:0,
                             ' / ', plusionfiltervoltage:4:0);
   until filterOK;
   fan_on := 1; setcontrols (15000);
   wait_for_starttime;
   end; {of failurepause}

Var i, j : integer;
Begin {measurement}
   if simulator then cycle := 1 else cycle := round (cycletime);
   textmode (259);
   textcolor (black);
   textbackground (lightgray);
   clrscr;
   for i := 1 to 64 do decimals [i] := 0;
   for i := 4 to  6 do decimals [i] := 1;
   decimals [ 3] := 3;
   decimals [60] := 1;
   if developer then begin
      decimals [ 7] := 2;
      decimals [60] := 2;
      decimals [62] := 2;
   end;
   storagecounter := 0; oldcounter := 0;
   for i := 0 to nfraction do {0.0316..3.16}
      fractionhighmobility [i] := exp ((i - 12) * ln (10) / 8);
   for i := 1 to nscanmax do for j := 1 to nfraction do
      fractiontable [i, j] := 0;
   column := 0;
   oldvoltage [1] := 0;
   datarecording := (initialcontrols [1] = '1');
   diagramtables := (initialcontrols [2] = '1');
   scandetails   := (initialcontrols [3] = '1');
   extrastorage  := (initialcontrols [4] = '1');
   run := true;
   showdata;
   if not simulator then init_flowrate_meteo_balance
   else begin flowpressure_zero := 0; flowrate := standardflowrate end;
   {flowrate will be corrected after the every scan}
   plusiongate_open := 0;
   minusiongate_open := 0;
   setcontrols (100);
   cyclenumber := 0;
   repeat
      failuremessage := '';
      run_cycle;
      if failuremessage <> '' then failurepause;
   until not run;
   savedata (true);
   startposition;
   clrscr;
   clearkeyboard;
End {of measurement};

Procedure Extratest; {hidden test operation available only for developer}
Type sc = array [1..150] of integer;
Const n = 200;
Var i, j, k, x : integer;
             f : file of sc;
          name : string;
             t : text;
             a : array [1..3000] of integer;
Begin
clrscr;
textmode (259);
startposition;

{writeln ('pressure noise test');
   fan_on := 1; setcontrols (15000);
   fan_on := 1; setcontrols (15000);
   for i := 1 to 3000 do a [i] := adc (4) - 3171;
   assign (t, 'x.txt'); rewrite (t);
   for i := 1 to 3000 do writeln (t, a [i]);
   close (t);
exit;}

{fan_on := 1; setcontrols (15000);
for k := 1 to 3 do begin
   for i := 1 to n do begin
      writeln ('k = ', k, '   i = ', i);
      HV_on := 1; setcontrols (0);
      for j := 1 to 10 do begin
         x := adc (2);
         write (x:8);
      end;
      HV_on := 0; setcontrols (0);
      for j := 1 to 150 do begin
         x := adc (2);
         write (x:8);
         y [i, j] := x;
      end;
      writeln;
   end;
   str (k, name);
   assign (f, name + 'm.int'); rewrite (f);
   for i := 1 to n do write (f, y [i]);
   close (f);
end;}

End; {of extratest}

Procedure Welcome;
   var y, m, d, h, s : word;
   begin
   blankscreen;
   writeln;
   writeln (' Welcome to BSMA3 control and logging program ', version);
   writeln;
   writeln (' Requirements for the computer:');
   writeln ('    running under MS DOS or DOS-regime of Windows9#,');
   writeln ('    free space on disk C for writing of results,');
   writeln ('    BSMA3 connected to PICO ADC-16 and the computer LPT1 port,');
   writeln ('    ADC-16 connected to the computer COM1 port.');
   writeln ('    The program can be interrupted using Ctrl+Break');
   writeln ('    (consider Fn key when working with a laptop).');
   writeln;
   writeln (' Local winter time:');
   getdate (y, m, d, s);
   writeln ('    Year ', y, '   Month ', m:2, '   Day ', d:2);
   gettime (h, m, s, y);
   writeln ('    Hour   ', h:2, '  Minute ', m:2, '   Sec ', s:2);
   writeln;
   writeln (' Selective keys and corresponding tasks are:');
   writeln ('    C - Check and adjust the computer clock,');
   writeln ('    T - Test operations,');
   writeln ('    M - Measure charged particles and clusters,');
   writeln
      ('    N - Noise test (measurement with permanenetly closed inlet gates),');
   writeln ('    X - eXit the program.');
   writeln (' Please press a selective key!');
   end;

CONST          idlelimit = 3; {minutes before autorun of measurement regime}
VAR         taskselector : char;
          adcclose, idle : integer;
                   prime : boolean;
  hh, mm, ss, tt, m0, s0 : word;

BEGIN {main}
prime := true;
blankscreen;
writeln (version);
read_ini;
balance := 7;
startposition;
read_dlf;
repeat
   welcome;
   if prime then writeln (' (Measurement will automatically start after ',
                          'about ', idlelimit, ' idle minutes)');
   gettime (hh, m0, s0, tt);
   repeat
      gettime (hh, mm, ss, tt);
      idle := 60 + mm - m0;
      if idle >= 60 then idle := idle - 60;
      if prime and (idle > idlelimit) then taskselector := 'M'
      else if keypressed then taskselector := upcase (readkey)
      else taskselector := ' ';
      if ss <> s0 then begin
         gotoxy (1, 14);
         write ('    Hour   ', hh:2, '  Minute ', mm:2, '   Sec ', ss:2);
         s0 := ss;
      end;
   until taskselector in ['C', 'T', 'M', 'N', 'Z', 'X'];
   if taskselector = 'N' then begin
      blankscreen;
      gotoxy (12, 12);
      writeln ('NB: N is the noise test, please repeat the selection!');
      repeat until keypressed;
      while keypressed do taskselector := upcase (readkey);
   end;
   blankscreen;
   case taskselector of
      'C' : Checktime;
      'T' : Test;
      'M' : Measurement (true);  {normal}
      'N' : Measurement (false); {noise recording}
      'Z' : if developer then Extratest;
   end;
   prime := false;
until taskselector = 'X';
startposition;
if not simulator then adcclose := adc (-1);
clrscr;
END.