Documentation

This is machine translation

Translated by Microsoft
Mouseover text to see original. Click the button below to return to the English verison of the page.

Note: This page has been translated by MathWorks. Please click here
To view all translated materals including this page, select Japan from the country navigator on the bottom of this page.

Adding a Quadrature Mirror Filter

This example shows how to add an orthogonal quadrature mirror filter (QMF) pair to the Wavelet Toolbox™. While Wavelet Toolbox™ already contains many of the most widely used orthogonal QMF families, including the Daubechies' extremal-phase, the Daubechies' least-asymmetric phase, the coiflet, and the Fejer-Korovkin filters, you can easily add your own QMF filters and use the filter in any of the discrete wavelet or wavelet packet algorithms. This example adds the Beylkin(18) QMF filter to the toolbox and shows how to subsequently use the filter in discrete wavelet analysis. Finally, the example demonstrates how to verify the necessary and sufficient conditions for the QMF pair to constitute a scaling and wavelet filter.

First, you must have some way of obtaining the coefficients. In this case, here are the coefficients for the lowpass (scaling) Beylkin(18) filter. You only need a valid scaling filter, wfilters creates the corresponding wavelet filter for you.

beyl = [9.93057653743539270E-02
    4.24215360812961410E-01
    6.99825214056600590E-01
    4.49718251149468670E-01
    -1.10927598348234300E-01
    -2.64497231446384820E-01
    2.69003088036903200E-02
    1.55538731877093800E-01
    -1.75207462665296490E-02
    -8.85436306229248350E-02
    1.96798660443221200E-02
    4.29163872741922730E-02
    -1.74604086960288290E-02
    -1.43658079688526110E-02
    1.00404118446319900E-02
    1.48423478247234610E-03
    -2.73603162625860610E-03
    6.40485328521245350E-04];

Save the Beylkin(18) filter and add the new filter to the toolbox.

save beyl beyl

Use wavemngr to add the wavelet filter to the toolbox. Define the wavelet family name and the short name used to access the filter. Define the wavelet type to be 1. Type 1 wavelets are orthogonal wavelets in the toolbox. Because you are adding only one wavelet in this family, define the NUMS variable input to wavemngr to be an empty string.

familyName      = 'beylkin';
familyShortName = 'beyl';
familyWaveType  = 1;
familyNums      = '';
fileWaveName    = 'beyl.mat';

Add the wavelet using wavemngr.

wavemngr('add',familyName,familyShortName,familyWaveType, ...
    familyNums,fileWaveName)

Verify that the wavelet has been added to the toolbox.

wavemngr('read')
ans =

  19x35 char array

    '==================================='
    'Haar              		haar           '
    'Daubechies        		db             '
    'Symlets           		sym            '
    'Coiflets          		coif           '
    'BiorSplines       		bior           '
    'ReverseBior       		rbio           '
    'Meyer             		meyr           '
    'DMeyer            		dmey           '
    'Gaussian          		gaus           '
    'Mexican_hat       		mexh           '
    'Morlet            		morl           '
    'Complex Gaussian  		cgau           '
    'Shannon           		shan           '
    'Frequency B-Spline		fbsp           '
    'Complex Morlet    		cmor           '
    'Fejer-Korovkin    		fk             '
    'beylkin           		beyl           '
    '==================================='

You can now use the wavelet to analyze signals or images. For example, load an ECG signal and obtain the MODWT of the signal down to level four using the Beylkin(18) filter.

load wecg
wtecg = modwt(wecg,'beyl',4);

Load a box image, obtain the 2-D DWT using the Beylkin(18) filter. Show the level-one diagonal detail coefficients.

load xbox
[C,S] = wavedec2(xbox,1,'beyl');
[H,V,D] = detcoef2('all',C,S,1);
subplot(2,1,1)
imagesc(xbox)
axis off
title('Original Image')
subplot(2,1,2)
imagesc(D)
axis off
title('Level-One Diagonal Coefficients')

Finally, verify that the new filter satisfies the conditions for an orthogonal QMF pair. Obtain the scaling (lowpass) and wavelet (highpass) filters.

[Lo,Hi] = wfilters('beyl');

Sum the lowpass filter coefficients to verify that the sum equals . Sum the wavelet filter coefficients and verify that the sum is 0.

sum(Lo)
sum(Hi)
ans =

    1.4142


ans =

  -1.9873e-16

Verify that the autocorrelation of the scaling and wavelet filters at all even nonzero lags is 0. You have to have the Signal Processing Toolbox™ to use xcorr.

[Clow,lags] = xcorr(Lo,Lo,10);
Chigh = xcorr(Hi,Hi,10);
subplot(2,1,1)
stem(lags,Clow,'markerfacecolor',[0 0 1])
grid on;
title('Autocorrelation of Scaling Filter');
subplot(2,1,2)
stem(lags,Chigh,'markerfacecolor',[0 0 1])
grid on;
title('Autocorrelation of Wavelet Filter');

Note that the autocorrelation values in both plots is zero for nonzero even lags. Show that the cross-correlation of the scaling and wavelet filter is zero at all even lags.

[xcr,lags] = xcorr(Lo,Hi,10);
figure
stem(lags,xcr,'markerfacecolor',[0 0 1]);
title('Cross-correlation of QMF filters')

The final criterion states the sum of squared magnitudes of the Fourier transforms of scaling and wavelet filters at each frequency is equal to 2. In other words, let be the Fourier transform of the scaling filter and be the Fourier transform of the wavelet filter. The following holds for all : . The DFT version of this equality is: for any . Check this for the Beylkin(18) filter with .

N = numel(Lo);
LoDFT = fft(Lo);
HiDFT = fft(Hi);
k = 0:N-1;
m = 0;
sumDFTmags = abs(LoDFT(1+mod(2^m*k,N))).^2+abs(HiDFT(1+mod(2^m*k,N))).^2
sumDFTmags =

    2.0000
    2.0000
    2.0000
    2.0000
    2.0000
    2.0000
    2.0000
    2.0000
    2.0000
    2.0000
    2.0000
    2.0000
    2.0000
    2.0000
    2.0000
    2.0000
    2.0000
    2.0000

All the values are equal to 2 as expected. To understand why these filters are called quadrature mirror filters, visualize the squared-magnitude frequency responses of the scaling and wavelet filters.

nfft = 512;
F = 0:1/nfft:1/2;
LoDFT = fft(Lo,nfft);
HiDFT = fft(Hi,nfft);
figure
plot(F,abs(LoDFT(1:nfft/2+1)).^2,'DisplayName','Scaling Filter');
hold on
plot(F,abs(HiDFT(1:nfft/2+1)).^2,'r','DisplayName','Wavelet Filter');
xlabel('Frequency'); ylabel('Squared Magnitude')
title('Beylkin(18) QMF Filter Pair')
grid on
plot([1/4 1/4], [0 2],'k','HandleVisibility','off');
legend show;

Note the magnitude responses are symmetric, or mirror images, of each other around the quadrature frequency of 1/4.

The following code removes the Beylkin(18) wavelet filter.

wavemngr('del',familyShortName);
delete('beyl.mat')
Was this topic helpful?