% getCPFdata.m
% Extract science data from engineering file

% 10/10/16 fixed cp data parsing problem - only pars lines w 3 csv's

% fname    = '2016-1-29T18-22-32';
% fname    = '2015-8-11T15-31-52'; % newer format
%fname    = '2016-4-8T16-3-5';
fname    = '2016-10-6T16-26-47';
%fname    = '2015-9-11T17-10-19'; % older format
data_dir = '\\atlas\ProjectLibrary\CPF\Data\CPF-P2\Paragon Tests\2016Oct06\';

fig_switch =[1 0 0 0 0];

% ***********************************************************************
% FORMAT STRINGS


       
FLBB_format = '99/99/99,99:99:99,%*f,%f,%*f,%f,%f';
gps_format  = '%*s %*s %*s %f %s %f %s';

% ***********************************************************************
% ***********************************************************************
% BUILD FORMAT STRING BASED ON VARIABLE LIST FROM HEADER DESIRED VARIABLES
% DETERMINED FROM "keep_vars"
% HEADER
% $	   Date	  Time	State	Comment	Pressure	Bel Pos	Batt Volt
% CTD Temp	CTD Sal	 LPS331 Press	LPS311 Temp	  SHT21 Humid	 SHT21 Temp	
% FLBB	Topt	Tphase	Rphase	NO3	pH	LPF PVel	
% Diff PVel 	LPF PAcc	Press Avg	Press Var	PID SP	PID PV	
% PID Y1	PID Y2	PID I	PID v	PID u	PID uCPS	
% Batt Amps	 12V Amps	3.3V Amps	Debug	GPS	Irid Modem	
% ***********************************************************************
fid       = fopen([data_dir,fname]);

tline     = fgetl(fid); % 1st line of data is always header line
hdr_vars  = textscan(tline,'%s','Delimiter','\t');
hdr_vars  = hdr_vars{1,1};

keep_vars = {'Date', 'Time', 'State', 'Comment', 'Pressure', ...
             'CTD Sal', 'FLBB', 'Topt','NO3','pH','GPS'};         

ct = 0;
txt_format = '';
for i = 1:length(hdr_vars) % Finish building format string
    hdr_test = strcmp(hdr_vars{i,1},keep_vars);
    if sum(hdr_test) == 0
        txt_format = [txt_format,' %*s']; % Un-wanted variables
    else
        txt_format = [txt_format,' %s'];  % desired variables
        ct = ct+1;
        saved_vars(ct,1) = hdr_vars(i,1); %Header for discovered vars
    end
end

% *************************************************************************
% GET INDICES FROM  "saved_vars" IF THEY EXIST
state_col = find(strcmp('State',saved_vars));
comnt_col = find(strcmp('Comment',saved_vars));
pres_col  = find(strcmp('Pressure',saved_vars));
CTD_col   = find(strcmp('CTD Sal',saved_vars));
FLBB_col  = find(strcmp('FLBB',saved_vars));
Topt_col  = find(strcmp('Topt',saved_vars));
NO3_col   = find(strcmp('NO3',saved_vars));
pH_col    = find(strcmp('pH',saved_vars));    
GPS_col   = find(strcmp('GPS',saved_vars)); 

% *************************************************************************
% Bring text into cell arrays
txt_data = textscan(fid, txt_format, 'Delimiter','\t');
fclose(fid);
clear ct tline hdr_test i

% *************************************************************************
% CHECK IF OLDER OR NEWER VERSION OF CPF SOFTWARE
% *************************************************************************
states = {'Surface Ops', 'Ascend', 'Dump CP Data'};
CPF_version = 'new';
O2_format   = ['%*s%*s%*s%*s%f%*s%f%*s%f%*s%f%*s%f%*s%f%*s%*f', ...
               '%*s%*f%*s%*f%*s%*f'];

% state_ind = strcmp('Park',txt_data{1,state_col});
% 
% if sum(state_ind) > 0
%     CPF_version = 'new';
%     states = {'Surface Ops', 'Ascend', 'Dump CP Data'};
%     O2_format   = ['%*s%*s%*s%*s%f%*s%f%*s%f%*s%f%*s%f%*s%f%*s%*f', ...
%                    '%*s%*f%*s%*f%*s%*f'];
%     %O2_format   = '%*s %*s %*s %*s %f %*s %f %*s %f'; 
% else
%     CPF_version = 'old';
%     states = {'SurfaceOps      ', 'Ascend          ', 'Dumping CP Data  '};
%     O2_format   = '%*s %*s %*s %f %*s %f %*s %f'; 
%     disp('Older data format - still need to finish the code for it')
%     return
%     
% end


% *************************************************************************
% *************************************************************************
%                         NOW CONVERT TEXT TO NUMBERS
% *************************************************************************
% *************************************************************************
[rows,~] = size(txt_data{1,2}); % Will be sdn and decimal seconds
time_stamp  = ones(rows,2)*NaN; %predimmension

% *************************************************************************
% PARSE TIME STAMP
% *************************************************************************
time_stamp(:,1) =  datenum(txt_data{1,1}, 'dd-mmm-yyyy'); % Convert day to SDN
for i = 1:rows
    t_tmp = sscanf(txt_data{1,2}{i},'%f:%f:%f');
    time_stamp(i,2) = t_tmp(1)*3600 + t_tmp(2)*60 + t_tmp(3) ; %seconds in day
end
time_stamp(:,1) = time_stamp(:,1) + time_stamp(:,2) ./ 24 /3600;

% *************************************************************************
% NEXT GET GPS FIX
% gps data returns cell { 'lat' 'N or S' 'lon' 'E or W'}
% *************************************************************************
fix_ind = find(strncmp('UBX',txt_data{1,GPS_col},3)); % get gps string
if ~isempty(fix_ind)
    lat = ones(length(fix_ind))*NaN; % predim, the size should usually be 1
    lon = lat;
    GPS = ones(length(fix_ind),3) * NaN;
    for i = 1:length(lat) % usually only 1 fix per cycle file
        gps = textscan(txt_data{1,GPS_col}{fix_ind(i)}, gps_format, 1,...
            'Delimiter',','); % Only use format string once
        lat(i) = fix(gps{1}/100) + rem(gps{1}/100,fix(gps{1}/100))*100/60;
        lon(i) = fix(gps{3}/100) + rem(gps{3}/100,fix(gps{3}/100))*100/60;
        if strcmp(gps{2},'S')
            lat = lat * -1;
        end
        if strcmp(gps{4},'W')
            lon = -lon + 360;
        end
        GPS(i,:) = [time_stamp(fix_ind,1), lon, lat]; %sdn lon lat
    end
else
    disp('NO GPS FIX FOUND IN DATA FILE!')
end
GPS = mean(GPS,1);
clear fix_ind lat lon i gps_format fid

% *************************************************************************
% OK NOW STEP THROGH LISTED STATES AND LOOK FOR LISTED VARIABLES
% Build matrices of data for each state
% *************************************************************************

for i = 1:length(states)
    state_ind = strcmp(states{i},txt_data{1,state_col});
    tmp       = cellfun(@(x) x(state_ind), txt_data, 'UniformOutput',false);
    
    d = ones(sum(state_ind),18)*NaN; %predim matrix for data
    d(:,1:2)  = time_stamp(state_ind,:); % add time
    
    % ********************************************************************
    % PARSE CTD STRINGS
    % STRING & COMMENT DIFF FOR DIFF STATES KIND OF A PAIN
    SBE_ct = 1; % 1 or 3 p vs pst
    if strcmp('Surface Ops', states{i}) || ....
        strcmp('SurfaceOps      ', states{i})
    
        SBE_ind = find(strncmp('SBE41',tmp{1,comnt_col},5));
        sbe_format = '%f';
    elseif strncmp('Dump', states{i},4)
        SBE_ind = find(strncmp('Unparsed',tmp{1,comnt_col},8));
        sbe_format = '%f,%f,%f';
        SBE_ct = 3;
    else % ALL OTHER STATES
        SBE_ind = find(strncmp('SBE41',tmp{1,comnt_col},5));
        sbe_format = '%f,%f,%f'; 
        SBE_ct = 3;
    end
    
    if ~isempty(SBE_ind)
        for j = 1:length(SBE_ind) % PUT NUMBERS IN MATRIX
            str = tmp{1,CTD_col}{SBE_ind(j)};
            
            if strncmp('Dump', states{i},4)
                if size(regexp(str,','),2) == 3
                    d(SBE_ind(j),4:4+SBE_ct-1) = sscanf(str,sbe_format);
                end
            elseif ~isempty(sscanf(str,'%f')) % test for numeric string
                d(SBE_ind(j),4:4+SBE_ct-1) = sscanf(str,sbe_format);
            end
            
            str = tmp{1,pres_col}{SBE_ind(j)}; 
            if ~isempty(sscanf(str,'%f'))
                d(SBE_ind(j),3) = sscanf(str,'%f');
            end
        end
    end
    % ********************************************************************
    % PARSE O2 STRINGS - 3 types of strings!!!
    O2_ind1 = find(strncmp('MEAS',tmp{1,Topt_col},4)); % New ver
    O2_ind2 = find(strncmp('%!M',tmp{1,Topt_col},3));  % new ver
    O2_ind3 = find(strncmp('%#',tmp{1,Topt_col},3));  % new ver
    O2_ind = [O2_ind1; O2_ind2; O2_ind3];
    if ~isempty(O2_ind)   
        for j = 1:length(O2_ind) % PUT NUMBERS IN MATRIX
            d(O2_ind(j),7:12) = cell2mat(textscan(tmp{1,Topt_col}{O2_ind(j)}, ...
                        O2_format,1,'Delimiter',',', 'CollectOutput',1));
            d(O2_ind(j),3) = sscanf(tmp{1,pres_col}{O2_ind(j)},'%f');
        end
    end
    
    % ********************************************************************
    % PARSE FLBB STRINGS
    FLBB_ind = find(strncmp('99/',tmp{1,FLBB_col},3));
    if ~isempty(FLBB_ind)   
        for j = 1:length(FLBB_ind) % PUT NUMBERS IN MATRIX
            d(FLBB_ind(j),13:15) = (sscanf(tmp{1,FLBB_col}{FLBB_ind(j)}, ...
                                  FLBB_format))';
            d(FLBB_ind(j),3) = sscanf(tmp{1,pres_col}{FLBB_ind(j)},'%f');
        end
    end
    
    % *******************************************************************
    
    %             *** ADD NO3 & pH LOOPS WHEN NEEDED ***
    
    % *******************************************************************
   
%    SORT DATA BY TIME BEFORE ANY TABLE LOOKUP TYPE OPPERATION 
%    DATA IN TXT FILE IS ASYNC WITH TIME DUE TO ALL THE ASYN TASKS
     [~,IX] = sort(d(:,1)); % Sort data by time due to async dataflow
     d = d(IX,:); % RAW DATA SORTED BY TIMESTAMP

     % *******************************************************************
     % EXTRACTION COMPLETED NOW CLEAN UP AND CONVERT #'S BY STATE
     % *******************************************************************
     
     % *******************************************************************
     if strcmp('Ascend', strtrim(states{i}))
         % Match data - closet neighboor    
         % MERGE ASCENDING DATA - 'Sampling Instruments' COMMENT
         % FOR SAMPLE COUNT, CTD FOR TIME STAMP
         % FOLLOW SIMPLE RULE:
         %  1)  GET 1st SAMPLE IN TIME FOR EACH INSTRUMENT AFTER
         %     "Sampling Instruments" COMMENT
         
         % GET ROW INDICES FOR START OF INSTUMENT SAMPLING, AS WELL AS
         % INDICIES FOR INDIVIDUAL RETURNED CTD, O2 & FLBB SAMPLES
             samp_ind = find(strcmp('Sampling Instruments',tmp{1,comnt_col}));
             t_ctd    = find(~isnan(d(:,4))); % Get CTD indices
             t_o2     = find(~isnan(d(:,7))); % Get O2 indices
             t_flbb   = find(~isnan(d(:,13))); % Get FLBB indices
         
         lr_data = ones(size(samp_ind,1),size(d,2))*NaN; % predim merged
         % CTD WILL BE MASTER TIME STAMP
         %if strcmp('new',CPF_version)
         for j = 1: length(samp_ind)
             % FIND INDICES IN d FROM INDIVIDUAL SAMPLE INDICE VECTORS
             ind_ctd = find(t_ctd >= samp_ind(j),1,'first'); % 1st CTD after samp
             ind_o2 = find(t_o2 >= samp_ind(j),1,'first'); % 1st o2 after samp
             ind_flbb = find(t_flbb >= samp_ind(j),2,'first'); %first 2 FLBB after samp
             
             % MERGE DATA
             lr_data(j,1:6) = d(t_ctd(ind_ctd),1:6); % time and ctd
             lr_data(j,7:12) = d(t_o2(ind_o2),7:12); % O2 added
             lr_data(j,13:15) = mean(d(t_flbb(ind_flbb),13:15),1); % Added avged 2 FLBB
         end
%          else % OLD VERSION NO "SAMPLING INSTRUMENT" COMMENT
% %              for j = 1: length(t_o2)
%          end   
         
         % CORRECT O2 FOR SALINITY & DEPTH
         % Benson & Krause (cm^3/dm^3) *** USE THESE COEFFICIENTS!!!! ***
         pA = [3.88767 -0.256847 4.94457 4.05010 3.22014 2.00907]; % Temperature Coeff
         pB = [-8.17083e-3 -1.03410e-2 -7.37614e-3 -6.24523e-3]; % Salinity coff
         Co = -4.88682e-7;
         
         O2_volume =  22.3916;
         
         % SALINITY CORRECTION
         Ts     = log((298.15 - lr_data(:,5))./(273.15 + lr_data(:,5)));
         part1  = polyval(pB,Ts);
         S_corr = lr_data(:,5).*part1 + lr_data(:,5).^2 * Co;
         L      = polyval(pA,Ts) + S_corr;
         O2sol  = (1000/O2_volume) * exp(L); 
         
         O2c = lr_data(:,7) .* exp(S_corr); % Aanderaa S correction to O2 concentration
         
         % DEPTH COMPENSATION - 3.2% is an emperical value from Uchida 2008
         O2c = O2c .* ((lr_data(:,4) * 0.032/1000) + 1);
     end
     
     % *******************************************************************
     if strncmp('Dump', states{i},4) 
         % NO MATCHING NEEDED - Just QUICK CLEAN UP OF DATA MATRIX
                           % 
          tNaN = sum(~isnan(d(:,4:end)),2) == 0; % rows with no CTD data
          cp_data = d(~tNaN,:);
     end   
end



% ************************************************************************
% ************************************************************************
%                             MAKE SOME PLOTS
% ************************************************************************
% ************************************************************************

if fig_switch(1) == 1
    
figure(1)
set(gcf,'Units', 'Normalized');
set(gcf,'Position',[0.2, 0.2, 0.6, 0.7]);
set(gcf,'NumberTitle','Off');
set(gcf,'Name',' CPF RAW DATA PLOTS');

h_bnds  = [0.1 1.0];
v_bnds  = [0.1 1.0];
h_plots = 3; % horiz plots
v_plots = 2; % vert plots

pos = make_SubplotPos(h_bnds,v_bnds, h_plots, v_plots); % [L B W H]
pos(1:3,2) = pos(1:3,2) - 0.025;
pos(:,2) = pos(:,2) - 0.025;

a1(1) = axes('Position',pos(1,:)); 
plot(lr_data(:,5),lr_data(:,4),'bo')
hold on
plot(cp_data(:,5),cp_data(:,4),'k-','LineWidth',2)
hold off
set(gca,'Ydir','Reverse')
ylabel('Depth, m')
xlabel('Temperature, C')
%title(fname)
legend('low res', 'cp', 'Location','SouthEast')

a1(2) = axes('Position',pos(2,:)); 
plot(lr_data(:,6),lr_data(:,4),'bo')
hold on
plot(cp_data(:,6),cp_data(:,4),'k-','LineWidth',2)
hold off
set(gca,'Ydir','Reverse')
%ylabel('Depth, m')
xlabel('Salinity, psu')
%title(datestr(GPS(1),'mm/dd/yy'))

a1(3) = axes('Position',pos(3,:)); 
plot(lr_data(:,7),lr_data(:,4),'ro')
hold on
plot(O2c, lr_data(:,4),'bo')
plot(O2sol, lr_data(:,4),'ko')
hold off

set(gca,'Ydir','Reverse')
%ylabel('Depth, m')
xlabel('Oxygen, % sat')
%s1 =[num2str(GPS(2),'%3.3f'),' N & ',num2str(GPS(3),'%3.3f'), ' E'];
%title(s1)
set(gca, 'Units', 'normalized')
leg3 = legend('raw','S cor', ' O2sol');
ax_pos = get(gca,'Position');
leg_pos = [ax_pos(1)+ax_pos(3)/2, ax_pos(2)+ax_pos(4)+0.005, 0.1, 0.075];

set(leg3, 'Position',leg_pos)


a1(4) = axes('Position',pos(4,:)); 
plot(lr_data(:,13),lr_data(:,4),'bo')
set(gca,'Ydir','Reverse')
ylabel('Depth, m')
xlabel('FLBB channel 1')

a1(5) = axes('Position',pos(5,:)); 
plot(lr_data(:,14),lr_data(:,4),'bo')
set(gca,'Ydir','Reverse')
%ylabel('Depth, m')
xlabel('FLBB channel 2')

a1(6) = axes('Position',pos(6,:)); 
plot(lr_data(:,1),lr_data(:,4),'b-','LineWidth',2)
set(gca,'Ydir','Reverse')
datetick('x','HH:MM')
%ylabel('Depth, m')
xlabel('Time')

txt_ax = axes('units','normalized','Position',[0 0 1 1],...
              'Visible', 'off');
txt1_str = {['Coastal Profiling float file : ',fname], ...
          datestr(GPS(1),'mmm dd yyyy'), ...
          [num2str(GPS(2),'%3.3f'),' E & ',num2str(GPS(3),'%3.3f'), ' N']};

txt1 = text(0.5, 0.93, txt1_str, 'fontsize', 16,...
    'HorizontalAlignment', 'center');



end
 
 