corTsgLinear.m

function [error] = corTsgLinear(hMainFig, PARA, dateMin, dateMax)
% Correct the TSG salinity time series with the Water sample.
% Use a linear fit to the water sample/tsg difference
% 
% Input
% hMainFig ..... Handle to the main GUI
% PARA ......... Parameter
% 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 correction is computed with
% less than 4 samples
%
% The error minimum cannot be lower than 0.01
%
% 2009/03/23 Y. Gouriou
% 1) Compute linear fit witha minimum of 3 points.
% 2) If only 2 points: no linear fit, but compute the average difference
%    with the 2 points.

% Get application data
% --------------------
tsg    = getappdata( hMainFig, 'tsg_data');
SAMPLE = tsg.plot.sample;

% -------------------------------------------------------------------------
% 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 = keys(tsg.qc.hash);

% TODO: define size of keptCode
% -----------------------------
%keptCode = zeros(numel(qc_list), 1);

% iterate (loop) on each key store inside hastable
% ------------------------------------------------
keptCode = [];

nKeptCode = 0;
for key = qc_list

  % get handle of checkbox
  % ----------------------
  hCb = findobj(hMainFig, 'tag', ['TAG_CHECK_CORRECTION_' char(key)]);
    
  if get( hCb, 'value' )
    nKeptCode = nKeptCode + 1;
    keptCode(nKeptCode) = tsg.qc.hash.(key).code;
  end
end

% Get PROBABLY_GOOD, PROBABLY_BAD and VALUE_CHANGED codes
% -------------------------------------------------------
PROBABLY_GOOD = tsg.qc.hash.PROBABLY_GOOD.code;
PROBABLY_BAD  = tsg.qc.hash.PROBABLY_BAD.code;
VALUE_CHANGED = tsg.qc.hash.VALUE_CHANGED.code;

% Intialisation
% 01/09/2009 : intialisation to NaN
% ---------------------------------
if isempty( tsg.([PARA '_ADJUSTED']) )  
  tsg.([PARA '_ADJUSTED'])       = NaN*ones(size(tsg.(PARA)));
  tsg.([PARA '_ADJUSTED_QC'])    = NaN*ones(size(tsg.([PARA '_QC'])));
  tsg.([PARA '_ADJUSTED_ERROR']) = NaN*ones(size(tsg.(PARA))); 
end

if dateMax > dateMin

  % Find samples within TIME_WINDOWS with Good, probably Good, QC
  % -------------------------------------------------------------
  ind = find( tsg.DAYD_EXT    >= dateMin &  tsg.DAYD_EXT    <= dateMax &...
    tsg.([SAMPLE '_EXT_QC']) <= PROBABLY_GOOD);

  if ~isempty(ind)

    % detect NaN in sample.SSPS_DIF due to bad QC code for tsg.SSPS
    % -------------------------------------------------------------
    ind2 = find(~isnan(tsg.EXT_DIF(ind)));

    % Compute linear fit of the TSG/SAMPLE difference
    % -----------------------------------------------
    if ~isempty(ind2) && length(ind2) > 2
      %       if ~isempty(tsg.EXT_DIF(ind(ind2)))    Is this line useful?

      % Linear fit applied to the difference tsg-sample
      % -----------------------------------------------
      X = tsg.DAYD_EXT(ind(ind2));
      Y = tsg.EXT_DIF(ind(ind2));

      [p, S, mu] = polyfit( X, Y, 1);

      % The correction is applied to the TSG between dateMin and dateMax using
      % a linear interpolation only on measurements with keptCode Quality Codes
      % ------------------------------------------------------------------------
      for icode = 1 : length( keptCode )
        dtTsg = find( tsg.DAYD    >= dateMin  & tsg.DAYD <= dateMax &...
          tsg.([PARA '_QC']) == keptCode( icode ));

        if ~isempty( dtTsg )

          % Compute the correction + the error
          % ----------------------------------
          [tsg.([PARA '_ADJUSTED'])(dtTsg),...
            tsg.([PARA '_ADJUSTED_ERROR'])(dtTsg)] =...
                                            polyval( p, tsg.DAYD(dtTsg), S, mu);

          % Compute the corrected value : orignal value + correction
          % --------------------------------------------------------
          tsg.([PARA '_ADJUSTED'])(dtTsg) =...
                            tsg.(PARA)(dtTsg) + tsg.([PARA '_ADJUSTED'])(dtTsg);

          % Transfer the QC
          % ---------------
          tsg.([PARA '_ADJUSTED_QC'])(dtTsg) = tsg.([PARA '_QC'])(dtTsg);
        end
      end

      % If there is only 2 points, compute the mean difference. Apply an error max.
      % --------------------------------------------------------------------------
    elseif ~isempty(ind2) && length(ind2) == 2

      meanDif = mean( tsg.EXT_DIF(ind(ind2)) );

      % The correction is applied to the TSG between dateMin and dateMax using
      % a linear interpolation only on measurements with keptCode Quality Codes
      % ------------------------------------------------------------------------
      for icode = 1 : length( keptCode )
        dtTsg = find( tsg.DAYD    >= dateMin  & tsg.DAYD <= dateMax &...
                      tsg.([PARA '_QC']) == keptCode( icode ));

        if ~isempty( dtTsg )
          tsg.([PARA '_ADJUSTED'])(dtTsg)       = tsg.(PARA)(dtTsg) + meanDif;
          tsg.([PARA '_ADJUSTED_ERROR'])(dtTsg) = meanDif / 2;
          tsg.([PARA '_ADJUSTED_QC'])(dtTsg)    = tsg.([PARA '_QC'])(dtTsg);
        end
      end

    end

    % The error minimum cannot be lower than 0.01 or equal to Inf
    % -----------------------------------------------------------
    if ~isempty(ind2) && length(ind2) >= 2
      tsg.([PARA '_ADJUSTED_ERROR'])...
                               ( tsg.([PARA '_ADJUSTED_ERROR']) < 0.01 ) = 0.01;

      tsg.([PARA '_ADJUSTED_ERROR'])...
                            (isinf(tsg.([PARA '_ADJUSTED_ERROR'])) == 1) = 0.01;
    end
    %       end
  end

  % Update tsg application data
  % ---------------------------
  setappdata( hMainFig, 'tsg_data', tsg);

  % everything OK
  % -------------
  error = 1;

else

  % DateMax <= DateMin
  % ------------------
  error = -1;

end