Цифровая фильтрация кардиосигнала на основе БПФ

В данной заметке мы расскажем о классическом построении цифрового фильтра на основе быстрого преобразования Фурье. Похожий алгоритм цифровой фильтрации кардиограммы используется в ПО ECG Control . Для примера отфильтруем фрагмент реальной кардиограммы в программном пакете Mathcad.

В двух словах о том, что такое цифровой фильтр. Цифровой фильтр - это некий программный алгоритм, направленный на подавление помех в исходном сигнале за счет коррекции его (сигнала) спектрального состава. Работает это довольно просто.

С любовью вспомним студенческие годы и среди множества полезных и не очень 🙂 занятий найдем в глубинах памяти основы теории радиоэлектронных цепей. А гласила она о следующем: любой периодический сигнал можно представить как сумму синусоидальных сигналов разной амплитуды и частоты. Эти синусоидальные сигналы называются гармоники или, в рамках радиотехники, спектральные составляющие (а сам результат разложения сигнала на синусоидальные функции называется спектром). Если совсем по-простому (да простят меня теоретики, которые считают себя физиками-практиками) спектр сигнала представляет собой зависимость амплитуды гармоники сигнала от ее частоты. Представление функции суммой ее гармоник в матанализе называется преобразованием Фурье. Представляем дальше: коль скоро любой сигнал раскладывается на гармоники, что и кардиограмма тоже, ведь она еще какая периодическая (напомню, что периодический сигнал - тот, который повторяется через равные промежутки времени). Разложить мы кардиограмму можем, но назревает вопрос - а зачем? Затем, что в реальных условиях, как бы мы ни давили помеху осветительной сети 50 Гц аналоговыми фильтрами, частотно-зависимыми отрицательными обратными связями усилителей, цепями активной нейтрализации, в конечном сигнале, который придет на АЦП, она все равно будет. Причем, помехи этой будет много. Очень много, чтобы там не говорили теоретики-разработчики аналоговых схем, у которых разработка заканчивается программными пакетами Matlab/Simulink или NI Multisim. Как от нее, гадины, избавится, полностью? С использованием вычислительных мощностей вашего ПК изгнание бесов помех становится задачей решаемой: раскладываем сигнал на спектральные составляющие (получаем его спектр), т.е. осуществляем прямое преобразование Фурье. Находим в спектре кардиосигнала помеху 50 Гц и ее амплитуду устанавливаем в ноль (так её!). После данной процедуры проводим обратное преобразование Фурье - то есть возвращаем наш сигнал с задавленной помехой обратно во временную область - получаем кардиосигнал без 50 Гц. Вот собственно и вся хитрость цифровой фильтрации. Дополню, что процедура дискретного преобразования Фурье, хоть прямого, хоть обратного, в классическом ее виде, муторная и медленная. Но нашлись умники, которые придумали быстрое преобразование Фурье, уменьшив, тем самым, количество требуемых операций с N^2 до N * log(N), где N - количество отсчетов в фильтруемом фрагемнте. Останавливаться на алгоритме БПФ не будем - кому интересно, читаем здесь. А мы переходим к практике.

Цифровую фильтрацию произведем с фрагментом кардиосигнала, длительность которого составляет 5 секунд, частота дискретизации – 2 кГц (т.е. размер выборки составляет 10 000 отсчетов). На рисунке слева представлен исследуемый фрагмент сигнала во временной области, в котором присутствует сигнал помехи 50 Гц и высшие гармоники. Произведем уменьшение тестовой выборки до размера 8 192 точки для корректного применения алгоритмов быстрого преобразования Фурье (алгоритм работает лишь с таким количеством выборок, которое представляют собой степень двойки - в нашем случае, ближайшее число к 10 000 - это два в тринцадцатой степени = 8192 ). На рисунке справа показан спектр мощности исследуемого сигнала, полученный с помощью алгоритма БПФ.

Для коррекции спектрального состава сигнала реализуем функцию полосовой цифровой фильтрации, в качестве аргумента которой передаются коэффициенты Фурье преобразования k, частота среза  f₀  и число высших гармоник помехи numH, которые необходимо подавить.

В представленном выражении fd  определяет ширину частотной полосы подавления, fi связывает частоту спектральной составляющей с ее индексом k в общей последовательности, res – переменная, которая возвращает результат вычислений. Подавление помехи 50 Гц и четырех ее высших гармоник  производится путем умножения спектра сигнала на реализованную корректирующую функцию полосовой цифровой фильтрации.

На рисунке слева представлен спектр кардиосигнала после применения полосовой фильтрации, на рисунке справа – ЭКГ сигнал после коррекции его спектрального состава с использованием цифрового фильтра на основе БПФ.

Работу алгоритма можно проверить в программе Mathcad; запись кардиограммы, зарегистрированной без цифровой фильтрации, в текстовом файле можно скачать здесь.

Обработка ЭКГ, цифровой фильтр на основе БПФ, фильтрация кардиограммы, как убрать помеху, зашумленная ЭКГ, спектр кардиограммы