% Calc for acid pump control 
% 
% This is based on Ed Peltzers calcs for the HCL concentration and
% predicted currents
% 
% copyright (c) 2005 by W.J.Kirkwood of MBARI
% February 2, 2005 
% Last Modified 3/6/05 
% Revision: 1.0 
%
% ///////////////////////////////////\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\
% 
% clean up from last activity 
clear all ; clf;
%
pH_SetPoint = 7.3;                % pH target value for control volume 
pH_Valid_Max = 8.0;               % Maximum allowable value for inclusion in average
pH_Valid_Min = 6.0;               % Minimum allowable value for inclusion in average
PumpH20 = 720;                    % Sets the water pump at 720 mL/min
PumpAcid = 720;                   % Sets the acid pump at 720 mL/min
PumpFixedFlow = 1440;             % Sets the total pump flow at maximum current

% Main Control Algorithm for Pumps
% Concept: We would like the pump control to be exponential as the pH
% wanders from the set point. This makes our control take the form 
% y = f(x) + K = a*2^(b*(x-h))+k 
% y is the acid pump commanded output

% 
% Test matrices used to create pH change (bounce) 
PP(1,:) = 1:1000;
PP(2,:) = sin(rand (1,1000));
PP(3,:) = .01:.01:10;

I = 1;
while I < 1000;
    J = PP(2,I);
    if J > .40
        PP(4,I) = 0.0;
    else 
        PP(4,I) = J;
    end
    I = I + 1; 
end
        
figure (1)                          % Plot for visual purposes during development 
plot (PP(1,:), PP(2,:)) 
hold on
plot (PP(1,:), PP(4,:), 'r')


%
% a = scalar multiplier and moves the zero crossover on the x axis, also
% change of sign of a flips the control curve around the x axis
a = .1;
% b = exponential multiplier and sets the rate of change of the slope of
% the control curve, also a sign change flips the control curve around the y axis
b = 1.5;
% h is the x offset constant and sets the x crossover at y = 0
h = 1.0;
% K is the y offset constant and sets the y crossover at f(x) = 0
K = 0;                  % scaled in ml/min 
% 
MaxCurrent = 30; 

% ////////////////////// Pump Flow Control \\\\\\\\\\\\\\\\\\\\\\\
% Two pumps run as an inverse ratio to maintain a set over flow volume
% uses the estimated current relationship. 
% ............ Note: Pumps ordered are 4.8 ml / rev ............... 
Current = input (' Enter current between 10 and 30 cm/sec : ');  

    if Current < 10
        Current = 10;    
    elseif Current > 30
        Current = 30;
    else
        Current = Current;
    end
Current   % Prints current for check for input by user

% Bounds the input to allowable current ranges 
   if Current == 10
          Total_Flow = 480;             % Minimum ml/min from combined pumps
          PumpAcid = Total_Flow / 2;
          % PumpFixedFlow =  Total_Flow - PumpAcid;  % Default total output from pumps in ml/sec
   elseif Current == 30
          Total_Flow = 1440;             % Maximum ml/min from combined pumps
          PumpAcid = Total_Flow / 2;
          % PumpFixedFlow =  Total_Flow - PumpAcid;  % Default total output from pumps in ml/sec
   else 
          Total_Flow = 1440 * (Current / MaxCurrent); % Minimum ml/min from combined pumps
          PumpAcid = Total_Flow / 2;
          % PumpFixedFlow =  Total_Flow - PumpAcid;  % Default total output from pumps in ml/sec
   end

% Quick Check on Pump Numbers
Total_Flow
PumpAcid

% Lets user input pH setpoint within allowable bounds
pH_SetPoint = input (' Entet pH desired set point within 6.0 and 8.0 : '); 
    if pH_SetPoint < 7.0
        pH_SetPoint = 7.0;    
    elseif pH_SetPoint > 8.0
        pH_SetPoint = 8.0;
    else
        pH_SetPoint = pH_SetPoint;
    end

% Sets dead band amounts - can be indepentantly adjusted
dbLow = .10;
dbHigh = .10;
    
% pH deadband set: lower is - .1 pH unit and upper is + .1 pH unit 
pH_BoundLow = pH_SetPoint - dbLow;
pH_BoundHigh = pH_SetPoint + dbHigh; 

% Scales bounds for later use 
pHdbLow = - 10 * dbLow; 
pHdbHigh = 10 * dbHigh;
    
% Initiate indices for counters 
i = 1;

% Set scaling for pump adjust
while i < 1000                 % iterates for testing = 1000 delta pH readings
    
  % create too high and too low pH conditions  
    if i < 200 | i > 600              
        % fake high pH measurement for testing
        pH_Measured = pH_SetPoint + PP(4,i);   
        % fake low pH measurement for testing
    else
        pH_Measured = pH_SetPoint - PP(4,i);        
    end
    
% takes measurements and finds delta pH - basically undoes above for later purposes
        pH_Delta = pH_Measured - pH_SetPoint;
        
% boolean for sign of pH_Delta (needed later for up or down corrections)       
        Y = sign (pH_Delta);
        
% lets one equation calc for up or down (uses symmetry of correction curve)      
        if Y == -1                                  
            pH_Deltb = -1 * pH_Delta; 
        else
            pH_Deltb = pH_Delta; 
        end
        
% scaling makes pH change correlate to correction curve  
        pH_Control = 10 * pH_Deltb;
        
% limits max corection to 100% for pH meausrements well away from set point       
            if pH_Control > 3.3                         
                pH_Control = 3.3;
            end 
            
% just stuff to make following math easier for me            
        pH_expo = b*(pH_Control - h);                   
        pH_X = a*(2.0^pH_expo);
        
% sets the percentage to increase acid flow - I left in the constant in case it is needed later        
        PumpAdjNumUp  = pH_X + K;                       
        PumpAdjNum = PumpAdjNumUp;
        
% sets limit of 100% pump correction
        if PumpAdjNumUp > 1                             
            PumpAdjNum = 1;
        end
        
% checks and sets up or down pump correction
        if Y == -1                  
            % sets the correction direction
            PumpAdjNum = -1* PumpAdjNum;
            % for plotting purposes, shows acid reduction more clearly
            pH_Control = -1* pH_Control;                
        end

    % ..... inserts model Data for plotting purposes ...... 
        pH_Data (1,i) = PumpAdjNum;
        pH_Data (2,i) = PP(4,i);
        pH_Data (3,i) = pH_Measured
        pH_Data (4,i) = pH_Delta;    
        pH_Data (5,i) = pH_Control;
        % used for giving dead band line + side on graph
        pH_Data (6,i) = pH_BoundHigh - pH_SetPoint; 
        % used for giving dead band line - side on graph
        pH_Data (7,i) = -pH_Data (6,i);                       
        pH_Data (8,i) = i;
        
% Corretion to acid pump  
PumpMod = PumpAdjNum; 

% Stores PumpMod for later analysis 
pH_Data (9,i) = PumpMod; 

    % sets first adjustment as no change 
    if i == 1
        PumpAcidAdj = 1 * PumpAcid; 
        PumpAcidAdj
        NUMB = input (' Any num : ')
        
    % Checks and installs dead band 
    %elseif pHdbLow < pH_Control < pHdbHigh
            % makes it easy  to look at past PumpMod cells by subtract number
            % j = i - 1;    
            % keeps pump set at the same amount 
            % PumpAcid = PumpAcid * pH_Data (9,j);     
            
    else
        % outside of deadband multipier is applied to modify acid pump
        PumpAcidAdj = PumpAcid * PumpMod;   
        
    end 

% Calcs Sea Water pump 
PumpAcidNew = PumpAcid + PumpAcidAdj;
PumpH20New = Total_Flow - PumpAcidNew; 

% Stores Pump adjustments 
pH_Data (10,i) = PumpAcidNew;
pH_Data (11,i) = PumpH20New; 
     
% iterates counter
    i = i+1;
end

%  Plotting of figures for visual of performance
figure(2) 
plot (pH_Data(8,:,:),pH_Data(4,:,:))

figure(3)
scatter (pH_Data(5,:),pH_Data(1,:))
hold on
plot (pH_Data(4,:,:),pH_Data(6,:,:),'--r')
plot (pH_Data(4,:,:),pH_Data(7,:,:),'--r')

figure(4)
scatter (pH_Data(8,:,:),pH_Data(5,:,:))

figure(5) 
plot (pH_Data(8,:,:),pH_Data(10,:,:),'--g')
hold on
plot (pH_Data(8,:,:),pH_Data(11,:,:),'r')
