function [ wcs_out ] = insertMaskedWaveletCoeffs( wcs_in, wcv, level, index )
%insertWaveletCoeffs Insert a vector of wavelet coefficients in to the WCS
%at level 'level'
%   [ wcs ] = insertWaveletCoeffs( wcs, wcv, level ) inserts a vector of
%   wavelet coefficients in wcv into a Wacelet Coefficient Structure (WCS).
%   The wavelet coefficients in wcv may be those governing trends in x-, y-
%   or xy-directions (index = 1, 2, 3), or they may contain wavelet
%   coefficients in x, y and xy directions, in the same sequence that
%   extractWaveletCoeffs() does things.
% 
% INPUT:
%  wcs_in     : The WCS into which the vector of wavelet coefficients has to
%            be inserted.
%
%  wcv     : The vector of wavelet coefficients which are to be inserted.
%            On level 'l', there are a total of (2^l X 2^l) - 
%            (2^{l-1} X 2^{l-1}) wavelet coefficients
%
%  level   : The wavelet coefficients at at level 'level'
%
%  index   : The wavelet coefficients correspond to those that determine
%            variation in the 'index' direction. index = {1, 2, 3} = {x-,
%            y-, xy-} directions.  If you specify an index, we expect that
%            the length of the wcv vector is  [ (2^l X 2^l) - 
%            (2^{l-1} X 2^{l-1})] / 3.
%
% OUTPUT :
%  wcs_out      : The wavelet coefficient structure (WCS) with the wavelet
%             coefficients from wcv inserted into it at the right level.
%
% See also :
%
%     extractMaskedWaveletCoeffs
%
% Jaideep Ray, SNL, 2/25/211
% =========================================================================

% ---- make a copy
wcs_out = wcs_in ;
wcs_out.waveletType    = wcs_in.waveletType ;
wcs_out.waveletOrder   = wcs_in.waveletOrder ;
wcs_out.qmf            = wcs_in.qmf ;
wcs_out.origDataSize   = wcs_in.origDataSize ;
wcs_out.paddedDataSize = wcs_in.paddedDataSize ;
wcs_out.nlevels        = wcs_in.nlevels ;
wcs_out.baseVal        = wcs_in.baseVal ;
wcs_out.mask           = wcs_in.mask ;
wcs_out.maskBaseVal    = wcs_in.maskBaseVal ;
wcs_out.coeffs         = cell(wcs_out.nlevels, 3) ;
wcs_out.maskCoeffs     = cell(wcs_out.nlevels, 3) ;
wcs_out.maskNonZeroIndices = cell(wcs_out.nlevels, 3) ;

for l = 1 : wcs_in.nlevels
    for index = 1 : 3
        wcs_out.coeffs{l, index} = wcs_in.coeffs{l, index} ;
        wcs_out.maskCoeffs{l, index} = wcs_in.maskCoeffs{l, index} ;
        wcs_out.maskNonZeroIndices{l, index} = wcs_in.maskNonZeroIndices{l, index} ;
    end
end

% if (nargin == 3)
%     % Print out what came in
%     fprintf('insertmaskedWaveletCoeff() Level = %d, len(wcv) = %d, norm(wcv) = %f\n',...
%             level, numel(wcv), norm(wcv)) ;
% end

% ---- What is the size of the incoming vector?
if (nargin == 4) % somebody gave me an index
    reqdSize = numel(wcs_out.maskNonZeroIndices{level, index}) ;
else
    if (level == 0)
        reqdSize = 1 ; % just the baseVal for level 0
    else
        reqdSize = 0 ;    
        for i = 1:3
            reqdSize = reqdSize + numel(wcs_in.maskNonZeroIndices{level, i}) ;
        end
    end
end

if ( numel(wcv) ~= reqdSize )
    error('insertWaveletCoeffs(): Incoming wavelet vector is of the wrong size'); 
end

% ---- If the index is specified, simply stick this into the WCS
if (level == 0)
    wcs_out.baseVal = wcv(1) * wcs_in.maskBaseVal; % wcv has only 1 element if level == 0
    return ;
end
% fprintf('insertmaskedWaveletCoeff() wcs_out.baseVal = %f\n', wcs_out.baseVal) ;

% ---- If the index is specified, reshape the vector to 2D arrays and copy
%      into the right place in the WCS
if (nargin == 4)
    if ( level <=  wcs_out.nlevels )
        tmp = wcs_out.coeffs{level, index} ;
        tmp(wcs_out.maskNonZeroIndices{level, index}) = wcv ;
        wcs_out.coeffs{level, index} = tmp ;
    end
    return ;
end

% ---- No index specified, so the vector has all the wavelet coeffs on this
% level.
lb = 1 ;
for index = 1 : 3
    if ( level <= wcs_out.nlevels )
        numElements = numel(wcs_out.maskNonZeroIndices{level, index}) ;
        ub = lb + numElements -1 ;   % ub of the index in wcv related to this set of coeffs
        
        tmp = wcs_out.coeffs{level, index} ;
        tmp(wcs_out.maskNonZeroIndices{level, index}) = wcv(lb:ub) ;
        wcs_out.coeffs{level, index} = tmp ;
        lb = ub + 1 ; % update the indices for wcv
    end
end

end

