FINDPROFILES Identify individual profiles and compute vertical direction from depth sequence.
Syntax:
[PROFILE_INDEX, PROFILE_DIRECTION] = FINDPROFILES(STAMP, DEPTH)
[PROFILE_INDEX, PROFILE_DIRECTION] = FINDPROFILES(STAMP, DEPTH, OPTIONS)
[PROFILE_INDEX, PROFILE_DIRECTION] = FINDPROFILES(STAMP, DEPTH, OPT1, VAL1, ...)
[PROFILE_INDEX, PROFILE_DIRECTION] = FINDPROFILES(DEPTH, ...)
Description:
[PROFILE_INDEX, PROFILE_DIRECTION] = FINDPROFILES(DEPTH) and
[PROFILE_INDEX, PROFILE_DIRECTION] = FINDPROFILES(STAMP, DEPTH)
identify upcast and downcast profiles in depth (or pressure) vector DEPTH,
with optional timestamps in vector STAMP, and computes a vector of profile
indices PROFILE_INDEX and a vector of vertical direction PROFILE_DIRECTION.
STAMP, DEPTH, PROFILE_DIRECTION and PROFILE_INDEX are the same length N,
and if STAMP is not specified, it is assumed to be the sample index [1:N].
PROFILE_DIRECTION entries may be 1 (down), 0 (flat), -1 (up).
PROFILE_INDEX entries associate each sample with the number of the profile
it belongs to. Samples in the middle of a profile are flagged with a whole
number, starting at 1 and increased by 1 every time a new cast is detected,
while samples between profiles are flagged with an offset of 0.5.
See note on identification algorithm below.
[PROFILE_INDEX, PROFILE_DIRECTION] = FINDPROFILES(..., OPTIONS) and
[PROFILE_INDEX, PROFILE_DIRECTION] = FINDPROFILES(..., OPT1, VAL1) accept
the following options given in key-value pairs OPT1, VAL1... or in a struct
OPTIONS with field names as option keys and field values as option values:
STALL: maximum range of a stalled segment (in the same units as DEPTH).
Only intervals of constant vertical direction spanning a depth range
not less than the given value are considered valid cast segments.
Shorter intervals are considered stalled segments inside or between
casts.
Default value: 0 (all segments are valid cast segments)
SHAKE: maximum duration of a shake segment (in the same units as STAMP).
Only intervals of constant vertical direction with duration
not less than the given value are considered valid cast segments.
Briefer intervals are considered shake segments inside or between
casts.
Default value: 0 (all segments are valid cast segments)
INVERSION: maximum depth inversion between cast segments of a profile.
Consecutive valid cast segments with the same direction are joined
together in the same profile if the range of the introduced depth
inversion, if any, is less than the given value.
Default value: 0 (never join cast segments)
INTERRUPT: maximum time separation between cast segments of a profile.
Consecutive valid cast segments with the same direction are joined
together in the same profile if the duration of the lapse (sequence of
stalled segments or shakes between them) is less than the given value.
When STAMP is not specified, the duration will be the number of samples
between them.
Default value: 0 (never join cast segments)
LENGTH: minimum length of a profile.
A sequence of joined cast segments will be considered a valid profile
only if the total spanned depth is greater or equal than the given.
value.
Default value: 0 (all profiles are valid)
PERIOD: minimum duration of a profile.
A sequence of joined cast segments will be considered a valid profile
only if the total duration is greater or equal than the given value.
Default value: 0 (all profiles are valid)
The following options are deprecated and should not be used:
RANGE: minimum depth range.
Deprecated in v1.1.0:
Superseded by STALL.
JOIN: whether to join consecutive profiles with the same direction.
Deprecated in v1.1.0:
Superseded by INTERRUPT and INVERSION:
Equivalent to INVERSION = INF and INTERRUPT = INF, and
INVERSION = 0 and INTERRUPT = 0.
Notes:
Profiles are identified as sequences of cast segments with the same
vertical direction, allowing for stalled or shake segments in between.
Vertical segments are intervals of constant vertical direction,
and are delimited by the changes of vertical direction computed
as the sign of forward differences of the depth sequence.
A segment is considered stalled if it is to short in depth,
or a shake if it is to short in time. Otherwise it is a cast segment.
Consecutive cast segments with the same direction are joined together
if the introduced depth inversion and the lapse between the segments
are not significant according to the specified thresholds.
Invalid samples (NaN) in input are ignored. In output, they are marked as
belonging to the previous profile, and with the direction of the previous
sample.
Examples:
depth = [4 4 3 2 3 4 4 5 6 6 6 5 4 4 5 3 2 3 1 1 0 4 4]
[profile_index, profile_direction] = findProfiles(depth)
figure
subplot(3, 1, 1, 'XGrid','on','YGrid','on', 'NextPlot', 'add')
stairs(profile_direction, '-g')
subplot(3, 1, 2, 'XGrid','on','YGrid','on', 'NextPlot', 'add')
plot(depth, '-db')
subplot(3, 1, 3, 'XGrid','on','YGrid','on', 'NextPlot', 'add')
stairs(profile_index, '-r')
[profile_index, profile_direction] = findProfiles(depth, 'stall', 1.5)
stairs(profile_index, '-c')
[profile_index, profile_direction] = ...
findProfiles(depth, 'stall', 1.5, 'inversion', 1.5, 'interrupt', inf)
stairs(profile_index, '-m')
Authors:
Joan Pau Beltran <joanpau.beltran@socib.cat>0001 function [profile_index, profile_direction] = findProfiles(varargin) 0002 %FINDPROFILES Identify individual profiles and compute vertical direction from depth sequence. 0003 % 0004 % Syntax: 0005 % [PROFILE_INDEX, PROFILE_DIRECTION] = FINDPROFILES(STAMP, DEPTH) 0006 % [PROFILE_INDEX, PROFILE_DIRECTION] = FINDPROFILES(STAMP, DEPTH, OPTIONS) 0007 % [PROFILE_INDEX, PROFILE_DIRECTION] = FINDPROFILES(STAMP, DEPTH, OPT1, VAL1, ...) 0008 % [PROFILE_INDEX, PROFILE_DIRECTION] = FINDPROFILES(DEPTH, ...) 0009 % 0010 % Description: 0011 % [PROFILE_INDEX, PROFILE_DIRECTION] = FINDPROFILES(DEPTH) and 0012 % [PROFILE_INDEX, PROFILE_DIRECTION] = FINDPROFILES(STAMP, DEPTH) 0013 % identify upcast and downcast profiles in depth (or pressure) vector DEPTH, 0014 % with optional timestamps in vector STAMP, and computes a vector of profile 0015 % indices PROFILE_INDEX and a vector of vertical direction PROFILE_DIRECTION. 0016 % STAMP, DEPTH, PROFILE_DIRECTION and PROFILE_INDEX are the same length N, 0017 % and if STAMP is not specified, it is assumed to be the sample index [1:N]. 0018 % PROFILE_DIRECTION entries may be 1 (down), 0 (flat), -1 (up). 0019 % PROFILE_INDEX entries associate each sample with the number of the profile 0020 % it belongs to. Samples in the middle of a profile are flagged with a whole 0021 % number, starting at 1 and increased by 1 every time a new cast is detected, 0022 % while samples between profiles are flagged with an offset of 0.5. 0023 % See note on identification algorithm below. 0024 % 0025 % [PROFILE_INDEX, PROFILE_DIRECTION] = FINDPROFILES(..., OPTIONS) and 0026 % [PROFILE_INDEX, PROFILE_DIRECTION] = FINDPROFILES(..., OPT1, VAL1) accept 0027 % the following options given in key-value pairs OPT1, VAL1... or in a struct 0028 % OPTIONS with field names as option keys and field values as option values: 0029 % STALL: maximum range of a stalled segment (in the same units as DEPTH). 0030 % Only intervals of constant vertical direction spanning a depth range 0031 % not less than the given value are considered valid cast segments. 0032 % Shorter intervals are considered stalled segments inside or between 0033 % casts. 0034 % Default value: 0 (all segments are valid cast segments) 0035 % SHAKE: maximum duration of a shake segment (in the same units as STAMP). 0036 % Only intervals of constant vertical direction with duration 0037 % not less than the given value are considered valid cast segments. 0038 % Briefer intervals are considered shake segments inside or between 0039 % casts. 0040 % Default value: 0 (all segments are valid cast segments) 0041 % INVERSION: maximum depth inversion between cast segments of a profile. 0042 % Consecutive valid cast segments with the same direction are joined 0043 % together in the same profile if the range of the introduced depth 0044 % inversion, if any, is less than the given value. 0045 % Default value: 0 (never join cast segments) 0046 % INTERRUPT: maximum time separation between cast segments of a profile. 0047 % Consecutive valid cast segments with the same direction are joined 0048 % together in the same profile if the duration of the lapse (sequence of 0049 % stalled segments or shakes between them) is less than the given value. 0050 % When STAMP is not specified, the duration will be the number of samples 0051 % between them. 0052 % Default value: 0 (never join cast segments) 0053 % LENGTH: minimum length of a profile. 0054 % A sequence of joined cast segments will be considered a valid profile 0055 % only if the total spanned depth is greater or equal than the given. 0056 % value. 0057 % Default value: 0 (all profiles are valid) 0058 % PERIOD: minimum duration of a profile. 0059 % A sequence of joined cast segments will be considered a valid profile 0060 % only if the total duration is greater or equal than the given value. 0061 % Default value: 0 (all profiles are valid) 0062 % 0063 % The following options are deprecated and should not be used: 0064 % RANGE: minimum depth range. 0065 % Deprecated in v1.1.0: 0066 % Superseded by STALL. 0067 % JOIN: whether to join consecutive profiles with the same direction. 0068 % Deprecated in v1.1.0: 0069 % Superseded by INTERRUPT and INVERSION: 0070 % Equivalent to INVERSION = INF and INTERRUPT = INF, and 0071 % INVERSION = 0 and INTERRUPT = 0. 0072 % 0073 % Notes: 0074 % Profiles are identified as sequences of cast segments with the same 0075 % vertical direction, allowing for stalled or shake segments in between. 0076 % Vertical segments are intervals of constant vertical direction, 0077 % and are delimited by the changes of vertical direction computed 0078 % as the sign of forward differences of the depth sequence. 0079 % A segment is considered stalled if it is to short in depth, 0080 % or a shake if it is to short in time. Otherwise it is a cast segment. 0081 % Consecutive cast segments with the same direction are joined together 0082 % if the introduced depth inversion and the lapse between the segments 0083 % are not significant according to the specified thresholds. 0084 % 0085 % Invalid samples (NaN) in input are ignored. In output, they are marked as 0086 % belonging to the previous profile, and with the direction of the previous 0087 % sample. 0088 % 0089 % Examples: 0090 % depth = [4 4 3 2 3 4 4 5 6 6 6 5 4 4 5 3 2 3 1 1 0 4 4] 0091 % [profile_index, profile_direction] = findProfiles(depth) 0092 % figure 0093 % subplot(3, 1, 1, 'XGrid','on','YGrid','on', 'NextPlot', 'add') 0094 % stairs(profile_direction, '-g') 0095 % subplot(3, 1, 2, 'XGrid','on','YGrid','on', 'NextPlot', 'add') 0096 % plot(depth, '-db') 0097 % subplot(3, 1, 3, 'XGrid','on','YGrid','on', 'NextPlot', 'add') 0098 % stairs(profile_index, '-r') 0099 % [profile_index, profile_direction] = findProfiles(depth, 'stall', 1.5) 0100 % stairs(profile_index, '-c') 0101 % [profile_index, profile_direction] = ... 0102 % findProfiles(depth, 'stall', 1.5, 'inversion', 1.5, 'interrupt', inf) 0103 % stairs(profile_index, '-m') 0104 % 0105 % Authors: 0106 % Joan Pau Beltran <joanpau.beltran@socib.cat> 0107 0108 % Copyright (C) 2013-2016 0109 % ICTS SOCIB - Servei d'observacio i prediccio costaner de les Illes Balears 0110 % <http://www.socib.es> 0111 % 0112 % This program is free software: you can redistribute it and/or modify 0113 % it under the terms of the GNU General Public License as published by 0114 % the Free Software Foundation, either version 3 of the License, or 0115 % (at your option) any later version. 0116 % 0117 % This program is distributed in the hope that it will be useful, 0118 % but WITHOUT ANY WARRANTY; without even the implied warranty of 0119 % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 0120 % GNU General Public License for more details. 0121 % 0122 % You should have received a copy of the GNU General Public License 0123 % along with this program. If not, see <http://www.gnu.org/licenses/>. 0124 0125 error(nargchk(1, 14, nargin, 'struct')); 0126 0127 0128 %% Set options and default values. 0129 options.length = 0; 0130 options.period = 0; 0131 options.inversion = 0; 0132 options.interrupt = 0; 0133 options.stall = 0; 0134 options.shake = 0; 0135 0136 options.range = []; % Deprecated since v1.1.0. 0137 options.join = []; % Deprecated since v1.1.0. 0138 0139 0140 %% Parse basic input arguments. 0141 % Get numeric (non-option) arguments. 0142 nargnum = find(~cellfun(@isnumeric, varargin), 1, 'first') - 1; 0143 if isempty(nargnum) 0144 nargnum = nargin; 0145 end 0146 switch(nargnum) 0147 case 1 0148 depth = varargin{1}; 0149 stamp = reshape(1:numel(depth), size(depth)); 0150 case 2 0151 [stamp, depth] = varargin{1:2}; 0152 end 0153 0154 0155 %% Parse optional arguments. 0156 % Get option key-value pairs in any accepted call signature. 0157 argopts = varargin(nargnum+1:end); 0158 if isscalar(argopts) && isstruct(argopts{1}) 0159 % Options passed as a single option struct argument: 0160 % field names are option keys and field values are option values. 0161 opt_key_list = fieldnames(argopts{1}); 0162 opt_val_list = struct2cell(argopts{1}); 0163 elseif mod(numel(argopts), 2) == 0 0164 % Options passed as key-value argument pairs. 0165 opt_key_list = argopts(1:2:end); 0166 opt_val_list = argopts(2:2:end); 0167 else 0168 error('glider_toolbox:findProfiles:InvalidOptions', ... 0169 'Invalid optional arguments (neither key-value pairs nor struct).'); 0170 end 0171 % Overwrite default options with values given in extra arguments. 0172 for opt_idx = 1:numel(opt_key_list) 0173 opt = lower(opt_key_list{opt_idx}); 0174 val = opt_val_list{opt_idx}; 0175 if isfield(options, opt) 0176 options.(opt) = val; 0177 else 0178 error('glider_toolbox:findProfiles:InvalidOption', ... 0179 'Invalid option: %s.', opt); 0180 end 0181 end 0182 0183 0184 %% Handle deprecated options. 0185 if ~isempty(options.range) 0186 warning('glider_toolbox:findProfiles:DeprecatedOption', ... 0187 'Deprecated option: range. See option: stall.'); 0188 options.stall = options.range; 0189 end 0190 if ~isempty(options.join) 0191 warning('glider_toolbox:findProfiles:DeprecatedOption', ... 0192 'Deprecated option: join. See option: shake.'); 0193 options.interrupt = options.join * inf; 0194 options.inversion = options.join * inf; 0195 end 0196 0197 0198 %% Identify the profiles. 0199 valid_index = find(~(isnan(depth(:)) | isnan(stamp(:)))); 0200 sdy = sign(diff(depth(valid_index))); 0201 depth_peak = true(size(valid_index)); 0202 depth_peak(2:end-1) = diff(sdy) ~= 0; 0203 depth_peak_index = valid_index(depth_peak); 0204 sgmt_frst = stamp(depth_peak_index(1:end-1)); 0205 sgmt_last = stamp(depth_peak_index(2:end)); 0206 sgmt_strt = depth(depth_peak_index(1:end-1)); 0207 sgmt_fnsh = depth(depth_peak_index(2:end)); 0208 sgmt_sinc = sgmt_last - sgmt_frst; 0209 sgmt_vinc = sgmt_fnsh - sgmt_strt; 0210 sgmt_vdir = sign(sgmt_vinc); 0211 cast_sgmt_valid = ... 0212 ~(abs(sgmt_vinc(:)) <= options.stall | sgmt_sinc(:) <= options.shake); 0213 cast_sgmt_index = find(cast_sgmt_valid); 0214 cast_sgmt_lapse = ... 0215 (sgmt_frst(cast_sgmt_index(2:end)) - sgmt_last(cast_sgmt_index(1:end-1))); 0216 cast_sgmt_space = -sgmt_vdir(cast_sgmt_index(1:end-1)) .* ... 0217 (sgmt_strt(cast_sgmt_index(2:end)) - sgmt_fnsh(cast_sgmt_index(1:end-1))); 0218 cast_sgmt_dirch = diff(sgmt_vdir(cast_sgmt_index)); 0219 cast_sgmt_bound = ~(cast_sgmt_dirch(:) == 0 & ... 0220 cast_sgmt_lapse(:) <= options.interrupt & ... 0221 cast_sgmt_space(:) <= options.inversion); 0222 cast_sgmt_head_valid = true(size(cast_sgmt_index)); 0223 cast_sgmt_tail_valid = true(size(cast_sgmt_index)); 0224 cast_sgmt_head_valid(2:end) = cast_sgmt_bound; 0225 cast_sgmt_tail_valid(1:end-1) = cast_sgmt_bound; 0226 cast_head_index = depth_peak_index(cast_sgmt_index(cast_sgmt_head_valid)); 0227 cast_tail_index = depth_peak_index(cast_sgmt_index(cast_sgmt_tail_valid) + 1); 0228 cast_length = abs(depth(cast_tail_index) - depth(cast_head_index)); 0229 cast_period = stamp(cast_tail_index) - stamp(cast_head_index); 0230 cast_valid = ... 0231 ~(cast_length(:) <= options.length | cast_period(:) <= options.period); 0232 cast_head = zeros(size(depth)); 0233 cast_tail = zeros(size(depth)); 0234 cast_head(cast_head_index(cast_valid) + 1) = 0.5; 0235 cast_tail(cast_tail_index(cast_valid)) = 0.5; 0236 profile_index = 0.5 + cumsum(cast_head + cast_tail); 0237 profile_direction = nan(size(depth)); 0238 for i = 1:numel(valid_index)-1 0239 profile_direction(valid_index(i):valid_index(i+1)-1) = sdy(i); 0240 end 0241 0242 end