wilorth.m

linear time－frequency toolbox

function [gamma]=wilorth(p1,p2,p3)
%WILORTH  Wilson orthonormal window.
%   Usage:   gamma=wilorth(M,L);
%            gamma=wilorth(g,M);
%            gamma=wilorth(g,M,L);
%
%   Input parameters:
%         g     : Auxiliary window window function (optional).
%         M     : Number of modulations.
%         L     : Length of window (optional).
%   Output parameters:
%         gamma : Window generating an orthonormal Wilson basis.
%
%   gamma=WILORTH(M,L) computes a nice window of length L generating a
%   Wilson or MDCT orthonormal basis with M frequency bands for signal of
%   length L.
%
%   gamma=WILORTH(g,M) creates an orthonormal window from the window g, 
%   g must be whole-point even.
%
%   gamma=WILORTH(g,M,L) pads or truncates g to length L before
%   calculating the orthonormal window. gamma will also be of length L.
%
%   SEE ALSO:  DWILT, WILDUAL, ISEVEN

wasrow=0;

% Detect which parameters was entered, and do simple transformations.
if nargin == 3
  g=p1;
  M=p1;
  L=p2;

  if (rem(L,M)~=0)
    error('L must be dividable by M.');
  end;

  if size(g,2)>1
    if size(g,1)>1
      error('g must be a vector');
    else
      % g was a row vector.
      wasrow=1;
      g=g(:);
    end;
  end;

% This check has been removed!
%  if ~iseven(g)
%    error('g must be whole point even around middle point.See help on pgauss.');
%  end;

  % firextend is now safe.
  g=firextend(g,L);
  
else

  if nargin==2

    [x,y]=size(p1);
    if x*y>1
      % First parameter is a vector.
      g=p1;
      M=p2;

      if size(g,2)>1
	if size(g,1)>1
	  error('g must be a vector');
	else
	  % g was a row vector.
	  wasrow=1;
	  g=g(:);
	end;
      end;

      L=length(g);

      if rem(L,M)~=0
	error('The length of the window must be dividable by M.');
      end;

      % This check has been removed.
      %  if ~iseven(g)
      %    error('g must be whole point even around middle point.See help on pgauss.');
      %  end;	  

    else
      M=p1;
      L=p2;
      a=M;
      b=L/(2*M);

      if rem(L,M)~=0
	error('L must be dividable by M.');
      end;

      % Create default window, a Gaussian.
      g=comp_pgauss(L,a/b,0);

    end;
  else

    error('Incorrect number of calling parameters.')

  end;

end;

if size(M,1)>1 || size(M,2)>1
  error('M must be a scalar');
end;

a=M;
b=L/(2*M);

% Multiply by sqrt(2), because comp_cantight will return a normalized
% tight frame, i.e. the framebounds are A=B=1 instead of A=B=2. This means
% that the returned gamma only has norm=.701 and not norm=1.

gf=comp_wfac(g,a,2*M);
gammaf=comp_cantight_fac(gf,L,a,2*M);
gamma=sqrt(2)*comp_iwfac(gammaf,L,a,2*M);

if wasrow
  gamma=gamma.';
end;

% Wilson windows are always real.
gamma=real(gamma);