Index: gliderproc/trunk/extractDbd.m =================================================================== --- gliderproc/trunk/extractDbd.m (revision 501) +++ gliderproc/trunk/extractDbd.m (revision 503) @@ -6,5 +6,5 @@ % NOTE: At top of file, set glider index and deployment number % -% MATLAB folder - contains util, seawater, +% MATLAB folder - contains util % % Index: gliderproc/trunk/gliderOptode_Generate_L1_Data.m =================================================================== --- gliderproc/trunk/gliderOptode_Generate_L1_Data.m (revision 498) +++ gliderproc/trunk/gliderOptode_Generate_L1_Data.m (revision 503) @@ -4,393 +4,370 @@ % Purpose: Generate Level 1 data mat files that include optode DO data. % -% NOTE: At top of file, set glider index and deployment number -% % Requires: -% MATLAB folder - contains util, matutil, seawater, plots, strfun, opnml -% -% -% Author: William Stark -% Marine Sciences Department -% UNC-Chapel Hill -% -% Created: 21 May 2012 -% -%////////////////////////////////////////////////////////////////////////// +% MATLAB folder - contains util +% +% Authors: Adapted from William Stark by Harvey Seim and Chris Calloway +% Marine Sciences Department +% UNC-Chapel Hill +% +% Created: March 2013 +% +%/////////////////////////////////////////////////////////////////////////////// clear all; - - -%!!!!!! SET PRIOR TO RUNNING CODE !!!!!!!!!!!!!!!!!!!!!!!! -%!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! -% -% SET THE GLIDER INDEX (Pelagia = 1, Ramses = 2) ... -gliderIndex = 1; - -% SET THE DEPLOYMENT NUMBER (1, 2 or 3) ... -deploymentNumber = 1; - -% -%!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! -%!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - - % add paths for required files... addpath('MATLAB/util/'); -addpath('MATLAB/matutil/'); -addpath('MATLAB/seawater/'); -addpath('MATLAB/plots/'); -addpath('MATLAB/strfun/'); -addpath('MATLAB/opnml/'); -addpath('MATLAB/opnml/FEM/'); - - - -% glider name string... -if (gliderIndex==1) - strGliderName = 'Pelagia'; -else - strGliderName = 'Ramses'; -end % populate arrays for the deployment start and end dates... % ex. strStart(2, 3) is start date for Ramses, Deployment 3 -strStart = {'26-Jan-2012', '16-Feb-2012', '16-Mar-2012'; '26-Jan-2012', '16-Feb-2012', '16-Mar-2012'}; -strEnd = {'14-Feb-2012', '08-Mar-2012', '04-Apr-2012'; '14-Feb-2012', '12-Mar-2012', '03-Apr-2012'}; - -% deployment number string... -strDeploymentNumber = num2str(deploymentNumber); - -% deployment start date string... -strStartDate = strStart(gliderIndex, deploymentNumber); - -% deployment end date string... -strEndDate = strEnd(gliderIndex, deploymentNumber); - -% define the path to the glider ascii files... -%datadir = strcat('/Users/haloboy/Documents/MASC/MATLAB/CTD_data_correction/GLIDER_CTD_DATA_LEVEL0/',... -datadir = strcat('GLIDER_DATA_LEVEL0/', strGliderName, '_Deployment', strDeploymentNumber, '/'); - - - -% define default bounds for use in plots... -% switch gliderIndex -% case 1 % Pelagia -% switch deploymentNumber -% case 1 % Deployment 1 -% tempBounds = [17.0 24.0]; -% salinBounds = [36.0 36.4]; -% densBounds = [1025.0 1026.6]; -% chlorBounds = [0.0 4.0]; -% -% case 2 % Deployment 2 -% tempBounds = [17.0 24.0]; -% salinBounds = [36.0 36.5]; -% densBounds = [1024.5 1026.8]; -% chlorBounds = [0.0 4.0]; -% -% case 3 % Deployment 3 -% tempBounds = [17.0 24.0]; -% salinBounds = [35.9 36.7]; -% densBounds = [1024.4 1026.4]; -% chlorBounds = [0.0 4.0]; -% end -% case 2 % Ramses -% switch deploymentNumber -% case 1 % Deployment 1 -% tempBounds = [8.0 23.0]; -% salinBounds = [35.0 36.4]; -% densBounds = [1024.5 1027.5]; -% chlorBounds = [0.0 4.0]; -% -% case 2 % Deployment 2 -% tempBounds = [9.0 25.0]; -% salinBounds = [35.2 36.6]; -% densBounds = [1024.0 1027.5]; -% chlorBounds = [0.0 4.0]; -% -% case 3 % Deployment 3 -% tempBounds = [10.0 24.5]; -% salinBounds = [35.3 36.7]; -% densBounds = [1024.4 1027.4]; -% chlorBounds = [0.0 4.0]; -% end -% end - - -%########################################################################################## - - - - - -%*** READ IN EBD DATA ***************************************************** -% declare variables for storing data... -oxy_dphase=[]; -ptime_oxy=[]; - -% try to load all *.ebdasc files at once... -[files, Dstruct] = wilddir(datadir, '.ebdasc'); -nfile = size(files, 1); - -for i=1:nfile-1 - % protect against empty ebd file - if(Dstruct(i).bytes>0) - data = read_gliderasc2([datadir, files(i,:)]); - %data = read_gliderasc3([datadir, files(i,:)]); - - % if the number of values (in data.data) is less than the number of - % vars (in data.vars), this means that the data were not completely read - % in. To correct this, pad data.data with NaNs until its length - % equals that of data.vars... - if (length(data.data) < length(data.vars)) - data.data = padarray(data.data, [0 length(data.vars)-length(data.data)], NaN, 'post'); +strStart = {'26-Jan-2012', '16-Feb-2012', '16-Mar-2012'; ... + '26-Jan-2012', '16-Feb-2012', '16-Mar-2012'}; +strEnd = {'14-Feb-2012', '08-Mar-2012', '04-Apr-2012'; ... + '14-Feb-2012', '12-Mar-2012', '03-Apr-2012'}; + +% SET THE GLIDER INDEX (Pelagia = 1, Ramses = 2) ... +for gliderIndex=1:2 + + % SET THE DEPLOYMENT NUMBER (1, 2 or 3) ... + for deploymentNumber=1:3 + + clearvars -except gliderIndex deploymentNumber strStart strEnd; + + % glider name string... + if (gliderIndex==1) + strGliderName = 'Pelagia'; + else + strGliderName = 'Ramses'; end - % populate variables with data... - if(~isempty(data.data)) - oxy_dphase = [oxy_dphase;data.data(:,strmatch('sci_oxy3835_wphase_dphase', data.vars, 'exact'))]; % oxygen (dphase) - ptime_oxy = [ptime_oxy;data.data(:,strmatch('sci_m_present_time', data.vars, 'exact'))]; % present time - %should read in other oxy values for comparison with - %re-calculated values + % deployment number string... + strDeploymentNumber = num2str(deploymentNumber); + + % deployment start date string... + strStartDate = strStart(gliderIndex, deploymentNumber); + + % deployment end date string... + strEndDate = strEnd(gliderIndex, deploymentNumber); + + % define the path to the glider ascii files... + datadir = strcat('GLIDER_DATA_LEVEL0/', strGliderName, ... + '_Deployment', strDeploymentNumber, '/'); + + % define optode foil calibration coefficients + % by glider and deployment + switch gliderIndex + case 1 % Pelagia + switch deploymentNumber + case {1, 2} % Deployment 1 or 2; + % before foil replacement; + % 21 June 2004 calibration date; + % batch 2204 + C0coef=[3.12899E+03,-1.05764E+02,2.06640E+00,-1.68983E-02]; + C1coef=[-1.67671E+02,4.80582E+00,-8.73323E-02,6.61507E-04]; + C2coef=[3.73685E+00,-8.78300E-02,1.47560E-03,-9.96701E-06]; + C3coef=[-3.96052E-02,7.46930E-04,-1,17804E-05,6.677619E-08]; + C4coef=[1.61999E-04,-2.37870E-06,3.63223E-08,-1.62194E-10]; + case 3 % Deployment 3; + % replaced foil on 12 March 2012; + % 18 August 2010 calibration date; + % batch 1023 + C0coef=[4.27019336E+03,-1.32723585E+02,2.15629751E+00,-1.40275831E-02]; + C1coef=[-2.29729690E+02,5.74242078E+00,-6.85357898E-02,1.88612346E-04]; + C2coef=[5.06401550E+00,-9.62084932E-02,5.22180779E-04,7.70889717E-06]; + C3coef=[-5.26332308E-02,7.15467419E-04,3.31185072E-06,-1.86124024E-07]; + C4coef=[2.10916841E-04,-1.84087896E-06,-4.28645540E-08,1.11120317E-09]; + end + case 2 % Ramses all deployments; + % 2 June 2010 calibration date; + % batch 5009 + C0coef=[4.53793e3 -1.62595e2 3.29574 -2.79285e-2]; + C1coef=[-2.50953e2 8.02322 -1.58398e-1 1.31141e-3]; + C2coef=[5.66417 -1.59647e-1 3.07910e-3 -2.46265e-5]; + C3coef=[-5.99449e-2 1.48326e-3 -2.82110e-5 2.15156e-7]; + C4coef=[2.43614e-4 -5.26759e-6 1.00064e-7 -7.14320e-10]; end - data = []; - end + %*** READ IN EBD DATA **** + % declare arrays for accumulating data + oxyw_oxygen = []; + oxyw_saturation = []; + oxyw_temp = []; + oxyw_dphase = []; + oxyw_bphase = []; + oxyw_rphase = []; + oxyw_bamp = []; + oxyw_bpot = []; + oxyw_ramp = []; + oxyw_rawtemp = []; + oxyw_time = []; + oxyw_installed = []; + ptime_ebd = []; + + % try to load all *.ebdasc files at once... + [files, Dstruct] = wilddir(datadir, '.ebdasc'); + nfile = size(files, 1); + + for i=1:nfile-1 + % protect against empty ebd file + if(Dstruct(i).bytes>0) + data = read_gliderasc2([datadir, files(i,:)]); + + % if the number of values (in data.data) is less than the number + % of vars (in data.vars), this means that the data were not + % completely read in. To correct this, pad data.data with NaNs + % until its length equals that of data.vars... + if (length(data.data) < length(data.vars)) + data.data = padarray(data.data, ... + [0 length(data.vars)-length(data.data)], ... + NaN, 'post'); + end + + % concatenate variables with data... + if(~isempty(data.data)) + oxyw_oxygen = [oxyw_oxygen; ... + data.data(:,strmatch('sci_oxy3835_wphase_oxygen',... + data.vars, 'exact'))]; + oxyw_saturation = [oxyw_saturation; ... + data.data(:,strmatch('sci_oxy3835_wphase_saturation',... + data.vars, 'exact'))]; + oxyw_temp = [oxyw_temp; ... + data.data(:,strmatch('sci_oxy3835_wphase_temp',... + data.vars, 'exact'))]; + oxyw_dphase = [oxyw_dphase; ... + data.data(:,strmatch('sci_oxy3835_wphase_dphase',... + data.vars, 'exact'))]; + oxyw_bphase = [oxyw_bphase; ... + data.data(:,strmatch('sci_oxy3835_wphase_bphase',... + data.vars, 'exact'))]; + oxyw_rphase = [oxyw_rphase; ... + data.data(:,strmatch('sci_oxy3835_wphase_rphase',... + data.vars, 'exact'))]; + oxyw_bamp = [oxyw_bamp; ... + data.data(:,strmatch('sci_oxy3835_wphase_bamp',... + data.vars, 'exact'))]; + oxyw_bpot = [oxyw_bpot; ... + data.data(:,strmatch('sci_oxy3835_wphase_bpot',... + data.vars, 'exact'))]; + oxyw_ramp = [oxyw_ramp; ... + data.data(:,strmatch('sci_oxy3835_wphase_ramp',... + data.vars, 'exact'))]; + oxyw_rawtemp = [oxyw_rawtemp; ... + data.data(:,strmatch('sci_oxy3835_wphase_rawtemp',... + data.vars, 'exact'))]; + oxyw_time = [oxyw_time; ... + data.data(:,strmatch('sci_oxy3835_wphase_timestamp',... + data.vars, 'exact'))]; + oxyw_installed = [oxyw_installed; ... + data.data(:,strmatch('sci_oxy3835_wphase_is_installed',... + data.vars, 'exact'))]; + ptime_ebd = [ptime_ebd; ... + data.data(:,strmatch('sci_m_present_time',... + data.vars, 'exact'))]; + end + + data = []; + end + end + %*** END READ IN EBD DATA **** + + % remove nans from EBD data... + i = find(~isnan(oxyw_dphase)); + oxyw_dphase = oxyw_dphase(i); + ptime_ebd = ptime_ebd(i); + oxyw_oxygen = oxyw_oxygen(i); + oxyw_saturation = oxyw_saturation(i); + oxyw_temp = oxyw_temp(i); + oxyw_bphase = oxyw_bphase(i); + oxyw_rphase = oxyw_rphase(i); + oxyw_bamp = oxyw_bamp(i); + oxyw_bpot = oxyw_bpot(i); + oxyw_ramp = oxyw_ramp(i); + oxyw_rawtemp = oxyw_rawtemp(i); + oxyw_time = oxyw_time(i); + oxyw_installed = oxyw_installed(i); + + % apply the sort() function to ptime_oxy + % to make sure it increases monotonically... + [Y,I] = sort(ptime_ebd); + ptime_ebd = Y; + oxyw_oxygen = oxyw_oxygen(I); + oxyw_saturation = oxyw_saturation(I); + oxyw_temp = oxyw_temp(I); + oxyw_dphase = oxyw_dphase(I); + oxyw_bphase = oxyw_bphase(I); + oxyw_rphase = oxyw_rphase(I); + oxyw_bamp = oxyw_bamp(I); + oxyw_bpot = oxyw_bpot(I); + oxyw_ramp = oxyw_ramp(I); + oxyw_rawtemp = oxyw_rawtemp(I); + oxyw_time = oxyw_time(I); + oxyw_installed = oxyw_installed(I); + + % load CTD data from existing mat file... + strMatFilePath = strcat('GLIDER_DATA_LEVEL1/', strGliderName, ... + '_Deployment', strDeploymentNumber, '_CTD_L1.mat'); + load(strMatFilePath); + + % interpolate CTD to times of oxyw_dphase + tempi = interp1(ptime,temp,ptime_ebd); + salini = interp1(ptime,salinCorrected,ptime_ebd); + depthi = interp1(ptime,depth,ptime_ebd); + + % first, implement temperature-dependent correction to DO concentration + % that utilizes dphase as the input (manual page 30), coefficients from cal + % sheets (a 5x4 matrix of values) and will be glider (and foil) dependent. + C0 = C0coef(1) + (C0coef(2) .* tempi) + ... + (C0coef(3) .* (tempi.^2)) + (C0coef(4) .* (tempi.^3)); + C1 = C1coef(1) + (C1coef(2) .* tempi) + ... + (C1coef(3) .* (tempi.^2)) + (C1coef(4) .* (tempi.^3)); + C2 = C2coef(1) + (C2coef(2) .* tempi) + ... + (C2coef(3) .* (tempi.^2)) + (C2coef(4) .* (tempi.^3)); + C3 = C3coef(1) + (C3coef(2) .* tempi) + ... + (C3coef(3) .* (tempi.^2)) + (C3coef(4) .* (tempi.^3)); + C4 = C4coef(1) + (C4coef(2) .* tempi) + ... + (C4coef(3) .* (tempi.^2)) + (C4coef(4) .* (tempi.^3)); + + o2_tcorr = C0 + (C1 .* oxyw_dphase) + (C2 .* (oxyw_dphase.^2)) + ... + (C3 .* (oxyw_dphase.^3)) + (C4 .* (oxyw_dphase.^4)); + + % second, implement the salinity correction to DO concentration (page 31) + % ASSUMES DEFAULT SALINITY IN SENSOR WAS SET TO ZERO, OTHERWISE ANOHTER + % SCALING IS REQUIRED + B0=-6.24097e-3; + B1=-6.93498e-3; + B2=-6.90358e-3; + B3=-4.29155e-3; + C0=-3.11680e-7; + + Ts = log((298.15 - tempi) ./ (273.15 + tempi)); + + o2_tscorr = o2_tcorr .* ... + exp((salini .* (B0 + (B1 .* Ts) + (B2 .* (Ts.^2)) + ... + (B3 .* (Ts.^3)))) + ... + (C0 .* (salini.^2))); + + % third, implement the pressure correction to DO concentration (page 32) + os_tspcoor = o2_tscorr .* (1 + (0.04 .* depthi ./ 1000)); + + % use polynomial to calculate DO satuations using the measured temp and + % sal (manual page 30) + A0 = 2.00856; + A1 = 3.22400; + A2 = 3.99063; + A3 = 4.80299; + A4 = 9.78188e-1; + A5 = 1.71069; + B0 = -6.24097e-3; + B1 = -6.93498e-3; + B2 = -6.90358e-3; + B3 = -4.29155e-3; + C0 = -3.11680e-7; + + % need interpolated temp, salinity, pressure at times of DO obs + rslt = (A0 + (A1 .* Ts) + (A2 .* (Ts.^2)) + (A3 .* (Ts.^3)) + ... + (A4 .* (Ts.^4)) + (A5 .* (Ts.^5))) + ... + (salini .* (B0 + (B1 .* Ts) + (B2 .* (Ts.^2)) + (B3 .*(Ts.^3)))) + ... + (C0 .* (salini.^2)); + + o2_sol = exp(rslt); + + o2_sat = (os_tspcoor .* 2.2414) ./ o2_sol; + + % convert ptime into datenum style... + ptime_ebd_datenum = ptime_ebd/3600/24+datenum(1970, 1, 1, 0, 0, 0); + + % create configuration struct... + units = struct( ... + 'oxyw_oxygen', '10?6 mol/dm3', ... + 'oxyw_saturation', 'percent', ... + 'oxyw_temp', 'degrees C', ... + 'oxyw_dphase', 'degrees', ... + 'oxyw_bphase', 'degrees', ... + 'oxyw_rphase', 'degrees', ... + 'oxyw_bamp', 'mA', ... + 'oxyw_bpot', 'mV', ... + 'oxyw_ramp', 'mA', ... + 'oxyw_rawtemp', 'degrees C', ... + 'oxyw_time', 'timestamp', ... + 'oxyw_installed', 'bool', ... + 'ptime_ebd', 'seconds since 0000-01-01T00:00', ... + 'ptime_ebd_datenum', 'days since 1970-01-01T00:00', ... + 'tempi', 'degrees C', ... + 'salini', 'psu', ... + 'depthi', 'decibars', ... + 'o2_tcorr', '10?6 mol/dm3', ... + 'o2_tscorr', '10?6 mol/dm3', ... + 'o2_tspcoor', '10?6 mol/dm3', ... + 'o2_sol', 'cm3/liter at 1031 hPa', ... + 'os_sat', 'percent'); + + variable_description = struct( ... + 'oxyw_oxygen', 'dissolved oxygen', ... + 'oxyw_saturation', 'dissolved oxygen saturation', ... + 'oxyw_temp', 'water temperature', ... + 'oxyw_dphase', 'phase difference', ... + 'oxyw_bphase', 'blue phase', ... + 'oxyw_rphase', 'red phase', ... + 'oxyw_bamp', 'blue current bias', ... + 'oxyw_bpot', 'blue voltage bias', ... + 'oxyw_ramp', 'red current bias', ... + 'oxyw_rawtemp', 'raw water temperature', ... + 'oxyw_time', 'optode timestampe', ... + 'oxyw_installed', 'bool', ... + 'ptime_ebd', 'seconds since 0000-01-01T00:00', ... + 'ptime_ebd_datenum', 'days since 1970-01-01T00:00', ... + 'tempi', 'degrees C', ... + 'salini', 'psu', ... + 'depthi', 'decibars', ... + 'o2_tcorr', '10?6 mol/dm3', ... + 'o2_tscorr', '10?6 mol/dm3', ... + 'o2_tspcoor', '10?6 mol/dm3', ... + 'o2_sol', 'cm3/liter at 1031 hPa', ... + 'os_sat', 'percent'); + + correction_coefficients = struct('C0coef', C0coef, ... + 'C1coef', C1coef, ... + 'C2coef', C2coef, ... + 'C3coef', C3coef, ... + 'C4coef', C4coef); + + config = struct('glider_name', strGliderName,... + 'deployment_number', strDeploymentNumber,... + 'start_date', strStartDate,... + 'end_date', strEndDate,... + 'correction_coefficients', correction_parameters,... + 'var_descriptions', variable_description,... + 'var_units', units); + + % set Level 1 data mat file name... + strMatFileName = strcat(strGliderName, '_Deployment', strDeploymentNumber, '_DO_L1.mat'); + + % save flight data to mat file... + save(strMatFileName,... + 'config', ... + 'oxyw_oxygen', ... + 'oxyw_saturation', ... + 'oxyw_temp', ... + 'oxyw_dphase', ... + 'oxyw_bphase', ... + 'oxyw_rphase', ... + 'oxyw_bamp', ... + 'oxyw_bpot', ... + 'oxyw_ramp', ... + 'oxyw_rawtemp', ... + 'oxyw_time', ... + 'oxyw_installed', ... + 'ptime_ebd', ... + 'ptime_ebd_datenum', ... + 'tempi', ... + 'salini', ... + 'depthi', ... + 'o2_tcorr', ... + 'o2_tscorr', ... + 'o2_tspcoor', ... + 'o2_sol', ... + 'os_sat'); + + end end -%************************************************************************** - - -% load glider data from existing mat file... -strMatFilePath = strcat('GLIDER_DATA_LEVEL1/', strGliderName, '_Deployment', strDeploymentNumber, '_DataL1.mat'); -load(strMatFilePath); - - -% first, apply the sort() function to ptime_oxy to make sure it increases monotonically... -[Y,I] = sort(ptime_oxy); -ptime_oxy = Y; -oxy_dphase = oxy_dphase(I); - - -% remove nans from EBD data... -i = find(~isnan(oxy_dphase)); -oxy_dphase = oxy_dphase(i); ptime_oxy = ptime_oxy(i); - -% interpolate CTD to times of oxy_dphase -tempi = interp1(ptime,temp,ptime_oxy); -salini = interp1(ptime,salinCorrected,ptime_oxy); -depthi = interp1(ptime,depth,ptime_oxy); - -% first, implement temperature-dependent correction to DO concentration -% that utilizes dphase as the input (manual page 30), coefficients from cal -% sheets (a 5x4 matrix of values) and will be glider (and foil) dependent. - -% for ramses, 2 June 2010 calibration sheet - -C0coef=[4.53793e3 -1.62595e2 3.29574 -2.79285e-2]; -C1coef=[-2.50953e2 8.02322 -1.58398e-1 1.31141e-3]; -C2coef=[5.66417 -1.59647e-1 3.07910e-3 -2.46265e-5]; -C3coef=[-5.99449e-2 1.48326e-3 -2.82110e-5 2.15156e-7]; -C4coef=[2.43614e-4 -5.26759e-6 1.00064e-7 -7.14320e-10]; - -C0 = C0coef(1)+C0coef(2)*tempi+C0coef(3)*tempi.^2+C0coef(4)*tempi.^3; -C1 = C1coef(1)+C1coef(2)*tempi+C1coef(3)*tempi.^2+C1coef(4)*tempi.^3; -C2 = C2coef(1)+C2coef(2)*tempi+C2coef(3)*tempi.^2+C2coef(4)*tempi.^3; -C3 = C3coef(1)+C3coef(2)*tempi+C3coef(3)*tempi.^2+C3coef(4)*tempi.^3; -C4 = C4coef(1)+C4coef(2)*tempi+C4coef(3)*tempi.^2+C4coef(4)*tempi.^3; - -o2 = C0+C1.*oxy_dphase+C2.*oxy_dphase.^2+C3.*oxy_dphase.^3+C4.*oxy_dphase.^4; - -% second, implement the salinity correction to DO concentration (page 31) -% ASSUMES DEFAULT SALINITY IN SENSOR WAS SET TO ZERO, OTHERWISE ANOHTER -% SCALING IS REQUIRED - -B0=-6.24097e-3; -B1=-6.93498e-3; -B2=-6.90358e-3; -B3=-4.29155e-3; -C0=-3.11680e-7; - -Ts = log((298.15 - tempi) ./ (273.15 + tempi)); - -o2_salcorr = o2*exp(salini*(B0+B1*Ts+B2*Ts.*Ts+B3*Ts.*Ts.*Ts)+C0*salini.^2); - -% third, implement the pressure correction to DO concentration (page 32) - -o2_salprescorr = os_salcorr*(1 + 0.04*depthi/1000); - -% use polynomial to calculate DO satuations using the measured temp and -% sal (manual page 30) - -A0 = 2.00856; -A1 = 3.22400; -A2 = 3.99063; -A3 = 4.80299; -A4 = 0.0978188; -A5 = 1.71069; -B0 = -00624097; -B1 = -00693498; -B2 = -00690358; -B3 = -00429155; -C0 = -000000311680; - -% need interpolated temp, salinity, pressure at times of DO obs - -rslt = A0 + A1*Ts + A2*Ts.^2 + A3*Ts.^3 + A4*Ts.^4 + A5*Ts.^5 +... - salini.*(B0 + B1*Ts + B2*Ts.^2 + B3*Ts.^3) + C0*salini.^2; - -O2sol = exp(rslt); - -o2_sat = (o2_salprescorr * 2.2414) ./ O2sol; - - -% convert ptime into datenum style... -ptime_oxy_datenum = ptime_oxy/3600/24+datenum(1970, 1, 1, 0, 0, 0); -% ptimehrly = fix(ptime_datenum*24)/24; -% ptimedaily = fix(ptime_datenum); -% ptimedaily2 = unique(ptimedaily); -% ptimedaily2 = ptimedaily2(1:2:end); - - - - -lb = datenum('10-Feb-2012 12:22:00'); -ub = datenum('10-Feb-2012 12:23:00'); -%lb = datenum('22-Mar-2012 12:01:00'); -%ub = datenum('23-Mar-2012 12:02:00'); - -% figure('Position', [500,500,1000,1000]); -% subplot(311); -% ind = find(ptime_oxy_datenum >= lb & ptime_oxy_datenum <= ub); -% plot(ptime_oxy_datenum(ind), oxy_dphase(ind), 'b.'); -% datetick; -% xlabel('time', 'fontsize', 12); -% ylabel('DO (dphase)', 'fontsize', 12); -% subplot(312); -% ind2 = find(ptime_datenum >= lb & ptime_datenum <= ub); -% plot(ptime_datenum(ind2), temp(ind2), 'r.'); -% datetick; -% xlabel('time', 'fontsize', 12); -% ylabel('temp', 'fontsize', 12); -% subplot(313); -% ind = find(ptime_oxy_datenum >= lb & ptime_oxy_datenum <= ub); -% ind2 = find(ptime_datenum >= lb & ptime_datenum <= ub); -% plot(ptime_oxy_datenum(ind), ones(length(ind), 1), 'b.', ptime_datenum(ind2), ones(length(ind2), 1), 'r.'); -% datetick; - - -figure('Position', [500,500,1000,1000]); -subplot(211); -ind = find(ptime_datenum >= lb & ptime_datenum <= ub); -plot(ptime_datenum(ind), oxy_dphase(ind), 'b.', ptime_datenum(ind), oxy_dphase_adjusted(ind), 'r.'); -datetick; -xlabel('time', 'fontsize', 12); -ylabel('DO (dphase)', 'fontsize', 12); -subplot(212); -ind2 = find(ptime_datenum >= lb & ptime_datenum <= ub); -plot(ptime_datenum(ind2), temp(ind2), 'r.'); -datetick; -xlabel('time', 'fontsize', 12); -ylabel('temp', 'fontsize', 12); - - - - - - - - -% interpolate DBD data to align with EBD data... -% horizontalVelocity = interp1(ptime1_dbd, horizontalVelocity, ptime); -% depth = interp1(ptime1_dbd, depth, ptime); -% pitch = interp1(ptime1_dbd, pitch, ptime); -% avgDepthRate = interp1(ptime1_dbd, avgDepthRate, ptime); -% glideAngle = interp1(ptime1_dbd, glideAngle, ptime); -% lat = interp1(ptime1_dbd, lat, ptime); -% lon = interp1(ptime1_dbd, lon, ptime); - - - - - - - - - - - - -% create configuration struct... -% units = struct('chlor', 'micrograms liter-1',... -% 'dens', 'kg m-3',... -% 'densCorrected', 'kg m-3',... -% 'depth', 'm',... -% 'gpsLat', 'decimal degrees',... -% 'gpsLon', 'decimal degrees',... -% 'pres', 'decibars',... -% 'ptime', 'seconds since 0000-01-01T00:00',... -% 'ptime_datenum', 'seconds since 0000-01-01T00:00',... -% 'salin', 'psu',... -% 'salinCorrected', 'psu',... -% 'temp', 'deg C'); - -% variable_description = struct('chlor', 'chlorophyll measured by glider',... -% 'chlorBounds', 'default chlorophyll bounds for plots',... -% 'dens', 'density measured by glider',... -% 'densBounds', 'default density bounds for plots',... -% 'densCorrected', 'density corrected for thermal lag',... -% 'depth', 'depth calculated as function of pressure and position latitude',... -% 'gpsLat', 'position latitude measured by glider GPS',... -% 'gpsLon', 'position longitude measured by glider GPS',... -% 'pres', 'pressure measured by glider',... -% 'ptime', 'time vector reported by glider',... -% 'ptime_datenum', 'Serial Date Number string',... -% 'salin', 'salinity measured by glider',... -% 'salinBounds', 'default salinity bounds for plots',... -% 'salinCorrected', 'salinity corrected for thermal lag',... -% 'temp', 'temperature measured by glider',... -% 'tempBounds', 'default temperature bounds for plots'); - -% correction_parameters = struct('alpha_offset', correctionParams(1),... -% 'alpha_slope', correctionParams(2),... -% 'tau_offset', correctionParams(3),... -% 'tau_slope', correctionParams(4)); - -% config = struct('glider_name', strGliderName,... -% 'deployment_number', strDeploymentNumber,... -% 'start_date', strStartDate,... -% 'end_date', strEndDate,... -% 'thermal_lag_correction_parameters', correction_parameters,... -% 'var_descriptions', variable_description,... -% 'var_units', units); - - - - -% set Level 1 data mat file name... -%strMatFileName = strcat(strGliderName, '_Deployment', strDeploymentNumber, '_WithDO_LEVEL1.mat'); - - - -% save glider/deployment data to mat file... -% save(strMatFileName,... -% 'config',... -% 'chlor',... -% 'chlorBounds',... -% 'dens',... -% 'densBounds',... -% 'densCorrected',... -% 'depth',... -% 'gpsLat',... -% 'gpsLon',... -% 'pres',... -% 'ptime',... -% 'ptime_datenum',... -% 'salin',... -% 'salinBounds',... -% 'salinCorrected',... -% 'temp',... -% 'tempBounds'); - - -