Пакет вейвлета преобразовывает и базовая перемена

В этом примере мы используем пакетный подход вейвлета [3], который имеет потенциал, чтобы быть свободным искажением. Смотрите Пакеты Вейвлета: Разложение Деталей для анализа пакета вейвлета преобразовывает.

Базис пакетного метода вейвлета для гармонической интерференции состоит из следующих шагов:

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

Интерференционный компонент на 60 Гц

Загрузите и постройте сигнал ECG, соответствующий, чтобы записать 200 из Базы данных Аритмии MIT-BIH [1], [2]. Данные производятся на уровне 360 Гц.

load mit200
plot(tm,ecgsig)
xlabel('Seconds')
ylabel('Amplitude')
axis tight
title('ECG Signal with 60-Hz Interference')

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

pspectrum(ecgsig,360)
xlim([45 80])
title('Power Spectrum of ECG Signal with 60-Hz Interference')

Существует также гармоника 60 Гц, видимых в данных на уровне 120 Гц. Однако уровень той интерференции на приблизительно 20 дБ ниже уровня интерференции на 60 Гц. Соответственно, мы фокусируемся здесь на компоненте на 60 Гц.

Во-первых, загрузите и анализируйте два БИХ-фильтра метки созданный с использованием iirnotch от DSP System Toolbox™. Фильтры были спроектированы с двумя факторами Q, 35 и 100, с полосой пропускания по умолчанию-3 дБ.

load Q100Filter3
load Q35Filter3

Исследуйте величину и фазовые отклики для двух фильтров с уровнями полосы пропускания-3 дБ.

helperFreqPhasePlot(Q100Filter3.num,Q100Filter3.den,...
    Q35Filter3.num,Q35Filter3.den);

Фильтр нулевой фазы данные о ECG с помощью метки фильтрует с Q-фактором 100. Постройте спектр мощности отфильтрованных данных.

y = filtfilt(Q100Filter3.num,Q100Filter3.den,ecgsig);
figure
pspectrum([ecgsig,y],360)
xlim([55 65])
title('Power Spectrum of ECG Signal and Notch-filtered Signal')
legend('ECG Signal', 'Notch Filtered')

Фильтр метки сделал хорошее задание удаления интерференции на 60 Гц, но существует ясно некоторое повреждение спектра сигнала вокруг интерференционной частоты.

Если вы уменьшаете масштаб, чтобы видеть целый спектр мощности от 0 до Найквиста, искажение спектра сигнала около 60 Гц ясно отображается, в то время как в целом спектр отфильтрованного сигнала совпадает со спектром исходного сигнала хорошо.

xlim([0 180])

Функция помощника helperWPHarmonicFilter реализует переключенный базовой линией пакетный метод вейвлета, описанный в [3]. Здесь мы используем Добечи меньше всего - асимметричный вейвлет с 8 исчезающими моментами (16 коэффициентов).

sigWP = helperWPHarmonicFilter(ecgsig,60,'sym8',360,'reflection');
pspectrum([ecgsig,sigWP],360)
xlim([55 65])
title({'Power Spectrum of ECG Signal and' ; 'Baseline-Shifted Wavelet Packet Signal'})
legend('ECG Signal', 'Wavelet Packet Cancellation','Location','SouthWest')

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

Уменьшите масштаб, чтобы видеть качество спектра мощности для пакетного метода вейвлета, сравненного с исходным сигналом от 0 до частоты Найквиста.

xlim([0 180])

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

Одно из преимуществ пакетного метода вейвлета - то, что он имеет способность в теории действовать в нескольких гармониках одновременно. Чтобы проиллюстрировать это, загрузите ecgHarmonics сигнала. Это - mit200 сигнал с сильными гармоническими интерференционными компонентами, введенными на уровне 60 и 120 Гц.

load ecgHarmonics
figure
pspectrum(ecgHarmonics,360)
title('Power Spectrum of ECG Signal with Artificial Harmonics')

Удалите эту гармоническую интерференцию с помощью пакетного метода вейвлета. Снова, мы используем наименьшее количество асимметричного вейвлета фазы Добечиса с 8 исчезающими моментами.

sigWP = helperWPHarmonicFilter(ecgHarmonics,60,'sym8',360,'reflection');

Постройте исходную и переключенную базовой линией пакетную реконструкцию вейвлета сигнала.

subplot(2,1,1)
plot(tm,ecgHarmonics)
ylabel('Amplitude')
axis tight
title('Signal with Interference Components')
subplot(2,1,2)
plot(tm,sigWP)
xlabel('Seconds')
ylabel('Amplitude')
title('Baseline-Shifted Wavelet Packet Signal')
axis tight

Исследуйте спектры мощности обоих сигналов.

figure
pspectrum([ecgHarmonics,sigWP],360)
xlim([30 150])
title({'Power Spectrum of ECG Signal and'; 'Wavelet Packet Baseline-Shifted Signal'})
legend('ECG Signal with Interference', 'Baseline-Shifted Wavelet Packet Signal','Location','NorthEast')

Результат довольно хорош с обеими гармониками, удаленными без спектрального искажения в окружении интерференционных гармоник.

Сравните это с БИХ-фильтрацией метки. Здесь мы должны отфильтровать данные дважды, однажды для фильтра метки, действующего на уровне 60 Гц и однажды для фильтра метки, действующего на уровне 120 Гц. В качестве альтернативы можно создать каскад фильтра. Здесь мы хотим использовать filtfilt для нулевой фазы, фильтрующей, таким образом, мы выбираем реализовывать фильтрацию дважды.

Загрузите фильтр метки на 120 Гц созданный с использованием iirnotch с фактором Q 100 и уровнем полосы пропускания по умолчанию-3 дБ. Отфильтруйте сигнал сначала с фильтром метки на 60 Гц и затем отфильтруйте выход той операции с фильтром метки на 120 Гц.

load Q100Filter3_120Hz
filt60Hz = filtfilt(Q100Filter3.num,Q100Filter3.den,ecgHarmonics);
filtCascade = filtfilt(Q100Filter3_120Hz.num,Q100Filter3_120Hz.den,filt60Hz);

Постройте исходный и отфильтрованный меткой сигнал сигнала.

subplot(2,1,1)
plot(tm,ecgHarmonics)
ylabel('Amplitude')
axis tight
title('Signal with Interference Components')
subplot(2,1,2)
plot(tm,filtCascade)
xlabel('Seconds')
ylabel('Amplitude')
title('Notch-filtered Signal')
axis tight

Исследуйте спектры мощности обоих сигналов.

figure
pspectrum([ecgHarmonics,filtCascade],360)
xlim([30 150])
title('Power Spectrum of ECG Signal and Notch-filtered Signal')
legend('ECG Signal', 'Notch-Filtered Signal')

Здесь снова мы видим искажение спектра в местоположениях фильтров метки. Базовая перемена пакетных коэффициентов вейвлета была более эффективной при удалении высоко-амплитудных синусоидальных интерференционных компонентов.

Сводные данные

Этот пример посмотрел на гармоническую интерференционную отмену с помощью переключенного базовой линией пакетного метода вейвлета, описанного в [3]. В данных о ECG, рассмотренных здесь, пакетный метод вейвлета был эффективным при удалении гармонической интерференции и в истинные и в симулированные условия. Пакетный метод вейвлета ввел меньше гармонического искажения в спектре сигнала около интерференционных частот, чем фильтры метки. Было некоторое искажение, введенное пакетным методом вейвлета около частоты Найквиста, которая гарантирует дальнейшее расследование.

Ссылки

[1] Голдбергер, Ари Л., Луис А. Н. Амараль, Леон Гласс, Джеффри М. Гаусдорф, Plamen Ch. Иванов, Роджер Г. Марк, Джозеф Э. Митус, Джордж Б. Муди, Чанг-Канг Пенг и Х. Юджин Стэнли. “PhysioBank, PhysioToolkit и PhysioNet: Компоненты Нового Ресурса Исследования для Комплексных Физиологических Сигналов”. Циркуляция 101, № 23 (13 июня 2000). https://doi.org/10.1161/01. CIR.101.23.e215.

[2] Капризный, G.B., и Р.Г. Марк. “Удар Базы данных Аритмии MIT-BIH”. Разработка IEEE в Журнале 20 Медицины и Биологии, № 3 (июнь 2001): 45–50. https://doi.org/10.1109/51.932724.

[3] Лицзюнь Сюй. “Отмена Гармонической Интерференции Базовой Переменой Пакетных Коэффициентов Разложения Вейвлета”. Транзакции IEEE на Обработке сигналов 53, № 1 (январь 2005): 222–30. https://doi.org/10.1109/TSP.2004.838954.

function helperFreqPhasePlot(num1,den1,num2,den2)
figure
[H100,W] = freqz(num1,den1,[],360);
phi100 = phasez(num1,den1);
H35 = freqz(num2,den2,[],360);
phi35 = phasez(num2,den2);
subplot(2,2,1)
plot(W,10*log10(abs(H100)))
title('Magnitude - Q-factor 100')
ylabel('dB')
grid on
axis tight
subplot(2,2,2)
plot(W,phi100*180/pi)
title('Phase  - Q-factor 100')
ylabel('Degrees')
grid on
axis tight
subplot(2,2,3)
plot(W,10*log10(abs(H35)))
title('Magnitude - Q-factor 35')
xlabel('Hz')
ylabel('dB')
grid on
axis tight
subplot(2,2,4)
plot(W,phi35*180/pi)
title('Phase - Q-factor 35')
xlabel('Hz')
ylabel('Degrees')
grid on
axis tight


end