Пусть и обозначают фильтры анализа нижних и верхних частот и и обозначают соответствующие фильтры синтеза нижних и верхних частот. Банк фильтров с двухканальной критической выборкой фильтрует входной сигнал с использованием фильтра нижних и верхних частот. Выходные сигналы поддиапазонов фильтров понижаются на два для сохранения общего числа выборок. Чтобы восстановить входные данные, выполните повторную выборку на два, а затем интерполяцию результатов с использованием фильтров синтеза нижних и верхних частот. Если фильтры удовлетворяют определенным свойствам, можно добиться идеальной реконструкции входных данных. Чтобы продемонстрировать это, отфильтруйте ЭКГ-сигнал, используя экстремальный фазовый вейвлет Daubechies с двумя моментами исчезновения. В примере объясняется понятие исчезающих моментов в более позднем разделе.

load wecg;
plot(wecg);
title('ECG Signal')

Получение фильтров низкочастотного (масштабирование) и высокочастотного (вейвлет) анализа и синтеза.

[gtilde,htilde,g,h] = wfilters('db2');

В этом примере установите режим заполнения для DWT на периодизацию. Это не расширяет сигнал.

origmodestatus = dwtmode('status','nodisplay');
dwtmode('per','nodisplay');

Получить дискретное вейвлет-преобразование (DWT) уровня один сигнала ЭКГ. Это эквивалентно ветви анализа (с понижающей выборкой) двухканального набора фильтров на рисунке.

[lowpass,highpass] = dwt(wecg,gtilde,htilde);

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

xrec = idwt(lowpass,highpass,g,h);
max(abs(wecg-xrec))
subplot(2,1,1)
plot(wecg); title('Original ECG Waveform')
subplot(2,1,2)
plot(xrec); title('Reconstructed ECG Waveform');

ans =

   1.3658e-12

Фильтры анализа и синтеза для вейвлета «db2» являются только реверсами времени друг друга. Это можно увидеть, сравнив следующее.

scalingFilters = [flip(gtilde); g]
waveletFilters = [flip(htilde); h]

scalingFilters =

    0.4830    0.8365    0.2241   -0.1294
    0.4830    0.8365    0.2241   -0.1294


waveletFilters =

   -0.1294   -0.2241    0.8365   -0.4830
   -0.1294   -0.2241    0.8365   -0.4830

Это относится ко всем ортогональным блокам вейвлет-фильтров. Ортогональными семействами вейвлетов, поддерживаемыми Wavelet Toolbox, являются «dbN», «fkN», «symN» и «cofN», где «N» является допустимым номером фильтра.

Вместо предоставления dwt с фильтрами в предыдущем примере вместо этого используется строка «db2». Используя короткое имя семейства вейвлетов и номер фильтра, нет необходимости правильно указывать фильтры анализа и синтеза.

[lowpass,highpass] = dwt(wecg,'db2');
xrec = idwt(lowpass,highpass,'db2');

Номер фильтра в экстремальной фазе Daubechies и наименьшей асимметричной фазе вейвлетов ('db' и 'sym') относится к числу исчезающих моментов. В основном, вейвлет с N моментами исчезновения удаляет полином порядка, N-1 в коэффициентах вейвлета. Для иллюстрации этого создайте сигнал, который состоит из линейного тренда с аддитивным шумом.

n = (0:511)/512;
x = 2*n+0.2*randn(size(n));
plot(n,x)

Линейный тренд - многочлен степени 1. Поэтому вейвлет с двумя моментами исчезновения удаляет этот многочлен. Линейный тренд сохраняется в коэффициентах масштабирования, и вейвлет-коэффициенты могут рассматриваться как состоящие только из шума. Получите DWT уровня один сигнала с вейвлетом «db2» (два момента исчезновения) и постройте график коэффициентов.

[A,D] = dwt(x,'db2');
subplot(2,1,1)
plot(A); title('Scaling Coefficients');
subplot(2,1,2)
plot(D); title('Wavelet Coefficients');

Вы можете использовать dwt и idwt для реализации двухканального ортогонального набора фильтров, но часто удобнее реализовать многоуровневый двухканальный набор фильтров, используя wavedec. Многоуровневый DWT выполняет итерацию на выходе фильтра нижних частот (масштабирования). Другими словами, вход на второй уровень набора фильтров является выходом фильтра нижних частот на уровне 1. Двухуровневый набор вейвлет-фильтров показан на следующем рисунке.

На каждом последовательном уровне количество коэффициентов масштабирования и вейвлет уменьшается на два, так что общее количество коэффициентов сохраняется. Получение DWT третьего уровня сигнала ЭКГ с помощью набора ортогональных фильтров «sym4».

[C,L] = wavedec(wecg,3,'sym4');

Число коэффициентов по уровню содержится в векторе, L. Первый элемент L равен 256, что представляет количество коэффициентов масштабирования на уровне 3 (конечный уровень). Второй элемент L - это количество вейвлет-коэффициентов на уровне 3. Последующие элементы дают количество вейвлет-коэффициентов на более высоких уровнях до достижения конечного элемента L. Конечный элемент L равен количеству выборок в исходном сигнале. Масштабные и вейвлет-коэффициенты сохраняются в векторе С в том же порядке. Чтобы извлечь коэффициенты масштабирования или вейвлета, используйте appcoef или detcoef. Извлеките все вейвлет-коэффициенты в массиве ячеек и коэффициенты масштабирования конечного уровня.

wavcoefs = detcoef(C,L,'dcells');
a3 = appcoef(C,L,'sym4');

Можно построить график коэффициентов вейвлета и масштабирования в их приближенных положениях.

cfsmatrix = zeros(numel(wecg),4);
cfsmatrix(1:2:end,1) = wavcoefs{1};
cfsmatrix(1:4:end,2) = wavcoefs{2};
cfsmatrix(1:8:end,3) = wavcoefs{3};
cfsmatrix(1:8:end,4) = a3;
subplot(5,1,1)
plot(wecg); title('Original Signal');
axis tight;
for kk = 2:4
    subplot(5,1,kk)
    stem(cfsmatrix(:,kk-1),'marker','none','ShowBaseLine','off');
    ylabel(['D' num2str(kk-1)]);
    axis tight;
end
subplot(5,1,5);
stem(cfsmatrix(:,end),'marker','none','ShowBaseLine','off');
ylabel('A3'); xlabel('Sample');
axis tight;

Поскольку банк вейвлет-фильтров с критической выборкой понижает выборку данных на каждом уровне, анализ должен быть остановлен, когда остается только один коэффициент. В случае сигнала ЭКГ с 2048 выборками, это должно произойти, когда.

[C,L] = wavedec(wecg,log2(numel(wecg)),'sym4');
fprintf('The number of coefficients at the final level is %d. \n',L(1));

The number of coefficients at the final level is 1. 

Если вы хотите реализовать банк ортогональных вейвлет-фильтров без понижающей дискретизации, вы можете использовать modwt.

ecgmodwt = modwt(wecg,'sym4',3);
ecgmra = modwtmra(ecgmodwt,'sym4');
subplot(5,1,1);
plot(wecg); title('Original Signal');

title('MODWT-Based Multiresolution Analysis');
for kk = 2:4
    subplot(5,1,kk)
    plot(ecgmra(kk-1,:));
    ylabel(['D' num2str(kk-1)]);
end
subplot(5,1,5);
plot(ecgmra(end,:));
ylabel('A3'); xlabel('Sample');

В банке биорогональных фильтров фильтры синтеза являются не просто обратимыми во времени версиями фильтров анализа. Семейство биортогональных сплайновых вейвлет-фильтров является примером таких банков фильтров.

[LoD,HiD,LoR,HiR] = wfilters('bior3.5');

При изучении фильтров анализа (LoD, HiD) и фильтров синтеза (LoR, HiR) видно, что они сильно отличаются. Эти банки фильтров все еще обеспечивают идеальную реконструкцию.

[A,D] = dwt(wecg,LoD,HiD);
xrec = idwt(A,D,LoR,HiR);
max(abs(wecg-xrec))

ans =

   6.6613e-16

Биоргональные фильтры полезны, когда линейная фаза является требованием для вашего банка фильтров. Ортогональные фильтры не могут иметь линейную фазу, за исключением вейвлет-фильтра Хаара. При наличии Toolbox™ обработки сигналов можно просмотреть фазовые отклики для ортогональной и биоргональной пары вейвлет-фильтров.

[Lodb6,Hidb6] = wfilters('db6');
[PHIdb6,W] = phasez(Hidb6,1,512);
PHIbior35 = phasez(HiD,1,512);
figure;
subplot(2,1,1)
plot(W./(2*pi),PHIdb6); title('Phase Response for db6 Wavelet');
grid on;
xlabel('Cycles/Sample'); ylabel('Radians');
subplot(2,1,2)
plot(W./(2*pi),PHIbior35); title('Phase Response for bior3.5 Wavelet');
grid on;
xlabel('Cycles/Sample'); ylabel('Radians');

Установите в dwtmode исходное значение.

dwtmode(origmodestatus,'nodisplay');