% program runs data before strutton calibration and regresses data 

% (once) run C:\Documents and Settings\Lisa Taylor\My Documents\Data\COASOSU\chl_night.m
  % to get chl_night.mat (Strutton Night Chl)
% (once) run C:\Documents and Settings\Lisa Taylor\My Documents\Data\COASOSU\compile-noref.m
  % to get ps-compile.mat (Strutton data in addition to chl)
% (daily) run plot_apex_data.m to get apex_data_comparison.mat (Float 0005 profile data)
  
close all
clear all
warning('off')

%% APEX DATA, FLOAT 0005
% load 'C:\Documents and Settings\Lisa Taylor\My Documents\Data\APEX\0005\apex_data_comparison'

% apex_data_comparison.mat contains:
  % pressure_apex
  % temp_apex
  % sal_apex
  % oxygen_apex
  % blueRef_apex
  % blueSig_apex counts
  % LSSRef_apex
  % LSSSig_apex counts
  % theta_apex
  % sigma_apex
  % O2sat_apex
  % O2_apex
  % Chl_apex
  % B_LSS_apex
  % date_apex (y m d h m s)
  % tt_apex, which is datenum(2004,6,17)+data_apex(:,16)
  % t_apex, which is 168+data_apex(:,16);
  % B_LSS_b_apex=(LSSSig_apex-51.8)*0.00004637; % instrument calibration ONLY
  % cp_apex, which is B_LSS_apex*63;
  % cp_b_apex, which is B_LSS_b_apex*63;
  % chl_b_apex, which is (((blueSig_apex-25)*0.01079772) instrument-calibrated ONLY
  % POC_b_apex, which is 10000* the corrected B_LSS_b_apex;

%   allApex=[date_apex, pressure_apex, temp_apex, sal_apex, oxygen_apex, chl_b_apex, B_LSS_b_apex]; % create one array to work with

load 'C:\Documents and Settings\Lisa Taylor\My Documents\data\0005\APEX\apex_data'

for i=1:N % for all files (pos, dat)
    file_i=int2str(i); % change i to string and store as file_i
    n=length(file_i); % check length of file_i
        if n==1 %if 1, 2, 3...9
            name=['pos00',file_i]; % concatenates pos00 and i (pos001) and assigns to array name
            name2=['dat00',file_i]; % concatenates dat00 and i (dat001) and assigns to array name2
            name3=['date00',file_i];
        elseif n==2 %if 10, 11, 12...99 
             name=['pos0',file_i]; % concatenates pos0 and i (pos010) and assigns to array name
             name2=['dat0',file_i]; % concatenates dat0 and i (dat010) and assigns to array name2
             name3=['date0',file_i];
        else %if 100, 101, 102...999
             name=['pos',file_i]; % concatenates pos and i (pos100) and assigns to array name
             name2=['dat',file_i]; % concatenates dat and i (dat100) and assigns to array name2
             name3=['date',file_i];
        end           
    A=eval(name); % eval(pos001) - takes strings in pos001 (lat, lon, minutes) and changes to numbers, then assigns to array A
    latlong=A(1,1:2); 
    long=[latlong(:,1)-60]; % longitude
    lat=latlong(:,2); % latitude
    B=eval(name2); % data
    C=eval(name3); % date
    [mm nn]=size(B);

    % open each file (pos, dat and date) 000, 001, 000+ in turn
    if i==1 % if first time through (pos001, dat001 and date001)
        dat=[B,zeros(mm,1)]; % create array dat from data in dat file with last column (16) of zeros
        alldat=[dat,zeros(mm,1:2)]; % create array alldat with dat plus columns 17 and 18 of zeros
        alldat(:,17)=long(1,1); % replace column 17 with longitude
        alldat(:,18)=lat(1,1); % replace column 18 with latitude
         
    else %if not the first time around
        AA=[B,(i-1)*5*ones(mm,1)]; % create array AA from file and replace/increment last column (16) to get profile stages 
        dat=[dat;AA]; % append AA to dat
        CC=[AA,zeros(mm,1:2)]; % create array CC from AA, and add two columns of zeros
        CC(:,17)=long(1,1); % replace column 17 with longitude
        CC(:,18)=lat(1,1); % replace column 18 with latitude
        alldat=[alldat;CC]; % append CC to alldat array
    end
end

% remove negative values (no data collected)
I=find(alldat<-100);
alldat(I)=NaN;

%variables are:
pressure_apex=alldat(:,1); % depth (dB)
temp_apex=alldat(:,2); % temperature (degrees C)
sal_apex=alldat(:,3); % salinity (psu)
oxygen_apex=alldat(:,4);
blueRef_apex=alldat(:,5);

blueSig_apex=alldat(:,6);
LSSRef_apex=alldat(:,7);

LSSSig_apex=alldat(:,8);
theta_apex=alldat(:,9); % temperature (degrees C)
sigma_apex=alldat(:,10); % density-1000 (kg/m3)

O2sat_apex=alldat(:,11); % oxygen (% saturation)
O2_apex=alldat(:,12); % oxygen (um/kg)

Chl_apex=alldat(:,13); % chlorophyll (mg/m3, not calibrated)
B_LSS_apex=alldat(:,14); % bb (m-1, not calibrated)
B_LSS_b_apex=(LSSSig_apex-51.8)*0.00004637; % instrument calibration ONLY

date_apex=alldat(:,15); % date (datenum format)
% tt_apex=datenum(2004,6,17)+data_apex(:,16)
% t_apex=168+date_apex(:,16);

% cp_apex=B_LSS_apex*63;
% cp_b_apex=B_LSS_b_apex*63;
chl_b_apex=(((blueSig_apex-25)*0.01079772)); % instrument-calibrated ONLY
% POC_b_apex=10000*B_LSS_b_apex;
% profile=alldat(:,16); % profile (in days from start)
% long=alldat(:,17); % longitude (degrees)
% lat=alldat(:,18); % latitude (degrees)

allApex=[date_apex, pressure_apex, temp_apex, sal_apex, oxygen_apex, chl_b_apex, B_LSS_b_apex]; % create one array to work with

%clear unnecessary variables from workspace
clear dat0*;
clear dat1*;
clear date0*; 
clear date1*;
clear pos*; 

%% STRUTTON CHL NIGHT ONLY
load 'C:\Documents and Settings\Lisa Taylor\My Documents\Data\0005\COASOSU\chl_night'
% chl_night.mat contains:
  % date_night at chlNight(:,1) (needs correction, see below)
  % chl_2m_night at chlNight(:,2) 
  % chl_5m_night at chlNight(:,3)
  % chl_10m_night at chlNight(:,4)
  % chl_15m_night at chlNight(:,5)
  % chl_20m_night at chlNight(:,6)
time_strut=datenum(2004,1,1)+chlNight(:,1)-1; % Strutton date changed to 2004
% chlNightAll=[time_strut,chl_2m_night,chl_5m_night,chl_10m_night,chl_15m_night,chl_20m_night,];

%% STRUTTON DATA (other than chl night)
load 'C:\Documents and Settings\Lisa Taylor\My Documents\Data\0005\COASOSU\ps-compile'

% ps-compile.mat contains:
  % PsCompile (6977x38), containing:
  % 1 date_pscompile (all arrays 6977x1)
  % 2 exceldate_pscompile
  % 3 year_day_pscompile
  % 4 lon_pscompile
  % 5 lat_pscompile
  % 6 flr_2m_pscompile
  % 7 chl_2m_pscompile
  % 8 flr_10m_pscompile
  % 9 chl_10m_pscompile
  % 10 flr_15m_pscompile
  % 11 chl_15m_pscompile
  % 12 b100_2m_pscompile
  % 13 b125_2m_pscompile
  % 14 b150_2m_pscompile
  % 15 b_back_2m_pscompile
  % 16 b100_10m_pscompile
  % 17 b125_10m_pscompile
  % 18 b150_10m_pscompile
  % 19 b_back_10m_pscompile
  % 20 b100_15m_pscompile
  % 21 b125_15m_pscompile
  % 22 b150_15m_pscompile
  % 23 b_back_15m_pscompile
  % 24 par_1m_pscompile
  % 25 unknown_1m_pscompile
  % 26 par_5m_pscompile
  % 27 flr_5m_pscompile
  % 28 par_10m_pscompile
  % 29 c_10m_pscompile
  % 30 par_20m_pscompile
  % 31 flr_20m_pscompile
  % 32 par_35m_pscompile
  % 33 c_35m_pscompile
  % 34 poc_bb_2m_pscompile
  % 35 poc_bb_10m_pscompile
  % 36 poc_beamc_10m_pscompile
  % 37 poc_bb_15m_pscompile
  % 38 poc_beamc_35m_pscompile
  
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% CHLOROPHYLL PROGRAMS

% set apex depth range
I1=find(allApex(:,2)<11.1 & allApex(:,2)>4.9); % apex depths from 5 to 11m
I2=find(allApex(:,2)<17.1 & allApex(:,2)>11.9); % apex depths from 12 to 17m
I3=find(allApex(:,2)<22.1 & allApex(:,2)>17.9); % apex depths from 18 to 22m

% match strutton 24-hr night data with apex chlorophyll for 10, 15, and 20m
% get apex measurements taken between 5 and 11m
apex10m(:,1)=date_apex(I1);
apex10m(:,2)=pressure_apex(I1);
apex10m(:,3)=temp_apex(I1);
apex10m(:,4)=sal_apex(I1);
apex10m(:,5)=oxygen_apex(I1);
apex10m(:,6)=chl_b_apex(I1);
apex10m(:,7)=B_LSS_b_apex(I1);

% get apex measurements taken between 12 and 17m
apex15m(:,1)=date_apex(I2);
apex15m(:,2)=pressure_apex(I2);
apex15m(:,3)=temp_apex(I2);
apex15m(:,4)=sal_apex(I2);
apex15m(:,5)=oxygen_apex(I2);
apex15m(:,6)=chl_b_apex(I2);
apex15m(:,7)=B_LSS_b_apex(I2);

% get apex measurements taken between 18 and 22m
apex20m(:,1)=date_apex(I3);
apex20m(:,2)=pressure_apex(I3);
apex20m(:,3)=temp_apex(I3);
apex20m(:,4)=sal_apex(I3);
apex20m(:,5)=oxygen_apex(I3);
apex20m(:,6)=chl_b_apex(I3);
apex20m(:,7)=B_LSS_b_apex(I3);
 
% following code to obtain first measurement of the day not necessary for
% apex data since only 1 measurement taken every 5 days
% j=1;
% for i=1:155 %number of different profiles to date NEED TO EDIT THIS NEXT TIME NEW APEX DATA FILES ARE RUN
%     I4=find(datenum(datestr(apexSurface(:,1)))-datenum(datestr(apexSurface(j,1)))<1 & datenum(datestr(apexSurface(:,1)))-datenum(datestr(apexSurface(j,1)))>=0);
%     % compares apexSurface(1,1) to apexSurface(1,1)
%     % compares apexSurface(2,1) to apexSurface(1,1)
%     % compares apexSurface(3,1) to apexSurface(1,1) and so on...all <1 are found (all dates from same day)
%     % then stores first measurement (min) of the day in array firstApex:
%     firstApex10(i,1)=datenum(datestr(apexSurface(min(I4),1))); % date
%     firstApex10(i,2)=apexSurface(min(I4),2); % pressure
%     firstApex10(i,3)=apexSurface(min(I4),3); % temperature
%     firstApex10(i,4)=apexSurface(min(I4),4); % salinity
%     firstApex10(i,5)=apexSurface(min(I4),5); % oxygen
%     firstApex10(i,6)=apexSurface(min(I4),6); % chl_b
%     firstApex10(i,7)=apexSurface(min(I4),7); % LSS_b
%     j=max(I4)+1;
% end

%% match strutton data to within 24 hours of apex measurements by depth
%% 10m
strutbin_mean_chl_10m = nan(size(apex10m,1),size(chlNight,2)); % initialize meaned array (entire data set) to size you want
strutbin_std_chl_10m = nan(size(apex10m,1),size(chlNight,2)); % initialize std deviation arrray
hours12a=datenum(0,0,0,12,0,0); % set range of hours to include in meaned values, 12 hrs either side of apex data

for mn=1:size(apex10m,1)
   thisApexDateChl_10m=apex10m(mn,1); 
   idx10=find((time_strut>=(thisApexDateChl_10m-hours12a)) & (time_strut<=thisApexDateChl_10m+hours12a)); % get all strutton within 24hrs of apex data
   strutbin_mean_chl_10m(mn,1)=mean(chlNight(idx10,4)); % get mean of 24hr strutton at 10m
   strutbin_std_chl_10m(mn,1)=std(chlNight(idx10,4)); % get std dev of 24hr strutton at 10m
end

% put into one array for regression analysis
chlReg10m=[apex10m(:,1), apex10m(:,6), strutbin_mean_chl_10m(:,1)];

%% 15m
strutbin_mean_chl_15m = nan(size(apex15m,1),size(chlNight,2)); % initialize meaned array (entire data set) to size you want
strutbin_std_chl_15m = nan(size(apex15m,1),size(chlNight,2)); % initialize std deviation arrray

for mn=1:size(apex15m,1)
   thisApexDateChl_15m=apex15m(mn,1); 
   idx15=find((time_strut>=(thisApexDateChl_15m-hours12a)) & (time_strut<=thisApexDateChl_15m+hours12a)); % get all strutton within 24hrs of apex data
   strutbin_mean_chl_15m(mn,1)=mean(chlNight(idx15,5)); % get mean of 24hr strutton at 15m
   strutbin_std_chl_15m(mn,1)=std(chlNight(idx15,5)); % get std dev of 24hr strutton at 15m
end

% put into one array for regression analysis
chlReg15m=[apex15m(:,1), apex15m(:,6), strutbin_mean_chl_15m(:,1)];

%% 20m
strutbin_mean_chl_20m = nan(size(apex20m,1),size(chlNight,2)); % initialize meaned array (entire data set) to size you want
strutbin_std_chl_20m = nan(size(apex20m,1),size(chlNight,2)); % initialize std deviation arrray

for mn=1:size(apex20m,1)
   thisApexDateChl_20m=apex20m(mn,1); 
   idx20=find((time_strut>=(thisApexDateChl_20m-hours12a)) & (time_strut<=thisApexDateChl_20m+hours12a)); % get all strutton within 24hrs of apex data
   strutbin_mean_chl_20m(mn,1)=mean(chlNight(idx20,6)); % get mean of 24hr strutton at 20m
   strutbin_std_chl_20m(mn,1)=std(chlNight(idx20,6)); % get std dev of 24hr strutton at 20m
end

% put into one array for regression analysis
chlReg20m=[apex20m(:,1), apex20m(:,6), strutbin_mean_chl_20m(:,1)];

% chlorophyll variables to graph:
  % strutbin_mean_chl_10m and apex10m(:,6)
  % strutbin_mean_chl_15m and apex15m(:,6)
  % strutbin_mean_chl_20m and apex20m(:,6)

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% BACKSCATTER PROGRAMS
% match Strutton data to within 24 hours of Apex
%% 2m
strutbin_mean_b_2m = nan(size(apex10m,1),size(PsCompile,2)); % initialize meaned array (entire data set) to size you want
strutbin_std_b_2m = nan(size(apex10m,1),size(PsCompile,2)); % initialize std deviation array
hours12b=datenum(0,0,0,12,0,0); % set range of hours to include in meaned values, 12 hrs either side of apex data

for la=1:size(apex10m,1)
   thisApexDateB_2m=apex10m(la,1); 
   idx2=find((date_pscompile>=(thisApexDateB_2m-hours12b)) & (date_pscompile<=thisApexDateB_2m+hours12b)); % get all strutton within 24hrs of apex data
   strutbin_mean_b_2m(la,1)=mean(PsCompile(idx2,15)); % get mean of 24hr strutton at 2m
   strutbin_std_b_2m(la,1)=std(PsCompile(idx2,15)); % get std dev of 24hr strutton at 2m
end

% put into one array for regression analysis
bReg2m=[apex10m(:,1), apex10m(:,7), strutbin_mean_b_2m(:,1)];

%% 10m % strutton sensor may be off, try using 2m instead with apex 10m (see above)
strutbin_mean_b_10m = nan(size(apex10m,1),size(PsCompile,2)); % initialize meaned array (entire data set) to size you want
strutbin_std_b_10m = nan(size(apex10m,1),size(PsCompile,2)); % initialize std deviation arrray

for la=1:size(apex10m,1)
   thisApexDateB_10m=apex10m(la,1); 
   idx10=find((date_pscompile>=(thisApexDateB_10m-hours12b)) & (date_pscompile<=thisApexDateB_10m+hours12b)); % get all strutton within 24hrs of apex data
   strutbin_mean_b_10m(la,1)=mean(PsCompile(idx10,19)); % get mean of 24hr strutton at 10m
   strutbin_std_b_10m(la,1)=std(PsCompile(idx10,19)); % get std dev of 24hr strutton at 10m
end

% put into one array for regression analysis
bReg10m=[apex10m(:,1), apex10m(:,7), strutbin_mean_b_10m(:,1)];

%% 15m
strutbin_mean_b_15m = nan(size(apex15m,1),size(PsCompile,2)); % initialize meaned array (entire data set) to size you want
strutbin_std_b_15m = nan(size(apex15m,1),size(PsCompile,2)); % initialize std deviation arrray

for la=1:size(apex15m,1)
   thisApexDateB_15m=apex15m(la,1); 
   idx15=find((date_pscompile>=(thisApexDateB_15m-hours12b)) & (date_pscompile<=thisApexDateB_15m+hours12b)); % get all strutton within 24hrs of apex data
   strutbin_mean_b_15m(la,1)=mean(PsCompile(idx15,23)); % get mean of 24hr strutton at 10m
   strutbin_std_b_15m(la,1)=std(PsCompile(idx15,23)); % get std dev of 24hr strutton at 10m
end

% put into one array for regression analysis
bReg15m=[apex15m(:,1), apex15m(:,7), strutbin_mean_b_15m(:,1)];

% backscattering variables to graph:
% strutbin_mean_b_2m and apex10m(:,7)
% strutbin_mean_b_10m and apex10m(:,7) ? sensor may be off
% strutbin_mean_b_15m and apex15m(:,7)

% match Apex data to nearest Strutton data
b_back_2m_pscompile(isnan(b_back_2m_pscompile)) = []; % initialize strutton b2m array
b_back_10m_pscompile(isnan(b_back_10m_pscompile)) = []; % initialize strutton b10m array
b_back_15m_pscompile(isnan(b_back_15m_pscompile)) = []; % initialize strutton b15m array
date_pscompile(isnan(date_pscompile)) = []; % initialize strutton date array
B_LSS_b_apex(isnan(B_LSS_b_apex)) = []; % initialize apex b_lss_b array
date_apex(isnan(date_apex)) = []; % initialize apex date array
 
strut_back2m_AP = interp1(date_pscompile,b_back_2m_pscompile,date_apex, 'nearest');
strut_back10m_AP = interp1(date_pscompile,b_back_10m_pscompile,date_apex, 'nearest');
strut_back15m_AP = interp1(date_pscompile,b_back_15m_pscompile,date_apex, 'nearest');

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% CHLOROPHYLL PLOTS
%%%%%%%%%%%%%%%%% all Apex Chl_b vs date
% plot(allApex(:,1), allApex(:,6), 'b.') %all apex data
% datetick('x', 'mmmyy')

%%%%%%%%%%%%%%%%% apex vs date by depth
% figure
% subplot(2,2,1)
% plot(date_apex(I1),chl_b_apex(I1),'b.') % all apex data from 5 to 11m
% set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
% set(gca,'XTick',[732099 732129 732160 732191])
% datetick('x', 'mmmyy','keeplimits','keepticks')
% xlabel('Date')
% ylabel('Chlorophyll')
% title({'Chlorophyll','Float 0005 Chl b (5-11m)'})
% 
% subplot(2,2,2) % all apex data from 12 to 17m
% plot(date_apex(I2),chl_b_apex(I2),'b.') %all apex data
% set(gca,'XTick',[732099 732129 732160 732191])
% set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
% datetick('x', 'mmmyy','keeplimits','keepticks')
% xlabel('Date')
% ylabel('Chlorophyll')
% title({'Chlorophyll','Float 0005 Chl b (12-17m)'})
% 
% subplot(2,2,3) % all apex data from 18 to 22m
% plot(date_apex(I3),chl_b_apex(I3),'b.') %all apex data
% set(gca,'XTick',[732099 732129 732160 732191])
% set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
% datetick('x', 'mmmyy','keeplimits','keepticks')
% xlabel('Date')
% ylabel('Chlorophyll')
% title({'Chlorophyll','Float 0005 Chl b (18-22m)'})

%%%%%%%%%%%%%%%%% Strutton night vs date by depth
% figure
% subplot(3,2,1)
% plot(time_strut, chl_2m_night,'b.')
% set(gca,'XTick',[732099 732129 732160 732191])
% set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
% datetick('x', 'mmmyy','keeplimits','keepticks')
% xlabel('Date')
% ylabel('Chlorophyll')
% title({'Chlorophyll','COASOSU (chl 2m night)'})
% 
% subplot(3,2,2)
% plot(time_strut, chl_5m_night,'b.')
% set(gca,'XTick',[732099 732129 732160 732191])
% set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
% datetick('x', 'mmmyy','keeplimits','keepticks')
% xlabel('Date')
% ylabel('Chlorophyll')
% title({'Chlorophyll','COASOSU (chl 5m night)'})
% 
% subplot(3,2,3)
% plot(time_strut, chl_10m_night,'b.')
% set(gca,'XTick',[732099 732129 732160 732191])
% set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
% datetick('x', 'mmmyy','keeplimits','keepticks')
% xlabel('Date')
% ylabel('Chlorophyll')
% title({'Chlorophyll','COASOSU (chl 10m night)'})
% 
% subplot(3,2,4)
% plot(time_strut, chl_15m_night,'b.')
% set(gca,'XTick',[732099 732129 732160 732191])
% set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
% datetick('x', 'mmmyy','keeplimits','keepticks')
% xlabel('Date')
% ylabel('Chlorophyll')
% title({'Chlorophyll','COASOSU (chl 15m night)'})
% 
% subplot(3,2,5)
% plot(time_strut, chl_20m_night,'b.')
% set(gca,'XTick',[732099 732129 732160 732191])
% set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
% datetick('x', 'mmmyy','keeplimits','keepticks')
% xlabel('Date')
% ylabel('Chlorophyll')
% title({'Chlorophyll','COASOSU (chl 20m night)'})

%%%%%%%%%%%%%%%%% Strutton and Apex vs date
figure(1)
subplot(2,2,1) % all apex data from 5 to 11m and strutton data at 10m 
plot(date_apex(I1),chl_b_apex(I1),'b.')
hold on
plot(time_strut,chl_10m_night,'rx')
set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
set(gca,'XTick',[732099 732129 732160 732191])
datetick('x', 'mmmyy','keeplimits','keepticks')
title({'ASIS Buoy (10m)','Float 0005 (5-11m)'})
hold off

subplot(2,2,2) % all apex data from 12 to 17m and strutton data at 15m
plot(date_apex(I2),chl_b_apex(I2),'b.')
hold on
plot(time_strut,chl_15m_night,'rx')
set(gca,'XTick',[732099 732129 732160 732191])
set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
datetick('x', 'mmmyy','keeplimits','keepticks')
title({'ASIS Buoy (15m)','Float 0005 (12-17m)'})
hold off

subplot(2,2,3) % all apex data from 18 to 22m and strutton data at 20m
plot(date_apex(I3),chl_b_apex(I3),'b.')
hold on
plot(time_strut,chl_20m_night,'rx')
set(gca,'XTick',[732099 732129 732160 732191])
set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
datetick('x', 'mmmyy','keeplimits','keepticks')
title({'ASIS Buoy (20m)','Float 0005 (18-22m)'})
hold off
% legend('Float 0005','COASOSU')

%%%%%%%%%%%%%%%%% Strutton and Apex vs date -2 y axes
figure(2)
% subplot(2,2,1) % all apex data from 5 to 11m and strutton data at 10m 
% plot(date_apex(I1),chl_b_apex(I1),'b.')
% hold on
% plot(time_strut,chl_10m_night,'rx')
[AX,H1,H2] = plotyy(date_apex(I1),chl_b_apex(I1),time_strut,chl_10m_night,'plot');
set(get(AX(1),'Ylabel'),'String','Counts')
set(get(AX(2),'Ylabel'),'String','Chlorophyll')
set(AX(1),'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
set(AX(2),'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
set(AX(1),'XTick',[732099 732129 732160 732191])
set(AX(2),'XTick',[],'YColor','red')
datetick(AX(1),'x', 'mmmyy','keeplimits','keepticks')
set(H1,'LineStyle','none','Marker','.')
set(H2,'LineStyle','none','Marker','x','MarkerEdgeColor','red')
title({'ASIS Buoy (10m)','Float 0005 (5-11m)'})

figure(3)
% subplot(2,2,2) % all apex data from 12 to 17m and strutton data at 15m
plot(date_apex(I2),chl_b_apex(I2),'b.')
[AX,H1,H2] = plotyy(date_apex(I2),chl_b_apex(I2),time_strut,chl_15m_night,'plot');
set(get(AX(1),'Ylabel'),'String','Counts')
set(get(AX(2),'Ylabel'),'String','Chlorophyll')
set(AX(1),'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
set(AX(2),'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
set(AX(1),'XTick',[732099 732129 732160 732191])
set(AX(2),'XTick',[],'YColor','red')
datetick(AX(1),'x', 'mmmyy','keeplimits','keepticks')
set(H1,'LineStyle','none','Marker','.')
set(H2,'LineStyle','none','Marker','x','MarkerEdgeColor','red')
title({'ASIS Buoy (15m)','Float 0005 (12-17m)'})

figure(4)
% subplot(2,2,3) % all apex data from 18 to 22m and strutton data at 20m
[AX,H1,H2] = plotyy(date_apex(I3),chl_b_apex(I3),time_strut,chl_20m_night,'plot');
set(get(AX(1),'Ylabel'),'String','Counts')
set(get(AX(2),'Ylabel'),'String','Chlorophyll')
set(AX(1),'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
set(AX(2),'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
set(AX(1),'XTick',[732099 732129 732160 732191])
set(AX(2),'XTick',[],'YColor','red')
datetick(AX(1),'x', 'mmmyy','keeplimits','keepticks')
set(H1,'LineStyle','none','Marker','.')
set(H2,'LineStyle','none','Marker','x','MarkerEdgeColor','red')
title({'ASIS Buoy (20m)','Float 0005 (18-22m)'})
% legend('Float 0005','COASOSU')

%%%%%%%%%%%%%%%% 24 HRS Strutton chlorophyll and Apex data vs date
% figure
% subplot(2,2,1)
% plot(apex10m(:,1),(apex10m(:,6)/2.5), 'b.') 
% hold on
% plot(apex10m(:,1),strutbin_mean_chl_10m, 'rx')
% set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
% set(gca,'XTick',[732099 732129 732160 732191])
% datetick('x', 'mmmyy','keeplimits','keepticks')
% xlabel('Date')
% ylabel('Chlorophyll')
% title({'Chlorophyll','24-Hr Mean COASOSU (chl 10m night)','Apex 5-11m'})
% legend('Float 0005','COASOSU')
% hold off
% 
% subplot(2,2,2)
% plot(apex15m(:,1),(apex15m(:,6)/2.5), 'b.') 
% hold on
% plot(apex15m(:,1),strutbin_mean_chl_15m, 'rx')
% set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
% set(gca,'XTick',[732099 732129 732160 732191])
% datetick('x', 'mmmyy','keeplimits','keepticks')
% xlabel('Date')
% ylabel('Chlorophyll')
% title({'Chlorophyll','24-Hr Mean COASOSU (chl 15m night)','Apex 12-17m'})
% legend('Float 0005','COASOSU')
% hold off
% 
% subplot(2,2,3)
% plot(apex20m(:,1),(apex20m(:,6)/2.5), 'b.') 
% hold on
% plot(apex20m(:,1),strutbin_mean_chl_20m, 'rx')
% set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
% set(gca,'XTick',[732099 732129 732160 732191])
% datetick('x', 'mmmyy','keeplimits','keepticks')
% xlabel('Date')
% ylabel('Chlorophyll')
% title({'Chlorophyll','24-Hr Mean COASOSU (chl 20m night)','Apex 18-22m'})
% legend('Float 0005','COASOSU')
% hold off

%%%%%%%%%%%%%%%%% 24 HRS Strutton chlorophyll and Apex correlations
figure(5) %% 6
subplot(2,2,1)
plot(strutbin_mean_chl_10m,apex10m(:,6), 'r.')
xlabel('COASOSU chl 10m night')
ylabel('Float 0005 chl b 5-11m')
title({'Mean COASOSU Data Within','24 Hours of Float 0005 chl b'})

subplot(2,2,2)
plot(strutbin_mean_chl_15m,apex15m(:,6), 'b*')
xlabel('COASOSU chl 15m night')
ylabel('Float 0005 chl b 12-17m')
title({'Mean COASOSU Data Within','24 Hours of Float 0005 chl b'})

subplot(2,2,3)
plot(strutbin_mean_chl_20m,apex20m(:,6), 'gx')
xlabel('COASOSU chl 20m night')
ylabel('Float 0005 chl b 18-22m')
title({'Mean COASOSU Data Within','24 Hours of Float 0005 chl b'})

subplot(2,2,4)
plot(strutbin_mean_chl_10m,apex10m(:,6), 'r.')
hold on
plot(strutbin_mean_chl_15m,apex15m(:,6), 'b*')
plot(strutbin_mean_chl_20m,apex20m(:,6), 'gx')
plot([0 10],[0 30])
xlabel('COASOSU chl night')
ylabel('Float 0005')
title({'Mean COASOSU Data Within','24 Hours of Float 0005 chl b'})
% legend('10m','15m','20m')
hold off

corrcoef(strutbin_mean_chl_10m(1:13,1),apex10m(1:13,6)); % 0.9222
corrcoef(strutbin_mean_chl_15m(1:13,1),apex15m(1:13,6)); % 0.6660
corrcoef(strutbin_mean_chl_20m(1:13,1),apex20m(1:13,6)); % 0.3764

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% BACKSCATTER PLOTS
%%%%%%%%%%%%%%%%% all Apex B_LSS_b vs date
% plot(allApex(:,1), allApex(:,7), 'b.') %all apex data
% datetick('x', 'mmmyy')

%%%%%%%%%%%%%%%%% Strutton vs date by depth
% figure %% 7
% subplot(2,2,1)
% plot(date_pscompile, b_back_2m_pscompile,'r.')
% set(gca,'XTick',[732099 732129 732160 732191])
% set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
% datetick('x', 'mmmyy','keeplimits','keepticks')
% xlabel('Date')
% ylabel('Backscatter')
% title({'Backscatter','COASOSU (b back 2m)'})
% 
% subplot(2,2,2)
% plot(date_pscompile, b_back_10m_pscompile,'r.')
% set(gca,'XTick',[732099 732129 732160 732191])
% set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
% datetick('x', 'mmmyy','keeplimits','keepticks')
% xlabel('Date')
% ylabel('Backscatter')
% title({'Backscatter','COASOSU (b back 10m)'})
% 
% subplot(2,2,3)
% plot(date_pscompile, b_back_15m_pscompile,'r.')
% set(gca,'XTick',[732099 732129 732160 732191])
% set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
% datetick('x', 'mmmyy','keeplimits','keepticks')
% xlabel('Date')
% ylabel('Backscatter')
% title({'Backscatter','COASOSU (b back 15m)'})

%%%%%%%%%%%%%%%%% Strutton and Apex vs date by depth
figure(6) %% 8
subplot(2,2,1)% all apex data from 5 to 11m and strutton data at 2m
plot(date_apex(I1),(B_LSS_b_apex(I1)/2.5),'b.')
hold on
plot(date_pscompile,b_back_2m_pscompile,'rx')
set(gca,'XTick',[732099 732129 732160 732191])
set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
datetick('x', 'mmmyy','keeplimits')
% xlabel('Date')
ylabel('Backscatter')
title({'ASIS Buoy (2m)','Float 0005 (5-11m)'})
% legend('Float 0005','Strutton')
hold off

% figure(7)
subplot(2,2,2)% all apex data from 5 to 11m and strutton data at 10m
plot(date_apex(I1),(B_LSS_b_apex(I1)/2.5),'b.')
hold on
plot(date_pscompile,b_back_10m_pscompile,'rx')
set(gca,'XTick',[732099 732129 732160 732191])
set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
datetick('x', 'mmmyy','keeplimits')
% xlabel('Date')
ylabel('Backscatter')
title({'ASIS Buoy (10m)','Float 0005 (5-11m)'})
% legend('Float 0005','Strutton')
hold off

% figure(8)
subplot(2,2,3)% all apex data from 12 to 17m and strutton data at 15m
plot(date_apex(I2),(B_LSS_b_apex(I2)/2.5),'b.')
hold on
plot(date_pscompile,b_back_15m_pscompile,'rx')
set(gca,'XTick',[732099 732129 732160 732191])
set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
datetick('x', 'mmmyy','keeplimits')
% xlabel('Date')
ylabel('Backscatter')
title({'ASIS Buoy (15m)','Float 0005 (12-17m)'})
legend('Float 0005','ASIS Buoy')
hold off

%%%%%%%%%%%%%%%% 24 HRS Strutton backscatter and Apex data vs date
% figure %% 9
% subplot(2,2,1)
% plot(apex10m(:,1),(apex10m(:,7)/2.5), 'b.') 
% hold on
% plot(apex10m(:,1),strutbin_mean_b_2m, 'rx')
% set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
% set(gca,'XTick',[732099 732129 732160 732191])
% datetick('x', 'mmmyy','keeplimits','keepticks')
% xlabel('Date')
% ylabel('Backscatter')
% title({'Backscatter','24-Hr Mean COASOSU (b back 2m)','Float 0005 (B LSS b, 5-11m)'})
% legend('Float 0005','COASOSU')
% hold off
% 
% subplot(2,2,2)
% plot(apex10m(:,1),(apex10m(:,7)/2.5), 'b.') 
% hold on
% plot(apex10m(:,1),strutbin_mean_b_10m, 'rx')
% set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
% set(gca,'XTick',[732099 732129 732160 732191])
% datetick('x', 'mmmyy','keeplimits','keepticks')
% xlabel('Date')
% ylabel('Backscatter')
% title({'Backscatter','24-Hr Mean COASOSU (b back 10m)','Float 0005 (B LSS b, 5-11m)'})
% legend('Float 0005','COASOSU')
% hold off
% 
% subplot(2,2,3)
% plot(apex15m(:,1),(apex15m(:,7)/2.5), 'b.') 
% hold on
% plot(apex15m(:,1),strutbin_mean_b_15m, 'rx')
% set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
% set(gca,'XTick',[732099 732129 732160 732191])
% datetick('x', 'mmmyy','keeplimits','keepticks')
% xlabel('Date')
% ylabel('Backscatter')
% title({'Backscatter','24-Hr Mean COASOSU (b back 15m)','Float 0005 (B LSS b, 12-17m)'})
% legend('Float 0005','COASOSU')
% hold off

%%%%%%%%%%%%%%%%% 24 HRS Strutton backscatter and Apex correlations
figure(9) %% 10
subplot(2,2,1)
plot(strutbin_mean_b_2m,apex10m(:,7), 'b.')
xlabel('COASOSU b back 2m')
ylabel('Float 0005 B LSS b (5-11m)')
title({'Mean COASOSU Data Within','24 Hours of Float 0005'})

subplot(2,2,2)
plot(strutbin_mean_b_10m,apex10m(:,7), 'rx')
xlabel('COASOSU b back 10m')
ylabel('Float 0005 B LSS b (5-11m)')
title({'Mean COASOSU Data Within','24 Hours of Float 0005'})

subplot(2,2,3)
plot(strutbin_mean_b_15m,apex15m(:,7), 'g*')
xlabel('COASOSU b back 15m')
ylabel('Float 0005 B LSS b (12-17m)')
title({'Mean COASOSU Data Within','24 Hours of Float 0005'})

subplot(2,2,4)
plot(strutbin_mean_b_2m,apex10m(:,7), 'b.')
hold on
plot(strutbin_mean_b_10m,apex10m(:,7), 'rx')
plot(strutbin_mean_b_15m,apex15m(:,7), 'g*')
plot([0 0.012],[0 0.04])
xlabel('COASOSU b back')
ylabel('Float 0005 B LSS b')
title({'Mean COASOSU Data Within','24 Hours of Float 0005'})
% legend('10m','15m')
hold off

corrcoef(strutbin_mean_b_2m(1:14),apex10m(1:14,7)) % 0.8720 
corrcoef(strutbin_mean_b_10m(1:14),apex10m(1:14,7)) % 0.8727 
corrcoef(strutbin_mean_b_15m(1:14),apex15m(1:14,7)) % 0.5540 


%%%%%%%%%%%%%%%%% Regression analysis (before Strutton correction)

%% bb
%%%%%%%%%% bb, strutton 2m and apex 10m
[poly_bb1] = polyfit(strutbin_mean_b_2m(~isnan(strutbin_mean_b_2m)), ...
                    apex10m(~isnan(strutbin_mean_b_2m),7), 1);
% for bb values after instrument correction: 3.8591 and 0.0055

[poly_bb2] = polyfit(apex10m(~isnan(strutbin_mean_b_2m),7), ...
                     strutbin_mean_b_2m(~isnan(strutbin_mean_b_2m)),1);
% for bb values after instrument correction: 0.1970 and -0.0004

slope_bb=realsqrt(3.8591*(1/0.1970)); % 4.4260

%%%%%%%%%% bb, strutton 10m and apex 10m
[poly_bb1] = polyfit(strutbin_mean_b_10m(~isnan(strutbin_mean_b_10m)), ...
                    apex10m(~isnan(strutbin_mean_b_10m),7), 1);
% for bb values after instrument correction: 3.0624 and -0.0017

[poly_bb2] = polyfit(apex10m(~isnan(strutbin_mean_b_10m),7), ...
                     strutbin_mean_b_10m(~isnan(strutbin_mean_b_10m)),1);
% for bb values after instrument correction: 0.2487 and 0.0018

slope_bb=realsqrt(3.0624*(1/0.2487)); % 3.5091

%%%%%%%%%% all bb, strutton 2m, 10m, and 15m and apex 10m (to match 2m), 10m (to match 10m), and 15m
bb_strut=[strutbin_mean_b_2m(~isnan(strutbin_mean_b_2m));strutbin_mean_b_10m(~isnan(strutbin_mean_b_10m));strutbin_mean_b_15m(~isnan(strutbin_mean_b_15m))];
bb_apex=[apex10m(~isnan(strutbin_mean_b_2m),7);apex10m(~isnan(strutbin_mean_b_10m),7);apex15m(~isnan(strutbin_mean_b_15m),7)];
[poly_bb1]=polyfit(bb_strut,bb_apex,1)
[poly_bb2]=polyfit(bb_apex,bb_strut,1)

figure;plot(bb_strut,bb_apex,'.')
slope_bb=sqrt(poly_bb1(1)*(1/poly_bb2(1))); % 3.1171
corrcoef(bb_strut,bb_apex); % correlation coefficient: 0.6417

%% chl
%%%%%%%%%% chl, strutton 10m and apex 10m
[poly_chl1] = polyfit(strutbin_mean_chl_10m(~isnan(strutbin_mean_chl_10m)), ...
                    apex10m(~isnan(strutbin_mean_chl_10m),6), 1);
% for chl values after instrument correction: 2.5223 and -1.6276

[poly_chl2] = polyfit(apex10m(~isnan(strutbin_mean_chl_10m),6), ...
                     strutbin_mean_chl_10m(~isnan(strutbin_mean_chl_10m)),1);
% for chl values after instrument correction: 0.3372 and 1.0579

slope_chl=realsqrt(2.5223*(1/0.3372)); % 2.7350

%%%%%%%%%% all chl, strutton 10m, 15m, and 20m and apex 10m, 15m, and 20m
chl_strut=[strutbin_mean_chl_10m(~isnan(strutbin_mean_chl_10m));strutbin_mean_chl_15m(~isnan(strutbin_mean_chl_15m));strutbin_mean_chl_20m(~isnan(strutbin_mean_chl_20m))];
chl_apex=[apex10m(~isnan(strutbin_mean_chl_10m),6);apex15m(~isnan(strutbin_mean_chl_15m),6);apex15m(~isnan(strutbin_mean_chl_20m),6)];
[poly_chl1]=polyfit(chl_strut,chl_apex,1)
[poly_chl2]=polyfit(chl_apex,chl_strut,1)

figure;plot(chl_strut,chl_apex,'.')
slope_chl=sqrt(poly_chl1(1)*(1/poly_chl2(1))); % 2.5562
corrcoef(chl_strut,chl_apex)
% correlation coefficient: 0.7600

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% POC PLOTS
%%%%%%%%%%% Strutton POC vs Apex POC
figure(10)
subplot(2,2,1) % all apex data from 5 to 11m and strutton data at 2m
plot(date_apex(I1),POC_b_apex(I1),'b.')
hold on
plot(date_pscompile,poc_bb_2m_pscompile,'rx')
set(gca,'XTick',[732099 732129 732160 732191])
set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
datetick('x', 'mmmyy','keeplimits','keepticks')
% xlabel('Date')
ylabel('POC')
title({'ASIS Buoy (2m)','Float 0005 (5-11m)'})
legend('Float 0005','ASIS Buoy')
hold off

subplot(2,2,2) % all apex data from 5 to 11m and strutton data at 10m 
plot(date_apex(I1),POC_b_apex(I1),'b.')
hold on
plot(date_pscompile,poc_bb_10m_pscompile,'rx')
set(gca,'XTick',[732099 732129 732160 732191])
set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
datetick('x', 'mmmyy','keeplimits','keepticks')
% xlabel('Date')
ylabel('POC')
title({'ASIS Buoy (10m)','Float 0005 (5-11m)'})
legend('Float 0005','ASIS Buoy')
hold off

subplot(2,2,3) % all apex data from 12 to 17m and strutton data at 15m 
plot(date_apex(I2),POC_b_apex(I2),'b.')
hold on
plot(date_pscompile,poc_bb_15m_pscompile,'rx')
set(gca,'XTick',[732099 732129 732160 732191])
set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
datetick('x', 'mmmyy','keeplimits','keepticks')
% xlabel('Date')
ylabel('POC')
title({'ASIS Buoy (15m)','Float 0005 (12-17m)'})
legend('Float 0005','ASIS Buoy')
hold off

%%%%%%%%%%% Strutton POC vs Apex B_LSS_b_apex
figure(11)
subplot(2,2,1) % all apex data from 5 to 11m and strutton data at 2m
plot(date_apex(I1),B_LSS_b_apex(I1),'b.')
hold on
plot(date_pscompile,poc_bb_2m_pscompile,'rx')
set(gca,'XTick',[732099 732129 732160 732191])
set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
datetick('x', 'mmmyy','keeplimits','keepticks')
% xlabel('Date')
ylabel('')
title({'ASIS Buoy (2m)','Float 0005 (5-11m)'})
legend('Float 0005','ASIS Buoy')
hold off

subplot(2,2,2) % all apex data from 5 to 11m and strutton data at 10m 
plot(date_apex(I1),B_LSS_b_apex(I1),'b.')
hold on
plot(date_pscompile,poc_bb_10m_pscompile,'rx')
set(gca,'XTick',[732099 732129 732160 732191])
set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
datetick('x', 'mmmyy','keeplimits','keepticks')
% xlabel('Date')
ylabel('')
title({'ASIS Buoy (10m)','Float 0005 (5-11m)'})
legend('Float 0005','ASIS Buoy')
hold off

subplot(2,2,3) % all apex data from 12 to 17m and strutton data at 15m 
plot(date_apex(I2),B_LSS_b_apex(I2),'b.')
hold on
plot(date_pscompile,poc_bb_15m_pscompile,'rx')
set(gca,'XTick',[732099 732129 732160 732191])
set(gca,'xlim',[datenum('01 jun 2004') datenum('01 sep 2004')])
datetick('x', 'mmmyy','keeplimits','keepticks')
% xlabel('Date')
ylabel('')
title({'ASIS Buoy (15m)','Float 0005 (12-17m)'})
legend('Float 0005','ASIS Buoy')
hold off


