view lot.m @ 17:0da5eb32e49d tip

minor changes

function y = lot(x,wLength,orthoBasis,tailLength,tailFunc)
% LOT lapped orthogonal transform
%
% Example
% y = lappedorthotransform(x,wLength,orthoBasis,tailLength,tailFunc)
%
% Input
% - x: vector or matrix containing the input signal. 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 x)
% - 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 portion
% 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
% - y: vector or matrix of coefficients of the lapped orthogonal transform.
%
% References:
% S. Mallat, A Wavelet Tour of Signal Processing
%
% See also LAPPEDWINDOW, ilot

% 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(x)/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(x)<sum(wLength), error('you must pad the signal!'); end
% define function handle to the forward local transform
orthoFun = orthohandle(orthoBasis,'for');

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

%% Compute transform
nWindows  = length(wLength);
sigLength = sum(wLength);
y         = zeros(sigLength,1);
% length of the support of the p-th window
suppLen   = @(p) wLength(p)+2*tailLength;

% first frame
% project input signal onto orthogonal subspace associated with the first
% window
h = lappedprojector(x(1:wLength(1)+tailLength),wLength(1),tailLength,tailFunc,'first');
% transform projected signal 
y(1:wLength(1)) = orthoFun(h(1:wLength(1)));

% central frames
for p=2:nWindows-1
    % project input signal onto orthogonal subspace associated with the
    % p-th window
    h = lappedprojector(x(sum(wLength(1:p-1))-tailLength+(1:suppLen(p))),...
        wLength(p),tailLength,tailFunc);
    % transform projected signal
    y(sum(wLength(1:p-1))+(1:wLength(p))) = orthoFun(h(tailLength+(1:wLength(p))));
end
    
%last frame
% project input signal onto orthogonal subspace associated with the
% last window
h = lappedprojector(x(sum(wLength(1:end-1))-tailLength+(1:wLength(end)+tailLength)),...
    wLength(end),tailLength,tailFunc,'last');
% transform projected signal
if strcmpi(orthoBasis,'mdct') %id mdct use DCT I at last frame to avoid artefacts
    y(sum(wLength(1:end-1))+(1:wLength(end))) = dcti(h(tailLength+(1:wLength(end))),'I');
else
    y(sum(wLength(1:end-1))+(1:wLength(end))) = orthoFun(h(tailLength+(1:wLength(end))));
end

function h = lappedprojector(f,wLength,tailLength,tailFunc,special)
% Projects the input f into the space of functions with support contained
% in the interval defined by the lapped window.

% check inputs
error(nargchk(2, 5, nargin, 'struct'));
if ~exist('special','var') || isempty(special), special = []; end
if ~exist('tailFunc','var') || isempty(tailFunc), tailFunc = 'sin2'; end
if ~exist('tailLength','var') || isempty(tailLength)
    tailLength = floor(min(wLength)/2);
end

% generate the window and project
g = lappedwindow(wLength,tailLength,tailFunc,special);
h = f;
if ~isempty(special)
    if strcmpi(special,'first')
        % overlapping region
        O = wLength-tailLength+(1:2*tailLength);
        % project
        h(O) = g(O).*f(O) - g(fliplr(O)).*f(fliplr(O));
    elseif strcmpi(special,'last')
        % overlapping portion
        O = 1:2*tailLength;
        % project
        h(O) = g(O).*f(O) + g(fliplr(O)).*f(fliplr(O));
    end
else
    % start overlapping region
    Os = 1:2*tailLength;
    % end overlapping region
    Oe = wLength+(1:2*tailLength);
    % project
    h(Os) = g(Os).*f(Os) + g(fliplr(Os)).*f(fliplr(Os));
    h(Oe) = g(Oe).*f(Oe) - g(fliplr(Oe)).*f(fliplr(Oe));
end