unit Simhh;
{ ==================================================================
  WinWCP - Voltage-activated current simulation (c) J. Dempster 1996
  ================================================================== }
interface

uses
  SysUtils, WinTypes, WinProcs, Messages, Classes, Graphics, Controls,
  Forms, Dialogs, StdCtrls, TabNotBk, ExtCtrls, Global, Shared, FileIO, maths ;

type
  TVClampSim = class(TForm)
    pbDisplay: TPaintBox;
    lbTMin: TLabel;
    lbTMax: TLabel;
    QuantumGrp: TGroupBox;
    Label9: TLabel;
    edNoiseRMS: TEdit;
    GroupBox1: TGroupBox;
    bStart: TButton;
    bAbort: TButton;
    bClose: TButton;
    edProgress: TEdit;
    Timer: TTimer;
    TabbedNotebook1: TTabbedNotebook;
    edGMax: TEdit;
    Label1: TLabel;
    Label6: TLabel;
    edVRev: TEdit;
    Label10: TLabel;
    Label13: TLabel;
    Label14: TLabel;
    Label17: TLabel;
    edPower: TEdit;
    Label18: TLabel;
    GroupBox2: TGroupBox;
    Label8: TLabel;
    edVHold: TEdit;
    edVStep: TEdit;
    edNumRecords: TEdit;
    GroupBox7: TGroupBox;
    edLeakDivBy: TEdit;
    rbNoLeakSubtraction: TRadioButton;
    rbLeakSubtraction: TRadioButton;
    edGLeak: TEdit;
    Label2: TLabel;
    ckInactivation: TCheckBox;
    edHInfVHalf: TEdit;
    Label15: TLabel;
    edTauH: TEdit;
    edHInfVSlope: TEdit;
    edMInfVHalf: TEdit;
    Label12: TLabel;
    edTauM: TEdit;
    edMInfVSlope: TEdit;
    edGSeries: TEdit;
    Label3: TLabel;
    Label24: TLabel;
    Label4: TLabel;
    edCm: TEdit;
    procedure TimerTimer(Sender: TObject);
    procedure bStartClick(Sender: TObject);
    procedure bAbortClick(Sender: TObject);
    procedure bCloseClick(Sender: TObject);
    procedure pbDisplayPaint(Sender: TObject);
    procedure FormClose(Sender: TObject; var Action: TCloseAction);
    procedure FormShow(Sender: TObject);
    procedure FormHide(Sender: TObject);
    procedure FormCreate(Sender: TObject);
  private
    { Private declarations }
    procedure HeapBuffers( Operation : THeapBufferOp ) ;
  public
    { Public declarations }
  end;

function Boltz( V, VHalf, VSlope : Single ) : Single ;
function AlphaFunction( V, VHalf, VSlope, Alpha0 : Single ) : Single ;
function HAlphaFunction( V, Alpha0 : Single ) : Single ;
Function GetParameter( var ed : TEdit ;
                       Default, Min, Max : Single ;
                       const Units : string ) : Single ;
var
  VClampSim: TVClampSim;

implementation

uses mdiform ;

{$R *.DFM}
const
     ChIm = 0 ;
     ChVm = 1 ;
    VMax = 0.2 ; { Note in volts}
type
    TState = (Idle, DoSimulation, ClearDisplay) ;
    THHParameter = record
                 value : single ;
                 Inf : single ;
                 Zero : single ;
                 Tau : single ;
                 p : single ;
                 VSlope : Single ;
                 VHalf : Single ;
                 Alpha : Single ;
                 Alpha0 : Single ;
                 end ;
var
   State : TState ;
   NumRecordsDone,NumRecordsRequired,GroupNumber,NumPulses : LongInt ;
   Rec : LongInt ;
   FirstRecord : Boolean ;
   RH : ^TRecHeader ; { Record header }
   ADC : ^TIntArray ; { Digitised signal buffer }
   GMax,GLeak,GSeries,Cm,TauC,Ic,ITot : Single ;
   VRev,VHold,VStep,VOld,IMax,IScale,VScale,VTest,VLeak,Vp,NoiseRMS: Single ;

   m,h : THHParameter ;
   TimerBusy : boolean ;
   BuffersAllocated : boolean ;{ Indicates if memory buffers have been allocated }


procedure TVClampSim.HeapBuffers( Operation : THeapBufferOp ) ;
{ -----------------------------------------------
  Allocate/deallocation dynamic buffers from heap
  -----------------------------------------------}
begin
     case Operation of
          Allocate : begin
             if not BuffersAllocated then begin
                New(ADC) ;
                New(RH) ;
                BuffersAllocated := True ;
                end ;
             end ;
          Deallocate : begin
             if BuffersAllocated then begin
                Dispose(ADC) ;
                Dispose(RH) ;
                BuffersAllocated := False ;
                end ;
             end ;
          end ;
     end ;


procedure TVClampSim.FormShow(Sender: TObject);
var
   dt : single ;
begin
     { Create buffers }
     HeapBuffers( Allocate ) ;

     Timer.Enabled := True ;
     TimerBusy := False ;
     State := ClearDisplay ;
     bStart.Enabled := True ;
     bAbort.Enabled := False ;
     { dt := (RawFH.dt*RawFH.NumSamples)*0.02*Settings.TScale ;
     edTauM.text := format( '%.3g %s', [dt,Settings.TUnits]);
     edTauH.text := format( '%.3g %s', [dt*10.,Settings.TUnits]); }
     end;


procedure TVClampSim.TimerTimer(Sender: TObject);
var
   i,j,iStart,iEnd,xPix,yPix,ch,ChOffset,MaxCh, nOffScreen : LongInt ;
  x,y,dx,NoiseTemp,Vm,Im,Gm : Single ;
begin
     if not TimerBusy then begin
        TimerBusy := True ;
        case State of
          ClearDisplay : Begin
                       { Erase display }
                       lbTMin.Left := pbDisplay.Left ;
                       lbTMin.Top := pbDisplay.Top + pbDisplay.Height + 2 ;
                       lbTMax.Left := pbDisplay.Left + pbDisplay.Width
                                      - canvas.TextWidth(lbTMax.caption) ;
                       lbTMax.Top := lbTMin.Top ;
                       lbTMin.caption := '' ;
                       lbTMax.caption := '' ;
                       EraseDisplay(pbDisplay) ;
                       State := Idle
                       end ;

          DoSimulation : Begin

              { Initialisations when first record is done }
              if FirstRecord then begin
                 WriteToLogFile( 'Hodgkin-Huxley Voltage Clamp Simulation ' ) ;
                 GMax := GetFromEditBox(edGMax, 1., 0., 1E4, '%.1f','nS' ) * nStoSeimens ;
                 WriteToLogFile( 'GMax = ' + edGMax.text ) ;
                 GLeak := GetFromEditBox( edGLeak, 1., 0., 1E4,'%.1f', 'nS' ) * nStoSeimens ;
                 WriteToLogFile( 'GLeak = ' + edGLeak.text ) ;
                 GSeries := GetFromEditBox( edGSeries, 100., 0., 1E4,'%.1f', 'nS' ) * nStoSeimens ;
                 WriteToLogFile( 'GSeries = ' + edGSeries.text ) ;
                 Cm := GetFromEditBox( edCm, 1., 0., 1E4,'%.1f', 'pF' ) * pFToFarad ;
                 WriteToLogFile( 'Cm = ' + edCm.text ) ;
                 VRev := GetFromEditBox(edVRev, 40., -200., 200.,'%.1f', 'mV' ) * mVToVolts ;
                 WriteToLogFile( 'VRev = ' + edVRev.text ) ;
                 VHold := GetFromEditBox(edVHold, -90., -200., 200.,'%.1f', 'mV' ) * mVToVolts ;
                 VStep := GetFromEditBox(edVStep, 10., -200., 200.,'%.1f', 'mV' ) * mVToVolts ;
                 { Activation gate (m) parameters }
                 m.VHalf := GetFromEditBox(edMInfVHalf,-25.,-200.,200.,'%.1f','mV')*mVToVolts ;
                 WriteToLogFile( 'm V.half = ' + edMInfVHalf.text ) ;
                 m.VSlope := GetFromEditBox(edMInfVSlope,10.,1.,200.,'%.1f','mV')*mVToVolts ;
                 WriteToLogFile( 'm V.slope = ' + edMInfVSlope.text ) ;
                 m.Tau := GetFromEditBox(edTauM,
                          1E-3*Settings.TScale, 1E-4*Settings.TScale,Settings.TScale,
                          '%.3g', Settings.TUnits ) * Settings.TUnScale ;
                 WriteToLogFile( 'm Tau = ' + edTauM.text ) ;
                 m.Alpha := 0.5 / m.Tau ;
                 m.Alpha0 := AlphaFunction( m.VHalf, m.VHalf, m.VSlope, 1. ) ;
                 m.Alpha0 :=  m.Alpha / m.Alpha0 ;
                 m.P := GetFromEditBox( edPower, 3., 1., 10., '%.1f', '' ) ;
                 WriteToLogFile( 'm Power factor = ' + edPower.text ) ;
                 { Sampling interval }
                 RawFH.dt := (m.Tau*6. ) / RawFH.NumSamples ;
                 if ckInactivation.checked then begin
                    { Inactivation gate (h) parameters }
                    h.VHalf := GetFromEditBox(edHInfVHalf,-25., -200., 200.,
                                               '%.1f','mV')*mVToVolts ;
                    WriteToLogFile( 'h V.Half = ' + edHInfVHalf.text ) ;
                    h.VSlope := GetFromEditBox(edHInfVSlope,10., 1., 1E3 ,
                                               '%.1f','mV' ) *mVToVolts ;
                    WriteToLogFile( 'h V.Slope = ' + edHInfVSlope.text ) ;
                    h.Tau := GetFromEditBox(edTauH,
                             1E-2*Settings.TScale, 1E-3*Settings.TScale,Settings.TScale,
                             '%.3g', Settings.TUnits) * Settings.TUnScale ;
                    WriteToLogFile( 'h Tau = ' + edTauH.text ) ;
                    h.Alpha := 0.5 / h.Tau ;
                    h.Alpha0 := HAlphaFunction( h.VHalf, 1. ) ;
                    h.Alpha0 :=  h.Alpha / h.Alpha0 ;

                    h.P := 1. ;
                    RawFH.dt := MaxFlt([(h.Tau*6. )/RawFH.NumSamples,RawFH.dt]) ;
                    end ;

                 { Capacity current time constant }
                 TauC := Cm / GSeries ;
                 { Background noise RMS }
                 NoiseRMS := ExtractFloat( edNoiseRMS.text, 0. ) * nAToAmps ;


                 { Set channel parameters }

                 { If this is the first record in the file create
                   an appropriate scaling factor }
                 if RawfH.NumRecords = 0 then begin

                    RawfH.NumChannels := 2 ;
                    Channel[0].ChannelOffset := 0 ;
                    Channel[1].ChannelOffset := 1 ;

                    { Set current channel scaling factor }
                    IMax := 0. ;
                    for i := 1 to NumRecordsRequired do begin
                        Vm := VHold + VStep*i ;
                        m.Inf := Boltz( Vm, m.VHalf, m.VSlope ) ;
                        Im := Abs( GMax*Power( m.Inf, m.p )*(Vm-VRev)
                                   + GLeak*Vm ) ;
                        if IMax < Im then IMax := Im ;
                        end ;
                    IMax := IMax*1.1 ;
                    if IMax = 0. then IMax := 1E-9 ;
                    Channel[ChIm].ADCScale := AmpsTonA*IMax / (MaxADCValue) ;
                    rH^.ADCVoltageRange[ChIm] := 1. ;
                    Channel[ChIm].ADCCalibrationFactor := rH^.ADCVoltageRange[ChIm] /
                                  ( Channel[ChIm].ADCScale * (MaxADCValue+1) ) ;

                    { Set voltage channel scaling factor }
                    Channel[ChVm].ADCScale := (VMax*VoltsTomV) / (MaxADCValue) ;
                    rH^.ADCVoltageRange[ChVm] := 1. ;
                    Channel[ChVm].ADCCalibrationFactor := rH^.ADCVoltageRange[ChVm] /
                                  ( Channel[ChVm].ADCScale * (MaxADCValue+1) ) ;
                    end
                 else begin
                      GetRecordHeaderOnly( RawfH, rH^, RawfH.NumRecords ) ;
                      GroupNumber := Trunc(rH^.Number) ;
                      IMax := Channel[ChIm].ADCScale*MaxADCValue*nAToAmps ;
                      end ;

                 { Set channel magnification }
                 for ch := 0 to RawfH.NumChannels-1 do begin
                     Channel[Ch].ADCZero := 0 ;
                     Channel[Ch].xMin := 0 ;
                     Channel[Ch].xMax := RawfH.NumSamples-1 ;
                     Channel[Ch].yMin := MinADCValue ;
                     Channel[Ch].yMax := MaxADCValue ;
                     Channel[Ch].InUse := True ;
                     end ;

                 Channel[CHVm].ADCName := 'Vm' ;
                 Channel[CHVm].ADCUnits := 'mV' ;
                 Channel[CHIm].ADCName := 'Im' ;
                 Channel[CHIm].ADCUnits := 'nA' ;

                 GroupNumber := 0 ;

                 { Set voltage channel to fixed zero level }
                 Channel[ChVm].ADCZeroAt := -1 ;
                 CHannel[ChVm].ADCZero := 0 ;
                 { Current channel zero level calculated frm samples 0..19}
                 Channel[ChIm].ADCZeroAt := 0 ;
                 SaveHeader( RawFH ) ;
                 FirstRecord := False ;
                 end ;

              { Volts --> bits scale factor }
              VScale := Channel[ChVm].ADCScale / VoltsTomV ;
              { Amps --> bits scale factor }
              IScale := Channel[ChIm].ADCScale / AmpsTonA ;

              VTest := (NumRecordsDone+1)*VStep ;
              iStart := RawFH.NumSamples div 10 ;
              iEnd := RawFH.NumSamples - iStart ;
              NoiseTemp := 1. ;
              Inc(GroupNumber) ;

              { ** Create leak subtraction pulse ** }

              if rbLeakSubtraction.checked then begin
                 { ** P/N leak subtraction ** }
                 NumPulses := ExtractInt(edLeakDivBy.text)  ;
                 If (NumPulses = 0) or (Abs(NumPulses)>200) then NumPulses := -4 ;
                 VLeak := VTest / NumPulses ;

                 {Create leak pulse record }
                 VOld := VHold ;
                 x := 0. ;
                 Ic := 0. ;
                 for i := 0 to RawFH.NumSamples-1 do begin

                     if (i >= iStart) and (i <= iEnd) then Vp := VHold + VLeak
                                                      else Vp := VHold ;

                     {Voltage channel }
                     j := i*RawFH.NumChannels + ChVm ;
                     ADC^[j] := Trunc( Vp/VScale ) ;

                     { Find membrane potential after voltage drop along
                     pipette series resistance }
                     Vm := Vp / ( 1. + GLeak/GSeries ) ;

                     { Current channel }
                     Im := (Vm * GLeak) +
                           (GaussianRandom(NoiseTemp)*NoiseRMS) /
                            sqrt(MaxFlt([Abs(NumPulses),1. ])) ;

                     { Add capacity current }
                     if Vp <> VOld then begin
                        Ic := (Vp - Vold) * GSeries ;
                        x := 0. ;
                        end ;
                     VOld := Vp ;
                     Im := Im + Ic*exp( -x/TauC ) ;

                     { Keep within display limits }
                     Im := MaxFlt([MinFlt([Im,IMax]),-IMax]) ;
                     x := x + RawFH.dt ;

                     j := i*RawFH.NumChannels + ChIm ;
                     ADC^[j] := Trunc( Im/IScale ) ;
                     end ;

                 { Save leak record to file }
                 Inc(RawFH.NumRecords) ;
                 rH^.Status := 'ACCEPTED' ;
                 rH^.RecType := 'LEAK' ;
                 rH^.Number := GroupNumber ;
                 rH^.Time := rH^.Number ;
                 rH^.dt := RawfH.dt ;
                 rH^.Equation.Available := False ;
                 rH^.Analysis.Available := False ;
                 rH^.Ident := ' ' ;
                 PutRecord( RawfH, rH^, RawfH.NumRecords, ADC^ ) ;
                 end ;

              { ** Create test pulse ** }

              x := 0. ;
              Gm := GLeak ;
              Vm := VHold / ( 1. + Gm/GSeries ) ;
              m.value := Boltz( Vm, m.VHalf, m.VSlope ) ;
              if ckInactivation.checked then h.value := Boltz(Vm,h.VHalf,-h.VSlope) ;
              Ic := 0. ;
              VOld := VHold ;
              for i := 0 to RawFH.NumSamples-1 do begin

                  if (i > iStart) and (i < iEnd) then Vp := VHold + VTest
                                                 else Vp := VHold ;
                  {Voltage channel }
                  j := i*RawFH.NumChannels + ChVm ;

                  { Find membrane potential after voltage drop along
                    pipette series resistance }

                  Vm := Vp / ( 1. + Gm/GSeries ) ;
                  ADC^[j] := Trunc( Vp/VScale ) ;
                  { Current channel }
                  { ACTIVATION KINETICS }
                  {Calculate m.zero & m.infinity from Boltzmann function }
                  m.Inf := Boltz( Vm, m.VHalf, m.VSlope ) ;
                  { Calculate activation rate constant }
                  m.Alpha := AlphaFunction( Vm, m.VHalf, m.VSlope, m.Alpha0 ) ;
                  { Calculate time constant }
                  m.Tau := m.Inf / m.Alpha ;
                  { Calculate value of activation parameter }
                  m.value := m.Inf - (m.Inf - m.value)*exp( -RawFH.dt/m.Tau ) ;
                  { Calculate conductance }
                  Gm := GMax*Power(m.value,m.p) ;

                  { INACTIVATION KINETICS }
                  if ckInactivation.checked then begin
                        h.Inf := Boltz( Vm, h.VHalf, -h.Vslope ) ;
                        h.Zero := Boltz( VHold, h.VHalf, -h.Vslope ) ;
                        h.Alpha := HAlphaFunction( Vm, h.Alpha0 ) ;
                        h.Tau := h.Inf / h.Alpha ;
                        h.value := h.Inf - (h.Inf - h.value)*exp( -RawFH.dt/h.Tau ) ;
                        Gm := Gm*h.value ;
                        end ;
                  { Add leak conductance }
                  Gm := Gm + GLeak ;
                  { Calculate current }
                  Im := (Vm - VRev)*Gm + GaussianRandom(NoiseTemp)*NoiseRMS ;
                  { Add capacity current }
                  if Vp <> VOld then begin
                     Ic := (Vp - Vold) * GSeries ;
                     x := 0. ;
                     end ;
                  VOld := Vp ;
                  x := x + RawFH.dt ;

                  Im := Im + Ic*exp( -x/TauC ) ;
                  { Keep within display limits }
                  Im := MaxFlt([MinFlt([Im,IMax]),-IMax]) ;

                  j := i*RawFH.NumChannels + ChIm ;
                  ADC^[j] :=  Trunc( Im/IScale ) ;
                  end ;

              { Display simulated record }

              InitializeDisplay( Channel, RawFH.NumChannels, rH^,
                                 lbTMin,lbTMax, pbDisplay) ;

              { Erase display }
              EraseDisplay(pbDisplay) ;

              MaxCh := RawFH.NumChannels - 1;
              for ch := 0 to RawFH.NumChannels-1 do begin
                  if Channel[ch].InUse then begin
                     dx := Channel[ch].xScale ;
                     x := -Channel[ch].xScale*Channel[ch].xMin + Channel[ch].Left ;
                     nOffScreen := 0 ;
                     for i := 0 to RawFH.NumSamples-1 do begin

                         y := ADC^[(i*RawFH.NumChannels) + ch ] ;
                         xPix := Trunc(x) ;
                         yPix := Trunc(Channel[ch].Bottom
                                 - Channel[ch].yScale*(y - Channel[ch].yMin));

                         if (xPix < Channel[ch].Left)
                            or (Channel[ch].Right < xPix )
                            or (yPix < Channel[ch].Top)
                            or (Channel[ch].Bottom < yPix ) then begin
                            xPix := MinInt( [ MaxInt( [xPix,Channel[ch].Left]),
                                              Channel[ch].Right] ) ;
                            yPix := MinInt( [ MaxInt( [yPix,Channel[ch].Top]),
                                              Channel[ch].Bottom] ) ;
                            nOffScreen := nOffScreen + 1 ;
                            end
                         else nOffScreen := 0 ;

                         if (nOffScreen > 1) or (i=0) then
                              pbDisplay.canvas.moveto(xPix,yPix)
                         else pbDisplay.canvas.lineto(xPix,yPix);

                         x := x + dx ;
                         end ;

                     { Display Channel Name }
                     pbDisplay.Canvas.TextOut( Channel[ch].Left,
                              (Channel[ch].Top + Channel[ch].Bottom) div 2,
                              ' ' + Channel[ch].ADCName ) ;

                     { Draw cursors }
                     DrawHorizontalCursor(pbDisplay,Channel[ch],
                                          Channel[ch].ADCZero) ;
                     end ;
                  end ;

              { Save Record to file }
              Inc(RawFH.NumRecords) ;
              rH^.Status := 'ACCEPTED' ;
              rH^.RecType := 'TEST' ;
              rH^.Number := GroupNumber ;
              rH^.Time := rH^.Number ;
              rH^.dt := RawfH.dt ;
              rH^.Equation.Available := False ;
              rH^.Ident := ' ' ;
              PutRecord( RawfH, rH^, RawfH.NumRecords, ADC^ ) ;

              { Terminate when all records done }
              Inc(NumRecordsDone) ;
              edProgress.text := format( 'Rec. %d/%d Vm= %f mV', [NumRecordsDone,
                                          NumRecordsRequired,
                                          (VHold+VStep*NumRecordsDone)*VoltsTomV ] ) ;
              if NumRecordsDone >= NumRecordsRequired then bABort.click ;
              end ;
          end ;
        TimerBusy := False ;
        end ;
     end;


procedure TVClampSim.bStartClick(Sender: TObject);
begin
     State := DoSimulation ;
     bStart.Enabled := False ;
     bAbort.Enabled := True ;
     NumRecordsRequired := Trunc(GetFromEditBox(edNumRecords,16.,1.,1E6,'%.0f','')) ;
     NumRecordsDone := 0 ;
     FirstRecord := True ;
     end;

procedure TVClampSim.bAbortClick(Sender: TObject);
begin
     State := Idle ;
     bStart.Enabled := True ;
     bAbort.Enabled := False ;
     if RawFH.NumRecords > 0 then SaveHeader( RawFH ) ;
     WriteToLogFile( format('(%d records)',[RawFH.NumRecords]));
     end;


procedure TVClampSim.bCloseClick(Sender: TObject);
begin
     close ;
     end;

function Boltz( V, VHalf, VSlope : Single ) : Single ;
begin
     Boltz := 1. / ( 1. + Exp( -(V - VHalf)/VSlope ) ) ;
     end ;

function AlphaFunction( V, VHalf, VSlope, Alpha0 : Single ) : Single ;
var
   A0,A1,V0,V1 : Single ;
begin
     if Abs(VHalf-V) > 1E-5 then
        AlphaFunction := (Alpha0*(VHalf-V)) / (exp((VHalf-V)/VSlope) - 1. )
     else begin
        V0 := -1E-5 ;
        A0 := (Alpha0*(VHalf-V0)) / (exp((VHalf-V0)/VSlope) - 1. ) ;
        V1 := 1E-5 ;
        A1 := (Alpha0*(VHalf-V1)) / (exp((VHalf-V1)/VSlope) - 1. ) ;
        AlphaFunction := (A0+A1)*0.5 ;
        end ;
     end ;

function HAlphaFunction( V, Alpha0 : Single ) : Single ;
begin
     HAlphaFunction := Alpha0*exp( -(V-0.06) / 0.02 ) ;
     end ;

Function GetParameter( var ed : TEdit ;
                       Default, Min, Max : Single ;
                       const Units : string ) : Single ;
var
   Value : single ;
begin
     Value := ExtractFloat( ed.text, Default ) ;
     if Value < Min then Value := Abs(Value) ;
     if Value < Min then Value := Min ;
     if Value > Max then Value := Max ;
     ed.text := format( '%f %s', [Value,Units] ) ;
     GetParameter := Value ;
     end ;


procedure TVClampSim.pbDisplayPaint(Sender: TObject);
begin
     State := ClearDisplay ;
     end;

procedure TVClampSim.FormClose(Sender: TObject; var Action: TCloseAction);
begin

     HeapBuffers( Deallocate ) ;

     if RawFH.NumRecords > 0 then begin
        Main.ShowRaw.Enabled := True ;
        Main.ShowRaw.checked := True ;
        Main.ShowAveraged.checked :=  False ;
        Main.ShowAveraged.Enabled :=  False ;
        Main.ShowLeakSubtracted.checked := False ;
        Main.ShowLeakSubtracted.Enabled := False ;
        FH := RawFH ;
        Main.UpdateWindows ;
        end ;

     { Enable "Hodgkin-Huxley" item in Simulations" menu}
     Main.HHSimulation.enabled := true ;
     Action := caFree ;
     end;



procedure TVClampSim.FormHide(Sender: TObject);
begin
     HeapBuffers( Deallocate ) ;
     Timer.Enabled := False ;
     end;

procedure TVClampSim.FormCreate(Sender: TObject);
begin
     { Disable "Hodgkin-Huxley" item in Simulations" menu}
     Main.HHSimulation.enabled := false ;
     end;

end.
