Добавление квадратурных зеркальных и биортогональных вейвлет-фильтров

Открыть сценарий в реальном времени

В этом примере показано, как добавить пару ортогонального квадратурного зеркального фильтра (QMF) и четырехкратный биоргональный вейвлет-фильтр к вейвлет- Toolbox™. В то время как Vavelet Toolbox™ уже содержит множество наиболее широко используемых ортогональных и биоргональных семейств вейвлетов, включая экстремальную фазу Daubechies, наименее асимметричную фазу Daubechies, коифлет, фильтры Фейера-Коровкина и биортогональные сплайновые вейвлеты, можно легко добавить собственные фильтры и использовать фильтр в любом из дискретных вейвлет.

Этот пример добавляет пару фильтров Бейлкина (18) QMF к панели инструментов и показывает, как впоследствии использовать фильтр в дискретном вейвлет-анализе. Затем пример демонстрирует, как проверить необходимые и достаточные условия для пары QMF, чтобы составить масштабный и вейвлет-фильтр. После добавления пары QMF в примере добавляется почти ортогональный четверной биоргональный вейвлет на основе схемы пирамиды Лапласа Берта и Адельсона (таблица 8.4 на стр. 283 в [1]).

Добавление QMF

Во-первых, у вас должен быть какой-то способ получения коэффициентов. В этом случае здесь приведены коэффициенты для фильтра Бейлкина (18) нижних частот (масштабирования). Вам нужен только допустимый фильтр масштабирования ,wfilters создает для вас соответствующий вейвлет-фильтр.

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];

Сохраните фильтр Бейлкина (18) и добавьте новый фильтр к панели инструментов.

save beyl beyl

Использовать wavemngr для добавления вейвлет-фильтра к панели инструментов. Определите имя семейства вейвлетов и краткое имя, используемое для доступа к фильтру. Определите тип вейвлета как 1. Вейвлеты типа 1 являются ортогональными вейвлетами на панели инструментов. Поскольку в этом семействе добавляется только один вейвлет, определите ввод переменной NUMS в wavemngr быть пустой строкой.

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

Добавление вейвлета с помощью wavemngr.

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

Убедитесь, что вейвлет добавлен в панель инструментов.

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           '
    '==================================='

Теперь можно использовать вейвлет для анализа сигналов или изображений. Например, загрузите ЭКГ-сигнал и получите MODWT сигнала до четвертого уровня с помощью фильтра Бейлкина (18).

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

Загрузите изображение коробки, получите 2-D DWT с помощью фильтра Бейлкина (18). Отображение коэффициентов детализации диагонали первого уровня.

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')

Figure contains 2 axes. Axes 1 with title Original Image contains an object of type image. Axes 2 with title Level-One Diagonal Coefficients contains an object of type image.

Наконец, убедитесь, что новый фильтр удовлетворяет условиям для ортогональной пары QMF. Получение фильтров масштабирования (lowpass) и вейвлет (high pass).

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

Суммирование коэффициентов фильтра нижних частот для проверки того, что сумма равна $\sqrt{2}$ . Суммируйте коэффициенты вейвлет-фильтра и убедитесь, что сумма равна 0.

sum(Lo)

ans = 1.4142

sum(Hi)

ans = -1.9873e-16

Убедитесь, что автокорреляция масштабного и вейвлет-фильтров на всех даже ненулевых лагах равна 0. Для использования необходимо иметь Toolbox™ обработки сигналов 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');

Figure contains 2 axes. Axes 1 with title Autocorrelation of Scaling Filter contains an object of type stem. Axes 2 with title Autocorrelation of Wavelet Filter contains an object of type stem.

Обратите внимание, что значения автокорреляции на обоих графиках равны нулю для ненулевых четных лагов. Показать, что взаимная корреляция масштабирования и вейвлет фильтра равна нулю на всех четных лагах.

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

Figure contains an axes. The axes with title Cross-correlation of QMF filters contains an object of type stem.

Конечный критерий определяет, что сумма квадратов преобразований Фурье масштабирования и вейвлет-фильтров на каждой частоте равна 2. Другими словами, пусть $G (f$ ) - преобразование Фурье масштабного фильтра и $H ($ f) - преобразование Фурье вейвлет-фильтра. Для $всех$ f: $| H (^{f}) | 2 +^{} |$ G (f) | 2 = 2. DFT версия этого $_{^{} равенства:}^{} |_{^{}}^{G2mkmodN} |$ 2 + | H2mkmodN | 2 = 2 для любого m. Проверьте это для фильтра Бейлкина (18) с m = 0.

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 = 18×1

    2.0000
    2.0000
    2.0000
    2.0000
    2.0000
    2.0000
    2.0000
    2.0000
    2.0000
    2.0000
      ⋮

Все значения равны 2, как и ожидалось. Чтобы понять, почему эти фильтры называются квадратурными зеркальными фильтрами, визуализируйте частотные характеристики масштабирования и вейвлет-фильтров в квадрате.

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;

Figure contains an axes. The axes with title Beylkin(18) QMF Filter Pair contains 2 objects of type line. These objects represent Scaling Filter, Wavelet Filter.

Обратите внимание, что амплитудные характеристики являются симметричными или зеркальными изображениями друг друга вокруг квадратурной частоты 1/4.

Следующий код удаляет вейвлет-фильтр Бейлкина (18).

wavemngr('del',familyShortName);
delete('beyl.mat')

Добавление биорогонального вейвлета

Добавление биорогонального импульса к панели инструментов аналогично добавлению QMF. Предоставляется действительная пара фильтров нижних частот (масштабирование), используемая при анализе и синтезе. wfilters функция генерирует фильтры верхних частот.

Подлежит признанию wfilters, фильтр масштабирования анализа должен быть назначен переменной Df, и фильтр масштабирования синтеза должен быть назначен переменной Rf. Фильтры биорогонального масштабирования не обязательно должны иметь даже одинаковую длину. Создаваемые выходные пары биоргональных фильтров будут иметь даже равные длины. Вот пары масштабных функций почти ортогонального квадрупла биоргональных вейвлетов, основанные на схеме лапласианской пирамиды Берта и Адельсона.

Df = [-1 5 12 5 -1]/20*sqrt(2);
Rf = [-3 -15 73 170 73 -15 -3]/280*sqrt(2);

Сохранение фильтров в .mat файл.

save burt Df Rf

Использовать wavemngr для добавления биорогональных вейвлет-фильтров к панели инструментов. Определите имя семейства вейвлетов и краткое имя, используемое для доступа к фильтру. Поскольку вейвлеты являются биортогональными, установите тип вейвлета равным 2. Поскольку в этом семействе добавляется только один вейвлет, определите ввод переменной NUMS в wavemngr быть пустой строкой.

familyName      = 'burtAdelson';
familyShortName = 'burt';
familyWaveType  = 2;
familyNums      = '';
fileWaveName    = 'burt.mat';
wavemngr('add',familyName,familyShortName,familyWaveType,...
    familyNums,fileWaveName)

Убедитесь, что биорогональный вейвлет добавлен в панель инструментов.

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             '
    'burtAdelson       ->->burt           '
    '==================================='

Теперь можно использовать вейвлет на панели инструментов. Создание банка фильтров DWT анализа с помощью burt вейвлет. Подтвердите, что банк фильтров DWT является биорогональным. Постройте график частотно-амплитудных характеристик полосовых вейвлет-фильтров и функции масштабирования с наибольшим разрешением.

fb = dwtfilterbank('Wavelet','burt');
isBiorthogonal(fb)

ans = logical
   1

freqz(fb)

Figure contains an axes. The axes with title DWT Filter Bank burt contains 8 objects of type line. These objects represent D 1, D 2, D 3, D 4, D 5, D 6, D 7, A 7.

Получение функций вейвлета и масштабирования набора фильтров. Постройте график функций вейвлета и масштабирования в самом крупном масштабе.

[fb_phi,t] = scalingfunctions(fb);
[fb_psi,~] = wavelets(fb);
subplot(2,1,1)
plot(t,fb_phi(end,:))
axis tight
grid on
title('Analysis - Scaling')
subplot(2,1,2)
plot(t,fb_psi(end,:))
axis tight
grid on
title('Analysis - Wavelet')

Figure contains 2 axes. Axes 1 with title Analysis - Scaling contains an object of type line. Axes 2 with title Analysis - Wavelet contains an object of type line.

Создание банка фильтров синтеза DWT с помощью burt вейвлет. Вычислите границы кадров.

fb2 = dwtfilterbank('Wavelet','burt','FilterType','Synthesis','Level',4);
[synthesisLowerBound,synthesisUpperBound] = framebounds(fb2)

synthesisLowerBound = 0.9800

synthesisUpperBound = 1.0509

Получение фильтров анализа и синтеза нижних и верхних частот, связанных с burt. Обратите внимание, что выходные фильтры имеют одинаковую длину. Подтвердите сумму коэффициентов фильтра нижних частот sqrt(2) и коэффициенты фильтра верхних частот суммируются до 0.

[LoD,HiD,LoR,HiR] = wfilters('burt');
[LoD' HiD' LoR' HiR']

ans = 8×4

         0    0.0152   -0.0152         0
         0   -0.0758   -0.0758         0
   -0.0707   -0.3687    0.3687   -0.0707
    0.3536    0.8586    0.8586   -0.3536
    0.8485   -0.3687    0.3687    0.8485
    0.3536   -0.0758   -0.0758   -0.3536
   -0.0707    0.0152   -0.0152   -0.0707
         0         0         0         0

sum([(LoD'/sqrt(2)) HiD' (LoR'/sqrt(2)) HiR'])

ans = 1×4

    1.0000   -0.0000    1.0000         0

Удалите фильтр Burt-Adelson из панели инструментов.

wavemngr('del',familyShortName);
delete('burt.mat')

Ссылки

[1] Daubechies, I. Десять лекций по вейвлетам. Серия региональных конференций CBMS-NSF по прикладной математике. Филадельфия, Пенсильвания: Общество промышленной и прикладной математики, 1992.

Документация

Добавление квадратурных зеркальных и биортогональных вейвлет-фильтров

Добавление QMF

Добавление биорогонального вейвлета

Ссылки

Документация инструментария вейвлет

Поддержка