function TelluricTemplate=spec_telluric_template(VecWave,Ratio,FlagTelluric,varargin)
%--------------------------------------------------------------------------
% spec_telluric_template function                                   ImSpec
% Description: Generate a Telluric template. Give a spectrum of the
%              ratio between the calibrated observed spectrum and
%              the standard star theoretical spectrum, and a vector of
%              flags indicating if a wavelength is a Telluric line,
%              return a template of the Telluric lines.
% Input  : - A column vector of wavelengths corresponding to the Ratio
%            vector (second input argument).
%          - A column vector of the ratio between the calibrated observed
%            spectrum and the standard star theoretical spectrum,
%            generated by spec_response.m (output: Tran.Ratio).
%          - A column vector of flags (true/false) indicating if a
%            wavelength is a Telluric line. This is available from
%            spec_response.m (output: Tran.FlagTelluric).
%          * Arbitrary number of pairs of ...,key,val,... arguments.
%            The following keywords are available:
%            'TiltMethod' - Method by which to tilt each Telluric region
%                         such that it edges will intersect with 1.
%                         Options are:
%                         'no'    - No tilting.
%                         'lin'   - Linear multiplicative tilting
%                                   (default).
%            'SizeOR'   - Number of elements on each side of a Telluric
%                         region to use in the linear fitting of the
%                         tilt. Default is 5.
%            'RemoveAboveTh' - If the mean of a Telluric region is above
%                         a threshold then remove it {true|false}.
%                         Default is true.
%            'RemoveTh' - Threshold for removal of Telluric region.
%                         Default is 1.
%            'AdjustAbove1' - Adjust the mean of a Telluric region to one
%                         if its mean is above one {true|false}.
%                         Default is true.
%            'MedFiltSize' - Run a median filter on the final output.
%                         This is the filter size (in "pixels").
%                         If empty or 0, do not run median filter.
%                         Default is 7.
% Tested : Matlab R2011b
%     By : Eran O. Ofek                    Mar 2014
%    URL : http://weizmann.ac.il/home/eofek/matlab/
% Example:
% TelluricTemplate=spec_telluric_template(VecWave,Ratio,FlagTelluric);
% Reliable: 2
%--------------------------------------------------------------------------
import Util.array.*

DefV.TiltMethod    = 'lin';
DefV.SizeOR        = 5;
DefV.RemoveAboveTh = true;
DefV.RemoveTh      = 1;
DefV.AdjustAbove1  = true;
DefV.MedFiltSize   = 7;
InPar = InArg.populate_keyval(DefV,varargin,mfilename);


Regions          = flag2regions(FlagTelluric);
Nreg             = size(Regions,1);
Nwave            = numel(Ratio);
TelluricTemplate = ones(Nwave,1);
for Ireg=1:1:Nreg
    IndReg = (Regions(Ireg,1):Regions(Ireg,2));
    
    % tilting
    switch lower(InPar.TiltMethod)
        case {'no','none'}
            % no tilt - use Ratio as is
            TelluricTemplate(IndReg) = Ratio(IndReg);
        case {'lin','linear'}
            % apply linear tilt to region according to its edge points
            IXY  = [ ((Regions(Ireg,1)-InPar.SizeOR-1):1:Regions(Ireg,1)-1), ( Regions(Ireg,2)+1:1:(Regions(Ireg,2)+InPar.SizeOR+1) ) ];
            IXY  = Util.array.index_outofbound(IXY,Nwave);
            Par  = polyfit(VecWave(IXY),Ratio(IXY),1);
            Line = polyval(Par,VecWave(IndReg));
            TelluricTemplate(IndReg) = Ratio(IndReg)./Line;
        otherwise
            error('Unknown TiltMethod option');
    end
    
    % if template mean height is above Threshold - remove
    if (InPar.RemoveAboveTh)
        if (mean(TelluricTemplate(IndReg))>InPar.RemoveTh)
            TelluricTemplate(IndReg) = 1;
        end
    end
    
    % adjust height if template region above 1
    if (InPar.AdjustAbove1)
        if (mean(TelluricTemplate(IndReg))>1)
            TelluricTemplate(IndReg) = TelluricTemplate(IndReg)./mean(TelluricTemplate(IndReg));
        end
    end
    
    % smooth template
    if (isempty(InPar.MedFiltSize) || InPar.MedFiltSize==0)
        % do nothing
    else
       TelluricTemplate = medfilt1(TelluricTemplate,InPar.MedFiltSize);
       TelluricTemplate(FlagTelluric==0) = 1;
    end
    
end

            
            
