## [ff23f8]: wavelets / wfilt_db.m  Maximize  Restore  History

### 96 lines (77 with data), 2.3 kB

 ``` 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88``` ```function [h, g, a, info] = wfilt_db(N) %WFILT_DB Daubechies FIR filterbank % Usage: [h,g] = wfilt_db(N); % % Input parameters: % N : Order of Daubechies filters. % Output parameters: % H : cell array of analysing filters impulse reponses % G : cell array of synthetizing filters impulse reponses % a : array of subsampling (or hop) factors accociated with % corresponding filters % % `[H,G] = dbfilt(N)` computes a two-channel Daubechies FIR filterbank % from prototype maximum-phase analysing lowpass filter obtained by % spectral factorization of the Lagrange interpolator filter. *N* also % denotes the number of zeros at \$z=-1\$ of the lowpass filters of length % \$2N\$. The prototype lowpass filter has the following form (all roots of % \$R(z)\$ are outside of the unit circle): % % .. H_l(z)=(1+z^-1)^N*R(z), % % .. math:: H_l(z)=\left(1+z^{-1}\right)^NR(z), % % where \$R(z)\$ is a spectral factor of the Legrange interpolator \$P(z)=2R(z)*R(z^{-1})\$ % All subsequent filters of the two-channel filterbank are derived as % follows: % % .. H_h(z)=H_l((-z)^-1) % .. G_l(z)=H_l(z^-1) % .. G_h(z)=-H_l(-z) % % .. math:: H_h(z)=H_l((-z)^-1) % .. math:: G_l(z)=H_l(z^-1) % .. math:: G_h(z)=-H_l(-z) % % making them an orthogonal causal perfect-reconstruction QMF. % % Examples: % --------- % % wfiltinfo('db8'); % % References: daub98tenlectures if(nargin<1) error('%s: Too few input parameters.',upper(mfilename)); end if(N<1) error('Minimum N is 1.'); end if(N>20) warning('Instability may occur when N is too large.'); end h = cell(2,1); flen = 2*N; % Calculating Lagrange interpolator coefficients sup = [-N+1:N]; a = zeros(1,N); for n = 1:N non = sup(sup ~= n); a(n) = prod(0.5-non)/prod(n-non); end P = zeros(1,2*N-1); P(1:2:end) = a; P = [P(end:-1:1),1,P]; R = roots(P); R = R(abs(R)<1 & real(R)>0); % roots of the 2*conv(lo_a,lo_r) filter hroots = [R(:);-ones(N,1)]; % building synthetizing low-pass filter from roots h{1}= real(poly(hroots)); % normalize h{1}= h{1}/norm(h{1}); % QMF modulation low-pass -> highpass h{2}= (-1).^(1:flen).*h{1}(end:-1:1); g=h; a = [2;2]; info.istight=1; ```