For completeness, a direct matlab implementation of the built-in
`filter` function (Eq.(3.3)) is given in Fig.3.2.
While this code is useful for study, it is far slower than the
built-in `filter` function. As a specific example, filtering
samples of data using an order 100 filter on a 900MHz Athlon
PC required 0.01 seconds for `filter` and 10.4 seconds for
`filterslow`. Thus, `filter` was over a thousand times
faster than `filterslow` in this case. The complete test is
given in the following matlab listing:

x = rand(10000,1); % random input signal B = rand(101,1); % random coefficients A = [1;0.001*rand(100,1)]; % random but probably stable tic; yf=filter(B,A,x); ft=toc tic; yfs=filterslow(B,A,x); fst=tocThe execution times differ greatly for two reasons:

- recursive feedback cannot be ``vectorized'' in general, and
- built-in functions such as
`filter`are written in C, precompiled, and linked with the main program.

function [y] = filterslow(B,A,x) % FILTERSLOW: Filter x to produce y = (B/A) x . % Equivalent to 'y = filter(B,A,x)' using % a slow (but tutorial) method. NB = length(B); NA = length(A); Nx = length(x); xv = x(:); % ensure column vector % do the FIR part using vector processing: v = B(1)*xv; if NB>1 for i=2:min(NB,Nx) xdelayed = [zeros(i-1,1); xv(1:Nx-i+1)]; v = v + B(i)*xdelayed; end; end; % fir part done, sitting in v % The feedback part is intrinsically scalar, % so this loop is where we spend a lot of time. y = zeros(length(x),1); % pre-allocate y ac = - A(2:NA); for i=1:Nx, % loop over input samples t=v(i); % initialize accumulator if NA>1, for j=1:NA-1 if i>j, t=t+ac(j)*y(i-j); %else y(i-j) = 0 end; end; end; y(i)=t; end; y = reshape(y,size(x)); % in case x was a row vector |

[How to cite this work] [Order a printed hardcopy] [Comment on this page via email]

Copyright ©

Center for Computer Research in Music and Acoustics (CCRMA), Stanford University