CAL = struct;

% Extract lines matching pattern
[indices, selectedLines] = indexLinesWithPatternO2("\\thalassa\ProjectLibrary\901805_Coastal_Biogeochemical_Sensing\Spray_Data\26702901\calibration_data\20260720_SN029_Aircal", 'err=   0');

% Convert selected lines to cell array of strings
data = cellstr(selectedLines);

% Define substring to remove
substringToRemove = 'err=   0';

% Initialize variables
numRows = length(data);
cleanedData = NaN(numRows, 3);

% Clean data and split into components
for i = 1:numRows
    % Remove substring and trim the line
    cleanLine = strtrim(strrep(data{i}, substringToRemove, ''));
    % Split cleaned line into components
    try
        cleanedData(i, 1:3) = str2double(strsplit(cleanLine));
    catch
        disp(['Check data for inconsistent array sizes or missing lines.'...
            ' Missing values are handled as NaNs'])
        continue
    end
end

% Reshape cleanedData for vectorized operations
numIterations = numRows/4;
aircal = NaN(4, 3, numIterations);
for i = 1:numIterations   
    startIdx = (i - 1) * 4 + 1;
    endIdx = startIdx + 4 - 1;
    aircal(:, :, i) = cleanedData(startIdx:endIdx, :);
end

% Calculate means and standard deviations
aircalMeans = mean(aircal, 3, 'omitnan');
aircalStds = std(aircal, 0, 3, 'omitnan');

% Ctd p, ctd t, ctd Sal
% Raw counts for those 3
% O2calc ml/l, O2temp, rawphase err = 0,
% Raw counts for 3 line items
%    1        2         3
% 1 CTD_p   | CTD t  | CTD s
% 2 raw     | raw    | raw
% 3 O2 mL/L | O2 t   | rawphase
% 4 raw     | raw    | raw

% Extract specific values and their uncertainties
CAL.Raw_Phase = aircalMeans(3, 3);
CAL.Raw_Thermistor = aircalMeans(3, 2);
CAL.O2_mL_L_aircal = aircalMeans(3, 1);
CAL.Temp_aircal = aircalMeans(1, 2);

CAL.Raw_Phase_std = aircalStds(3, 3);
CAL.Raw_Thermistor_std = aircalStds(3, 2);
CAL.O2_mL_L_std = aircalStds(3, 1);
CAL.Temp_aircal_std = aircalStds(1, 2);


CAL.RelHumidityEditField.Value = 1;
CAL.AtmpressureEditField.Value = 1014;
CAL.TempEditField.Value = 20.1;
%%
CAL.pH2O_100_hum = exp(52.57 - (6690.9 / (CAL.Temp_aircal + 273.15)) - 4.681 * log(CAL.Temp_aircal + 273.15));
CAL.pH2O_real = CAL.pH2O_100_hum * CAL.RelHumidityEditField.Value;
CAL.T_O2soly = log((298.15 - CAL.Temp_aircal) / (273.15 + CAL.Temp_aircal));
CAL.O2_soly_umol_kg = exp(5.80818 + 3.20684 * CAL.T_O2soly + 4.1189 * CAL.T_O2soly^2 + 4.93845 * CAL.T_O2soly^3 + 1.01567 * CAL.T_O2soly^4 + 1.41575 * CAL.T_O2soly^5);
CAL.Density_waterkg_L = (999.84847 + 0.06337563 * CAL.Temp_aircal - 0.008523829 * CAL.Temp_aircal^2 + 0.00006943248 * CAL.Temp_aircal^3 - 0.0000003821216 * CAL.Temp_aircal^4) / 1000;
CAL.O2_soly_uM = CAL.O2_soly_umol_kg * CAL.Density_waterkg_L;
CAL.O2raw_air_uM = CAL.O2_mL_L_aircal * 44.6596; 
CAL.Gain = (CAL.O2_soly_uM * (CAL.AtmpressureEditField.Value - CAL.pH2O_real)) / (1013.25 - CAL.pH2O_100_hum) / CAL.O2raw_air_uM;

%%
CAL = struct2table(CAL);
writetable(CAL,'CAL','FileType','spreadsheet')