view ilot.m @ 14:2e42f5fb764d

Updated documentation, deleted lotbasis (use ilot(eye) instead)

function x = ilot(y,wLength,orthoBasis,tailLength,tailFunc)
% ILOT Inverse lapped orthogonal transform
%
% Example
% x = ilot(y,wLength,orthoBasis,tailLength,tailFunc) returns the inverse
% lapped orthogonal transform of the signal y
%
% 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 lot

% Author(s): Daniele Barchiesi
% Copyright 2011-2011

%% 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
invOrthoFun = orthohandle(orthoBasis,'inv');


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

%% Compute transform
nWindows  = length(wLength);
sigLength = sum(wLength);
x         = zeros(sigLength,1);
p = 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)
if strcmpi(orthoBasis,'mdct') && p==nWindows-1 % if mdct use idct I
    hp = idcti(y(wLength(1)+(1:wLength(2))),'I');
else
    hp = invOrthoFun(y(wLength(1)+(1:wLength(2))));
end
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:p))+(1:wLength(p+1))),'I');
    else
        hp = invOrthoFun(y(sum(wLength(1:p))+(1:wLength(p+1))));
    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);