diff --git a/AutomaticQC/imosTiltVelocitySetQC.m b/AutomaticQC/imosTiltVelocitySetQC.m index aed2337fb..40e686ae2 100644 --- a/AutomaticQC/imosTiltVelocitySetQC.m +++ b/AutomaticQC/imosTiltVelocitySetQC.m @@ -196,7 +196,12 @@ tilt = acos(sqrt(1 - sin(roll*pi/180).^2 - sin(pitch*pi/180).^2))*180/pi; % initially everything is failing the tests -sizeCur = size(sample_data.variables{idWcur}.flags); +if idUcur + idVar = idUcur; +else + idVar = idCspd; +end +sizeCur = size(sample_data.variables{idVar}.flags); flags = ones(sizeCur, 'int8')*secondFlagThreshold; % tilt test diff --git a/Graph/scatterMooring1DVarAgainstDepth.m b/Graph/scatterMooring1DVarAgainstDepth.m index 89f936ed1..cf749b68a 100644 --- a/Graph/scatterMooring1DVarAgainstDepth.m +++ b/Graph/scatterMooring1DVarAgainstDepth.m @@ -108,7 +108,7 @@ function scatterMooring1DVarAgainstDepth(sample_data, varName, isQC, saveToFile, xMin = min(xMin); xMax = max(xMax); -markerStyle = {'+', 'o', '*', '.', 'x', 's', 'd', '^', 'v', '>', '<', 'p', 'h'}; +markerStyle = {'+', 'o', '*', 's', 'd', '^', 'v', '>', '<', 'p', 'h'}; lenMarkerStyle = length(markerStyle); instrumentDesc = cell(lenSampleData + 1, 1); @@ -128,17 +128,20 @@ function scatterMooring1DVarAgainstDepth(sample_data, varName, isQC, saveToFile, iDepth = getVar(sample_data{iSort(i)}.variables, 'DEPTH'); izVar = getVar(sample_data{iSort(i)}.variables, varName); - if izVar > 0 && iDepth > 0 && ... - size(sample_data{iSort(i)}.variables{izVar}.data, 2) == 1 && ... % we're only plotting 1D variables with depth variable but no current + if izVar > 0 && size(sample_data{iSort(i)}.variables{izVar}.data, 2) == 1 && ... % we're only plotting 1D variables with depth variable but no current all(~strcmpi(sample_data{iSort(i)}.variables{izVar}.name, {'UCUR', 'VCUR', 'WCUR', 'CDIR', 'CSPD', 'VEL1', 'VEL2', 'VEL3'})) iGood = true(size(sample_data{iSort(i)}.variables{izVar}.data)); if isQC %get time, depth and var QC information timeFlags = sample_data{iSort(i)}.dimensions{iTime}.flags; - depthFlags = sample_data{iSort(i)}.variables{iDepth}.flags; varFlags = sample_data{iSort(i)}.variables{izVar}.flags; - iGood = (timeFlags == 1 | timeFlags == 2) & (varFlags == 1 | varFlags == 2) & (depthFlags == 1 | depthFlags == 2); + iGood = (timeFlags == 1 | timeFlags == 2) & (varFlags == 1 | varFlags == 2); + + if iDepth + depthFlags = sample_data{iSort(i)}.variables{iDepth}.flags; + iGood = iGood & (depthFlags == 1 | depthFlags == 2); + end end if any(iGood) @@ -235,14 +238,17 @@ function scatterMooring1DVarAgainstDepth(sample_data, varName, isQC, saveToFile, if isQC %get time, depth and var QC information timeFlags = sample_data{iSort(i)}.dimensions{iTime}.flags; - depthFlags = sample_data{iSort(i)}.variables{iDepth}.flags; varFlags = sample_data{iSort(i)}.variables{izVar}.flags; varValues = sample_data{iSort(i)}.variables{izVar}.data; iGood = (timeFlags == 1 | timeFlags == 2) & ... (varFlags == 1 | varFlags == 2) & ... ~isnan(varValues); - iGoodDepth = (depthFlags == 1 | depthFlags == 2); + + if iDepth + depthFlags = sample_data{iSort(i)}.variables{iDepth}.flags; + iGoodDepth = (depthFlags == 1 | depthFlags == 2); + end end if all(~iGood) && isQC @@ -250,7 +256,17 @@ function scatterMooring1DVarAgainstDepth(sample_data, varName, isQC, saveToFile, ', there is not any ' varName ' data with good flags.']); continue; else - depth = sample_data{iSort(i)}.variables{iDepth}.data; + if iDepth + depth = sample_data{iSort(i)}.variables{iDepth}.data; + else + if isfield(sample_data{iSort(i)}, 'instrument_nominal_depth') + depth = sample_data{iSort(i)}.instrument_nominal_depth*ones(size(iGood)); + else + fprintf('%s\n', ['Error : in ' sample_data{iSort(i)}.toolbox_input_file ... + ', global attribute instrument_nominal_depth is not documented.']); + continue; + end + end depth(~iGoodDepth) = metaDepth(i); depth = depth(iGood); diff --git a/Graph/scatterMooring2DVarAgainstDepth.m b/Graph/scatterMooring2DVarAgainstDepth.m index 9bb368173..38841e132 100644 --- a/Graph/scatterMooring2DVarAgainstDepth.m +++ b/Graph/scatterMooring2DVarAgainstDepth.m @@ -120,7 +120,7 @@ function scatterMooring2DVarAgainstDepth(sample_data, varName, isQC, saveToFile, return; end -markerStyle = {'+', 'o', '*', '.', 'x', 's', 'd', '^', 'v', '>', '<', 'p', 'h'}; +markerStyle = {'+', 'o', '*', 's', 'd', '^', 'v', '>', '<', 'p', 'h'}; lenMarkerStyle = length(markerStyle); instrumentDesc = cell(lenSampleData + 1, 1); @@ -139,10 +139,9 @@ function scatterMooring2DVarAgainstDepth(sample_data, varName, isQC, saveToFile, iTime = getVar(sample_data{iSort(i)}.dimensions, 'TIME'); iHeight = getVar(sample_data{iSort(i)}.dimensions, 'HEIGHT_ABOVE_SENSOR'); if iHeight == 0, iHeight = getVar(sample_data{iSort(i)}.dimensions, 'DIST_ALONG_BEAMS'); end % is equivalent when tilt is negligeable - iDepth = getVar(sample_data{iSort(i)}.variables, 'DEPTH'); iVar = getVar(sample_data{iSort(i)}.variables, varName); - if iVar > 0 && iHeight > 0 && iDepth > 0 && ... + if iVar > 0 && iHeight > 0 && ... size(sample_data{iSort(i)}.variables{iVar}.data, 2) > 1 && ... size(sample_data{iSort(i)}.variables{iVar}.data, 3) == 1 % we're plotting ADCP 2D variables with DEPTH variable. isPlottable(i) = true; @@ -163,7 +162,7 @@ function scatterMooring2DVarAgainstDepth(sample_data, varName, isQC, saveToFile, yLimMax = max(yLimMax, max(max(sample_data{iSort(i)}.variables{iVar}.data(iGood)))); end - elseif iVar > 0 && iDepth > 0 && ... + elseif iVar > 0 && ... any(strcmpi(sample_data{iSort(i)}.variables{iVar}.name, {'UCUR', 'VCUR', 'WCUR', 'CDIR', 'CSPD', 'VEL1', 'VEL2', 'VEL3'})) && ... size(sample_data{iSort(i)}.variables{iVar}.data, 2) == 1 % we're plotting current metre 1D variables with DEPTH variable. iGood = true(size(sample_data{iSort(i)}.variables{iVar}.data)); @@ -307,11 +306,14 @@ function scatterMooring2DVarAgainstDepth(sample_data, varName, isQC, saveToFile, if isQC %get time and var QC information timeFlags = sample_data{iSort(i)}.dimensions{iTime}.flags; - depthFlags = sample_data{iSort(i)}.variables{iDepth}.flags; varFlags = sample_data{iSort(i)}.variables{iVar}.flags; varValues = sample_data{iSort(i)}.variables{iVar}.data; - iGoodDepth = (depthFlags == 1 | depthFlags == 2); + if iDepth + depthFlags = sample_data{iSort(i)}.variables{iDepth}.flags; + iGoodDepth = (depthFlags == 1 | depthFlags == 2); + end + iGoodTime = (timeFlags == 1 | timeFlags == 2); iGood = repmat(iGoodTime, [1, size(sample_data{iSort(i)}.variables{iVar}.data, 2)]); @@ -337,7 +339,17 @@ function scatterMooring2DVarAgainstDepth(sample_data, varName, isQC, saveToFile, yScatter = 0; end - dataDepth = sample_data{iSort(i)}.variables{iDepth}.data; + if iDepth + dataDepth = sample_data{iSort(i)}.variables{iDepth}.data; + else + if isfield(sample_data{iSort(i)}, 'instrument_nominal_depth') + dataDepth = sample_data{iSort(i)}.instrument_nominal_depth*ones(size(iGoodTime)); + else + fprintf('%s\n', ['Error : in ' sample_data{iSort(i)}.toolbox_input_file ... + ', global attribute instrument_nominal_depth is not documented.']); + continue; + end + end dataDepth(~iGoodTime) = NaN; dataDepth(~iGoodDepth) = metaDepth(i); diff --git a/NetCDF/template/global_attributes_timeSeries.txt b/NetCDF/template/global_attributes_timeSeries.txt index 46896f80e..b73adedb8 100644 --- a/NetCDF/template/global_attributes_timeSeries.txt +++ b/NetCDF/template/global_attributes_timeSeries.txt @@ -23,6 +23,7 @@ S, metadata = S, sensorML = S, instrument_serial_number = [mat updateIfEmpty('[ddb InstrumentID Instruments SerialNumber]', sample_data.meta.instrument_serial_no)] N, instrument_sample_interval = [mat sample_data.meta.instrument_sample_interval] +N, instrument_average_interval = [mat sample_data.meta.instrument_average_interval] N, instrument_beam_angle = [mat sample_data.meta.beam_angle] N, instrument_burst_interval = [mat sample_data.meta.instrument_burst_interval] N, instrument_burst_duration = [mat sample_data.meta.instrument_burst_duration] diff --git a/Parser/aquadoppProfilerParse.m b/Parser/aquadoppProfilerParse.m index f69a0d4a6..b4293ea73 100644 --- a/Parser/aquadoppProfilerParse.m +++ b/Parser/aquadoppProfilerParse.m @@ -241,6 +241,7 @@ sample_data.meta.instrument_serial_no = hardware.SerialNo; sample_data.meta.instrument_firmware = hardware.FWversion; sample_data.meta.instrument_sample_interval = median(diff(time*24*3600)); +sample_data.meta.instrument_average_interval= user.AvgInterval; sample_data.meta.beam_angle = 25; % http://www.hydro-international.com/files/productsurvey_v_pdfdocument_19.pdf sample_data.meta.beam_to_xyz_transform = head.TransformationMatrix; @@ -255,7 +256,7 @@ end iWellOriented = adcpOrientations == adcpOrientation; % we'll only keep data collected when ADCP is oriented as expected dims = { - 'TIME', time(iWellOriented), ''; ... + 'TIME', time(iWellOriented), ['Time stamp corresponds to the start of the measurement which lasts ' num2str(user.AvgInterval) ' seconds.']; ... 'DIST_ALONG_BEAMS', distance, 'Nortek instrument data is not vertically bin-mapped (no tilt correction applied). Cells are lying parallel to the beams, at heights above sensor that vary with tilt.' }; clear time distance; @@ -277,6 +278,9 @@ end clear dims; +% add information about the middle of the measurement period +sample_data.dimensions{1}.seconds_to_middle_of_measurement = user.AvgInterval/2; + % add variables with their dimensions and data mapped. % we assume no correction for magnetic declination has been applied if velocityProcessed diff --git a/Parser/aquadoppVelocityParse.m b/Parser/aquadoppVelocityParse.m index 879edb026..5bf7d500c 100644 --- a/Parser/aquadoppVelocityParse.m +++ b/Parser/aquadoppVelocityParse.m @@ -158,12 +158,13 @@ sample_data.meta.instrument_serial_no = hardware.SerialNo; sample_data.meta.instrument_firmware = hardware.FWversion; sample_data.meta.instrument_sample_interval = median(diff(time*24*3600)); +sample_data.meta.instrument_average_interval= user.AvgInterval; sample_data.meta.beam_angle = 45; % http://wiki.neptunecanada.ca/download/attachments/18022846/Nortek+Aquadopp+Current+Meter+User+Manual+-+Rev+C.pdf sample_data.meta.featureType = mode; % add dimensions with their data mapped dims = { - 'TIME', time + 'TIME', time, ['Time stamp corresponds to the start of the measurement which lasts ' num2str(user.AvgInterval) ' seconds.'] }; clear time; @@ -173,9 +174,13 @@ sample_data.dimensions{i}.name = dims{i, 1}; sample_data.dimensions{i}.typeCastFunc = str2func(netcdf3ToMatlabType(imosParameters(dims{i, 1}, 'type'))); sample_data.dimensions{i}.data = sample_data.dimensions{i}.typeCastFunc(dims{i, 2}); + sample_data.dimensions{i}.comment = dims{i, 3}; end clear dims; +% add information about the middle of the measurement period +sample_data.dimensions{1}.seconds_to_middle_of_measurement = user.AvgInterval/2; + % add variables with their dimensions and data mapped. % we assume no correction for magnetic declination has been applied vars = { diff --git a/Parser/awacParse.m b/Parser/awacParse.m index 49ffb190b..09ed00e23 100644 --- a/Parser/awacParse.m +++ b/Parser/awacParse.m @@ -241,6 +241,7 @@ sample_data.meta.instrument_model = 'AWAC'; sample_data.meta.instrument_serial_no = hardware.SerialNo; sample_data.meta.instrument_sample_interval = median(diff(time*24*3600)); +sample_data.meta.instrument_average_interval= user.AvgInterval; sample_data.meta.instrument_firmware = hardware.FWversion; sample_data.meta.beam_angle = 25; % http://www.hydro-international.com/files/productsurvey_v_pdfdocument_19.pdf sample_data.meta.beam_to_xyz_transform = head.TransformationMatrix; @@ -256,7 +257,7 @@ end iWellOriented = adcpOrientations == adcpOrientation; % we'll only keep data collected when ADCP is oriented as expected dims = { - 'TIME', time(iWellOriented), ''; ... + 'TIME', time(iWellOriented), ['Time stamp corresponds to the start of the measurement which lasts ' num2str(user.AvgInterval) ' seconds.']; ... 'DIST_ALONG_BEAMS', distance, 'Nortek instrument data is not vertically bin-mapped (no tilt correction applied). Cells are lying parallel to the beams, at heights above sensor that vary with tilt.' }; clear time distance; @@ -278,6 +279,9 @@ end clear dims; +% add information about the middle of the measurement period +sample_data.dimensions{1}.seconds_to_middle_of_measurement = user.AvgInterval/2; + % add variables with their dimensions and data mapped. % we assume no correction for magnetic declination has been applied if velocityProcessed @@ -370,20 +374,36 @@ [filePath, fileRadName, ~] = fileparts(filename); filename = fullfile(filePath, [fileRadName '.wap']); -sample_data{2}.toolbox_input_file = filename; -sample_data{2}.meta.head = []; -sample_data{2}.meta.hardware = []; -sample_data{2}.meta.user = []; -sample_data{2}.meta.instrument_sample_interval = median(diff(waveData.Time*24*3600)); +sample_data{2}.toolbox_input_file = filename; +sample_data{2}.meta.head = []; +sample_data{2}.meta.hardware = []; +sample_data{2}.meta.user = []; +sample_data{2}.meta.instrument_sample_interval = median(diff(waveData.Time*24*3600)); + +avgInterval = []; +if isfield(waveData, 'summary') + iMatch = ~cellfun(@isempty, regexp(waveData.summary, 'Wave - Number of samples [0-9]')); + if any(iMatch) + nSamples = textscan(waveData.summary{iMatch}, 'Wave - Number of samples %f'); + + iMatch = ~cellfun(@isempty, regexp(waveData.summary, 'Wave - Sampling rate [0-9\.] Hz')); + if any(iMatch) + samplingRate = textscan(waveData.summary{iMatch}, 'Wave - Sampling rate %f Hz'); + avgInterval = nSamples{1}/samplingRate{1}; + end + end +end +sample_data{2}.meta.instrument_average_interval = avgInterval; +if isempty(avgInterval), avgInterval = '?'; end % we assume no correction for magnetic declination has been applied % add dimensions with their data mapped dims = { - 'TIME', waveData.Time; ... - 'FREQUENCY_1', waveData.pwrFrequency; ... - 'FREQUENCY_2', waveData.dirFrequency; ... - 'DIR_MAG', waveData.Direction + 'TIME', waveData.Time, ['Time stamp corresponds to the start of the measurement which lasts ' num2str(avgInterval) ' seconds.']; ... + 'FREQUENCY_1', waveData.pwrFrequency, ''; ... + 'FREQUENCY_2', waveData.dirFrequency, ''; ... + 'DIR_MAG', waveData.Direction, '' }; nDims = size(dims, 1); @@ -395,6 +415,9 @@ end clear dims; +% add information about the middle of the measurement period +sample_data{2}.dimensions{1}.seconds_to_middle_of_measurement = sample_data{2}.meta.instrument_average_interval/2; + % add variables with their dimensions and data mapped vars = { 'TIMESERIES', [], 1; ... diff --git a/Parser/continentalParse.m b/Parser/continentalParse.m index efb73f779..79c37b91e 100644 --- a/Parser/continentalParse.m +++ b/Parser/continentalParse.m @@ -243,6 +243,7 @@ sample_data.meta.instrument_model = 'Continental'; sample_data.meta.instrument_serial_no = hardware.SerialNo; sample_data.meta.instrument_sample_interval = median(diff(time*24*3600)); +sample_data.meta.instrument_average_interval= user.AvgInterval; sample_data.meta.instrument_firmware = hardware.FWversion; sample_data.meta.beam_angle = 25; % http://www.hydro-international.com/files/productsurvey_v_pdfdocument_19.pdf sample_data.meta.beam_to_xyz_transform = head.TransformationMatrix; @@ -258,7 +259,7 @@ end iWellOriented = adcpOrientations == adcpOrientation; % we'll only keep data collected when ADCP is oriented as expected dims = { - 'TIME', time(iWellOriented), ''; ... + 'TIME', time(iWellOriented), ['Time stamp corresponds to the start of the measurement which lasts ' num2str(user.AvgInterval) ' seconds.']; ... 'DIST_ALONG_BEAMS', distance, 'Nortek instrument data is not vertically bin-mapped (no tilt correction applied). Cells are lying parallel to the beams, at heights above sensor that vary with tilt.' }; clear time distance; @@ -280,6 +281,9 @@ end clear dims; +% add information about the middle of the measurement period +sample_data.dimensions{1}.seconds_to_middle_of_measurement = user.AvgInterval/2; + % add variables with their dimensions and data mapped. % we assume no correction for magnetic declination has been applied if velocityProcessed diff --git a/Parser/readAWACWaveAscii.m b/Parser/readAWACWaveAscii.m index feb96bc2b..337a4913b 100644 --- a/Parser/readAWACWaveAscii.m +++ b/Parser/readAWACWaveAscii.m @@ -153,6 +153,7 @@ % transform the filename into processed wave data filenames [path, name] = fileparts(filename); +summaryFile = fullfile(path, [name '.hdr']); headerFile = fullfile(path, [name '.whd']); waveFile = fullfile(path, [name '.wap']); dirFreqFile = fullfile(path, [name '.wdr']); @@ -170,6 +171,15 @@ end try + waveData = struct; + + if exist(summaryFile, 'file') + summaryFileID = fopen(summaryFile); + summary = textscan(summaryFileID, '%s', 'Delimiter', ''); + waveData.summary = summary{1}; + fclose(summaryFileID); + end + header = importdata(headerFile); wave = importdata(waveFile); dirFreq = importdata(dirFreqFile); @@ -198,8 +208,6 @@ iWave = ismember(totalTime, waveTime); clear headerTime waveTime; - waveData = struct; - % copy data over to struct waveData.Time = totalTime; nTime = length(totalTime); diff --git a/Parser/readWorkhorseWaveAscii.m b/Parser/readWorkhorseWaveAscii.m index b6ccf31f5..2572d8c3d 100644 --- a/Parser/readWorkhorseWaveAscii.m +++ b/Parser/readWorkhorseWaveAscii.m @@ -66,7 +66,17 @@ waveData = struct; % transform the filename into a path -filePath = fileparts(filename); +[filePath, name, ext] = fileparts(filename); + +% read the summary file +summaryFile = fullfile(filePath, [name '.txt']); + +if exist(summaryFile, 'file') + summaryFileID = fopen(summaryFile); + summary = textscan(summaryFileID, '%s', 'Delimiter', ''); + waveData.summary = summary{1}; + fclose(summaryFileID); +end % Load the *_LOG9.TXT file logFile = dir(fullfile(filePath, '*_LOG9.TXT')); diff --git a/Parser/readXR420.m b/Parser/readXR420.m index aeb7b3841..3a68bf295 100644 --- a/Parser/readXR420.m +++ b/Parser/readXR420.m @@ -76,7 +76,7 @@ sample_data.meta.featureType = mode; if isfield(header, 'correction'), sample_data.meta.correction = header.correction; end - if isfield(header, 'averaging_time_period'), sample_data.meta.averaging_time_period = header.averaging_time_period; end + if isfield(header, 'averaging_time_period'), sample_data.meta.instrument_average_interval = header.averaging_time_period; end if isfield(header, 'burst_sample_rate'), sample_data.meta.burst_sample_rate = header.burst_sample_rate; end sample_data.dimensions = {}; @@ -287,6 +287,10 @@ sample_data.dimensions{1}.name = 'TIME'; sample_data.dimensions{1}.typeCastFunc = str2func(netcdf3ToMatlabType(imosParameters(sample_data.dimensions{1}.name, 'type'))); sample_data.dimensions{1}.data = sample_data.dimensions{1}.typeCastFunc(data.time); + if isfield(header, 'averaging_time_period') + sample_data.dimensions{1}.comment = ['Time stamp corresponds to the start of the measurement which lasts ' num2str(header.averaging_time_period) ' seconds.']; + sample_data.dimensions{1}.seconds_to_middle_of_measurement = header.averaging_time_period/2; %assume averaging time is always in seconds + end sample_data.variables{end+1}.name = 'TIMESERIES'; sample_data.variables{end}.typeCastFunc = str2func(netcdf3ToMatlabType(imosParameters(sample_data.variables{end}.name, 'type'))); @@ -465,7 +469,7 @@ data.time = header.start:header.interval:header.end; data.time = data.time(1:length(samples{1}))'; - if isfield(header, 'averaging_time_period') - data.time = data.time + (header.averaging_time_period/(24*3600))/2; %assume averaging time is always in seconds - end +% if isfield(header, 'averaging_time_period') +% data.time = data.time + (header.averaging_time_period/(24*3600))/2; %assume averaging time is always in seconds +% end end diff --git a/Parser/signatureParse.m b/Parser/signatureParse.m index 32e369f21..48947e817 100644 --- a/Parser/signatureParse.m +++ b/Parser/signatureParse.m @@ -49,101 +49,226 @@ % only one file supported filename = filename{1}; +[path, name, ext] = fileparts(filename); -% read in all of the structures in the raw file -structures = readAD2CPBinary(filename); +isMagBias = false; +magDec = 0; -instrumentConfig = ''; -if isfield(structures, 'IdA0') - % this is a string data record - if structures.IdA0.Data.Id == 16 - instrumentConfig = structures.IdA0.Data.String; - end -end - -if isfield(structures, 'Id15Version3') - % this is a burst data record - dataRecordType = 'Id15Version3'; -elseif isfield(structures, 'Id16Version3') - % this is an average data record - dataRecordType = 'Id16Version3'; -else - % string/bottom track/interleaved data record is not supported - error('Only burst/average data record version 3 format is supported'); -end - -nSamples = length(structures.(dataRecordType).Header); -nCells = structures.(dataRecordType).Data(1).nCells; +switch lower(ext) + case '.ad2cp' + % read in all of the structures in the raw file + structures = readAD2CPBinary(filename); + + instrumentModel = ''; + if isfield(structures, 'IdA0') + % this is a string data record + if structures.IdA0.Data.Id == 16 + % looking for instrument model + instrumentModel = regexp(structures.IdA0.Data.String, 'Signature[0-9]*', 'match', 'once'); + + % check for magnetic declination + stringCell = textscan(structures.IdA0.Data.String, '%s', 'Delimiter', ','); + iMagDec = strncmp('DECL=', stringCell{1}, 5); + if any(iMagDec) + magDec = stringCell{1}{iMagDec}; + magDec = textscan(magDec, 'DECL=%f'); + magDec = magDec{1}; + end + end + end + + if isfield(structures, 'Id15Version3') + % this is a burst data record + dataRecordType = 'Id15Version3'; + elseif isfield(structures, 'Id16Version3') + % this is an average data record + dataRecordType = 'Id16Version3'; + else + % string/bottom track/interleaved data record is not supported + error('Only burst/average data record version 3 format is supported'); + end + + nSamples = length(structures.(dataRecordType).Header); + nCells = structures.(dataRecordType).Data(1).nCells; + + serialNumber = num2str(unique(vertcat(structures.(dataRecordType).Data.SerialNumber))); + + cellSize = unique(vertcat(structures.(dataRecordType).Data.CellSize))*0.001; % m + blankDist = unique(vertcat(structures.(dataRecordType).Data.Blanking))*0.001; % m + if length(cellSize) > 1, error('Multiple cell sizes/blanking distance not supported'); end + + % retrieve sample data + Time = vertcat(structures.(dataRecordType).Data.Time); + speedOfSound = vertcat(structures.(dataRecordType).Data.SpeedOfSound)*0.1; % m/s + Temperature = vertcat(structures.(dataRecordType).Data.Temperature)*0.01; % degree Celsius + Pressure = vertcat(structures.(dataRecordType).Data.Pressure)*0.001; % dBar + Heading = vertcat(structures.(dataRecordType).Data.Heading)*0.01; % deg + Pitch = vertcat(structures.(dataRecordType).Data.Pitch)*0.01; % deg + Roll = vertcat(structures.(dataRecordType).Data.Roll)*0.01; % deg + Battery = vertcat(structures.(dataRecordType).Data.BatteryVoltage)*0.1; % Volt + Status = vertcat(structures.(dataRecordType).Data.Status); + velocity = cat(3, structures.(dataRecordType).Data.VelocityData); + velocityScaling = repmat(10.^vertcat(structures.(dataRecordType).Data.VelocityScaling), 1, nCells); + Velocity_E = squeeze(velocity(1,:,:))'.*velocityScaling; % m/s + Velocity_N = squeeze(velocity(2,:,:))'.*velocityScaling; + Velocity_U = squeeze(velocity(3,:,:))'.*velocityScaling; + Velocity_U2 = squeeze(velocity(4,:,:))'.*velocityScaling; + amplitude = cat(3, structures.(dataRecordType).Data.AmplitudeData); + Backscatter1 = squeeze(amplitude(1,:,:))'; % count + Backscatter2 = squeeze(amplitude(2,:,:))'; + Backscatter3 = squeeze(amplitude(3,:,:))'; + Backscatter4 = squeeze(amplitude(4,:,:))'; + correlation = cat(3, structures.(dataRecordType).Data.CorrelationData); + Correlation1 = squeeze(correlation(1,:,:))'; % percent [0 - 100] + Correlation2 = squeeze(correlation(2,:,:))'; + Correlation3 = squeeze(correlation(3,:,:))'; + Correlation4 = squeeze(correlation(4,:,:))'; + clear structures velocity velocityScaling amplitude correlation; + + case '.mat' + % read in all of the structures in the file + structures = load(filename); + dataFields = fieldnames(structures.Data); + + % look for potentially more data in other .mat files + [path, name, ext] = fileparts(filename); + iMatFile = 1; + nextFilename = fullfile(path, [name '_' num2str(iMatFile) ext]); + while exist(nextFilename, 'file') + nextStructures = load(nextFilename); + % for each field of the structure we append new data to the + % existing ones + for i=1:length(dataFields) + structures.Data.(dataFields{i}) = [structures.Data.(dataFields{i}); nextStructures.Data.(dataFields{i})]; + end + + iMatFile = iMatFile + 1; + nextFilename = fullfile(path, [name '_' num2str(iMatFile) ext]); + clear nextStructures; + end + + % investigate which mode has been set for acquisition + if strcmpi(structures.Config.Plan_BurstEnabled, 'TRUE') + acquisitionMode = 'Burst'; + end + + if strcmpi(structures.Config.Plan_AverageEnabled, 'TRUE') + acquisitionMode = 'Average'; + end -serialNumber = num2str(unique(vertcat(structures.(dataRecordType).Data.SerialNumber))); + % we need to sort the data by time just in case + [Time, iSort] = sort(structures.Data.([acquisitionMode '_Time'])); + for i=1:length(dataFields) + structures.Data.(dataFields{i}) = structures.Data.(dataFields{i})(iSort, :); + end + + nSamples = length(Time); + nCells = double(structures.Config.([acquisitionMode '_NCells'])); + + serialNumber = num2str(structures.Config.SerialNo); + instrumentModel = structures.Config.InstrumentName; + + cellSize = structures.Config.([acquisitionMode '_CellSize']); + blankDist = structures.Config.([acquisitionMode '_BlankingDistance']); + + % check for magnetic declination + magDec = structures.Config.Declination; + + [~, ~, cpuEndianness] = computer; + + Status = dec2bin(bytecast(structures.Data.([acquisitionMode '_Status']), 'L', 'uint32', cpuEndianness), 32); + % Error = structures.Data.Average_Error; % error codes for each cell of a velocity profile inferred from the beams. 0=good; otherwise error. See http://www.nortek-as.com/en/knowledge-center/forum/waves/20001875?b_start=0#769595815 + % AmbiguityVel = structures.Data.Average_AmbiguityVel; + + if isfield(structures.Data, [acquisitionMode '_Velocity_ENU']) + Velocity_E = squeeze(structures.Data.([acquisitionMode '_ENU'])(:, 1, :)); + Velocity_N = squeeze(structures.Data.([acquisitionMode '_ENU'])(:, 2, :)); + Velocity_U = squeeze(structures.Data.([acquisitionMode '_ENU'])(:, 3, :)); + Velocity_U2 = squeeze(structures.Data.([acquisitionMode '_ENU'])(:, 4, :)); + else + Velocity_E = structures.Data.([acquisitionMode '_VelEast']); + Velocity_N = structures.Data.([acquisitionMode '_VelNorth']); + Velocity_U = structures.Data.([acquisitionMode '_VelUp']); + Velocity_U2 = structures.Data.([acquisitionMode '_VelUp2']); + end + if isfield(structures.Data, [acquisitionMode '_Amplitude_Beam']) + Backscatter1 = squeeze(structures.Data.([acquisitionMode '_Amplitude_Beam'])(:, 1, :))*2; % looks like the .mat format is giving raw counts * 0.5 by default + Backscatter2 = squeeze(structures.Data.([acquisitionMode '_Amplitude_Beam'])(:, 2, :))*2; + Backscatter3 = squeeze(structures.Data.([acquisitionMode '_Amplitude_Beam'])(:, 3, :))*2; + Backscatter4 = squeeze(structures.Data.([acquisitionMode '_Amplitude_Beam'])(:, 4, :))*2; + else + Backscatter1 = structures.Data.([acquisitionMode '_AmpBeam1'])*2; + Backscatter2 = structures.Data.([acquisitionMode '_AmpBeam2'])*2; + Backscatter3 = structures.Data.([acquisitionMode '_AmpBeam3'])*2; + Backscatter4 = structures.Data.([acquisitionMode '_AmpBeam4'])*2; + end + if isfield(structures.Data, [acquisitionMode '_Correlation_Beam']) + Correlation1 = squeeze(structures.Data.([acquisitionMode '_Correlation_Beam'])(:, 1, :)); + Correlation2 = squeeze(structures.Data.([acquisitionMode '_Correlation_Beam'])(:, 2, :)); + Correlation3 = squeeze(structures.Data.([acquisitionMode '_Correlation_Beam'])(:, 3, :)); + Correlation4 = squeeze(structures.Data.([acquisitionMode '_Correlation_Beam'])(:, 4, :)); + else + Correlation1 = structures.Data.([acquisitionMode '_CorBeam1']); + Correlation2 = structures.Data.([acquisitionMode '_CorBeam2']); + Correlation3 = structures.Data.([acquisitionMode '_CorBeam3']); + Correlation4 = structures.Data.([acquisitionMode '_CorBeam4']); + end + Battery = structures.Data.([acquisitionMode '_Battery']); + Heading = structures.Data.([acquisitionMode '_Heading']); + Pitch = structures.Data.([acquisitionMode '_Pitch']); + Roll = structures.Data.([acquisitionMode '_Roll']); + Temperature = structures.Data.([acquisitionMode '_Temperature']); + speedOfSound = structures.Data.([acquisitionMode '_Soundspeed']); + Pressure = structures.Data.([acquisitionMode '_Pressure']); + clear structures; -% generate distance values -cellSize = unique(vertcat(structures.(dataRecordType).Data.CellSize))*0.001; % m -blankDist = unique(vertcat(structures.(dataRecordType).Data.Blanking))*0.001; % m -if length(cellSize) > 1, error('Multiple cell sizes/blanking distance not supported'); end -distance = nan(nCells, 1); -distance(:) = (blankDist): ... - (cellSize): ... - (blankDist + (nCells-1) * cellSize); - -% Note this is actually the distance between the ADCP's transducers and the -% middle of each cell -% See http://www.nortek-bv.nl/en/knowledge-center/forum/current-profilers-and-current-meters/579860330 -distance = (distance + cellSize); - -% retrieve sample data -time = vertcat(structures.(dataRecordType).Data.Time); -speedOfSound = vertcat(structures.(dataRecordType).Data.SpeedOfSound)*0.1; % m/s -temperature = vertcat(structures.(dataRecordType).Data.Temperature)*0.01; % degree Celsius -pressure = vertcat(structures.(dataRecordType).Data.Pressure)*0.001; % dBar -heading = vertcat(structures.(dataRecordType).Data.Heading)*0.01; % deg -pitch = vertcat(structures.(dataRecordType).Data.Pitch)*0.01; % deg -roll = vertcat(structures.(dataRecordType).Data.Roll)*0.01; % deg -battery = vertcat(structures.(dataRecordType).Data.BatteryVoltage)*0.1; % Volt -status = vertcat(structures.(dataRecordType).Data.Status); -velocity = cat(3, structures.(dataRecordType).Data.VelocityData); -velocityScaling = repmat(10.^vertcat(structures.(dataRecordType).Data.VelocityScaling), 1, nCells); -velocity1 = squeeze(velocity(1,:,:))'.*velocityScaling; % m/s -velocity2 = squeeze(velocity(2,:,:))'.*velocityScaling; -velocity3 = squeeze(velocity(3,:,:))'.*velocityScaling; -velocity4 = squeeze(velocity(4,:,:))'.*velocityScaling; -amplitude = cat(3, structures.(dataRecordType).Data.AmplitudeData); -backscatter1 = squeeze(amplitude(1,:,:))'; % count -backscatter2 = squeeze(amplitude(2,:,:))'; -backscatter3 = squeeze(amplitude(3,:,:))'; -backscatter4 = squeeze(amplitude(4,:,:))'; -correlation = cat(3, structures.(dataRecordType).Data.CorrelationData); -correlation1 = squeeze(correlation(1,:,:))'; % percent [0 - 100] -correlation2 = squeeze(correlation(2,:,:))'; -correlation3 = squeeze(correlation(3,:,:))'; -correlation4 = squeeze(correlation(4,:,:))'; -clear structures velocity velocityScaling amplitude correlation; + otherwise + error('Data format not supported'); + +end sample_data = struct; +if magDec ~= 0 + isMagBias = true; + magBiasComment = ['A compass correction of ' num2str(magDec) ... + 'degrees has been applied to the data by a technician using Nortek''s software ' ... + '(usually to account for magnetic declination).']; +end + sample_data.toolbox_input_file = filename; sample_data.meta.featureType = ''; % strictly this dataset cannot be described as timeSeriesProfile since it also includes timeSeries data like TEMP sample_data.meta.binSize = cellSize; sample_data.meta.instrument_make = 'Nortek'; -if isempty(instrumentConfig) +if isempty(instrumentModel) sample_data.meta.instrument_model = 'Signature'; else - sample_data.meta.instrument_model = regexp(instrumentConfig, 'Signature[0-9]*', 'match', 'once'); + sample_data.meta.instrument_model = instrumentModel; end sample_data.meta.instrument_serial_no = serialNumber; -sample_data.meta.instrument_sample_interval = median(diff(time*24*3600)); +sample_data.meta.instrument_sample_interval = median(diff(Time*24*3600)); switch sample_data.meta.instrument_model - case {'Signature500', 'Signature1000'} - sample_data.meta.beam_angle = 25; - otherwise % Signature250 + case 'Signature250' sample_data.meta.beam_angle = 20; + otherwise % Signature500 and Signature1000 + sample_data.meta.beam_angle = 25; end +% generate distance values +distance = nan(nCells, 1); +distance(:) = (blankDist): ... + (cellSize): ... + (blankDist + (nCells-1) * cellSize); + +% distance between the ADCP's transducers and the middle of each cell +% See http://www.nortek-bv.nl/en/knowledge-center/forum/current-profilers-and-current-meters/579860330 +distance = (distance + cellSize); + % add dimensions with their data mapped iStartOrientation = 26; iEndOrientation = 28; -adcpOrientations = bin2dec(status(:, end-iEndOrientation+1:end-iStartOrientation+1)); +adcpOrientations = bin2dec(Status(:, end-iEndOrientation+1:end-iStartOrientation+1)); adcpOrientation = mode(adcpOrientations); % hopefully the most frequent value reflects the orientation when deployed % we assume adcpOrientation == 4 by default "ZUP" if adcpOrientation == 5 @@ -152,10 +277,10 @@ end iWellOriented = all(adcpOrientations == repmat(adcpOrientation, nSamples, 1), 2); % we'll only keep data collected when ADCP is oriented as expected dims = { - 'TIME', time(iWellOriented), ''; ... + 'TIME', Time(iWellOriented), ''; ... 'DIST_ALONG_BEAMS', distance, 'Nortek instrument data is not vertically bin-mapped (no tilt correction applied). Cells are lying parallel to the beams, at heights above sensor that vary with tilt.' }; -clear time distance; +clear Time distance; nDims = size(dims, 1); sample_data.dimensions = cell(nDims, 1); @@ -169,6 +294,11 @@ % add variables with their dimensions and data mapped. % we assume no correction for magnetic declination has been applied +if isMagBias + magExt = ''; +else + magExt = '_MAG'; +end iDimVel = nDims; iDimDiag = nDims; vars = { @@ -176,30 +306,30 @@ 'LATITUDE', [], NaN; ... 'LONGITUDE', [], NaN; ... 'NOMINAL_DEPTH', [], NaN; ... - 'VCUR_MAG', [1 iDimVel], velocity2(iWellOriented, :); ... % V - 'UCUR_MAG', [1 iDimVel], velocity1(iWellOriented, :); ... % U - 'WCUR', [1 iDimVel], velocity3(iWellOriented, :); ... - 'WCUR_2', [1 iDimVel], velocity4(iWellOriented, :); ... - 'ABSIC1', [1 iDimDiag], backscatter1(iWellOriented, :); ... - 'ABSIC2', [1 iDimDiag], backscatter2(iWellOriented, :); ... - 'ABSIC3', [1 iDimDiag], backscatter3(iWellOriented, :); ... - 'ABSIC4', [1 iDimDiag], backscatter4(iWellOriented, :); ... - 'CMAG1', [1 iDimDiag], correlation1(iWellOriented, :); ... - 'CMAG2', [1 iDimDiag], correlation2(iWellOriented, :); ... - 'CMAG3', [1 iDimDiag], correlation3(iWellOriented, :); ... - 'CMAG4', [1 iDimDiag], correlation4(iWellOriented, :); ... - 'TEMP', 1, temperature(iWellOriented); ... - 'PRES_REL', 1, pressure(iWellOriented); ... + ['VCUR' magExt], [1 iDimVel], Velocity_N(iWellOriented, :); ... % V + ['UCUR' magExt], [1 iDimVel], Velocity_E(iWellOriented, :); ... % U + 'WCUR', [1 iDimVel], Velocity_U(iWellOriented, :); ... + 'WCUR_2', [1 iDimVel], Velocity_U2(iWellOriented, :); ... + 'ABSIC1', [1 iDimDiag], Backscatter1(iWellOriented, :); ... + 'ABSIC2', [1 iDimDiag], Backscatter2(iWellOriented, :); ... + 'ABSIC3', [1 iDimDiag], Backscatter3(iWellOriented, :); ... + 'ABSIC4', [1 iDimDiag], Backscatter4(iWellOriented, :); ... + 'CMAG1', [1 iDimDiag], Correlation1(iWellOriented, :); ... + 'CMAG2', [1 iDimDiag], Correlation2(iWellOriented, :); ... + 'CMAG3', [1 iDimDiag], Correlation3(iWellOriented, :); ... + 'CMAG4', [1 iDimDiag], Correlation4(iWellOriented, :); ... + 'TEMP', 1, Temperature(iWellOriented); ... + 'PRES_REL', 1, Pressure(iWellOriented); ... 'SSPD', 1, speedOfSound(iWellOriented); ... - 'VOLT', 1, battery(iWellOriented); ... - 'PITCH', 1, pitch(iWellOriented); ... - 'ROLL', 1, roll(iWellOriented); ... - 'HEADING_MAG', 1, heading(iWellOriented) + 'VOLT', 1, Battery(iWellOriented); ... + 'PITCH', 1, Pitch(iWellOriented); ... + 'ROLL', 1, Roll(iWellOriented); ... + ['HEADING' magExt], 1, Heading(iWellOriented) }; -clear analn1 analn2 time distance velocity1 velocity2 velocity3 velocity4 ... - backscatter1 backscatter2 backscatter3 backscatter4 ... - correlation1 correlation2 correlation3 correlation4 ... - temperature pressure speedOfSound battery pitch roll heading status; +clear Velocity_N Velocity_E Velocity_U Velocity_U2 ... + Backscatter1 Backscatter2 Backscatter3 Backscatter4 ... + Correlation1 Correlation2 Correlation3 Correlation4 ... + Temperature Pressure speedOfSound Battery Pitch Roll Heading Status; nVars = size(vars, 1); sample_data.variables = cell(nVars, 1); @@ -215,5 +345,10 @@ end end sample_data.variables{i}.data = sample_data.variables{i}.typeCastFunc(vars{i, 3}); + + if any(strcmpi(vars{i, 1}, {'VCUR', 'UCUR', 'HEADING'})) + sample_data.variables{i}.compass_correction_applied = magDec; + sample_data.variables{i}.comment = magBiasComment; + end end clear vars; diff --git a/Parser/workhorseParse.m b/Parser/workhorseParse.m index 0edcc8ebd..43f3f9e2a 100644 --- a/Parser/workhorseParse.m +++ b/Parser/workhorseParse.m @@ -138,8 +138,10 @@ variable.y2kMinute,... variable.y2kSecond + variable.y2kHundredth/100.0]); - % shift the timestamp to the middle of the burst - time = time + (fixed.pingsPerEnsemble .* (fixed.tppMinutes*60 + fixed.tppSeconds + fixed.tppHundredths/100) / (3600 * 24))/2; + timePerPing = fixed.tppMinutes*60 + fixed.tppSeconds + fixed.tppHundredths/100; + timePerEnsemble = fixed.pingsPerEnsemble .* timePerPing; +% % shift the timestamp to the middle of the burst +% time = time + (timePerEnsemble / (3600 * 24))/2; % % auxillary data @@ -253,6 +255,7 @@ sample_data.meta.instrument_model = [model ' Workhorse ADCP']; sample_data.meta.instrument_serial_no = serial; sample_data.meta.instrument_sample_interval = median(diff(time*24*3600)); + sample_data.meta.instrument_average_interval = mode(timePerEnsemble); sample_data.meta.instrument_firmware = ... strcat(num2str(fixed.cpuFirmwareVersion(1)), '.', num2str(fixed.cpuFirmwareRevision(1))); % we assume the first value is correct for the rest of the dataset if all(isnan(fixed.beamAngle)) @@ -272,7 +275,7 @@ end iWellOriented = adcpOrientations == adcpOrientation; % we'll only keep data collected when ADCP is oriented as expected dims = { - 'TIME', time(iWellOriented), ''; ... + 'TIME', time(iWellOriented), ['Time stamp corresponds to the start of the measurement which lasts ' num2str(sample_data.meta.instrument_average_interval) ' seconds.']; ... 'HEIGHT_ABOVE_SENSOR', height(:), 'Data has been vertically bin-mapped using tilt information so that the cells have consistant heights above sensor in time.'; ... 'DIST_ALONG_BEAMS', distance(:), 'Data is not vertically bin-mapped (no tilt correction applied). Cells are lying parallel to the beams, at heights above sensor that vary with tilt.' }; @@ -288,6 +291,9 @@ end clear dims; + % add information about the middle of the measurement period + sample_data.dimensions{1}.seconds_to_middle_of_measurement = sample_data.meta.instrument_average_interval/2; + % add variables with their dimensions and data mapped if isMagBias magExt = ''; @@ -415,14 +421,25 @@ sample_data{2}.meta.user = []; sample_data{2}.meta.instrument_sample_interval = median(diff(waveData.param.time*24*3600)); + avgInterval = []; + if isfield(waveData, 'summary') + iMatch = ~cellfun(@isempty, regexp(waveData.summary, 'Each Burst Contains [0-9]* Samples, Taken at [0-9\.]* Hz.')); + if any(iMatch) + avgInterval = textscan(waveData.summary{iMatch}, 'Each Burst Contains %f Samples, Taken at %f Hz.'); + avgInterval = avgInterval{1}/avgInterval{2}; + end + end + sample_data{2}.meta.instrument_average_interval = avgInterval; + if isempty(avgInterval), avgInterval = '?'; end + sample_data{2}.dimensions = {}; sample_data{2}.variables = {}; % add dimensions with their data mapped dims = { - 'TIME', waveData.param.time; ... - 'FREQUENCY', waveData.Dspec.freq; ... - ['DIR' magExt], waveData.Dspec.dir + 'TIME', waveData.param.time, ['Time stamp corresponds to the start of the measurement which lasts ' num2str(avgInterval) ' seconds.']; ... + 'FREQUENCY', waveData.Dspec.freq, ''; ... + ['DIR' magExt], waveData.Dspec.dir, '' }; nDims = size(dims, 1); @@ -438,6 +455,9 @@ end clear dims; + % add information about the middle of the measurement period + sample_data{2}.dimensions{1}.seconds_to_middle_of_measurement = sample_data{2}.meta.instrument_average_interval/2; + % add variables with their dimensions and data mapped vars = { 'TIMESERIES', [], 1; ... diff --git a/Preprocessing/depthPP.m b/Preprocessing/depthPP.m index cf2b3a3fd..f6fdb26e8 100644 --- a/Preprocessing/depthPP.m +++ b/Preprocessing/depthPP.m @@ -137,10 +137,10 @@ % samples without pressure information are excluded if (presCurIdx == 0 && presRelCurIdx == 0), continue; end - if isSensorHeight - samSensorZ = sam.instrument_nominal_height; - else + if isSensorTargetDepth samSensorZ = sam.instrument_nominal_depth; + else + samSensorZ = sam.instrument_nominal_height; end % current sample or samples without vertical nominal @@ -199,10 +199,12 @@ iFirst = []; iSecond = []; for l = 1:m - if isSensorHeight - diffWithOthers(l) = curSam.instrument_nominal_height - otherSam{l}.instrument_nominal_height; - else + if isSensorTargetDepth diffWithOthers(l) = curSam.instrument_nominal_depth - otherSam{l}.instrument_nominal_depth; + else + % below is reversed so that sign convention is the + % same + diffWithOthers(l) = otherSam{l}.instrument_nominal_height - curSam.instrument_nominal_height; end end @@ -229,7 +231,9 @@ % the calculated depth could be too far off the truth distMin = 10; while distance < distMin && ~all(isnan(newDiffWithOthers)) - iNextBelow = find(diffWithOthers == max(newDiffWithOthers(newDiffWithOthers < 0)), 1); + iNextBelow = diffWithOthers == max(newDiffWithOthers(newDiffWithOthers < 0)); + iNextBelow(isnan(newDiffWithOthers)) = 0; % deals with the case of same depth instrument previously found + iNextBelow = find(iNextBelow, 1, 'first'); distance = abs(diffWithOthers(iNextBelow) - diffWithOthers(iBelow)); if distance >= distMin iSecond = iNextBelow; @@ -239,6 +243,27 @@ end elseif isempty(iBelow) && ~isempty(iAbove) iFirst = iAbove; + + % extending reseach to further nearest above didn't + % lead to better results + +% % let's find the second nearest above +% newDiffWithOthers = diffWithOthers; +% newDiffWithOthers(iFirst) = NaN; +% distance = 0; +% +% % if those two sensors are too close to each other then +% % the calculated depth could be too far off the truth +% distMin = 10; +% while distance < distMin && ~all(isnan(newDiffWithOthers)) +% iNextAbove = find(diffWithOthers == min(newDiffWithOthers(newDiffWithOthers > 0)), 1); +% distance = abs(diffWithOthers(iNextAbove) - diffWithOthers(iAbove)); +% if distance >= distMin +% iSecond = iNextAbove; +% break; +% end +% newDiffWithOthers(iNextAbove) = NaN; +% end else iFirst = iAbove; iSecond = iBelow; @@ -248,7 +273,7 @@ fprintf('%s\n', ['Warning : ' curSam.toolbox_input_file ... ' computing actual depth from only one pressure sensor '... 'on mooring']); - % we found only one sensor (it is above) + % we found only one sensor otherSam = otherSam{iFirst}; presIdxOther = getVar(otherSam.variables, 'PRES'); presRelIdxOther = getVar(otherSam.variables, 'PRES_REL'); @@ -285,10 +310,19 @@ % vertical axis is positive down when talking about % depth % - if isSensorHeight - distOtherCurSensor = otherSam.instrument_nominal_height - curSam.instrument_nominal_height; - else + if isSensorTargetDepth distOtherCurSensor = curSam.instrument_nominal_depth - otherSam.instrument_nominal_depth; + signOtherCurSensor = sign(distOtherCurSensor); + + distOtherCurSensor = abs(distOtherCurSensor); + else + distOtherCurSensor = otherSam.instrument_nominal_height - curSam.instrument_nominal_height; + signOtherCurSensor = -sign(distOtherCurSensor); + % 0 => two sensors at the same depth + % 1 => current sensor is deeper than other sensor + % -1 => current sensor is lower than other sensor + + distOtherCurSensor = abs(distOtherCurSensor); end if ~isempty(curSam.geospatial_lat_min) && ~isempty(curSam.geospatial_lat_max) @@ -332,7 +366,7 @@ zOther = interp1(tOther, zOther, tCur); clear tOther tCur; - computedDepth = zOther + distOtherCurSensor; + computedDepth = zOther + signOtherCurSensor*distOtherCurSensor; clear zOther; else samFirst = otherSam{iFirst}; @@ -382,12 +416,12 @@ % compute pressure at current sensor using trigonometry and % assuming sensors repartition on a line between the two % nearest pressure sensors - if isSensorHeight + if isSensorTargetDepth + distFirstSecond = samSecond.instrument_nominal_depth - samFirst.instrument_nominal_depth; + distFirstCurSensor = curSam.instrument_nominal_depth - samFirst.instrument_nominal_depth; + else distFirstSecond = samFirst.instrument_nominal_height - samSecond.instrument_nominal_height; distFirstCurSensor = samFirst.instrument_nominal_height - curSam.instrument_nominal_height; - else - distFirstSecond = samFirst.instrument_nominal_depth - samSecond.instrument_nominal_depth; - distFirstCurSensor = samFirst.instrument_nominal_depth - curSam.instrument_nominal_depth; end % theta is the angle between the vertical and line diff --git a/imosToolbox.m b/imosToolbox.m index 9404f8b3d..adfccdd47 100644 --- a/imosToolbox.m +++ b/imosToolbox.m @@ -73,7 +73,7 @@ function imosToolbox(auto, varargin) end % Set current toolbox version -toolboxVersion = ['2.5.6 - ' computer]; +toolboxVersion = ['2.5.7 - ' computer]; switch auto case 'auto', autoIMOSToolbox(toolboxVersion, varargin{:}); diff --git a/imosToolbox_Linux64.bin b/imosToolbox_Linux64.bin index bac0fc3d4..e5e3e7d87 100755 Binary files a/imosToolbox_Linux64.bin and b/imosToolbox_Linux64.bin differ diff --git a/imosToolbox_Win32.exe b/imosToolbox_Win32.exe index 8701f239c..9fbc86f76 100644 Binary files a/imosToolbox_Win32.exe and b/imosToolbox_Win32.exe differ diff --git a/imosToolbox_Win64.exe b/imosToolbox_Win64.exe index fd75130a8..a18b04cea 100644 Binary files a/imosToolbox_Win64.exe and b/imosToolbox_Win64.exe differ