candual.m

Matlab时频分析工具箱,希望能对大家有所帮助啊

function gamma=candual(g,a,M,L)
%CANDUAL  Canonical dual window.
%   Usage:  gamma=candual(g,a,M);
%           gamma=candual(g,a,M,L);
%
%   Input parameters:
%         g     : Gabor window.
%         a     : Length of time shift.
%         M     : Number of channels.
%         L     : Length of window. (optional)
%   Output parameters:
%         gamma : Canonical dual window.
%
%   CANDUAL(g,a,M) computes the canonical dual window of the discrete Gabor
%   frame with window g and parameters a, M. It is assumed that g is an IIR
%   window, and the length of the returned tight window is the same as the
%   length of g.
%
%   If the length of g is equal to M, then the input window is assumed to
%   be a FIR window. In this case, the canonical dual window also has
%   length of g.
%
%   CANDUAL(g,a,M,L) pads or truncates g to a length L IIR window before the
%   canonical dual window is calculated. gamma will also have length L.
%
%   If a>M then the dual window of the Gabor Riesz sequence with window g
%   and parameters a and M will be calculated.
%
%   SEE ALSO:  CANTIGHT, FIR2IIR, DGT

% This program is free software: you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
% 
% This program is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
% GNU General Public License for more details.
% 
% You should have received a copy of the GNU General Public License
% along with this program.  If not, see <http://www.gnu.org/licenses/>.

error(nargchk(3,4,nargin));

if nargin<4
  L=[];
end;

wasrow=0;
if isnumeric(g)
  if size(g,2)>1
    if size(g,1)>1
      error('g must be a vector');
    else
      % g was a row vector.
      g=g(:);
      wasrow=1;
    end;
  end;
  Lwindow=size(g,1);
else
  if isempty(L)
    error(['You must specify a length L if a window is represented as a ' ...
           'text string or cell array.']);
  end;
  Lwindow=0;    
end;

assert_squarelat(a,M,1,'CANDUAL',0);

dofirdual=0;

[b,N,L]=assert_L(Lwindow,Lwindow,L,a,M,'CANDUAL');

g=comp_window(g,a,M,L,0,'CANDUAL');

if nargin<4
% XXL: Why ??????

  if Lwindow<L
    % We are dealing with a FIR window.
    if Lwindow<=M
      % The dual will also be FIR
      dofirdual=1;
      g=fir2iir(g,L);
%
% XXL: do it directly ??
%
%       N_win = ceil(Lwindow/a);
%       Lwin_new = N_win*a;
%       if Lwin_new ~= Lwindow
%         g_new = fir2iir(g,Lwin_new);
%       else
%         g_new = g;
%       end
%       weight = sum((reshape(abs(g_new).^2,a,N_win)).').'
%       gamma = g./repmat(weight,N_win,1);
%       if Lwin_new ~= Lwindow
%           gamma = iir2fir(gamma,Lwindow);
%       end
%       if wasrow
%          gamma=gamma.';
%       end;
%       return;
    else
      error(['The computed dual window is an IIR window. Please specify ', ...
             'the exact length L of the window: h=candual(g,a,M,L)']);
    end;
  end;
else
   g=fir2iir(g,L);
end;
  
% If input is real, output must be real as well.
inputwasreal = isreal(g);

scale=1;
if a>M
  % Handle the Riesz basis (dual lattice) case.
  % Swap a and M, and scale differently.
  scale=a/M;
  tmp=a;
  a=M;
  M=a;
end;

% Do the actual computation
gf=comp_wfac(g,a,M);
gammaf=comp_candual_fac(gf,L,a,M)*scale;
gamma=comp_iwfac(gammaf,L,a,M);

if dofirdual
  gamma=iir2fir(gamma,Lwindow);
end;

if inputwasreal
  gamma=real(gamma);
end;
      
if wasrow
  gamma=gamma.';
end;