function [error] = corTsgMedian(hMainFig, PARA, dateMin, dateMax) % % Correct the TSG salinity time series with the Water sample. % Use the median value of TIME_WINDOWS water sample to compute the % correction. see the documentation % % Input % hMainFig ..... Handle to the main GUI % PARA ......... Parameter (SSPS, SSJT, SSTP) % dateMin ...... the correction is applied between dateMin and date Max % dateMax ...... the correction is applied between dateMin and date Max % % Output % Error ........ 1 everything OK % ........ -1 dateMax <= date Min % % The error is maximum, equal to 1, if the median is computed with less % than 4 samples % The error minimum cannot be lower than 0.01 % Get application data % -------------------- tsg = getappdata( hMainFig, 'tsg_data'); SAMPLE = tsg.preference.sample; % Shorten the variable name % ------------------------- TIME_WINDOWS = tsg.cst.COR_TIME_WINDOWS; % ------------------------------------------------------------------------- % Get from the checkbox the QC code on which the correction will be applied % ------------------------------------------------------------------------- % get list of keys from hashtable tsg.qc.hash, defined inside % tsg_initialisation.m % ----------------------------------------------------------- qc_list = get(tsg.qc.hash); % iterate (loop) on each key store inside hastable % ------------------------------------------------ keptCode = []; nKeptCode = 0; for i=1:numel(qc_list) % get key and some values in hashtable % ------------------------------------ key = qc_list{i}; % get handle of checkbox % ---------------------- hCb = findobj(hMainFig, 'tag', ['TAG_CHECK_CORRECTION_' key]); if get( hCb, 'value' ) nKeptCode = nKeptCode + 1; keptCode(nKeptCode) = get(tsg.qc.hash, key, 'code'); end end % Get PROBABLY_GOOD, PROBABLY_BAD and VALUE_CHANGED codes % ------------------------------------------------------- PROBABLY_GOOD = get(tsg.qc.hash, 'PROBABLY_GOOD', 'code'); PROBABLY_BAD = get(tsg.qc.hash, 'PROBABLY_BAD', 'code'); BAD = get(tsg.qc.hash, 'BAD', 'code'); VALUE_CHANGED = get(tsg.qc.hash, 'VALUE_CHANGED', 'code'); % intialisation % ------------- if isempty( tsg.([PARA '_ADJUSTED']) ) tsg.([PARA '_ADJUSTED']) = tsg.(PARA); tsg.([PARA '_ADJUSTED_QC']) = tsg.([PARA '_QC']); end % If only calibration have been applied ERROR is empty % --------------------------------------------------- if isempty( tsg.([PARA '_ADJUSTED_ERROR']) ) tsg.([PARA '_ADJUSTED_ERROR']) = NaN* ones(size(tsg.(PARA))); end % Create a structure with an NaN % No other solution, as I can't add a structure to an empty one % ------------------------------------------------------------- cor = struct('DAYD', NaN, 'DIFF', NaN, 'ERROR', NaN, 'NVALUE', NaN); if dateMax > dateMin % Find the indices of samples within the time limits. % -------------------------------------------------- indSample = find(tsg.DAYD_EXT >= dateMin & tsg.DAYD_EXT <= dateMax); indCor = 0; for i = 1:length(indSample) % Find samples within TIME_WINDOWS with Good, probably Good, probably Bad QC % --------------------------------------------------------------- ind = find( tsg.DAYD_EXT >= (tsg.DAYD_EXT(indSample(i)) - TIME_WINDOWS/2) &... tsg.DAYD_EXT <= (tsg.DAYD_EXT(indSample(i)) + TIME_WINDOWS/2) &... tsg.([SAMPLE '_EXT_QC']) <= PROBABLY_GOOD); if ~isempty(ind) % detect NaN in sample.SSPS_DIF due to bad QC code in tsg.SSPS % ------------------------------------------------------------ ind2 = find(~isnan(tsg.EXT_DIF(ind))); % Compute the median difference and error within TIME_WINDOWS % ----------------------------------------------------------- if ~isempty(ind2) if ~isempty(tsg.EXT_DIF(ind(ind2))) A = tsg.EXT_DIF(ind(ind2)); meanA = mean(A); stdA = std(A); % Standard deviation test: keep these values % ------------------------------------------ ind3 = find( A >= meanA-3*stdA & A <= meanA+3*stdA); B = tsg.EXT_DIF(ind(ind2(ind3))); if ~isempty( B ) indCor = indCor + 1; cor.DAYD(indCor) = tsg.DAYD_EXT((indSample(i))); cor.DIFF(indCor) = median(B); cor.ERROR(indCor) = nanstd(B)/sqrt(length(B)); cor.NVALUE(indCor) = length(B); end % Standard deviation test: don't keep these values % ------------------------------------------------------ ind4 = find( A < meanA-3*stdA | A > meanA+3*stdA); if ~isempty( ind4 ) tsg.([SAMPLE '_EXT_QC'])(ind(ind2(ind4))) = BAD; end end end end end % Eliminate the first element if NaN % ---------------------------------- if isnan(cor.DAYD(1)) cor.DAYD(1) = []; cor.DIFF(1) = []; cor.ERROR(1) = []; cor.NVALUE(1) = []; end if ~isempty( cor.DAYD ) % The error is maximum, equal to 1, if the median is computed with less % than 4 samples % ----------------------------------------------------------------------- cor.ERROR( cor.NVALUE < 4 ) = 1; % The error minimum cannot be lower than 0.01 % ------------------------------------------- cor.ERROR( cor.ERROR < 0.01 ) = 0.01; % The correction is applied between dateMin and dateMax % We attribute to dateMin the first correction computed % and to dateMax the last one % % Find the tsg date in the interval dateMin-dateMax % ------------------------------------------------- dtTsg = find(tsg.DAYD >= dateMin & tsg.DAYD <= dateMax); if cor.DAYD(1) ~= dateMin cor.DAYD = [tsg.DAYD(dtTsg(1)) cor.DAYD]; cor.DIFF = [cor.DIFF(1) cor.DIFF]; cor.ERROR = [cor.ERROR(1) cor.ERROR]; cor.NVALUE = [cor.NVALUE(1) cor.NVALUE]; end if cor.DAYD(end) ~= dateMax cor.DAYD = [cor.DAYD tsg.DAYD(dtTsg(end))]; cor.DIFF = [cor.DIFF cor.DIFF(end)]; cor.ERROR = [cor.ERROR cor.ERROR(end)]; cor.NVALUE = [cor.NVALUE cor.NVALUE(end)]; end % Test if the date are strictly increasing otherwise interp1 % cannot be used % ---------------------------------------------------------- if ~isempty( find(diff(cor.DAYD) == 0) ) msgbox( {'Function corTsgMedian'; ' '; ... 'At least 2 samples have the same date';... 'The soft cannot make the interpolation'},... 'Correction Method', 'warn'); error = -1; return end % The correction is applied to the TSG between dateMin and dateMax using % a linear interpolation only on measurements with keptCode % ---------------------------------------------------------------------- for icode = 1 : length( keptCode ) dtTsg = find( tsg.DAYD >= dateMin & tsg.DAYD <= dateMax &... tsg.([PARA '_QC']) == keptCode( icode )); tsg.([PARA '_ADJUSTED'])(dtTsg) = tsg.(PARA)(dtTsg) + ... interp1(cor.DAYD, cor.DIFF, tsg.DAYD(dtTsg)); tsg.([PARA '_ADJUSTED_ERROR'])(dtTsg) = ... interp1(cor.DAYD, cor.ERROR, tsg.DAYD(dtTsg)); end end % Update tsg application data % --------------------------- setappdata( hMainFig, 'tsg_data', tsg); % everything OK % ------------- error = 1; else % DateMax <= DateMin % ------------------ msgbox( {'Function corTsgMedian'; ' '; ... 'Date limits are not correct' },... 'Correction Method', 'warn', 'modal'); error = -1; end