Mercurial > hg > lots

changeset 4:e23a23349e31

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Added inverse transform
author danieleb@code.soundsoftware.ac.uk
date Thu, 16 Jun 2011 13:16:29 +0100
parents ee2a86d7ec07
children 084000643ad1
files dcti.m idcti.m ilot.m invlappedorthotransform.m l1overl2.m lappedorthobasis.m lappedorthotransform.m lappedwindow.m lot.m
diffstat 9 files changed, 342 insertions(+), 227 deletions(-) [+]
line wrap: on
line diff
--- a/dcti.m	Wed Jun 15 11:26:03 2011 +0100
+++ b/dcti.m	Thu Jun 16 13:16:29 2011 +0100
@@ -14,14 +14,14 @@
 % c: vector or matrix contaning the dct coefficients
 %
 % REFERENCES
-%   1)
+% 1)
 %
 % SEE ALSO
 % DCT, FFT, IDCT, IFFT
 
-%   Author(s): Daniele Barchiesi
-%   Copyright QMUL
-%   $Revision: 1
+% Author(s): Daniele Barchiesi
+% Copyright QMUL
+% $Revision: 1
 
 
 %% Check inputs & defaults
--- a/idcti.m	Wed Jun 15 11:26:03 2011 +0100
+++ b/idcti.m	Thu Jun 16 13:16:29 2011 +0100
@@ -14,14 +14,14 @@
 % c: vector or matrix contaning the inverse dct coefficients
 %
 % REFERENCES
-%   1)
+% 1)
 %
 % SEE ALSO
 % DCTI, DCT, IDCT
 
-%   Author(s): Daniele Barchiesi
-%   Copyright QMUL
-%   $Revision: 1
+% Author(s): Daniele Barchiesi
+% Copyright QMUL
+% $Revision: 1
 
 
 %% Check inputs & defaults
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/ilot.m	Thu Jun 16 13:16:29 2011 +0100
@@ -0,0 +1,131 @@
+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);
--- a/invlappedorthotransform.m	Wed Jun 15 11:26:03 2011 +0100
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,59 +0,0 @@
-function x = invlappedorthotransform(y,wLength,orthoBasis,tailLength,tailFunc)
-%>INVLAPPEDORTHOTRANSFORM Inverse lapped orthogonal transform (iLOT)
-%>
-%> x = lappedorthotransform(y,wLength,orthoBasis,tailLength,tailFunc)
-%>
-%> returns the inverse lapped orthogonal transform of the signal y 
-%>
-%> INPUT:
-%> - y: vector containing the LOT
-%> - wLength: vector containing the lengths of the overlapping windows (the
-%> sum of wLength must equal the length of x)
-%> - orthobasis: either a string or a function handle corresponding to an
-%> orhogonal transform to be used in each overlapping window. The default
-%> is DCT-IV (which is globally referred as MDCT)
-%> - tailLength: 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: either a string or a function handle used to define the tails
-%> of the overlapping windows. The default is 'sin2' (twicely differentiable
-%> sinusoidal).
-%>
-%> OUTPUT:
-%> - x: inverse LOT
-%>
-%> REFERENCE:
-%> S. Mallat, A Wavelet Tour of Signal Processing
-%>
-%> SEE ALSO
-%> lappedorthobasis, lappedwindow
-%>
-%> ----------------------------------------------------------- %
-%>  Daniele Barchiesi, daniele.barchiesi@eecs.qmul.ac.uk       %
-%>  Centre for Digital Music, Queen Mary University of London  %
-%>  May 2011                                                  %
-%> ----------------------------------------------------------- %
-%%> Check input & defaults
-error(nargchk(2, 5, nargin, 'struct'));
-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
-if length(x)<sum(wLength), error('you must pad the signal!'); end
-
-if ischar(orthoBasis)
-    switch orthoBasis
-        case 'mdct'
-            orthoFun = @(x) dcti(x,'IV');
-        otherwise
-            error('invalid basis function');
-    end
-elseif isa(orthoBasis,'function_handle')
-    orthoFun = orthoBasis;
-else
-    error('invalid basis function');
-end
-
-%%> Compute transform
--- a/l1overl2.m	Wed Jun 15 11:26:03 2011 +0100
+++ b/l1overl2.m	Thu Jun 16 13:16:29 2011 +0100
@@ -1,33 +1,27 @@
 function z = l1overl2(X)
-%L1OVERL2 ratio between the l1 and the l2 norms.
+%L1OVERL2 ratio between l1 and l2 norm
 %
-% z = l1overl2(X)
+% USAGE
+% z = l1overl2(x)
 %
-% computes the ratio between the l1 and the l2 norm of the input X. If
-% input is a matrix, the norms are calculated along the columns and a mean
-% value is returned. If the input is a cell, the norms are calculated
-% independently for each cell entry and a mean value is returned.
+% INPUT
+% x: input vector or matrix. If matrix, the function acts columnwise
+% returning one value for each column.
 %
-% ----------------------------------------------------------- %
-%  Daniele Barchiesi, daniele.barchiesi@eecs.qmul.ac.uk       %
-%  Centre for Digital Music, Queen Mary University of London  %
-%  Apr. 2011                                                  %
-% ----------------------------------------------------------- %
-if iscell(X)
-    nCols = length(X);
-    z = zeros(nCols,1);
-    for i=1:nCols
-        z(i) = norm(X{:,i},1)/norm(X{:,i},2);
-    end
-else
-    nCols = size(X,2);
-    z = zeros(nCols,1);
-    for i=1:nCols
-        z(i) = norm(X(:,i),1)/norm(X(:,i),2);
-    end
+% OUTPUT
+% z: the ratio ||x||_1/||x||_2
+%
+% SEE ALSO
+% NORM
+
+% Author(s): Daniele Barchiesi
+% Copyright QMUL
+% $Revision: 1
+
+nCols = size(X,2);
+z = zeros(nCols,1);
+for i=1:nCols
+    z(i) = norm(X(:,i),1)/norm(X(:,i),2);
 end
 
-z(isnan(z)) = [];
-z = mean(z);
-
 %EOF
\ No newline at end of file
--- a/lappedorthobasis.m	Wed Jun 15 11:26:03 2011 +0100
+++ b/lappedorthobasis.m	Thu Jun 16 13:16:29 2011 +0100
@@ -60,41 +60,56 @@
 suppLen   = @(p) wLength(p)+2*tailLength;
 
 % first frame
+% generate window
 g = lappedwindow(wLength(1),tailLength,tailFunc,'first');
+% generate local orthogonal matrix
 orthoMat = basisFun(eye(wLength(1)));
+% periodically extend length of basis functions to overlapping regions
 orthoMatExtended = orthomatextend(orthoMat,tailLength,'first');
+% multiply local orthogonal matrix by windowing function
 phi(1:wLength(1)+tailLength,1:wLength(1)) = diag(sparse(g))*orthoMatExtended;
 
 % central frames
 for p=2:nWindows-1
+    % generate window
     g = lappedwindow(wLength(p),tailLength,tailFunc);
+    % generate local orthogonal matrix
     orthoMat = basisFun(eye(wLength(p)));
+    % periodically extend length of basis functions to overlapping regions
     orthoMatExtended = orthomatextend(orthoMat,tailLength);
+    % multiply local orthogonal matrix by windowing function
     phi(sum(wLength(1:p-1))-tailLength+(1:suppLen(p)),sum(wLength(1:p-1))+(1:wLength(p)))...
         = diag(sparse(g))*orthoMatExtended;
 end
 
 %last frame
+% generate window
 g = lappedwindow(wLength(end),tailLength,tailFunc,'last');
-if strcmpi(orthobasis,'mdct')
+% generate local orthogonal matrix
+if strcmpi(orthobasis,'mdct') %if mdct use dct I at last block to avoid artefacts
     orthoMat = dcti(eye(wLength(end)),'I');
 else
     orthoMat = basisFun(eye(wLength(end)));
 end
+% periodically extend length of basis functions to overlapping regions
 orthoMatExtended = orthomatextend(orthoMat,tailLength,'last');
+% multiply local orthogonal matrix by windowing function
 phi(sum(wLength(1:end-1))-tailLength+(1:wLength(end)+tailLength),sum(wLength(1:end-1))+(1:wLength(end)))...
     = diag(sparse(g))*orthoMatExtended;
 
 function orthomatextended = orthomatextend(orthomat,tailLength,special)
 % Periodic extension of local orthonormal basis
-
-if ~exist('special','var') || isempty(special), special = []; end 
+if ~exist('special','var') || isempty(special), special = []; end
 
 if strcmpi(special,'first')
+    % odd symmetry at the end of basis functions
     orthomatextended = [orthomat; -orthomat(end:-1:end-tailLength+1,:)];
 elseif strcmpi(special,'last')
+    % even symmetry at the beginning of basis functions
     orthomatextended = [flipud(orthomat(1:tailLength,:)); orthomat];
 else
+    % even symmetry at the beginning and odd symmetry at the end of basis
+    % functions
     orthomatextended = [flipud(orthomat(1:tailLength,:)); orthomat; ...
         -orthomat(end:-1:end-tailLength+1,:)];
 end
\ No newline at end of file
--- a/lappedorthotransform.m	Wed Jun 15 11:26:03 2011 +0100
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,117 +0,0 @@
-function y = lappedorthotransform(x,wLength,orthoBasis,tailLength,tailFunc)
-%LAPPEDORTHOTRANSFORM Lapped orthogonal transform
-%
-% y = lappedorthotransform(x,wLength,orthoBasis,tailLength,tailFunc)
-%
-% returns the lapped orthogonal transform of the signal x 
-%
-% INPUT:
-% - x: vector containing the input signal 
-% - wLength: vector containing the lengths of the overlapping windows (the
-% sum of wLength must equal the length of x)
-% - orthobasis: either a string or a function handle corresponding to an
-% orhogonal transform to be used in each overlapping window. The default
-% is DCT-IV (which is globally referred as MDCT)
-% - tailLength: 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: either a string or a function handle used to define the tails
-% of the overlapping windows. The default is 'sin2' (twicely differentiable
-% sinusoidal).
-%
-% OUTPUT:
-% - y: vector of coefficients of the lapped orthogonal transform
-%
-% REFERENCE:
-% S. Mallat, A Wavelet Tour of Signal Processing
-%
-% SEE ALSO
-% lappedorthobasis, lappedwindow
-%
-% ----------------------------------------------------------- %
-%  Daniele Barchiesi, daniele.barchiesi@eecs.qmul.ac.uk       %
-%  Centre for Digital Music, Queen Mary University of London  %
-%  Apr. 2011                                                  %
-% ----------------------------------------------------------- %
-%% Check input & defaults
-error(nargchk(2, 5, nargin, 'struct'));
-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
-if isscalar(wLength)
-    nWindows = floor(length(x)/wLength);
-    wLength = wLength*ones(nWindows,1);
-end
-if length(x)<sum(wLength), error('you must pad the signal!'); end
-if ischar(orthoBasis)
-    switch orthoBasis
-        case 'mdct'
-            orthoFun = @(x) dcti(x,'IV');
-        otherwise
-            error('invalid basis function');
-    end
-elseif isa(orthoBasis,'function_handle')
-    orthoFun = orthoBasis;
-else
-    error('invalid basis function');
-end
-
-%% Compute transform
-nWindows  = length(wLength);
-sigLength = sum(wLength);
-y         = zeros(sigLength,1);
-suppLen   = @(p) wLength(p)+2*tailLength;
-
-% first frame
-h = lappedprojector(x(1:wLength(1)+tailLength),wLength(1),tailLength,tailFunc,'first');
-y(1:wLength(1)) = orthoFun(h(1:wLength(1)));
-
-% central frames
-for p=2:nWindows-1
-    h = lappedprojector(x(sum(wLength(1:p-1))-tailLength+(1:suppLen(p))),...
-        wLength(p),tailLength,tailFunc);
-    y(sum(wLength(1:p-1))+(1:wLength(p))) = orthoFun(h(tailLength+(1:wLength(p))));
-end
-    
-%last frame
-h = lappedprojector(x(sum(wLength(1:end-1))-tailLength+(1:wLength(end)+tailLength)),...
-    wLength(end),tailLength,tailFunc,'last');
-if strcmpi(orthoBasis,'mdct')
-    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 defined according to the
-% input parameteres
-
-% check inputs
-error(nargchk(2, 5, nargin, 'struct'));
-if ~exist('special','var') || isempty(special), special = []; end
-if ~exist('tailFunx','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')
-        O = wLength-tailLength+(1:2*tailLength);
-        h(O) = g(O).*f(O) - g(fliplr(O)).*f(fliplr(O));
-    elseif strcmpi(special,'last')
-        O = 1:2*tailLength;
-        h(O) = g(O).*f(O) + g(fliplr(O)).*f(fliplr(O));
-    end
-else
-    Os = 1:2*tailLength;
-    Oe = wLength+(1:2*tailLength);
-    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
\ No newline at end of file
--- a/lappedwindow.m	Wed Jun 15 11:26:03 2011 +0100
+++ b/lappedwindow.m	Thu Jun 16 13:16:29 2011 +0100
@@ -1,13 +1,13 @@
 function g = lappedwindow(wLength,tailLength,tailFunc,special)
-%LAPPEDWINDOW Lapped window used for lapped orthogonal transform
+%LAPPEDWINDOW Lapped window for lapped orthogonal transform
 %
+% USAGE
 % g = lappedwindow(wLength,tailLength,tailFunc,special)
 %
-% returns the window g used in the lapped orthogonal transform
-%
 % INPUT:
 % - wLength: length of the window
-% - tailLength: length of the tail of the overlapping window. 
+% - tailLength: (optional, default=floor(wLength/2))length of the tail of 
+% the overlapping window. 
 %
 %     ......................
 %   |.                      .|
@@ -18,9 +18,9 @@
 %                         tailLength
 %
 % The maximum tailLength cannot exceed half the window length. 
-% - tailFunc: either a string or a function handle used to define the tails
-% of the overlapping windows. The default is 'sin2' (twicely differentiable
-% sinusoidal).
+% - tailFunc: (optional, default='sin2') either a string or a function 
+% handle used to define the tails of the overlapping windows. The default 
+% is 'sin2' (twicely differentiable sinusoidal).
 %
 % OUTPUT:
 % - g: vector containing the window
@@ -31,16 +31,16 @@
 % SEE ALSO
 % lappedorthobasis, lappedorthotransform
 %
-% ----------------------------------------------------------- %
-%  Daniele Barchiesi, daniele.barchiesi@eecs.qmul.ac.uk       %
-%  Centre for Digital Music, Queen Mary University of London  %
-%  Apr. 2011                                                  %
-% ----------------------------------------------------------- %
+
+% Author(s): Daniele Barchiesi
+% Copyright QMUL
+% $Revision: 1
 
 %% Check inputs
-error(nargchk(2, 4, nargin, 'struct'));
+error(nargchk(1, 4, nargin, 'struct'));
 if ~exist('special','var') || isempty(special), special = []; end
-if ~exist('tailFunx','var') || isempty(tailFunc), tailFunc = 'sin2'; end
+if ~exist('tailFunc','var') || isempty(tailFunc), tailFunc = 'sin2'; end
+if ~exist('tailLength','var') || isempty(tailLength), tailLength = floor(wLength/2); end
 
 %check that the window is not too short
 if wLength<2*tailLength
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/lot.m	Thu Jun 16 13:16:29 2011 +0100
@@ -0,0 +1,151 @@
+function y = lot(x,wLength,orthoBasis,tailLength,tailFunc)
+%LAPPEDORTHOTRANSFORM Lapped orthogonal transform
+%
+% USAGE
+% 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
+% LAPPEDORTHOBASIS, LAPPEDWINDOW, INVLAPPEDORTHOTRANSFORM
+
+% 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(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
+if ischar(orthoBasis)
+    switch orthoBasis
+        case 'mdct'
+            orthoFun = @(x) dcti(x,'IV');
+        otherwise
+            error('invalid basis function');
+    end
+elseif isa(orthoBasis,'function_handle')
+    orthoFun = orthoBasis;
+else
+    error('invalid basis function');
+end
+
+%% Matrix case: act columnwise
+if size(x,2)>1
+    y = zeros(size(x));
+    for iCol=1:size(x,2)
+        y(:,iCol) = lappedorthotransform(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
\ No newline at end of file