c_ana.m

Sparse Signal Representation using Overlapping Frames (matlab toolbox)

function Y = C_ana(X,a0,a1)
% C_ana     The analysis part of an IIR filter bank based on allpass filters
%           Each column of X is filtered independently of the other columns
% to give two columns in Y.
% Block diagram, x is a column of X. y0 and y1 are columns of Y.
%                      
%         even x       +-------+  v0             y0
%      +-------------->| A0(z) |------------>(+)----->
%  x   |               +-------+      \     /
% -----+                                x 
%      |  odd x        +-------+  v1  /     \    y1
%      +-------------->| A1(z) |-------->--->(+)----->
%                      +-------+       -1 
%
% Y = C_ana(X,a0,a1)
% ---------------------------------------------------------------------------
% arguments:
%  X        the signal to be filtered
%  a0,a1    the values of a0 and a1 used in the IIR filterbank.
% ---------------------------------------------------------------------------

%----------------------------------------------------------------------
% Copyright (c) 1998.  Karl Skretting.  All rights reserved.
% Hogskolen in Stavanger (Stavanger University), Signal Processing Group
% Mail:  [e-mail]   Homepage:  http://www.ux.his.no/~karlsk/
% 
% HISTORY:
% Ver. 1.0  27.03.1998  Karl Skretting, used in Master thesis 1998, HIS: 
%                       "Kompresjon av seismiske data"
% Ver. 1.1  20.07.2000  KS: added optional arguments a0, a1
% Ver. 1.2  02.12.2002  KS: moved from ..\Frames to ..\FrameTools
%----------------------------------------------------------------------

Mfile='C_ana';

% check input and output arguments, and assign values to arguments
if (nargin < 1); 
   error([Mfile,': function must have one input argument, see help.']); 
end
if (nargin < 3)
   % the default values
   a0=0.1576056086;
   a1=0.6148404766;
end
% periodic extension of the signal in the ends by "ekstra" samples 
ekstra=50;		

N=size(X,1);
if (rem(N,2) ~= 0)
   error([Mfile,': length of vectors in X must have a factor of 2.']); 
end
N2=N/2;

for kol=1:size(X,2)
   x=X(:,kol);
   v0=zeros(N2,1);
   v1=v0;
   % filtering using A0
   k=2-2*ekstra;
   i=mod(k,N);
   if (i==0); i=N; end
   v=0;
   xl=0;					% last x
   for k=1:ekstra
      v=a0*(x(i)-v)+xl;
      xl=x(i);
      i=i+2; if (i>N); i=i-N; end;
   end
   if (i~=2)
      error([Mfile,': logical error, i~=2']); 
   end
   for k=1:N2
      v=a0*(x(i)-v)+xl;
      xl=x(i);
      i=i+2;  % if (i>N); i=i-N; end;
      v0(k)=v;
   end
   % filtering using A1
   k=1-2*ekstra;
   i=mod(k,N);
   v=0;
   xl=0;					% last x
   for k=1:ekstra
      v=a1*(x(i)-v)+xl;
      xl=x(i);
      i=i+2; if (i>N); i=i-N; end;
   end
   if (i~=1)
      error([Mfile,': logical error, i~=1']); 
   end
   for k=1:N2
      v=a1*(x(i)-v)+xl;
      xl=x(i);
      i=i+2; % if (i>N); i=i-N; end;
      v1(k)=v;
   end
   y0=v0+v1;
   y1=v0-v1;
   if (kol==1)
      Y=[y0,y1];
   else
      Y=[Y,y0,y1];
   end
end

return