function x = ilot(y,wLength,orthoBasis,tailLength,tailFunc)
%INVLAPPEDORTHOTRANSFORM Inverse lapped orthogonal transform
%
% USAGE
% x = lappedorthotransform(y,wLength,orthoBasis,tailLength,tailFunc)
%
% INPUT:
% y: vector or matrix containing the LOT coefficients. If matrix,
% the function acts columnwise.
%
% wLength: scalar or vector containing the lengths of the overlapping
% windows (if vector, the sum of wLength must equal the length of y)
%
% orthobasis: (optional, default='mdct') either a string or a function
% handle corresponding to an orhogonal transform to be used in each
% overlapping window.
%
% tailLength: (optional, default=floor(min(wLength)/2)) length of the tail
% of the overlapping windows. Two consecutive windows overlap in a interval
% of dimension 2*tailLength. The maximum tailLength cannot exceed half the
% length of the smallest window.
%
% tailFunc: (optional, default='sin2') either a string or a function handle
% used to define the tails of the overlapping windows (see LAPPEDWINDOW for
% more details).
%
% OUTPUT:
% x: vector or matrix of coefficients of the inverse LOT.
%
% REFERENCES:
% S. Mallat, A Wavelet Tour of Signal Processing
%
% SEE ALSO
% LAPPEDORTHOBASIS, LAPPEDWINDOW, LAPPEDORTHOTRANSFORM

% Author(s): Daniele Barchiesi
% Copyright QMUL
% $Revision: 1

%% Check input & defaults
error(nargchk(2, 5, nargin, 'struct'));

if isscalar(wLength) %create vector of fixed window lengths
    nWindows = floor(length(y)/wLength);
    wLength = wLength*ones(nWindows,1);
end

if ~exist('orthoBasis','var') || isempty(orthoBasis), orthoBasis = 'mdct'; end
if ~exist('tailFunc','var')   || isempty(tailFunc), tailFunc = 'sin2'; end
if ~exist('tailLength','var') || isempty(tailLength)
    tailLength = floor(min(wLength)/2);
end

% check wLength
if length(y)<sum(wLength), error('you must pad the signal!'); end
% define function handle to the forward local transform
if ischar(orthoBasis)
    switch orthoBasis
        case 'mdct'
            invOrthoFun = @(y) idcti(y,'IV');
        otherwise
            error('invalid basis function');
    end
elseif isa(orthoBasis,'function_handle')
    invOrthoFun = orthoBasis;
else
    error('invalid basis function');
end

%% Matrix case: act columnwise
if size(y,2)>1
    x = zeros(size(y));
    for iCol=1:size(y,2)
        x(:,iCol) = invlappedorthotransform(y(:,iCol),wLength,orthoBasis,...
            tailLength,tailFunc);
    end
    return
end

%% Compute transform
nWindows  = length(wLength);
sigLength = sum(wLength);
x         = zeros(sigLength,1);

% first frame
g = lappedwindow(wLength(1),tailLength,tailFunc,'first');
h = invOrthoFun(y(1:wLength(1))); % folded ilot
C = 1:wLength(1)-tailLength;      % central interval
Om = wLength(1)-tailLength+(1:tailLength);  %overlapping interval (end)
hp = invOrthoFun(y(wLength(1)+(1:wLength(2)))); % next folded ilot
x(C) = h(C);
% unfold
x(Om) = g(Om).*h(Om) + g(end:-1:end-tailLength+1).*hp(tailLength:-1:1);

% central frames
for p=2:nWindows-1
    hm = h; % previous folded ilot
    h = hp; % current folded ilot
    % next folded ilot
    if strcmpi(orthoBasis,'mdct') && p==nWindows-1 % if mdct use idct I
        hp = idcti(y(sum(wLength(1:end-1))+(1:wLength(end))),'I');
    else
        hp = invOrthoFun(y(sum(wLength(1:end-1))+(1:wLength(end))));
    end
    g = lappedwindow(wLength(p),tailLength,tailFunc);
    % overlapping interval (start)
    Op = sum(wLength(1:p-1))+(1:tailLength);
    % central interval
    C  = sum(wLength(1:p-1))+(tailLength+1:wLength(p)-tailLength);
    % overlapping interval (end)
    Om = sum(wLength(1:p-1))+ wLength(p)-tailLength+(1:tailLength);
    % unfold
    x(Op) = g(tailLength+(1:tailLength)).*h(1:tailLength) - ...
        g(tailLength:-1:1).*hm(end:-1:end-tailLength+1);
    x(C) = h(tailLength+1:wLength(p)-tailLength);
    x(Om) = g(wLength(p)+(1:tailLength)).*h(wLength(p)-tailLength+(1:tailLength)) + ...
        g(end:-1:end-tailLength+1).*hp(tailLength:-1:1);
end

%last frame
hm = h; % previous folded ilot
h = hp; % current folded ilot
g = lappedwindow(wLength(end),tailLength,tailFunc);
% overlapping interval (start)
Op = sum(wLength(1:end-1))+(1:tailLength);
% central interval
C  = sum(wLength(1:end-1))+(tailLength+1:wLength(end));
% unfold
x(Op) = g(tailLength+(1:tailLength)).*h(1:tailLength) - ...
        g(tailLength:-1:1).*hm(end:-1:end-tailLength+1);
x(C) = h(tailLength+1:end);
