Как получить цифровой БИХ фильтр из аналогового или как самому «изобрести» цифровой БИХ фильтр
on 19.07.2014 at 17:27Несмотря на свою столь широкую распространенность, разработка цифровых БИХ фильтров (фильтры с Бесконечной Импульсной Характеристикой), даже в наше время, время высоких технологий, для многих остается загадкой. Время от времени мне попадаются статьи или же видео о том, каков принцип разработки БИХ фильтров. Становится печально, когда просмотрев очередное из них я осознаю, что благодаря именно такой подаче материала большинству учащихся действительно не остается ничего другого, как опустить руки и придти в унылое состояние. У многих моих знакомых имеются практические навыки в разработке БИХ фильтров, но они заключаются в банальном использовании, например пакета Matlab, для получения коэффициентов того или иного БИХ фильтра, иногда даже без понимания сути его работы. Лично я, пока полностью не разберусь с принципом работы того или иного участка кода, не могу с уверенностью заявлять о его надёжности. Возобновив свои университетские знания я решил поделиться с Вами своими познаниями в разработке цифрового БИХ фильтра «с нуля»! Дальнейшее изложение материала будет происходить в контексте разработки реального цифрового ФНЧ (фильтра низких частот) на основе аналогового прототипа RC — фильтра. Мы не будем вдаваться во все подробности преобразования Лапласа и Z-преобразования. Я поставил себе за задачу показать читателю сам скелет разработки и мне кажется для старта этого будет вполне достаточно. Все недостающие знания при желании можно без труда почерпнуть из Интернета. «Что нового в этой статье?» — спросите вы. Новшество заключается в том, что в данном относительно коротком посте будут собраны воедино такие вещи как:
- история возникновения аналоговой, а позже и цифровой фильтрации;
- описание схемы с помощью математической модели;
- преобразование полученной модели в другую, более удобную, с помощью преобразования Лапласа;
- будет объяснено, зачем нужен предыдущий пункт;
- Z преобразование;
- получение разностного уравнения (результат работы).
Итак, начнём! Откуда же взялись БИХ фильтры? Как часто бывает, они взяты с уже годами проверенной практикой схемотехнических решений. Исторически так сложилось, что такой электронный компонент как конденсатор был разработан задолго (1745 год) до разработки, скажем радиоприёмника (1895 год). Так вот, ещё на заре развития аудиотехники люди разрабатывали радиоприёмники. В любом из таких устройств в наличии конечно же имелись наушники или динамик. Разработчики заметили, что если подключать конденсатор к динамику (параллельно или же последовательно), то в звуке что-то меняется. Однозначно можно было подобрать конденсатор таким образом, что звук становился чище, т.е. конденсатор уменьшал количество шума (при параллельном подключении к наушникам). Время шло, наука стремительно развивалась и вскоре выросла целая физико-математическая машина, с помощью которой стало возможным объяснить подобные явления. Как мы знаем, цифровая фильтрация до недавнего времени не сильно использовалась из-за невозможности её физической реализации, ведь процессоры обладали низкими значениями тактовой частоты. В виду этого развивалась индустрия фильтрации в аналоговом виде. К слову, развивалась она довольно хорошо и что самое главное качественно. Прошло время и процессоры, с частотами выше 100 МГц, стали реальностью. Стало реальным производить фильтрацию используя механизм БИХ фильтров. Но для того, чтобы использовать модель проверенных годами БИХ фильтров в цифровой фильтрации их необходимо каким-то образом описать и преобразовать математической моделью с целью использования в цифровом виде. Дальше, мы с Вами пройдем путь от аналоговой RC цепочки до конечной реализации цифрового БИХ фильтра конечно же с примером программы в конце статьи, так как статья ориентированна на практиков, а не только теоретиков. Итак, имеем следующую схему:
Видим, в составе резистор и конденсатор. На вход подается напряжение, которое меняется во времени, на выходе будем иметь отличное от входного напряжение, которое, в свою очередь, тоже меняется во времени. Ток, который течет в цепи можно описать, используя следующее выражение для конденсатора:
В таком случае можно записать:
Только ради удобства сделаем замену:
Получаем следующее выражение:
На данном шаге нами было получено выражение, которые в полной мере описывает поведение нашей RC цепи. Но вместо простой зависимости выхода от входа, мы получили что-то, на первый взгляд страшное, ведь мы имеем производную выходного сигнала по времени. В общем главной нашей задачей является ответ на следующий вопрос: «Что нужно сделать с входным сигналом, чтобы получить выходной? Умножить, суммировать, сдвинуть во времени, что?«. Каким бы наивным не казался этот вопрос, но он является ключевым. Все последующие действия, связанные с математическими преобразованиями, будут нацелены именно на то, чтобы можно было получить чёткую запись, которая покажет зависимость выхода от входа. Иными словами, нам необходимо решить полученное дифференциальное уравнение. Конечно же, решать его есть возможным множественным количеством методов, но одним из мне известных, универсальным, является метод с использованием преобразования Лапласа. Называется этот метод в математике методом операционного исчисления. Что это такое? Я специально не буду писать формул, а постараюсь объяснить «на пальцах». Допустим мы имеем функцию (в дискретном виде ее можно представить, скажем, последовательностью чисел, осциллограммой). Так вот, под преобразованием Лапласа понимают разложение начальной функции в сумму затухающих, незатухающих и возрастающих осциллирующих экспонент. Вот как они выглядят:
Некоторые из читающих могут провести некую параллель с преобразованием Фурье. Действительно, преобразование Фурье является частным случаем для преобразования Лапласа при константной (незатухающей) осциллограмме.
Частота, фаза, и скорость затухания, возрастания, постоянство этих составляющих функций может быть разной у каждой из них. Преобразованием Лапласа пользоваться очень просто (во всяком случае в этом конкретном случае). Существуют таблицы преобразования Лапласа и нам необходимо просто посмотреть на них и подставить нужные. В результате такой подстановки мы получим:
Ну, уже совсем другое дело. Теперь не составит труда показать зависимость образа выходного сигнала от входного:
где:
Полученное отношение образов можно рассматривать как отношение выходного спектра (экспоненциальных осциллограмм) к входному спектру (экспоненциальных осциллограмм). Что же такое переменная S? Тут посложнее, хотя это всего лишь комплексная переменная. Откуда она тут взялась? Из преобразования Лапласа. Тут можно оценить достижения математики, так как благодаря операциям с комплексными числами стало возможным получать решение в более элегантном представлении. В общем я обещал не вдаваться в детали. Только скелет (ну, с минимальными отступлениями)! Должен признать, дальше будет ещё сложнее. Но если не вдаваться в подробности, то всё довольно таки просто. Следующим шагом будет переход от преобразования Лапласа к Z-преобразованию (оно же преобразование Фишера или же преобразование Лорана). Z-преобразование это то же самое преобразование Лапласа, но вот только для дискретного (а не непрерывного) сигнала. К счастью для нас, существует простой способ перехода от непрерывного преобразования Лапласа к дискретному Z-преобразованию и называется оно «Билинейное преобразование». Пользоваться им очень просто. Достаточно в преобразовании Лапласа произвести замену:
После подстановки, получим следующее:
Напомню, пока мы работаем в частотной области, не во временной. Наша задача получить зависимость выхода от входа именно во временной области, а не в частотной. Нужно что-то сделать для перехода от частотной во временную область. Этим «что-то» является обратное Z-преобразование. Оно также имеет свою таблицу преобразований и нашей задачей является приведение уравнения до вида «табличных составляющих». Если переписать полученное выражение в виде:
то получить обратное Z-преобразование очень просто:
Чтобы читатель понимал, мы уже находимся во временной области. Перепишем полученное выражение:
Полученное выражение можно считать конечным результатом, но если произвести замену:
то получим более элегантную запись:
Поздравляю нас с Вами! Мы получили выражение, пользуясь которым возможно фильтровать дискретный сигнал. Для тех кто не понимает как применить данное выражение на практике привожу пример кода, написанного на языке программирования Python:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 |
import numpy as np import matplotlib.pyplot as plot Fs = 8000 # Sampling frequency (Hz) Fc = 2000 # Cutoff frequency (Hz) N = Fs # Number of samples per 1 s x = range(0, N) y = range(0, N) f1 = 100 # Hz f2 = 3500 # Hz N_TO_DISPLAY = 2*Fs/f1 # Generating test signal for t in range(0, N): x[t] = (np.sin(2.0*np.pi*f1*t/Fs) + np.sin(2.0*np.pi*f2*t/Fs))/2.0 # Calculating necessary values for filter dt = 1.0/Fs tau = 1.0/(2.0*np.pi*Fc) alpha = dt / (dt + 2*tau) # Applying filter for i in range(1, N): y[i] = alpha*(x[i]+x[i-1]) + (1-2*alpha)*y[i-1] # Drawing original and filtered spectrums plot.clf() plot.xlabel("Time") plot.ylabel("Sample value") plot.plot(x[0:N_TO_DISPLAY], 'r') plot.plot(y[0:N_TO_DISPLAY], 'g') plot.draw() |
В этом коде начальным сигналом является сумма двух синусоид с частотой 100Гц и 3500Гц. Предположим, что синусоида с частотой 3500 Гц — это помеха и мы хотим, с помощью созданного нами цифрового фильтра, отфильтровать её. Частота сэмплирования в данном случае равна 8кГц, а частота среза была выбрана равной 2кГц. В результате отработки данного скрипта получим следующий рисунок:
На данном графике красным цветом показан начальный сигнал (полезный сигнал+помеха), зелёным цветом показан отфильтрованный сигнал (полезный сигнал). Как видим, фильтрация удалась!
И в завершение… Надеюсь, у меня получилось изложить тему доступным Вам, читателю, языком? В следующей статье мы рассмотрим амплитудно-частотную и фазо-частотную характеристики разработанного нами с Вами низкочастотного RC фильтра. Ждите новые посты касательно цифровой обработки сигналов в контексте их применения на практике. Я же буду ждать Ваши отзывы внизу страницы! Вы даже не представляете как приятно получать такого рода отдачу:) Всего хорошего!)
P.S. Обратите внимание, трансформация частотной шкалы «frequency-prewrapping» была преднамеренно упущена мною, дабы не отпугнуть читающего, тем более что эта версия фильтра имеет право на существование и без этой трансформации. Т. е. если сравнивать версию фильтра С использованием «frequency-prewrapping» и БЕЗ него, то АЧХ фильтра будут схожи в диапазоне частот от 0 до частоты дискретизации делённой на 8 (приблизительное значение взято из моей памяти!). Детальнее об «frequency-prewrapping» я собирался написать и напишу в следующих статьях!
P.P.S. А вот и долгожданная статья о трансформации частотной оси! В качестве примера выступает именно RC фильтр низких частот первого порядка. Это абсолютно полноценная версия!
P.P.P.S. И ещё. Дабы не возникало трудностей на этапе получения конечной формулы при подстановке , будет правильным показать данный переход:
Вопрос: а где трансформация частотной шкалы, которая ДОЛЖНА быть выполнена перед z- преобразованием?
Формула преобразования Wd=(2/T)tg(wT/2)
Кстати, и формулы получаются более элегантные — пропадает множитель 2/Т !
Для ФНЧ 1 порядка H(S)=Wp/(S+Wp), где Wp=tg(Fc/FN *Pi/2)
z — преобразование S=(1- Z^{-1})/(1+Z^{-1})
Результат: y[i]= Wp/(Wp+1)*(X[i]+X[i-1])- (Wp+1)/(Wp-1)*Y[i-1]
Спасибо огромное за материал. Перерыл весь интернет, такого подробного объяснения для чайников не нашел. Признаюсь, до этого не приходилось так плотно сталкиваться с теорией ЦОС, поэтому от преобразования Лапласса был в шоке. Сейчас необходимо реализовать на контроллере ФВЧ и ФНЧ. ФВЧ пытаюсь сделать по аналогии, только не понятно как у Вас из
с подстановкой
вышло
. Скорее всего, я упустил что-то очевидное.
Пожалуйста! Рад приносить пользу!)

А преобразования Лапласа не бойтесь. Смотрите на него как на «мудрённое» преобразование Фурье, если так будет проще. Недопонимание с формулой я расписал в «P.P.P.S.». Если собираетесь использовать данный тип фильтра во всей полосе частот и для вас критично важна именно та же самая АЧХ, что и у аналогового прототипа, особенно на высоких частотах, тогда внимательно перечитайте моё «P.P.S.» (frequency-prewrapping). Если же собираетесь просто выделять низкочастотные составляющие — можете не применять «frequency-prewrapping».
Не смог пройти мимо, что бы не выразить свою благодарность
Долгое время искал подобное объяснение, очень давно интересует этот вопрос, читал несколько книг, но результата не было.
Особенно интересно как строятся цифровые фильтры осциллографов.
А тут такое доступное объяснение, подобных сайтов по ЦОС не встречал, спасибо большое!!!
Baghear, спасибо вам за положительный отзыв и пожалуйста! Благодаря таким посетителям как вы, на этом ресурсе и появляются новые статьи
Добрый день!
Я новичок в цифровой обработке сигналов, и мне хотелось бы узнать, какая будет схема фильтра в цифровом виде?
Просто здесь я ее не смог найти
Здравствуйте lo_de00!
Считаю, что для начинающего будет не рыбка, а удочка для ее ловли.
Поэтому, могу посоветовать вам следующую замечательную статью по интересующему вас вопросу:
http://www.dsplib.ru/content/filters/ch10/ch10.html
Разобравшись с материалом вы без труда нарисуете схему полученного в этой статье фильтра.