function [smooth] = tsg_average(hMainFig, PARA, iTsg)
% Perform the average of a time series at the position of a WS sample
%
% The average is made for a period equal to 'tsg.cst.TSG_DT_SMOOTH'
% Data exceeding the average over that period by 'tsg.cst.TSG_STDMAX'
% are not taken into account.
%
% Input
% hMainFig ........ Handle to the main GUI : TSGQC
% PARA ............ Cell array
%                   PARA{1} contains the characters (SSP, SSJT, SSTP)
%                   PARA{2} contains either the cahracters (SSPS, SSJT, SSTP)
%                           or (SSPS_CAL, SSJT_CAL, SSTP_CAL)
% iTsg ............ Indice of the TSG measurement the closest to the sample
%                   taking into account the TSG quality code
%
% No computation : NaN
%

% Get the tsg structure from the application
% ------------------------------------------
tsg = getappdata( hMainFig, 'tsg_data');

% Get PROBABLY_GOOD code
% ----------------------
probablyGoodCode = tsg.qc.hash.PROBABLY_GOOD.code;

% Select the param data over 'tsg.cst.TSG_DT_SMOOTH' time interval
% taking into account the TSG quality code
% ind1 : indices of tsg.PARA in the 'tsg.cst.TSG_DT_SMOOTH' time interval
% ind2 : indices of tsg.PARA not rejected by the S.D. test
% -----------------------------------------------------------------------
ind1 = find( tsg.DAYD >= tsg.DAYD(iTsg) - tsg.cst.TSG_DT_SMOOTH/2 & ...
    tsg.DAYD <= tsg.DAYD(iTsg) + tsg.cst.TSG_DT_SMOOTH/2 & ...
    tsg.([PARA{1} '_QC']) <= probablyGoodCode);

ind2 = ind1;
if ~isempty(ind2)

  currentStd  = Inf;
  previousStd = 0;

  % Compare Standard Deviation to the MAX acceptable STD
  % ----------------------------------------------------
  while currentStd > tsg.([PARA{1} '_STDMAX']) && currentStd ~= previousStd

    previousStd = currentStd;

    % Standard deviation and average over timeInterval
    % ------------------------------------------------
    currentStd = nanstd(  tsg.(PARA{2})(ind2) );
    meanParam  = nanmean( tsg.(PARA{2})(ind2) );

    % Indices of 'good' values of Para
    % ---------------------------------
    ind2 = ind1( tsg.(PARA{2})(ind1) >= meanParam - currentStd & ...
                 tsg.(PARA{2})(ind1) <= meanParam + currentStd );
  end

  smooth = nanmean( tsg.(PARA{2})(ind2) );

else
  smooth = NaN;
end