Позвольте и обозначьте lowpass и highpass аналитические фильтры и и обозначьте соответствующий lowpass и highpass фильтры синтеза. Двухканальный критически произведенный набор фильтров фильтрует входной сигнал с помощью фильтра highpass и lowpass. Поддиапазон выходные параметры фильтров прорежен два, чтобы сохранить общее количество выборок. Чтобы восстановить вход, сверхдискретизируйте два и затем интерполируйте результаты с помощью lowpass и highpass фильтров синтеза. Если фильтры удовлетворяют определенным свойствам, можно достигнуть совершенной реконструкции входа. Чтобы продемонстрировать это, отфильтруйте сигнал ECG использование экстремального вейвлета фазы Добечиса с двумя исчезающими моментами. Пример объясняет понятие исчезающих моментов в более позднем разделе.

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

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

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

В данном примере установите дополнительный режим для DWT к periodization. Это не расширяет сигнал.

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

Получите уровень один дискретный вейвлет преобразовывает (DWT) сигнала ECG. Это эквивалентно аналитической ветви (с субдискретизацией) двухканального набора фильтров на рисунке.

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

Сверхдискретизируйте и интерполируйте lowpass (масштабные коэффициенты) и highpass (коэффициенты вейвлета) поддиапазоны с фильтрами синтеза и продемонстрируйте совершенную реконструкцию.

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', и 'coifN', где N является допустимым номером фильтра.

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

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

Номер фильтра в экстремальной фазе Добечиса и наименьшем количестве асимметричных вейвлетов фазы ('дб' и '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 выполняет итерации на выходе lowpass (масштабирование) фильтра. Другими словами, входом к второму уровню набора фильтров является выход фильтра lowpass на уровне 1. Двухуровневый набор фильтров вейвлета проиллюстрирован в следующем рисунке.

На каждом последовательном уровне количество масштабирования и коэффициентов вейвлета прорежено два, таким образом, общее количество коэффициентов сохраняется. Получите уровень три DWT сигнала ECG использование 'sym4' ортогонального набора фильтров.

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

Количество коэффициентов уровнем содержится в векторе, L. Первый элемент L равен 256, который представляет количество масштабных коэффициентов на уровне 3 (итоговый уровень). Вторым элементом L является количество коэффициентов вейвлета на уровне 3. Последующие элементы дают количество коэффициентов вейвлета в более высоких уровнях, пока вы не достигаете итогового элемента L. Итоговый элемент L равен количеству отсчетов в исходном сигнале. Масштабирование и коэффициенты вейвлета хранятся в векторе C в том же порядке. Чтобы извлечь масштабирование или коэффициенты вейвлета, используйте 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;

Поскольку критически произведенный набор фильтров вейвлета прореживает данные на каждом уровне, анализ должен остановиться, когда вы имеете только один коэффициент в запасе. В случае сигнала ECG с 2 048 выборками это должно произойти когда.

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

Если вы исследуете аналитические фильтры (ЛОД, HiD) и фильтры синтеза (LoR, HiR), вы видите, что они очень отличаются. Эти наборы фильтров все еще обеспечивают совершенную реконструкцию.

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

ans =

   4.4409e-16

Биоортогональные фильтры полезны, когда линейная фаза является требованием для вашего набора фильтров. Ортогональные фильтры не могут иметь линейной фазы за исключением фильтра вейвлета Хаара. Если у вас есть Signal Processing 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');