CORRECTSENSORLAG Correct data sequence from sensor lag effects. Syntax: COR = CORRECTSENSORLAG(TIME, RAW, PARAMS) COR = CORRECTSENSORLAG(TIME, RAW, FLOW, PARAMS) Description: COR = CORRECTSENSORLAG(TIME, RAW, PARAMS) corrects the sequence in vector RAW with timestamp in vector TIME and constant flow speed from sensor lag effects. The correction is done advancing the signal in time using linear interpolation. PARAMS should be a scalar with the number of time units the signal should be shifted. COR = CORRECTSENSORLAG(TIME, RAW, FLOW, PARAMS) performs the same correction but for a sequence with variable flow speed, given by FLOW. FLOW should be a vector with the same dimensions as TIME and RAW. PARAMS should be a two element vector with the offset and the slope of the time lag parameter with respect to the inverse flow speed. Notes: This function is a recoding of a previous function by Tomeu Garau, called CORRECTTIMERESPONSE. He is the true glider man. Main changes are: - Different correction using shift by interpolation instead of gradient. - Support for dynamic time lag parameters for sequences with variable flow speed. Invalid values (NaN) in input sequence are ignored but preserved in output. The first order approximation was done using forward differences with the function DIFF. Examples: % Example sequence: time = 0:25 tru = exp(-0.5*(time-10).^2/4^2) % Correct a sequence with a lag of 5 time units and constant flow speed: tau = 5 raw = exp(-0.5*(time-10-tau).^2/4^2) raw([5 6 10 19 20 21]) = nan cor1 = correctSensorLag(time, raw, 3) cor2 = correctSensorLag(time, raw, tau) subplot(2, 1, 1) plot(time, tru, 'k--', time, raw, 'k', time, cor1, 'b', time, cor2, 'r') legend('true', ... sprintf('raw (delay params %.2f)', tau), ... sprintf('corrected 1 (shift params %.2f)', 3), ... sprintf('corrected 2 (shift params %.2f)', tau)) title('constant flow speed') % Correct a sequence with a lag of 5 time units and variable flow speed: flow = 0.4 + 0.1 * randn(size(time)) tau_offset = 1.5 tau_scale = 0.2 tau = tau_offset + tau_scale ./ flow raw = exp(-0.5*((time-10-tau).^2)/4^2) raw([5 6 10 19 20 21]) = nan cor1 = correctSensorLag(time, raw, flow, [1.0 0.1]) cor2 = correctSensorLag(time, raw, flow, [tau_offset tau_scale]) subplot(2, 1, 2) plot(time, tru, 'k--', time, raw, 'k', time, cor1, 'b', time, cor2, 'r') legend('true', ... sprintf('raw (delay params %.2f %.2f)', tau_offset, tau_scale), ... sprintf('corrected 1 (shift params %.2f %.2f)', 1.0, 0.1), ... sprintf('corrected 2 (shift params %.2f %.2f)', tau_offset, tau_scale)) title('variable flow speed') See also: FINDSENSORLAGPARAMS INTERP1 Authors: Joan Pau Beltran <joanpau.beltran@socib.cat>
0001 function cor = correctSensorLag(varargin) 0002 %CORRECTSENSORLAG Correct data sequence from sensor lag effects. 0003 % 0004 % Syntax: 0005 % COR = CORRECTSENSORLAG(TIME, RAW, PARAMS) 0006 % COR = CORRECTSENSORLAG(TIME, RAW, FLOW, PARAMS) 0007 % 0008 % Description: 0009 % COR = CORRECTSENSORLAG(TIME, RAW, PARAMS) corrects the sequence in vector 0010 % RAW with timestamp in vector TIME and constant flow speed from sensor lag 0011 % effects. The correction is done advancing the signal in time using linear 0012 % interpolation. PARAMS should be a scalar with the number of time units the 0013 % signal should be shifted. 0014 % 0015 % COR = CORRECTSENSORLAG(TIME, RAW, FLOW, PARAMS) performs the same 0016 % correction but for a sequence with variable flow speed, given by FLOW. 0017 % FLOW should be a vector with the same dimensions as TIME and RAW. 0018 % PARAMS should be a two element vector with the offset and the slope 0019 % of the time lag parameter with respect to the inverse flow speed. 0020 % 0021 % Notes: 0022 % This function is a recoding of a previous function by Tomeu Garau, called 0023 % CORRECTTIMERESPONSE. He is the true glider man. Main changes are: 0024 % - Different correction using shift by interpolation instead of gradient. 0025 % - Support for dynamic time lag parameters for sequences with variable 0026 % flow speed. 0027 % 0028 % Invalid values (NaN) in input sequence are ignored but preserved in output. 0029 % 0030 % The first order approximation was done using forward differences with the 0031 % function DIFF. 0032 % 0033 % Examples: 0034 % % Example sequence: 0035 % time = 0:25 0036 % tru = exp(-0.5*(time-10).^2/4^2) 0037 % 0038 % % Correct a sequence with a lag of 5 time units and constant flow speed: 0039 % tau = 5 0040 % raw = exp(-0.5*(time-10-tau).^2/4^2) 0041 % raw([5 6 10 19 20 21]) = nan 0042 % cor1 = correctSensorLag(time, raw, 3) 0043 % cor2 = correctSensorLag(time, raw, tau) 0044 % subplot(2, 1, 1) 0045 % plot(time, tru, 'k--', time, raw, 'k', time, cor1, 'b', time, cor2, 'r') 0046 % legend('true', ... 0047 % sprintf('raw (delay params %.2f)', tau), ... 0048 % sprintf('corrected 1 (shift params %.2f)', 3), ... 0049 % sprintf('corrected 2 (shift params %.2f)', tau)) 0050 % title('constant flow speed') 0051 % 0052 % % Correct a sequence with a lag of 5 time units and variable flow speed: 0053 % flow = 0.4 + 0.1 * randn(size(time)) 0054 % tau_offset = 1.5 0055 % tau_scale = 0.2 0056 % tau = tau_offset + tau_scale ./ flow 0057 % raw = exp(-0.5*((time-10-tau).^2)/4^2) 0058 % raw([5 6 10 19 20 21]) = nan 0059 % cor1 = correctSensorLag(time, raw, flow, [1.0 0.1]) 0060 % cor2 = correctSensorLag(time, raw, flow, [tau_offset tau_scale]) 0061 % subplot(2, 1, 2) 0062 % plot(time, tru, 'k--', time, raw, 'k', time, cor1, 'b', time, cor2, 'r') 0063 % legend('true', ... 0064 % sprintf('raw (delay params %.2f %.2f)', tau_offset, tau_scale), ... 0065 % sprintf('corrected 1 (shift params %.2f %.2f)', 1.0, 0.1), ... 0066 % sprintf('corrected 2 (shift params %.2f %.2f)', tau_offset, tau_scale)) 0067 % title('variable flow speed') 0068 % 0069 % See also: 0070 % FINDSENSORLAGPARAMS 0071 % INTERP1 0072 % 0073 % Authors: 0074 % Joan Pau Beltran <joanpau.beltran@socib.cat> 0075 0076 % Copyright (C) 2013-2016 0077 % ICTS SOCIB - Servei d'observacio i prediccio costaner de les Illes Balears 0078 % <http://www.socib.es> 0079 % 0080 % This program is free software: you can redistribute it and/or modify 0081 % it under the terms of the GNU General Public License as published by 0082 % the Free Software Foundation, either version 3 of the License, or 0083 % (at your option) any later version. 0084 % 0085 % This program is distributed in the hope that it will be useful, 0086 % but WITHOUT ANY WARRANTY; without even the implied warranty of 0087 % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 0088 % GNU General Public License for more details. 0089 % 0090 % You should have received a copy of the GNU General Public License 0091 % along with this program. If not, see <http://www.gnu.org/licenses/>. 0092 0093 error(nargchk(3, 4, nargin, 'struct')); 0094 0095 % Parse input arguments. 0096 switch(nargin) 0097 case 3 0098 % Constant flow speed. 0099 [timestamp, raw, params] = varargin{:}; 0100 constant_flow = true; 0101 case 4 0102 % Variable flow speed. 0103 [timestamp, raw, flow, params] = varargin{:}; 0104 constant_flow = false; 0105 end 0106 0107 if constant_flow 0108 % Positive time check to filter out bad initial lines on Slocum data. 0109 valid = (timestamp > 0) & ~isnan(raw); 0110 % Time lag parameter is a constant scalar. 0111 tau = params; 0112 else 0113 % Positive time check to filter out bad initial lines on Slocum data. 0114 valid = (timestamp > 0) & ~(isnan(raw) | isnan(flow)); 0115 % Compute dynamic time lag parameter inversely proportional to flow speed. 0116 tau_offset = params(1); 0117 tau_slope = params(2); 0118 tau = tau_offset + tau_slope ./ flow(valid); 0119 end 0120 0121 cor = nan(size(raw)); 0122 timestamp_valid = timestamp(valid); 0123 raw_valid = raw(valid); 0124 [timestamp_unique, index_from, index_to] = unique(timestamp_valid); 0125 raw_unique = raw_valid(index_from); 0126 if any(raw_valid ~= raw_unique(index_to)) 0127 error('glider_toolbox:correctSensorLag:InconsistentData', ... 0128 'Inconsistent sensor data.'); 0129 end 0130 if numel(timestamp_unique) > 1 0131 % del = [0; diff(raw_valid) ./ diff(timestamp_valid)]; 0132 % cor(valid) = raw_val + tau .* del; 0133 cor(valid) = interp1(timestamp_unique, raw_unique, ... 0134 timestamp_valid + tau, 'linear', 'extrap'); 0135 end 0136 0137 end