Цифровые частотные фильтры высоких и низких частот

Цифровые частотные фильтры высоких и низких частот

Цифровые фильтры (Лекция)

По виду импульсной характеристики цифровые фильтры делятся на два больших класса:

· Фильтры с конечной импульсной характеристикой (КИХ – фильтры, трансверсальные фильтры, нерекурсивные фильтры). Знаменатель передаточной функции таких фильтров – некая константа.

КИХ – фильтры характеризуются выражением:

· Фильтры с бесконечной импульсной характеристикой (БИХ – фильтры, рекурсивные фильтры) используют один или более своих выходов в качестве входа, то есть образуют обратную связь. Основным свойством таких фильтров является то, что их импульсная переходная характеристика имеет бесконечную длину во временной области, а передаточная функция имеет дробно-рациональный вид.

БИХ – фильтры характеризуются выражением:

Отличие КИХ – фильтров от БИХ – фильтров заключается в том, что у КИХ – фильтров выходная реакция зависит от входных сигналов, а у БИХ – фильтров выходная реакция зависит от текущего значения.

Импульсная характеристика – это реакция схемы на единичный сигнал.

Е диничный сигнал определяется следующим образом:

Таким образом, единичный сигнал только в одной точке равен единице – в точке начала координат.

Задержанный е диничный сигнал определяется следующим образом:

Таким образом, задержанный единичный сигнал задерживает на k периодов дискретизации.

Сигналы и спектры

Дуальность (двойственность) представления сигналов.

Все сигналы можно представить во временной или частотной плоскости.

Причем, частотных плоскостей – несколько.

Для просмотра сигнала во временной плоскости существует прибор:

Представим, что здесь есть достаточно длинный синусоидальный сигнал (в 1 сек. 1000 раз повторилась синусоида):

Возьмем сигнал с частотой, в два раза больше:

Сложим эти сигналы. Получим не синусоиду, а искаженный сигнал:

Преобразования из временной плоскости в частотную плоскость производятся с помощью преобразований Фурье.

Для просмотра сигнала в частотной плоскости существует прибор:

Частота циклическая или круговая ( f ).

Частотная плоскость покажет засечку:

Величина засечки пропорциональна амплитуде синусоиды, а частота:

Для второго сигнала частотная область покажет другую засечку:

Во временной области суммарного сигнала появится 2 засечки:

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

Преобразования из временной плоскости в частотную плоскость может производиться различными путями. Например: с помощью преобразований Лапласа или с помощью преобразований Фурье.

Три формы записи рядов Фурье.

Существует три формы записи рядов Фурье:

· Синус – косинусная форма .

1.) В синус – косинусной форме ряд Фурье имеет вид:

Входящие в формулу кратные частоты 1 называются гармониками; гармоники нумеруются в соответствии с индексом k ; частота ωk = 1называется k-й гармоникой сигнала.

Это выражение говорит о следующем: что любую периодическую функцию можно представить в виде суммы гармоник, где:

где T – период повторений этой функции;

ω – круговая частота .

где t – текущее время;

При разложении по Фурье самое главное – это периодичность. За счет неё происходит дискретизация по частоте, начинается некоторое количество гармоник.

Для того, чтобы установить возможность тригонометрического разложения для заданной периодичной функции, нужно исходить из определенного набора коэффициентов. Прием для их определения придумал во второй половине XVIII века Эйлер и независимо от него в начале XIX века – Фурье.

Три формулы Эйлера для определения коэффициентов:

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

Мощность гармоник падает с увеличением их номера. Если добавить/отбросить некоторые гармонические составляющие, то перерасчет остальных членов (других гармоник) не требуется.

Практически все функции являются четными или нечетными:

Например , функция Cos :

Четная функция симметрична относительно

Если функция четная, то все синусные коэффициенты bk будут равны нулю и в формуле ряда Фурье будут присутствовать только косинусные слагаемые.

Например , функция Sin :

Нечетная функция симметрична относительно центра.

Если функция нечетная, то все косинусные коэффициенты ak будут равны нулю и в формуле ряда Фурье будут присутствовать только синусные слагаемые.

2.) Вещественная форма записи ряда Фурье.

Некоторое неудобство синусно-косинусной формы ряда Фурье состоит в том, что для каждого значения индекса суммирования k (т.е. для каждой гармоники с частотой 1) в формуле фигурирует два слагаемых – синус и косинус. Воспользовавшись формулами тригонометрических преобразований, сумму этих двух слагаемых можно трансформировать в косинус той же частоты с иной амплитудой и некоторой начальной фазой:

Если S ( t ) является четной функцией, фазы φ могут принимать только значения 0 и π , а если S ( t ) – функция нечетная, то возможные значения для фазы φ равны + π /2.

Если bk = 0, тогда tg φ = 0 и угол φ = 0

Если ak = 0, тогда tg φ – бесконечен и угол φ =

В этой формуле может быть и минус (смотря какое направление взято).

3.) Комплексная форма записи ряда Фурье.

Данная форма представления ряда Фурье является, пожалуй, наиболее употребимой в радиотехнике. Она получается из вещественной формы представлением косинуса в виде полусуммы комплексных экспонент (такое представление вытекает из формулы Эйлера e jθ = Cosθ + jSinθ ):

Применив данное преобразование к вещественной форме ряда Фурье, получим суммы комплексных экспонент с положительными и отрицательными показателями:

А теперь будем трактовать экспоненты со знаком «минус» в показателе как члены ряда с отрицательными номерами. В рамках этого же общего подхода постоянное слагаемое a /2 станет членом ряда с нулевым номером. В результате получится комплексная форма записи ряда Фурье :

Формула расчета коэффициентов Ck ряда Фурье:

Если S ( t ) является четной функцией, коэффициенты ряда Ck будут чисто вещественными, а если S ( t ) – функция нечетная, коэффициенты ряда окажутся чисто мнимыми.

Совокупность амплитуд гармоник ряда Фурье часто называют амплитудным спектром, а совокупность их фаз – фазовым спектром.

Спектром амплитуд является действительная часть коэффициентов Ck ряда Фурье:

Re ( Ck ) – спектр амплитуд.

Спектр прямоугольных сигналов.

Рассмотрим сигнал в виде последовательности прямоугольных импульсов с амплитудой A , длительностью τ и периодом повторения T . Начало отсчета времени примем расположенным в середине импульса.

Данный сигнал является четной функцией, поэтому для его представления удобнее использовать синусно-косинусную форму ряда Фурье – в ней будут присутствовать только косинусные слагаемые ak , равные:

Из формулы видно, что длительность импульсов и период их следования входят в нее не обособлено, а исключительно в виде отношения. Этот параметр – отношение периода к длительности импульсов – называют скважностью последовательности импульсов и обозначают буквой: g : g = T /τ. Введем этот параметр в полученную формулу для коэффициентов ряда Фурье, а затем приведем формулу к виду Sin ( x )/ x :

Читайте также:  Гаражный блок питания для ремонтных работ

Примечание: В зарубежной литературе вместо скважности используется обратная величина, называемая коэффициентом заполнения ( duty cycle ) и равная τ / T .

При такой форме записи становится хорошо видно, чему равно значение постоянного слагаемого ряда: поскольку при x → 0 Sin ( x )/ x →1, то

Теперь можно записать и само представление последовательности прямоугольных импульсов в виде ряда Фурье:

Амплитуды гармонических слагаемых ряда зависят от номера гармоники по закону Sin ( x )/ x .

График функции Sin ( x )/ x имеет лепестковый характер. Говоря о ширине этих лепестков, следует подчеркнуть, что для графиков дискретных спектров периодических сигналов возможны два варианта градуировки горизонтальной оси – в номерах гармоник и в частотах.

На рисунке градуировка оси соответствует номерам гармоник, а частотные параметры спектра нанесены на график с помощью размерных линий.

Итак, ширина лепестков, измеренная в количестве гармоник, равна скважности последовательности (при k = ng имеем Sin ( π k / g ) = 0, если n ≠ 0). Отсюда следует важное свойство спектра последовательности прямоугольных импульсов – в нем отсутствуют (имеют нулевые амплитуды) гармоники с номерами, кратными скважности.

Расстояние по частоте между соседними гармониками равно частоте следования импульсов – 2 π / T . Ширина лепестков спектра, измеренная в единицах частоты, равна 2 π / τ, т.е. обратно пропорциональна длительности импульсов. Это проявление общего закона – чем короче сигнал, тем шире его спектр.

Вывод: для любого сигнала известны его разложения в ряд Фурье. Зная τ и T можем посчитать сколько гармоник нужно, чтобы передать мощность.

Методы анализа линейных систем с постоянными коэффициентами.

Задача в постановке:

Имеется линейная система (не зависит от амплитуды сигнала):

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

Это типичная задача электротехники. Имеется мощный способ решения данной задачи во временной области.

Порядок уравнения зависит от числа реактивных элементов.

Может быть записано в виде системы уравнений первой степени.

Схема состоит из резистора и конденсатора

(интегрирующая цепь). На вход подали сигнал X ( t ). Определить Y вых.

RC + UC = X ( t )

UC – является Y выхода, поэтому: RC + U ВЫХ. = X ( t )

Дальнейшее решение сводится к решению сначала однородного уравнения, а затем неоднородного.

Это решение немного упрощается при переводе из временной плоскости в другую плоскость комплексной переменной. Перевод из временной плоскости в комплексную плоскость производится прямым преобразованием Лапласа.

RCY ‘ + Y = X ( t )

Вычисляется разностное уравнение.

Прямое преобразование Лапласа.

Преобразование Лапласа – интегральное преобразование, связывающее функцию S ( p ) комплексного переменного (изображение) с функцией s ( x ) действительного переменного ( оригинал).

Преобразования Лапласа играют очень важную роль при исследовании систем, описываемых линейными дифференциальными уравнениями. С помощью прямого преобразования Лапласа можно перейти от дифференциальных уравнений к алгебраическим, решить их в алгебраической форме, а затем с помощью обратного преобразования получить искомый результат. Аналогичный результат достигается при решении линейных разностных уравнений при использовании аппарата Z – преобразования.

Прямое преобразование Лапласа осуществляется по формуле: , где – комплексная переменная , где σ – затухание.

Честно простой цифровой фильтр

Вы работаете с АЦП. Получаете результаты преобразования, один за одним. И замечаете, что эти результаты «скачут». А хотелось бы, чтобы стояли, как… Ну, короче, чтобы стояли!
Есть много причин, почему отсчеты АЦП могут быть нестабильны. В своей заметке я не говорю об этих причинах. Я говорю о том, как успокоить показания, получая их AS IS. И как сделать это максимально просто. При этом, возможно, не имея ни малейшего понятия о науке под названием «цифровая обработка сигналов».
Эта заметка написана в качестве полной замены предыдущей заметки. Ту лучше не читать :)

ПОСТАНОВКА ЗАДАЧИ

Имеется последовательность кодов, дискретно во времени представляющих физический сигнал. Мы будем говорить о последовательности кодов с АЦП. Физический сигнал и полученная последовательность кодов имеют шумы, выбросы, «болтанку» — назвать можно как угодно. Наша задача — сгладить входную последовательность, то есть, выдать выходную последовательность такую, чтобы влияние шумов было уменьшено.
При этом мы стремимся выполнить задачу максимально простыми программными средствами обычного микроконтроллера.
Более того, у нас поставлено еще одно условие. Допустим, что нам неудобно накапливать данные, а потом обрабатывать их и выдавать результат один раз на N полученных кодов. То ли буфер для данных негде организовать, то ли темп выдачи результатов должен совпадать с темпом получения кодов от АЦП, то ли еще что. Но условие поставлено.
И тут нам на помощь приходит цифровой рекурсивный фильтр, самой простой реализацией которого является фильтр первого порядка. Вот его и будем делать.

ТЕРМИНОЛОГИЯ

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

В каждый из моментов времени t1-2Т , t1-T, t1, t1+T и так далее на входе появляется очередной отсчет последовательности Х, а цифровой фильтр выдает новый отсчет последовательности Y. Если на каждый входной отсчет фильтр выдает и выходной отсчет, то это и есть классический цифровой фильтр.
Иногда поступают не так: набирают некоторое количество входных отсчетов и затем обрабатывают их. Таким образом, на несколько входных получают один выходной отсчет. Это, как правило, уже задача получения того или иного интегрального значения (например, среднего). Она настолько распространена при фильтрации шумов, что часто воспринимается как «обычная». И поэтому классический цифровой фильтр, выдающий отсчеты с каждым входным, называют “скользящим усреднением“, как нечто не совсем обычное для усреднения. Что ж, значит мы рассматриваем скользящее усреднение. Но понимаем, что это обычный цифровой фильтр :)
Важной особенностью цифрового фильтра является то, что входные отсчеты приходят с постоянным темпом, в нашем случае — с периодом Т. Обратная величина называется частотой дискретизации и играет огромную роль во всех параметрах цифрового фильтра. Самое простое, что необходимо себе четко представлять: частотная характеристика цифрового фильтра симметрична относительно Fд/2 и полностью повторяется за пределами . Из чего следует, что фильтр нижних частот, который мы проектируем сейчас, не способен подавлять сигналы с частотой , 2Fд, 3Fд, и т.д. Этого понимания пока достаточно для нашей задачи.
Математически цифровой фильтр 1-го порядка описывается различными способами. Мы будем использовать такое представление:

Читайте также:  Преобразователь напряжения +U в -U на микросхеме CD4049, схема.

Y(n) = Alfa*Y(n-1) + Beta*X(n)

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

Alfa + Beta = 1

В рамках данной заметки задача расчета цифрового рекурсивного фильтра 1-го порядка состоит в нахождении коэффициентов Alfa иBeta с учетом удобства их использования в микроконтроллере (МК) для цифровой фильтрации отсчетов.

РАСЧЕТ ФИЛЬТРА

Цифровой фильтр (равно как и не цифровой) имеет как бы 2 лица: частотную характеристику и переходную (временнУю) характеристику. Задачу расчета можно ставить как получение требуемой частотной либо требуемой переходной характеристики. Я люблю исходить из заданного поведения во временной области. И объясню почему.
Дело в том, что частотная характеристика фильтра первого порядка… как бы это выразиться… простенькая. Возьмем 2 фильтра с частотой дискретизации 200 Гц. Первый с частотой среза по уровню -3 дБ, равной 5 Гц, второй — с частотой 1 Гц:

Как видим, в большей части частотной характеристики подавление шумов составляет всего 15. 30 дБ, а у второго — 20..40 дБ. Вспомним, я отмечал выше, что за пределами Fд/2 (у нас это 100 Гц) характеристика симметрична и в районе 200 Гц подавления шумов снова нет. Если нам нужно «взрослое зубастое подавление», то необходимо строить более серьезные фильтры.
Но все же, второй фильтр лучше давит помехи! А что, если частоту среза еще понизить?
И тут оказывается, что в погоне за парочкой децибел дополнительного подавления мы совсем забыли о временнОм «лице» фильтра. А давайте глянем, как реагируют эти два фильтра на ступенчатое входное воздействие (сигнала не было, а потом он вдруг стал равен 1 и таким и остался).

Что мы видим? Если говорить о времени переходного процесса, то первый фильтр (который похуже фильтрует) отрабатывает входной скачек примерно за 160 мс, а второй, наш передовик с подавлением 20. 40 дБ — почти за 800 мс. Вот так-то. Лучше фильтрация — хуже переходной процесс. Поэтому и нужно выбрать некий оптимум.
Вот, понимая это, я и предлагаю: исходить из требований по быстродействию. Задавшись временем реакции на ступенчатое входное воздействие, мы получим параметры фильтра (вот те самые Alfa иBeta), а параметры фильтрации примем уж какие получатся. Парочка децибел туды-сюды уже мало что изменят, а быстродействие будет известно.

Вот почему в качестве исходного требования я выбираю обеспечение заданного времени переходного процесса при подаче на вход ступенчатого воздействия. Рассмотрим чуть пристальнее, как же фильтр отрабатывает такой сигнал. Если дать скачек от нуля до некоторого значения, обозначенного 100%, то на выходе фильтра увидим:

Как видим, на 8-м отсчете фильтр уже отработал 90% входного воздействия, на 10-м — 95%, и так далее. Обычно принято говорить о «недоработанном», т.е. о тех 10 или 5 процентах, на которые фильтр еще «врет» к какому-то отсчету. Говорят еще, что погрешность установления выходного сигнала составляет столько-то процентов. Далее я привожу формулы для погрешности установления от 5 до 0,1%. В первом случае переходной процесс закончился «на глазок», а точности 0,1% обычно достаточно для того, чтобы считать процесс полностью законченным.
Разница между этими 5%-ым и 0,1%-ым фильтрами не столь уж велика. Я предполагаю, что все фильтры, которые будут разработаны по описываемой методике, находятся в континиуме между этими двумя крайними точками. В качестве характеристики фильтра по степени «законченности» переходного процесса введем такой параметр: Ntau — и вычислим его крайние значения:

LN(1/5%) = 2,996
LN(1/0,1%) = 6,908
Ntau = 2,996. 6,908

Смысл Ntau примерно таков: чем более жесткие требования мы предъявляем к завершенности переходного процесса, тем больше это число. Так что нижней границе значений Ntau соответствует «грубый» процесс, когда мы спешим считать отработку законченной, а верхней границе — «точный» переходной процесс.

Итак, разработчик задался временем переходного процесса Тпп, за которое фильтр отработает скачек на 95. 99,9%. Что еще нужно знать? Время выборки — тот период, с которым на вход фильтра поступают выборки и с таким же темпом с фильтра уходят результаты. Он у нас обозначен Т.
Ясно, что за время переходного процесса будет обработано много отсчетов. Сколько?

N = Тпп / Т (1)

И наша задача — выбрать значения Alfa иBeta так, чтобы уложиться в Тпп.
Оказывается, все очень просто. Должно выполняться условие:

Задавшись двумя граничными значениями Ntau, зная требуемое время переходного процесса и период выборки, мы находим 2 значения Alfa, соответствующие «грубому» и «точному» решениям. Выбираем значение Alfa (ниже я покажу как) в диапазоне:

где А потом легко вычисляем второй коэффициент:

Beta = 1 — Alfa (3)

Вот и все. Коэффициенты фильтра вычислены. Поставив их в утилиту, описанную в заметке коллеги mc2, можете полюбоваться частотной характеристикой. И убедиться, что за требуемое вам время Тпп переходной процесс заканчивается.
Если уж говорить о программе для расчета, то можно обойтись и без «моего» метода. Но тут уж вам решать. Расчет по формулам (1). (3) несложен, но гарантирует вам безытерационное решение для заданного переходного процесса.

ПРОГРАММНАЯ РЕАЛИЗАЦИЯ

Если мы работаем с МК, то очень хорошей практикой является использование целочисленной арифметики. Даже работа с 32-разрядными LONG переменными в 8-битном МК здесь выполняется намного быстрее, чем с FLOAT. Поэтому я рекомендую выбирать Alfa и Beta таким образом:

Alfa = Na / 2^k
Beta = Nb / 2^k
Na + Nb = 2^k

Иначе говоря, вместо дробных Alfa и Beta нужно взять целые числа, такие, чтобы их сумма была равна целой степени двойки: 2, 4, 8, 16, 32, 64, 128, 256, 512 и т.д.
Вот тогда процессору работы поменьше:

Y(n) = (Na*Y(n-1) + Nb*X(n)) >> k

Микроконтроллер умножает 2 раза целое на целое, складывает, производит сдвиг вправо на k разрядов. Несколько микросекунд — и у вас новый отсчет выходной величины.

ПРИМЕР ДЛЯ ЦЕЛОЧИСЛЕННОЙ АРИФМЕТИКИ

Рассмотрим пример, в котором я задался N=100. Тогда условие (2) будет выглядеть так:

0,9333
Для целочисленной арифметики выберем однобайтную базу коэффициентов. Иначе говоря, пусть сумма Na и Nb будет равна 256. Тогда k=8 и получаем такие границы целочисленного эквивалента Alfa:

239 239
Выбрав любое значение Na из этого диапазона, мы можем быть уверены, что за 100 выборок после подачи сигнала на вход выходной сигнал войдет в зону 0,1. 5 % от установившегося значения. Хотим именно 0,1% — берем меньшее значение Na из указанного диапазона, допускаем 5% — берем большее значение.
И уже от выбранного значения Na расчитываем второй коэффициент:

Читайте также:  Автоматизированное управление вентилятором на микросхемах К561ЛА7 и К561ИЕ9, схема.

Nb = 256 — Na

Теперь формула для работы микроконтроллера такова:

Y(n) = (Na*Y(n-1) + Nb*X(n)) >> 8 (4)

На этом расчет фильтра с целочисленными коэффициентами закончен.

Коллега _pv подсказал в обсуждении хорошую запись вместо (4), полностью эквивалентую по результату:

Y(n) = Y(n-1) + [Nb*(X(n) — Y(n-1)) >> 8]

Здесь имеем экономию на одно умножение (при произвольных соотношениях между коэффициентами) и формула приобретает вид приращения старого значения Y(n) на разность между входом и выходом, умноженную на коэффициент Beta. Чем меньше этот коэффициент, тем медленнее следит выход за входом — более сильная фильтрация.

Приведу сишный код, соответствующий формуле от _pv:

Здесь учтено, что при делении (сдвиге вправо) в целочисленной арифметике теряются разряды. Поэтому промежуточный результат сохраняется в переменной z в масштабе суммы (до деления), а точнее, в масштабе выходной величины, деленной на Beta. Кстати, здесь и ответ на вопрос о необходимом диапазоне представления z: по максимальному значению выходной величины, деленному на Beta. Например, однобайтные переменные и коэффициент Nb, равный 1 (при последующем сдвиге на 8 бит, как в показанной функции) требуют двухбайтного z.

Еще один фокус — и я откланяюсь :)
Если не стремиться к конкретной цифре погрешности установления (0,1% или 5% или что-то еще), то можно попытаться выбросить одно умножение, выбрав меньшее из чисел Na и Nb равным 1.
В нашем примере меньшее из чисел именно Nb и легко посчитать, что оно может быть в диапазоне 8. 17. Вычислим следующие границы:

256 / Nb = 32. 16

Здесь числа равны степени двойки случайно. Важно то, что из диапазона допустимых значений 256 / Nb мы выбираем одно из чисел, равное именно степени двойки. Например, в нашем случае, берем 16.
Тогда 16 – это сумма «удобных» целочисленных Na и Nb, а меньшее из них равно 1, то есть микроконтроллеру не нужно умножать:

Y(n) = (15*Y(n-1) + X(n)) >> 4 (5)

Видим, что расчет (5) имеет на одну операцию умножения меньше, чем расчет (4).

РЕЗЮМЕ

Расчет простого цифрового фильтра нижних частот сводится к таким действиям:

Первое. Зададимся временем переходного процесса, периодом выборок и вычислим, за сколько выборок процесс должен закончиться (1):

N = Тпп / Т

Второе. Вычислим границы возможных значений коэффициента при Y (2):

EXP(-6,908 / N)
Третье. Выбираем k (рекомендую 8, но можно и иное) и умножаем границы Alfa на 2^k. Выбираем любое из возможных значений в качестве первого коэффициента Na, а второй коэффициент находим как дополнение:

Nb = 2^k — Na

Результирующая формула обработки входной последовательности X(n):

Y(n) = (Na*Y(n-1) + Nb*X(n)) >> k

Как показано в примере, иногда можно подобрать Na или Nb, равными 1, что упрощает вычисления.

БЛАГОДАРНОСТИ
Хочу выразить искреннюю благодарность уважаемому коллеге angel5a за то, что он вернулся из детского садика (кто читал дискуссию, тот поймет) и очень помог мне с переделкой исходной заметки. Радикальной переделкой, хочу заметить.
Также в доработке статьи помогли дискуссии с коллегами _pv и известным в наших кругах товарищем avreal, выступившим на ненашем форуме.

    , , ,
  • +18
  • 21 августа 2012, 18:03

Комментарии ( 128 )

Первую статью не читал. Эту с первого прочтения понял. Досконально не разобрался, но я это сделаю при первой же необходимости. Ваша цель достигнута, спасибо )

drvlas, Спасибо за подробное толкование первого варианта.
Какую литературу посоветуете почитать человеку, несведущему в данной области, но желающему понять основные принципы?
(хотя Юкио Сато я прочел)

Так я ж… Эта… Не читатель. Я писатель. Современных книг не знаю. Учился лет 30 назад. ЩАС, погодь, коллеги обратят внимание, подскажут.
А вообще, основные принципы я понял не в институте (нам такого и не рассказывали), а только тогда, когда понадобилось практически. Брал книги, читал, зная, зачем мне это нужно. Тогда и обучение шло быстро. А иначе… Ой, как там всего много!

Введение

Цифровой фильтр – фильтр, обрабатывающий цифровой сигнал с целью выделения и (или) подавления определённых частот этого сигнала.

  • Разделяют два больших класса цифровых фильтров:
  • Фильтры с конечной импульсной характеристикой (КИХ-фильтры, FIR)
  • Фильтры с бесконечной импульсной характеристикой (БИХ-фильтры, IIR)

Зачем они нужны? Чем нас не устраивают аналговые фильтры? Давайте рассмотрим преимущества и недостатки цифровых фильтров.

Преимущества цифровых фильтров:

  • Высокая точность и воспроизводимость (у аналоговых фильтров это определяется разбросом номиналов компонентов)
  • Гибкость (можно изменять характеристику, не затрагивая аппаратную часть)
  • Компактность (например, аналоговые ФНЧ могут занимать много места, в то время как цифровые имеют одинаковые габариты при любой форме АЧХ в пределах одной частоты дискретизации)

Недостатки цифвровых фильтров:

  • Сложность работы в реальном времени (все вычисления должны пройти менее, чем за один период дискретизации)
  • Высокая стоимость

Начнём изучение фильтров с наиболее простого типа — КИХ-фильтров.

КИХ-фильтры

Давайте сразу рассмотрим пример. Сгенерируем в Matlab зашумлённый сигнал частотой 0.5 Гц, оцифрованный с частотй дискретизации 10 Гц. Для начала, объявим переменную fs , которой присвоим значение частоты дискретизации, создадим массив временных отсчётов ts , и переменную N , которой присвоим количество получившихся отсчётов:

Теперь создадим сам синусоидальный сигнал x с частотой 0.5 Гц, зашумлённый случайным сигналом, амплитуда которого изменяется в диапазоне от a=-0.01 до b=0.1 :

Построим график сигнала x :

И получим результат:

Сгенерированный сигнал x(n)

Возьмём первые пять отсчётов с x(1) по x(5) , найдём их среднее арифметическое, запишем его в новый массив в элемент y(1) , затем возьмём отсчёты с x(2) по x(6) , также рассчитаем их среднее значение и запишем его в элемент y(2) и так далее, пока не пройдём по всему сигналу:

А теперь построим сигналы x и y друг под другом в одном окне:

Результаты выполнения получившегося скрипта показаны ниже:

Результаты работы фильтра “скользящее среднее”

И что же мы видим? Сигнал y(n) (на рисунке снизу) похож на сигнал x(n) (сверху), однако имеет некоторую задержку и более гладкую форму. Получается, что мы только что применили к сигналу x(n) фильтр нижних частот (ФНЧ)! Спроектированный нами фильтр называется скользящее среднее. Уравнение данного фильтра имеет вид:

Читайте также:  Повышающий регулятор мощности для паяльника на микросхеме К176ЛА7 , схема.

begin{equation*} y[n]=sumlimits_{k=0}^{M-1} h(k) cdot x(n-k) text{, где} sumlimits_{k=0}^M h(k) = 1 end{equation*}

(1)

А его структурная схема выглядит так:

Структурная схема фильтра

Структурная схема фильтра “скользящее среднее состоящего из пяти отсчётов

В каждый момент времени берётся текущий отсчёт и 4 предыдущих, каждый из которых домножается на коэффициент 1/5, затем эти результаты складываются и полученная сумма поступает на выход. Процедура повторяется, пока весь сигнал не будет обработан.

Свёртка

А что, если вместо коэффициентов 1/5 взять что-то другое? На самом деле, уравнение (1) представляет собой частный случай операции свёртки:

begin{equation*} y[n] = x[n] ast h[n]=sumlimits_{k=0}^{M-1} h(k) cdot x(n-k), end{equation*}

(2)

где , — длина , — длина

Это уравнение и является уравнением КИХ-фильтра. Его структурная схема выглядит так:

Структурная схема КИХ-фильтра

Структурная схема КИХ-фильтра

Теорема о свёртке

Со свёрткой связана следующая теорема:

begin{equation*} x[n] * h[n] leftrightarrow X[m] cdot H[m] end{equation*}

Если два сигнала и имеют дискретные преобразования Фурье (ДПФ) и соответственно, то свёртка этих сигналов во временной области эквивалентна произведению их спектров в частотной области (и наоборот):

(3)

И действует это в обе стороны. Получается, что, когда мы делаем свёртку исходного сигнала с коэффициентами фильтра, мы перемножаем амплитудно-частотную характеристику (АЧХ) исходного сигнала и АЧХ фильтра, тем самым убираем ненужные частоты. А это и есть фильтрация.

Мы выше говорили: “фильтры с конечной импульсной характеристикой, фильтры с бесконечной импульсной характеристикой”. А что же такое “импульсная характеристика”? Давайте вспомним на примере фильтров.

Импульсная характеристика КИХ-фильтра

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

Возьмём предыдущий пример из Matlab и вместо исходного сигнала x подставим в него функцию Дирака:

Результат выполнения скрипта показан ниже:

Функция Дирака и реакция на неё скользящего среднего из пяти отсчётов

Из рисунка видно, что импульсная характеристика фильтра представляет собой 5 отсчётов амплитудой 0.2 (или 1/5). Получается, что мы с вами видим ни что иное, как коэффициенты фильтра. Поэтому коэффициенты КИХ-фильтра также называют его импульсной характеристикой.

Проектирование КИХ-фильтров

Теперь возникает вопрос: как рассчитать коэффициенты, чтобы получить требуемую АЧХ фильтра? Самый простой способ: “рисуем” необходимую АЧХ в частотной области, делаем обратное ДПФ от этой АЧХ и получаем набор коэффициентов во временной области. Наверняка, каждый из вас мечтает увидеть фильтр с идеально прямоугольной АЧХ. Возможно ли такое? Давайте разбираться. Что из себя представляет ДПФ от идеально прямоугольного сигнала? Правильно, функцию , которую мы рассматривали, когда изучали растекание спектра ДПФ:

Функция

Данная функция бесконечна во временной области, поэтому для реализации идеально прямоугольной АЧХ нам потребуется бесконечное количество коэффициентов фильтра. Получается, чем больше мы используем коэффициентов h[n] , тем больше АЧХ фильтра будет похожа на прямоугольную. Вроде бы логично — возьми побольше коэффициентов, получишь хороший фильтр. Но не всё так просто: в реальности мы не можем увеличивать количество коэффициентов фильтра, т.к. это приведёт к дополнительным операциям перемножения, каждая из которых вызывает задержку. Для того, чтобы система работала в реальном времени, нужно, чтобы все вычисления по формуле (2) производились за время, не превышающее один период дискретизации. Вот и простейшее ограничение.

Давайте попробуем на примере. Создадим сигнал x(n) , состоящий из N нулевых отсчётов, затем присвоим первым 200 отсчтётам значение 1. Это будет наша требуемая форма АЧХ:

Рассчитаем ОДПФ от сигнала x(n) , для удобства сдвинем его с помощью ifftshift , отбросим мнимую часть и присвоим полученный результат сигналу y(n) . Далее построим графики сигнала x(n) , и, чтобы лучше рассмотреть сигнал у(n) , график его 200 центральных отсчётов.

Результат выполнения скрипта показан ниже:

Требуемая форма АЧХ фильтра x(n) и 200 центральных отсчётов его импульсной характеристики y(n)

Теперь давайте анализировать реальную АЧХ от спроектированного фильтра. Для этого создадим сигнал a(n) , состоящий из N нулевых отсчётов, затем присвоим первому отсчёту 1. Таким образом, получим функцию Дирака, спектр которой является константой на всей частотной оси:

Теперь попробуем пропустить функцию Дирака через спроектированный ранее фильтр y(n) , при этом будем использовать разное количество отсчётов y(n) :

Результат выполнения скрипта показан ниже:

Форма АЧХ фильтра в зависимости от количества отсчётов

Из графиков видно, что, чем больше мы используем отсчётов нашего фильтра (а количество отсчётов – это порядок цифрового фильтра), тем более его характеристика становится похожа на идеальную. Однако, стоит обратить внимание, что на вершине АЧХ всех наших фильтров видны пульсации независимо от того, сколько используется отсчётов в импульсной характеристике. Эти пульсации называются пульсации Гиббса и возникают из-за медленной сходимости ряда Фурье, которая обусловлена наличием разрыва функции на частоте среза полосы пропускания фильтра. С увеличением числа отсчетов уменьшается длительность выбросов на вершине АЧХ, но их амплитуда не меняется и составляет примерно 9% от амплитуды АЧХ на частоте среза.

Окна при проектировании КИХ-фильтров

Запишем выражение (3) следующим образом:

begin{equation*} h^infty[k] cdot w[k] leftrightarrow H^infty [m] * W[m] end{equation*}

(4)

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

Результат выполнения скрипта показан ниже:

Результаты моделирования выражения (4)

Чтобы уменьшить амплитуду пульсаций, как и в случае с ДПФ, при проектировании КИХ-фильтров используют окна, отличные от прямоугольного. Давайте модифицируем наш листинг “Анализ АЧХ фильтра, часть 4” и добавим в строчках 26, 34, 42, 50 умножение на окно Хэмминга:

Полученные результаты показаны ниже:

Форма АЧХ фильтра, взвешенного окном, в зависимости от количества отсчётов

Другое дело! Пульсации ушли, однако АЧХ фильтра стала более “заваленной”. Выбор типа окна зависит от конкретной решаемой задачи. Самые распространённые: параметрическое окно Кайзера, Чебышёва, окно Блэкмана, Хэмминга и др.

БИХ-фильтры

Основное отличие в структуре БИХ-фильтра от КИХ-фильтра — наличие обратной связи. Что такое обратная связь? Это когда на вход устройства подаётся сигнал, пропорциональный сигналу на его выходе. Таким образом, значение выходного сигнала БИХ-фильтра зависит не только от текущего и предыдущего значения входного сигнала, но и от предыдущего значения выходного сигнала. Это одновременно и плюс, и минус: наличие обратной связи позволяет значительно сократить количество отсчётов (а значит, и порядок фильтра, а как следствие – количество умножений), однако может привести к неустойчивости его работы, или вообще превратить фильтр в генератор. Из-за того, что при отсутствии входного сигнала на выходе БИХ-фильтра может генерироваться бесконечное количество отсчётов, его и прозвали фильтром с бесконечной импульсной характеристикой.

Читайте также:  Зарядное Устройство для любого шуруповерта и не только

Давайте на примере сравним, сколько же потребуется коэффициентов КИХ- и БИХ-фильтра для реализации одинаковой АЧХ. Спроектируем с помощью filterDesigner два ФНЧ со следующими параметрами:

В результате получили два фильтра, АЧХ и ФЧХ которых представлены ниже:

АЧХ и ФЧХ КИХ-фильтра АЧХ и ФЧХ БИХ-фильтра

Из рисунка видно, что АЧХ фильтров действительно схожи, однако, для реализации КИХ-фильтра нам потребовалось 222 коэффициента, а для реализации БИХ-фильтра — всего 41. Получается, что в случае с БИХ-фильтром, операций умножения требуется в 5 раза меньше.

Следует также обратить внимание на график фазы. У КИХ-фильтров ФЧХ гарантированно линейная, у БИХ-фильтров она гарантированно нелинейная, о чём нужно помнить.

Рассмотрим структуру БИХ-фильтра:

Структурная схема БИХ-фильтра

Структурная схема БИХ-фильтра

Он состоит из двух частей: прямой и обратной связи. Прямая связь повторяет структурную схему КИХ фильтра, обратная связь напоминает зеркальную копию прямой связи.

Таким образом, разностное уравнение БИХ-фильтра имеет вид:

begin{equation*} y[n] = sumlimits_{i=0}^{N} b(i) cdot x(n-i) - sumlimits_{k=1}^{M} a(k) cdot y(n-k), end{equation*}

(5)

При проектировании БИХ-фильтров используют z-преобразование, корни которого идут от преобразования Лапласа. Второе преобразование вам точно знакомо, его должны были изучать на курсе математического анализа. Давайте вспомним, что это такое, а затем перейдём к изучению z-преобразования. Но это на следующей лекции.

Скачать конспект в pdf: Digital Filters Lecture – V.V. Leonidov.pdf

Симметричные НЧ-ВЧ фильтры

В задачах обработки сигналов часто возникает необходимость фильтрации сигналов, когда сигнал разбивается на узкополосные диапазоны. В бытовом плане мы с этим сталкиваемся при воспроизведении музыки через акустические системы, в которых каждый громкоговоритель (динамик) воспроизводит свою полосу частот, которых обычно три — низкие (НЧ), средние (СЧ) и высокие (ВЧ); для воспроизведения сверхнизких частот иногда выделяют отдельную акустическую систему под названием «сабвуфер». Конкретные границы частот зависят от реализации и ориентировочно находятся на границах 100 Гц, 1 кГц и 5 кГц. Для того, чтобы не было резких скачков громкости между динамиками, используют частичное перекрытие — когда амплитуда воспроизводимой полосы частот плавно спадает на одном, одновременно нарастая на другом.

Наиболее популярными фильтрами для такого разбиения являются фильтры Линквитца-Рейли 4-го порядка, представляющих из себя два последовательно соединённых фильтра Баттерворта, изображение АЧХ которых многим хорошо знакомо:

Популярность их обусловлена простотой схемотехнической реализации и гладкой АЧХ с отсутствием пульсаций. Порядок фильтра, обуславливающий скорость затухания, обычно выбирают 4-ым.

У всех аналоговых фильтров есть недостатки — нестабильность параметров и сложность в настройке. Другой немаловажным недостатком является неустранимый сдвиг фаз. Реализовать фазолинейный фильтр в аналоге невозможно, потому что это нарушает причинно-следственную связь выходного сигнала от входного.

Появление цифровой техники в целом и цифровых источников сигналов в частности позволило избавиться от недостатков аналоговых фильтров путём реализации их непосредственно программным образом, в цифре. Такой подход требует «мультиампинг» — когда используется отдельный усилитель мощности для каждой полосы частот — в отличие от классического подхода, когда сначала усиливается широкополосный сигнал, который разбивается на частотные полосы (обычно) пассивными фильтрами. Применение нескольких усилителей вместо одного очевидным образом удорожает технику, поэтому такое решение выбирают для особо качественных систем.

▍ Краткое введение в цифровые фильтры

Цифровой фильтр — это функция, которую можно рассматривать как во временной, так и в частотной области. Во временной области она определяет импульсную характеристику, получаемой при реакции некоторого устройства или его мат. модели на единичный импульс (на слух воспринимаемый как щелчок). В частотной области она определяет затухание или усиление амплитуды и сдвиг фазы на отдельно взятой частоте. Переход между этими областями (в англоязычной литературе используется слово domain) осуществляется через преобразование Фурье, которое между функциями во временной и в частотной области задаёт однозначное соответствие.

Разделяют два типа цифровых фильтров:

IIR (Infinite Impulse Response) — фильтры с бесконечной импульсной характеристикой. По сути, представляют собой всё те же аналоговые фильтры, но «моделируемые» в цифре — математический аппарат в обоих случаях идентичен (в основе которого лежит преобразование Лапласа). Отсюда следует бесконечность импульсной характеристики, которую можно представить в виде суммы экспоненциально затухающих синусоид.

Вот так выглядит (обрезанная справа) импульсная характеристика рекурсивного фильтра низких частот с высокой добротностью

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

Наиболее простым является фильтр низких частот первого порядка:

FIR (Finite Impulse Response) — фильтры с конечной импульсной характеристикой. Самый известный фильтр такого рода — это скользящее среднее, а сама фильтрация осуществляется посредством линейной свёртки отсчётов фильтра с отсчётами входного сигнала. Здесь сложность уже квадратичная — однако её можно уменьшить до , если использовать алгоритм быстрой свёртки с использованием быстрого преобразования Фурье (FFT). Конечность количества отсчётов накладывает свои ограничения на форму АЧХ, приводя к пульсациям.

На самом деле математически и IIR-фильтры, и FIR-фильтры описываются одинаково — через отношение полиномов. Разница в том, что у IIR-фильтров степень и у числителя, и знаменателя небольшая, а у FIR-фильтров высокая степень числителя компенсирует единичный знаменатель. Эта тема достаточно сложная и занимает немалую часть в теории обработки сигналов; для практических задач в неё углубляться не обязательно, интересующим же можно порекомендовать вспомнить про ряд Тейлора, затем перейти к производящей функции и уже после разбираться с z-преобразованием.

А так выглядит импульсная характеристика фазолинейного FIR-фильтра низких частот

Довольно популярной практикой при проектировании фильтров является линеаризация — когда берётся АЧХ известного IIR фильтра и формируется из неё FIR с линейной фазой.

▍ Постановка задачи

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

  1. Гладкая АЧХ без пульсаций;
  2. Фазолинейность (сдвиг фаз на всех частотах равен нулю);
  3. Симметричность АЧХ в логарифмическом масштабе частот.

пункт 1) означает, что мы не закладываем в фильтр пульсации так, как это делается, например, в фильтре Чебышёва. В теории, идеально гладкая АЧХ доступна только в IIR-фильтре бесконечной длины. На практике нам достаточно лишь, чтобы уровень пульсаций не превышал уровень шума входного сигнала, который для аудио-сигналов составляет -120 дБ в целом и -90 дБ для компакт-дисков и прочих 16-битных записей.

Читайте также:  Блок питания с регулировкой напряжения и тока

Cимметричность (зеркальность) АЧХ решает две задачи:
1) ВЧ-составляющую сигнала можно получить вычитанием из исходного НЧ-составляющей — как в частотной, так и во временной области;
2) избавляет от мучительного выбора, для какого фильтра — НЧ или ВЧ — предпочтительнее более пологий или крутой спад АЧХ.

Похожее требование есть у квадратурных зеркальных фильтров, в которых задаётся симметрия в линейном масштабе частот. Но нас интересует именно логарифмический масштаб, как более естественный для человеческого слуха.

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

▍ Общий алгоритм построения

Существует несколько подходов к проектированию КИХ-фильтров, из которых наиболее простым и интуитивно-понятным является следующий:

  1. В частотном домене задаётся желаемая АЧХ фильтра;
  2. Производится обратное преобразование Фурье:
  3. Накладывается оконная функция.

▍ Аналитическое решение

Озвученным выше требованиям гладкости и симметрии (но не фазолинейности) соответствует фильтр Линквитца-Рейли. Зная формулу АЧХ , где — порядок фильтра, формулу импульсной характеристики можно получить аналитически через интеграл Фурье, а конкретные отсчёты FIR фильтра считать непосредственным её вычислением.

График фильтра в логарифмическом масштабе:

Лучший способ для вычисления интеграла Фурье — это использовать какую-нибудь систему компьютерной алгебры. Например, в Wolfram Mathematica это будет выглядеть так:

InverseFourierTransform[1/(1 + w^2), w, x] // FullSimplify

График получившейся функции, она же импульсная характеристика:

Как видно из формулы и графика, импульсная характеристика бесконечна — стремится к нулю, но не достигает его; а также симметрична относительно нуля (за счёт фазолинейности). Ограничение её во времени осуществляется путём умножения на оконную функцию, за счёт чего значения функции за пределами окна приобретают нулевые значения, не требующих участия в расчётах. Побочным эффектом этого становится искажения исходного спектра фильтра — поскольку при умножении функций их спектры сворачиваются. Итоговый спектр можно увидеть также через преобразование Фурье. Например, для прямоугольного окна получим

в результате чего АЧХ фильтра изменится на

FourierTransform[E^-Abs[x] Sqrt[Pi/2] UnitBox[x/5], x,
w] // FullSimplify

Пульсации, которые мы видим, являются следствием от наличия «боковых лепестков» оконной функции. Их можно уменьшить, взяв более гладкое и широкое окно. Среди множества разработанных окон одним из наиболее оптимальных и удобным для аналитических вычислений является окно Нуттала, определяемое как

Помножив её на импульсную характеристику фильтра, получим:

а АЧХ примет вид

frnw = FourierTransform[E^-Abs[x] Sqrt[Pi/2] NuttallWindow[x/20],
x, w] // Simplify

Как видно, величина пульсаций значительно уменьшилась. По графику видно, что «полка» фильтра слегка опустилась — в идеале это тоже нужно учитывать и компенсировать. В данном случае она составила (на нулевой частоте):

FourierTransform[E^-Abs[x] Sqrt[Pi/2] NuttallWindow[x/20], x, w] /.w->0.0

Мнимая часть здесь получилась вследствие погрешности при численных вычислениях и её можно смело отбрасывать (об этом говорит значение , поскольку точность вычислений в формате double и представляет примерно 16 десятичных цифр).

Аналогичным образом можно посчитать импульсные характеристики и для более высоких (но только целых) порядков, например:

InverseFourierTransform[1/(1+w^4), w, x] // FullSimplify

Как видно, автоматического упрощения уже недостаточно и здесь требуется ручная работа по приведению формулы к удобочитаемому виду. В итоге получается следующая формула (в комплексных числах) для импульсной характеристики фильтра Линквитца-Рейли произвольного порядка:

При вычислении за счёт сложения комплексно-сопряжённых чисел мнимая часть обнуляется и в результате получится чисто действительная функция. Эту формулу можно выписать и непосредственно в действительных числах:

▍ Переход в дискретную область

Синтез фильтра в дискретном виде осуществляется сходным образом. Принципиальное отличие состоит в том, что:

  1. Частотный диапазон фильтра ограничен частотой дискретизации по определению;
  2. На выходе мы получаем периодический сигнал — в отличие от аналитического решения, в которой функция во времени затухает к бесконечности. Этот момент особенно важен, когда проектируются фазосдвигающие фильтры. Из этого также следует, что даже если мы не будем явно накладывать оконную функцию, по факту у нас в любом случае будет наложено прямоугольное окно.

Быстрые алгоритмы существуют и для степеней тройки, и для взаимно простых чисел, и для некоторых других частных случаев — но для степеней двойки, четвёрки и восьмёрки они наиболее быстрые.

Также имеются различия (но не результат) при работе с классическим (комплексным) FFT, Real-FFT, который работает только с половиной спектра, предполагая его симметрию, и FHT (быстрое преобразование Хартли), в котором изначально и данные, и спектр определены в поле действительных чисел. Здесь будет использоваться классический FFT.

Итак, определим формулу для фильтра — для примера 4-го порядка:

Выберем частоту дискретизации (в герцах), стандартную для современных ЦАП:

Выберем частоту среза, на которой подавления будет составлять 0.5 (-6 дБ) и соответствовать нормированной частоте w=1:

Выберем размер массива (маленький, для наглядности):

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

Заполняем первую половину (положительную часть спектра) массива:

Здесь фазу на каждой нечётной частоте мы повернули на 180° для того чтобы после БПФ максимум импульса был расположен по центру массива.

Вторую половину (отрицательную часть спектра) получаем реверсом первой, за исключением крайних частот (нулевой — постоянной составляющей и максимальной — частоты Найквиста, потому что для их описания достаточно одного числа).

spectrumdata =
Join[halfspectrumdata, Reverse[halfspectrumdata[[2;;size/2]]]];

Теперь делаем обратное преобразование Фурье:

И накладываем оконную функцию — Нутталла, как и в прошлый раз:

В качестве заключительного аккорда можно нормировать полку фильтра к единице. Для этого снова потребуется прямое FFT, чтобы узнать значение постоянной составляющей:

и затем просто поделить на него значение каждого отсчёта в импульсе

Сделав FFT ещё раз, можно убедиться, что постоянная составляющая стала равной единице:

Если бы мы проектировали фильтр высоких частот, но нормализацию нужно было бы делать по амплитуде не на нулевой частоте, а на частоте Найквиста.

▍ Кусочно-непрерывные фильтры

Исходя из изначально поставленной задачи, описывать АЧХ удобнее кусочно-непрерывной функцией, в которой полосы пропускания и подавления константны и равны 1 и 0 соответственно, заданной в логарифмическом масштабе. Выводу таких функций была посвящена отдельная статья, здесь для примера мы возьмём стандартную smooth-step функцию, построенной из кубического полинома 3-го порядка:

Читайте также:  Мощный, регулируемый блок питания на lm317

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

А вот в логарифмическом масштабе — уже нет, зато хорошо видно главное отличие от фильтра Линквитца-Рейли — спад не пологий, а имеет выраженную границу справа, за которой громкость (в децибелах) равна минус бесконечности. Дополнение до единицы даёт ВЧ-фильтр (жёлтым цветом), который и обладает желаемой симметрией:

На большем диапазоне громкости это свойство фильтров выглядит ещё более наглядным:

В отличие от классических фильтров, где параметризация крутизны спада АЧХ через порядок обусловлена их схемотехнической реализацией, здесь мы можем оперировать непосредственно шириной полосы сопряжения, задавая её в октавах от частоты раздела в -6 дБ. Для этого потребуется дополнительная функция для перевода логарифмического масштаба в линейный, значение которой будет передаваться в smooth-step функцию для определения амплитуды на частоте . Её можно определить несколькими путями:

1) через граничные частоты

2) через центральную частоту и ширину (в октавах)

3) через центральную и верхнюю граничную частоту

4) через центральную и нижнюю граничную частоту

Теперь, используя некоторую функцию фильтра, можно получить его импульсную характеристику по вышеописанному алгоритму, который можно определить отдельной функцией:

SymmetricFilterImpulseResponse[f0_, width_, size_, samplerate_,
ClipFunction_, WindowFunction_] :=
Table[ClipFunction[(2 Log[(index samplerate)/(f0 size)])/
(width Log[2]) ] (-1)^index, ] //
Join[#, Reverse[#[[2 ;; size/2]]]] & // InverseFourier[#]
Table[WindowFunction[(index-1)/size-1/2], ] & // #/
Fourier[#, FourierParameters -> ][[1]] & // Re

В разных реализациях свёртки через ДПФ существуют разные подходы к нормализации импульсной характеристики по амплитуде. Здесь она делается не отдельной операцией, а заданием опции FourierParameters при выполнении ДПФ.

Передав в эту функцию наш кубический фильтр и окно Нуттала получим

SymmetricFilterImpulseResponse[1000, 0.5, 512, 48000, FilterCubicClip, NuttallWindow]

▍ Несколько интересных фильтр-функций

Помимо уже рассмотренной кубической, можно придумать множество и других функций для фильтров, основанных на какой-нибудь математической идее. Также имеет значение, насколько легко может быть вычислена (и может ли быть вычислена вообще) обратная функция.

Параболическая, наименее вычислительно затратная. При необходимости, обратная к ней функция находится довольно элементарно:

Кубическая. Обратная к ней функция уже немного контринтуитивна (а если вам интересно, откуда здесь взялись тригонометрические функции — тогда сюда):

Используя полином 5-ой степени, с двумя нулевыми производными на границах сопряжения. Обратная функция здесь уже в элементарных функциях не выражается (потому и не приведена):

Используя полином 13-ой степени — с дополнительным «выпрямлением»:

Используя рациональный полином — даёт более гладкую характеристику и более простую обратную функцию, чем у кубической:

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

Линейно спадающая в логарифмическом масштабе:

Более плавный вариант предыдущего с возможностью регулировки «жёсткости». Здесь уже обратная функция для произвольного n в элементарных функциях не выражается:

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

С максимально быстро затухающей импульсной характеристикой — что при достаточно больших () позволяет обойтись без оконной функции вообще. Сама формула получена модуляцией аргумента кубической функции, т.е. :

Здесь — это специальная функция (функция ошибок), определяемая как интеграл от гауссианы. Вместо неё можно использовать функцию гиперболического тангенса — импульсная характеристика также будет достаточно быстро затухать, хоть и не так быстро, как в предыдущем случае. Причём при больших n визуально она будет походить на Линквитца-Рейли и это совсем не случайно —

▍ Заключение

Несмотря на очевидность постановки задачи, мне не встречалось ничего подобного ни в специализированных решениях, ни в мат. пакетах общего назначения. Возможно, что это связано с инерцией мышления в стиле «всё, что можно придумать, уже придумано». Однако это не так хотя бы потому, что сейчас доступны технологии и вычислительные мощности, о которых всего лишь 50 лет назад инженеры и математики не могли и мечтать.

В качестве proof of concept фильтры по такой методике были реализованы в плагине для популярного в аудиофильских кругах плеера foobar, извлекающего канал сабвуфера из стереосигнала (официальный репозиторий); и ещё в одном, просто для демонстрации звучания (неофициальный репозиторий). Для реализации таких фильтров в «мультиампинг»-системах можно использовать готовые решения, позволяющие загружать заранее посчитанные импульсные файлы. При отсутствии мат.пакетов их можно посчитать даже в excel-е (используя преобразование Фурье из пакета анализа).

Насколько такие фильтры «звучат» лучше классических решений — вопрос уже аудиофильский, но, как минимум, разницу в звучании обеспечивают вполне объективную и обеспечивают возможность для более точной настройки.

Цифровые фильтры без умножений

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

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

Введение

В операциях цифровой обработки сигналов особое внимание уделяется цифровой фильтрации, которая по объему вычислений в среднем занимает от 20 до 60%. В узком смысле цифровой фильтр — это частотно-избирательная цепь, обеспечивающая селекцию цифровых сигналов по частоте [2]. После выполнения цифровой фильтрации мы, как правило, получаем интересующий нас сигнал, то есть сигнал, несущий нужную нам информацию в виде, удобном для последующей обработки. Соответственно, к параметрам цифровых фильтров в современных системах цифровой обработки сигналов начинают предъявляться повышенные требования. Порядки фильтров нередко достигают тысячи и более. Это ведет к увеличению объема вычислений, а значит, и к резкому росту аппаратных затрат. При синтезе цифровых фильтров наибольшие затраты времени и оборудования приходятся на операции умножения [2]. Таким образом, задача минимизации времени вычислений и уменьшения аппаратных затрат сводится к минимизации количества умножений, требуемых для вычисления очередного отфильтрованного отсчета.

Читайте также:  Маломощный лабораторный источник питания на LM317

Структура фильтра

Как известно, классический цифровой фильтр описывается выражением [3]:

где y(n) — сигнал на входе фильтра; x(n) — сигнал на выходе фильтра; ai , bi — коэффициенты фильтра.

Все цифровые фильтры делятся на два обширных класса: нерекурсивные — фильтры с конечной импульсной характеристикой (КИХ) и рекурсивные — фильтры с бесконечной импульсной характеристикой (БИХ) [3]. Мы будем рассматривать только КИХ-фильтры, фазовая характеристика которых, в отличие от БИХ-фильтров, линейна. Для КИХ-фильтров выражение (1) принимает следующий вид [3]:

Таким образом, задача синтеза КИХ-фильтра сводится к вычислению коэффициентов bi такого фильтра.

Проблема минимизации количества операций умножения может быть решена различными путями. Одним из них является использование неклассических структур цифровых фильтров, которые бы описывались выражениями, отличными от (2), и в которых операция умножения исключена вообще. Такие структуры позволяют для некоторых классов цифровых фильтров обеспечить улучшенные характеристики при их реализации. Рассмотрим одну из упомянутых структур. В ее основе — субструктура, изображенная на рис. 1 [1].

Рис. 1. Структурная схема элементарного звена: K — количество задержек; W — весовой коэффициент (количество сумматоров)

Каскадным соединением субструктур, представленных на рис. 1, можно синтезировать фильтры, обладающие различными амплитудно-частотными характеристиками, которые бы удовлетворяли заданным требованиям по таким критериям, как ширина полосы пропускания фильтра, подавление боковых лепестков, прямоугольность и т. д.

Структуры и использование фильтров без умножений

Фильтры, синтезированные на основе структур, похожих на описанную выше, широко используются в микросхемах digital down converter (DDC) различных фирм. Например, таких как Analog Devices и Texas Instruments. В микросхеме AD6620 фирмы Analog Devices фильтры типа Cascaded Integrator Comb (CIC) CIC2 и CIC5 являются фильтрами без умножителей. CIC2 и CIC5 — фильтры второго и пятого порядков, состоящие, соответственно, из двух и пяти каскадов. Эти фильтры используются для программируемого понижения частоты входных отсчетов.

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

Коэффициенты фильтра

Попробуем получить аналитическое выражение для коэффициентов фильтра. Рассмотрим каскадное соединение описанных выше звеньев. Возьмем три звена, у каждого из которых W= 1, K = 1 и суммирование. Тогда на выходе первого звена получим:

На выходе второго звена получим:

На выходе третьего звена получим:

Значит, исходя из (3), (4), (5), сигнал на выходе очередного каскада определяется по формуле:

или, если произвести Z-преобразование и перейти в Z-область[2]:

где k — число задержек последнего звена; Hi–1(Z) — Z-преобразование предыдущего передаточного звена.

Для коэффициентов же bi данного фильтра будет справедливо выражение:

где I — порядковый номер каскада; k — число задержек в i-звене.

Как известно, bi представляют собой коэффициенты импульсной характеристики фильтра [2]. График импульсной характеристики для каскадного соединения четырех звеньев с K = 1 приведен на рис. 2. Очевидно, что при добавлении звеньев количество «ступенек» импульсной характеристики, соответствующих коэффициентам передаточной функции фильтра, будет расти.

Покажем, что при добавлении элементарного звена с вычитанием появляются отрицательные коэффициенты. На рис. 3 приведен график импульсной характеристики каскадного соединения четырех звеньев с K = 1 иодного звена с вычитанием, у которого K = 1.

Влияние структуры на параметры фильтра

Обратимся к амлитудно-частотным характеристикам (АЧХ) элементарных звеньев. В логарифмическом масштабе при каскадном соединении звеньев их АЧХ будут складываться. Значит, задача получения нужной АЧХ будет сводиться к сложению АЧХ типовых структур с различными параметрами.

Покажем, как при изменении каждого из трех параметров элементарного звена изменяется его АЧХ. На рис. 4 представлен график АЧХ элементарного звена, для которого K = 1, W = 1. При увеличении количества задержек происходит увеличение экстремумов АЧХ — пропорционально числу задержек. При этом частотный диапазон делится на равные части в соответствии с количеством задержек. Максимумы и минимумы АЧХ в данном случае чередуются через одинаковые частотные промежутки. На рис. 5 и 6 приведены примеры элементарных звеньев с K = 2, K = 4.

Графики, представленные на рис. 4–6, позволяют определить структуру фильтра нижних частот (ФНЧ). Сначала синтезируется звено с самым большим количеством задержек. Тогда задача синтеза ФНЧ сводится к уничтожению всех частот, лежащих правее полученной.

Изменение знака сумматора с плюса на минус приводит к появлению экстремума в области нулевой частоты и так называемому дифференцированию АЧХ. После точки перегиба АЧХ становится монотонно возрастающей. Это свойство позволяет синтезировать полосовые фильтры. На рис. 7 и 8 приведены примеры элементарного звена с вычитанием, у которого K = 1, а также каскадного соединения элементарных звеньев с K = 1, K = 2, K = 4 и вычитанием, у которого K = 1.

Рис. 8. АЧХ каскадного соединения элементарного звена с вычитанием, у которого K = 1, и звеньев с K = 1, K = 2, K = 4

Увеличение числа сумматоров, по сути дела, добавляет дополнительные элементарные звенья к заданному звену. Единственное отличие таких звеньев будет в том, что у них K = 0, то есть это звенья без задержек. АЧХ при увеличении числа сумматоров становится менее крутой, но форма АЧХ существенно не изменяется. На рис. 9 приведен пример элементарного звена с W = 2 и K = 4.

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

Синтез фильтров с оптимальными параметрами

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

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

Читайте также:  Бегущие огни на светодиодах

Другой путь — это изначальное задание точек АЧХ фильтра и получение по заданным точкам АЧХ фильтра его структуры, то есть комбинации элементарных звеньев с известными АЧХ. Естественно, что это не всегда возможно сделать точно, и структура фильтра может быть получена с некоторыми допущениями и приближениями. Кроме того, это применимо не для всех классов фильтров, в силу рассмотренных нами особенностей элементарных звеньев.

Рассмотрим фильтр, полученный каскадным соединением элементарных звеньев, подобранных так, чтобы получить максимальное подавление его боковых лепестков. График АЧХ этого фильтра представлен на рис. 10.

Как видно, мы получили фильтр нижних частот (ФНЧ) с конечной импульсной характеристикой (КИХ). У этого фильтра количество элементарных звеньев равно 44. Ширина полосы пропускания по критерию 3 дБ равна 0,062 ω. Подавление боковых лепестков — 109 дБ.

Синтезируем фильтр с вышеперечисленными параметрами с помощью утилиты Filter Design среды Matlab.

Получился фильтр нижних частот (ФНЧ) с конечной импульсной характеристикой (КИХ), порядок фильтра равен 74. Его параметры и АЧХ представлены на рис. 11.

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

Во-первых, фильтры этого класса имеют строго определенные частоты подавления, равные 1/2 ω, 1/3 ω, 1/4 ω, 1/6 ω и т. д.

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

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

На всех графиках (рис. 12–15) приведена нормированная частота ω/ωд.

Нерекурсивные частотные цифровые фильтры

Материал из Национальной библиотеки им. Н. Э. Баумана
Последнее изменение этой страницы: 18:22, 18 июня 2016.

Авторство
Чичварин Н. В.
Согласовано: 03.06.2016
Статья по учебной дисциплине
Название дисциплины:

Обнаружение и распознавание сигналов

2. Анализ регулярных сигналов

2.9 Цифровая обработка сигналов в электронном тракте. Квантование и дискретизация. Модель цифрового тракта на основе Z-преобразования.

Нерекурсивные фильтры всегда имеют конечную импульсную характеристику (КИХ-фильтры), при этом реализуют алгоритм свертки двух функций:

Содержание

Общие сведения

Свойства любого фильтра однозначно определяют его частотная (frequency response) и фазовая характеристики. Они показывают, какое влияние фильтр оказывает на амплитуду и фазу различных гармоник обрабатываемого сигнала.

К наиболее известным типам нерекурсивных цифровых фильтров (НЦФ) относятся частотные фильтры, алгоритм которых для симметричных НЦФ, не изменяющих фазу сигналов, имеет вид

Типы фильтров

В зависимости от вида частотной характеристики выделяют три основных группы частотных фильтров:

  • ФНЧ − фильтры низких частот (low-pass filters) – пропускание низких и подавление высоких частот во входном сигнале;
  • ФВЧ − фильтры высоких частот (high-pass filters) – пропускание высоких и подавление низких частот;
  • ПФ − полосовые фильтры, которые пропускают (band-pass filters) или подавляют (band-reject filters) сигнал в определенной частотной полосе.

Среди последних в отдельную группу иногда выделяют РФ – режекторные фильтры, понимая под ними фильтры с подавлением определенной гармоники во входном сигнале, и СФ – селекторные фильтры, обратные РФ.

Если речь идет о подавлении определенной полосы частот во входном сигнале, то такие фильтры называют заградительными. Ни теоретического, ни практического интереса к методам их расчета обычно не проявляется, так как их частотная характеристика обычно задается инверсией характеристики полосового фильтра ( 1 − H Π ( ω ) ) (omega)) ,!> и каких-либо дополнительных особенностей при своем проектировании не имеет.

Схематические частотные характеристики фильтров приведены на рисунке 1.1. Между частотными интервалами пропускания и подавления сигнала существует зона, которая называется переходной. Ширина переходной зоны определяет резкость характеристики фильтра. В этой зоне амплитудная характеристика монотонно уменьшается (или увеличивается) от полосы пропускания до полосы подавления (или наоборот).

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

Отсюда, условие инверсии симметричного низкочастотного фильтра в высокочастотный:

Пример обращения и спектры фильтров приведены на рис. 1.2 (в правой части главных диапазонов).

Последнее означает смену знака всех нечетных гармоник передаточной характеристики фильтра и, соответственно, всех нечетных членов фильтра:

Пример частотного реверса приведен на рис. 1.3. Физическую сущность такой операции инверсии спектра легко понять на постоянной составляющей сигнала. При изменении на противоположный знака каждого второго отсчета постоянной величины это постоянной значение превращается в “пилу”, частота которой равна частоте Найквиста главного частотного диапазона (отсчеты по амплитудным значениям этой частоты), равно как и наоборот, отсчеты гармоники сигнала на частоте Найквиста (знакочередующиеся в силу сдвига по интервалам дискретизации на π ) превращаются в постоянную составляющую.

Полосовой фильтр может реализоваться последовательным применением ФНЧ и ФВЧ с соответствующим перекрытием частот пропускания. В математическом представлении это означает последовательную свертку массива данных с массивами коэффициентов h H – низкочастотного, и h B – высокочастотного фильтров:

v k = h H ( n ) × s ( k − n ) , ) , y k = h B ( n ) × v k = h H ( n ) × h B ( n ) × s ( k − n )

Так как операция свертки коммутативна, то вместо отдельных массивов коэффициентов ФНЧ и ФВЧ их сверткой может быть определен непосредственно массив коэффициентов полосового фильтра:

h n = h H ( n ) × h B ( n )

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

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

Методика расчетов НЦФ

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

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

Наиболее простой является методика расчетов программных двусторонних симметричных фильтров без изменения фазы выходного сигнала относительно входного. В самом общем виде она включает:

  1. Задание идеальной амплитудно-частотной характеристики передаточной функции фильтра. Термин идеальности понимается здесь в том смысле, что на характеристике указываются полосы пропускания и подавления частот с коэффициентами передачи 1 и 0 соответственно без переходных зон.
  2. Расчет функции импульсного отклика идеального фильтра (обратное преобразование Фурье частотной характеристики фильтра). При наличии скачков функций на границах пропускания/подавления импульсный отклик содержит бесконечно большое количество членов.
  3. Ограничение функции отклика до определенного количества членов, при этом на передаточной характеристике фильтра возникает явление Гиббса – осцилляции частотной характеристики с центрами на скачках.
  4. Для нейтрализации явления Гиббса производится выбор весовой функции и расчет ее коэффициентов, на которые умножаются коэффициенты функции отклика фильтра. Результатом данной операции являются значения коэффициентов оператора фильтра (рабочий импульсный отклик фильтра). По существу, операции 3 и 4 представляют собой усечение ряда Фурье динамического (временного) представления передаточной функции фильтра определенной весовой функцией (умножение на весовую функцию).
  5. С использованием полученных значений коэффициентов оператора фильтра производится построение его частотной характеристики и проверяется ее соответствие поставленной задаче.

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

Фильтры с линейной фазовой характеристикой

Несколько сложнее расчет каузальных (односторонних) частотных фильтров, для которых требуется обеспечить линейность фазово-частотной характеристики для исключения изменения гармонии сочетания частотных составляющих сигнала на его выходе по отношению к входу. Чтобы фильтр имел линейную фазовую характеристику необходимо обеспечить выполнение условия:

Оно выполняется, если импульсная характеристика фильтра имеет положительную симметрию:

При этом фазовая характеристика будет определяться длиной фильтра:

Частотная характеристика фильтра:

Следует также учитывать, что частотную характеристику типа 2 нельзя использовать для проектирования фильтров верхних частот, т.к. она всегда равна нулю на частоте Найквиста.

Идеальные частотные фильтры

Импульсная реакция фильтра

Импульсная реакция фильтра (коэффициенты оператора) находится обратным преобразованием Фурье заданной передаточной функции H ( ω ) . В общем случае:

Введение

Во многих цифровых устройствах для преобразования аналоговых сигналов используется АЦП. Часто аналоговые сигналы содержат нежелательный высокочастотный шум.
Чтобы “очистить” сигнал от этих шумов применяются аналоговые RC фильтры низких частот, которые устанавливаются после источника сигнала. Такой подход не всегда идеален и практичен. Например, для больших постоянных времени требуются большие значения R и C.
В качестве альтернативы, можно “очистить” зашумленный сигнал с помощью цифрового эквивалента аналогового RC фильтра нижних частот.

Алгоритм цифрового фильтра

По сути, программа этого цифрового фильтра состоит всего из двух строчек на Си.:

Dacc = Dacc + Din – Dout
Dout = Dacc/K

где Dout – выходное значение фильтра, Din – входное значение фильтра, K – постоянный коэффициент, который рассчитывается по формуле:

K = T x SPS

где T – это постоянная времени фильтра, SPS – частота дискретизации АЦП.

Dacc и Dout должны сохранять свои значения, после выполнения алгоритма. Если алгоритм реализовать в виде функции, то эти переменные можно просто сделать статическими.

Для 8-ми разрядных входных данных алгоритм цифрового фильтра в Си коде может выглядеть так:

Пример реализации цифрового фильтра на mega16

Для наглядности рассмотрим реальный пример использования этого алгоритма в микроконтроллере AVR atmega16. Допустим мы хотим сымитировать RC фильтр с постоянной времени 0.001 c. Схема представлена на рисунке ниже.

аналоговый RC фильтр низкой частоты

R1 = 10 кОм
C1 = 0.1 мкФ
Trc = R1*C1 = 10000 Ом * (0.1/1000000) Ф = 0.001 с
F = 1/(2*Pi*R1*C1) = 1/(6.28 * Trc) = ~50 Гц

Тактовая частота модуля АЦП в микроконтроллерах AVR зависит от его тактовой частоты и внутреннего предделителя. Допустим, наш микроконтроллер тактируется от внутреннего генератора с частотой 8 МГц, а предделитель в модуле АЦП установлен равным 64. Тогда тактовая частота модуля АЦП будет равна:

Fadc = Fcpu/Pre = 8000000/64 = 125000 Гц

Из этой частоты можно рассчитать частоту дискретизации АЦП при работе в режиме непрерывного преобразования. Она равна отношению тактовой частоты АЦП к количеству тактов, которые требуются для выполнения одного преобразования. По даташиту можно узнать, что одно преобразование выполняется за 13 тактов (если это не первое преобразование).

Fs = Fadc/n = 125000/13 = ~9600 Гц

Итак, частота дискретизации равна 9600 Гц, а постоянная времени 0.001 с. Коэффициент фильтра будет равен:

K = SPS x T = 9600 * 0.001 = 9.6 = ~ 10

Теперь все данные известны и можно написать тестовую программу для проверки алгоритма.

Я сделал эту программу очень простой. АЦП работает в режиме непрерывного преобразования. В прерывании 8-ми разрядный результат преобразования обрабатывается алгоритмом и записывается в порт C. К порту C подключен схема R-2R ЦАПа на основе резисторов, чтобы можно было сравнить полученный сигнал с сигналом от аналогового RC фильтра. Тактовая частота микроконтроллера atmega16 – 8МГц, коэффициент предделителя АЦП – 64. Тестовая схема показана на рисунке ниже.

схема для тестирования цифрового фильтра


Код программы

Результат моделирования

Результат моделирования программы в Proteus`e показан на рисунке ниже. Красный сигнал – это входной меандр частотой 200 Гц, синий – сигнал на выходе RC фильтра с постоянной времени 0.001 с, а желтый – сигнал обработанный микроконтроллером. Он имеет ступенчатую форму, так как после ЦАПа не подвергался фильтрации.

Как видно из рисунка форма сигнала микроконтроллера достаточно точно повторяет сигнал с выхода RC фильтра.

результаты моделирования цифрового фильтра в Proteus`e

Заключение

Для наибольшего быстродействия коэффициент К лучше выбирать кратным степени 2 (например 2, 4, 8..), тогда компилятор будет заменять операцию деления сдвигами. В противном случае при высокой частоте дискретизации, микроконтроллер может не успевать рассчитывать следующее выходное значение фильтра. Я столкнулся с этим при моделировании.
Также необходимо учесть тот момент, что при больших значениях коэффициента К, переменная Dacc должна иметь достаточную разрядность, чтобы не наступало ее переполнение.

Рейтинг
( Пока оценок нет )
Понравилась статья? Поделиться с друзьями:
Добавить комментарий

;-) :| :x :twisted: :smile: :shock: :sad: :roll: :razz: :oops: :o :mrgreen: :lol: :idea: :grin: :evil: :cry: :cool: :arrow: :???: :?: :!: