Моделирование цифровой обработки сигналов ЦОС в MATLAB. Часть 2. Синтез оптимальных цифровых БИХ-фильтров программными средствами MATLAB

БИХ–фильтр описывается передаточной функцией общего вида:

и при (N−1)≤(M−1) (по умолчанию) имеет порядок R = M−1.

Сложность БИХ–фильтра определяется порядком R передаточной функции (1).

БИХ–фильтры характерны следующими особенностями:

• нелинейной ФЧХ;
• необходимостью проверки на устойчивость.

При синтезе частотно–избирательных БИХ–фильтров с существенно нелинейной ФЧХ последняя обычно не контролируется, и требования задаются к АЧХ. Они не отличаются от требований к АЧХ КИХ–фильтров, за тем исключением, что для рассматриваемого далее метода синтеза значение АЧХ в полосе пропускания не превышает единицы. Кроме того, для БИХ–фильтров требования задаются к характеристике затухания — АЧХ (дБ) (см. формулу (6) в [6]) и включают в себя:

• частоту дискретизации fд (Гц);
• граничные частоты полос пропускания (ПП) и полос задерживания (ПЗ), для которых введены условные обозначения:
  • fχ — граничная частота ПП для ФНЧ и ФВЧ;
  • fk — граничная частота ПЗ для ФНЧ и ФВЧ;
  • f−χ, fχ — левая и правая граничные частоты ПП для ПФ и РФ;
  • f−k, fk — левая и правая граничные частоты ПЗ для ПФ и РФ;
  • amax (дБ) — максимально допустимое затухание в ПП;
  • amin (дБ) — минимально допустимое затухание в ПЗ.

  В статье рассматривается синтез оптимальных БИХ–фильтров методом билинейного Z–преобразования на основе аналоговых фильтров–прототипов (АФП).

  Идея синтеза БИХ–фильтров на основе АФП возникла из желания воспользоваться давно известными и хорошо себя зарекомендовавшими методами синтеза аналоговых фильтров. Обоснование такой возможности вытекает из следующих положений:

  • передаточные функции АФП и БИХ–фильтров — дробно–рациональные;
  • импульсные характеристики АФП и БИХ–фильтров — бесконечные.

  Для того чтобы подчеркнуть контраст типа фильтра (аналоговый или цифровой), будем использовать аббревиатуры АФП и ЦФ, по умолчанию подразумевая под ЦФ БИХ–фильтр.

  Процедура синтеза БИХ–фильтра

  Процедура синтеза ЦФ на основе АФП включает в себя [4]:

  1. Задание требований к АЧХ ЦФ.
  2. Выбор метода синтеза.
  3. Формирование требований к АЧХ АФП.
  4. Выбор типа аппроксимирующей функции. Четырем типам аппроксимирующих функций соответствуют четыре разновидности аналоговых (и цифровых) фильтров:
  • Баттерворта (Butterwhorth) — с АЧХ, максимально плоской в ПП и монотонной в ПЗ;
  • Чебышева I рода (Chebyshov Type I) — с АЧХ, равноволновой в ПП и монотонной в ПЗ;
  • Чебышева II рода (Chebyshov Type II) — с АЧХ, максимально плоской в ПП и равноволновой в ПЗ;
  • Золотарева – Кауэра (эллиптические фильтры) (Eleptic) — с АЧХ, равноволновой в ПП и ПЗ.

  Для лучшего понимания синтеза в MATLAB ЦФ на основе АФП коротко познакомимся с синтезом АФП.

  Синтез аналоговых фильтров

  Синтез частотно–избирательных АФП Баттерворта, Чебышева I и II рода и Золотарева – Кауэра выполняется соответственно с помощью функций:

  Здесь R — порядок АФП; Wn — вектор частот среза в шкале ω = 2πf (рад/с), содержащий один элемент — для ФНЧ и ФВЧ и два — для ПФ и РФ (частотами среза называют частоты, на которых нормированная АЧХ АФП Â(f) равна 1/√2≈0,707, а Â(f) (дБ) — 3 дБ); rp, rs — соответственно максимально и минимально допустимые затухания amax (дБ) в ПП и amin (дБ) в ПЗ для Â(f) (дБ); ftype — параметр, указывающий тип избирательности и принимающий значения: ‘high’ — для ФВЧ; ‘stop’ — для РФ; по умолчанию (если значение параметра не задано явно) — для ФНЧ и ПФ; ‘s’ — признак аналогового фильтра, при его отсутствии по умолчанию подразумевается ЦФ; bs, as — векторы коэффициентов числителя и знаменателя передаточной функции АФП Ha(p) в порядке возрастания степеней p; as(1) = 1.

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

  Как правило, при синтезе АФП порядок фильтра (R) и частоты среза (Wn) заранее неизвестны. Их можно определить по требованиям к АЧХ с помощью следующих функций, соответственно для АФП Баттерворта, Чебышева I и II рода и Золотарева – Кауэра:

  Здесь R — минимальный порядок при заданных требованиях, соответствующий оптимальному АФП; Wp, Ws — соответственно векторы граничных частот ПП и ПЗ в порядке их следования слева направо в шкале частот ω = 2πf (рад/с). Остальные параметры были определены ранее.

  Пример 1

  Синтезировать оптимальные АФП Баттерворта, Чебышева I и II рода и Золотарева – Кауэра по заданным требованиям к АЧХ ФНЧ (см. табл. 3 в [6]). Значения amax = 0,4455 дБ и amin = 40 дБ (rp и rs) были вычислены в [6] в примере 1:

  Выведем значения порядков R1, R2, R3, и R4 соответственно оптимальных ФНЧ Баттерворта, Чебышева I и II рода и Золотарева – Кауэра:

  Наименьший порядок имеет ФНЧ Золотарева – Кауэра.

  Построим графики АЧХ аналоговых ФНЧ Баттерворта, Чебышева I и II рода и Золотарева – Кауэра на густой сетке частот (выберем 1000 точек) и выведем их в основной полосе частот [0; fд/2] ЦФ при частоте дискретизации 8000 Гц (для сравнения с ним впоследствии).

  Для построения графиков АЧХ АФП используем функцию:

  где bs, as — коэффициенты числителя и знаменателя передаточной функции АФП; W — вектор, задающий сетку частот в шкале ω = 2πf (рад/с).

  Выведем значения АЧХ всех АФП в одинаковом диапазоне [0;1] по оси ординат с помощью функции ylim([0 1]) (рис. 1):

  Синтез БИХ–фильтров методом билинейного Z–преобразования

  Отображение p–плоскости в z–плоскость выполняется в соответствии с формулой билинейного Z–преобразования:

  где γ = 1/T, получаемой из стандартного Z–преобразования z = e pT → p = (1/T)lnz, путем разложения lnz в ряд Тейлора:

  и сохранения первого члена.

  Формула (2) позволяет представить передаточную функцию ЦФ H(z) на основе передаточной функции АФП Ha(p):

  Используя (2), выражаем z через p:

  И подставляя z = re jˆω и p = jΩ, где jΩ — обозначение оси частот АФП (во избежание путаницы), получаем:

  Откуда имеем нелинейные зависимости между частотами АФП и ЦФ:

  Согласно (4), ось частот jΩ p–плоскости, как и при стандартном Z–преобразовании, отображается на z–плоскости в единичную окружность (радиус равен единице), однако каждому ее обороту (изменению нормированной частоты на Δˆω = 2π), а именно: …, −3π < ˆω < −π, −π < ˆω < π, π < ˆω < 3π, …, соответствует не отрезок оси, как при стандартном Z–преобразовании, а вся ось jΩ (так как зависимость между частотами определяется функцией arctg: …, −∞ < Ω < ∞, −∞ < Ω < ∞, −∞ < Ω < ∞…

  Связь между частотными характеристиками АФП и ЦФ, соответственно H(e jΩT ) и Ha(jΩ), имеет вид [1, 4]:

  При этом частотная характеристика АФП Ha(jΩ), бесконечная в шкале частот Ω, в шкале частот ω, согласно (6), сжимается в отрезок Δω = ωд, то есть становится финитной. Соответственно, частотная характеристика ЦФ H(e jΩT ), согласно (7), оказывается равной (с точностью до множителя 1/T) бесконечной сумме копий финитных частотных характеристик АФП длины Δω = ωд, сдвинутых друг относительно друга на частоту ωд. Вследствие этого элайсинг (Aliasing) принципиально отсутствует, и АЧХ ЦФ и АФП в основной полосе частот [0; ωд/2] совпадают (с учетом нелинейной зависимости между частотами).

  «Платой» за отсутствие элайсинга, помимо нелинейной зависимости между частотами АФП и ЦФ, будет полное несовпадение их импульсных характеристик и ФЧХ.

  Метод билинейного Z–преобразования реализуется по–разному, в зависимости от поставленной задачи, а именно:

  • Если ЦФ синтезируется на основе имеющегося АФП (копирует его АЧХ с учетом нелинейного соотношения между частотами), то в этом случае удобно воспользоваться функцией bilinear следующих основных форматов:

  Здесь bs, as — векторы коэффициентов числителя и знаменателя передаточной функции АФП Ha(p) в порядке возрастания степеней p; as(1)=1; Fs — частота дискретизации fд (Гц); Fp — необязательный параметр — частота f (Гц), на которой значения АЧХ АФП и ЦФ должны совпадать; b, a — векторы коэффициентов числителя и знаменателя передаточной функции ЦФ H(z) (1) в порядке возрастания отрицательных степеней z; a(1)=1; q, p, K — соответственно векторы нулей и полюсов и коэффициент усиления передаточной функции, представленной в виде произведения простейших множителей:

  где zok, z*k — соответственно k–е нуль и полюс передаточной функции (1), а b0 — коэффициент усиления; qs, ps, Ks — аналогичные параметры для передаточной функции АФП.

  • Если ЦФ синтезируется непосредственно по заданным требованиям к АЧХ, то в этом случае процедура синтеза подобна рассмотренной ранее для АФП, более того, используются те же функции, но без параметра ‘s’:

  Здесь R — порядок ЦФ; WDn — вектор нормированных частот среза в шкале ˆ f (частотами среза называют частоты, на которых нормированная АЧХ ЦФ Â(f) равна 1/√2≈0,707, а Â(f) (дБ) — 3 дБ), содержащий один элемент для ФНЧ и ФВЧ, равный:

  где f0 — абсолютная частота среза, и два — для ПФ и РФ, равные:

  где f01, f02 — абсолютные частоты среза; rp, rs — соответственно максимально и минимально допустимые затухания amax (дБ) в ПП и amin (дБ) в ПЗ для Â(f) (дБ) (они совпадают с допустимыми отклонениями для АФП); ftype — параметр, указывающий тип избирательности и принимающий значения: ‘high’ — для ФВЧ; ‘stop’ — для РФ; по умолчанию (если значение параметра не задано явно) — для ФНЧ и ПФ; b, a — векторы коэффициентов числителя и знаменателя передаточной функции ЦФ (1) в порядке возрастания отрицательных степеней z; где a(1)=1.

  Примечание. Здесь и далее в обозначениях частот ЦФ добавлена буква «D» от слова «Digital». Согласно (5)–(6) зависимость между частотами WDn и Wn нелинейная.

  При синтезе ЦФ Баттерворта, Чебышева I и II рода и Золотарева – Кауэра по методу билинейного Z–преобразования свойство оптимальности ЦФ сохраняется: синтезируемый ЦФ, подобно АФП, при заданных требованиях к АЧХ (характеристике затухания) имеет минимальный порядок.

  Расчет минимального порядка R и вектора частот среза WDn для ЦФ Баттерворта, Чебышева I и II рода и Золотарева – Кауэра выполняется с помощью тех же функций, что и для АФП, но без параметра ‘s’:

  Здесь WDp, WDs — соответственно векторы граничных нормированных частот ПП и ПЗ в порядке их следования слева направо в шкале ˆ f . Остальные параметры были определены ранее.

  При синтезе ЦФ с помощью данных функций в MATLAB реализуется алгоритм билинейного Z–преобразования, а именно:

  • по требованиям к АЧХ ЦФ автоматически формируются требования к АЧХ АФП путем пересчета граничных частот по формуле (5);
  • синтезируется АФП;
  • в соответствии с (3) передаточная функция АФП Ha(p) преобразуется в передаточную функцию ЦФ H(z) (1).

  Подобно функции bilinear, выходными параметрами функций butter, cheby1, cheby2 и ellip могут быть [q,p,K].

  Приведем примеры синтеза оптимальных БИХ–фильтров ФНЧ и ПФ непосредственно по заданным требованиям к АЧХ (дБ) ЦФ.

  Пример 2

  Заданы требования к АЧХ ФНЧ (табл. 3 и пример 2 в [6]). Значения amax = 0,4455 дБ и amin = 40 дБ (rp и rs) были вычислены в [6] в примере 1. Синтезировать оптимальные БИХ–фильтры Баттерворта, Чебышева I и II рода и Золотарева – Кауэра методом билинейного Z–преобразования:

  Выведем рассчитанные значения порядков R1, R2, R3, и R4 соответственно оптимальных ФНЧ Баттерворта, Чебышева I и II рода и Золотарева – Кауэра:

  Поскольку свойство оптимальности синтезируемых ЦФ сохраняется, их порядки совпадают с порядками соответствующих АФП (переменная R в примере 1).

  Построим графики АЧХ БИХ–фильтров ФНЧ Баттерворта, Чебышева I и II рода и Золотарева – Кауэра на густой сетке частот (выберем 1000 точек) в основной полосе [0; fд/2] и одинаковом диапазоне [0;1] по оси ординат, установленном с помощью функции ylim([0 1]). АЧХ рассчитывается с помощью функции freqz (рис. 2):

  Пример 3

  Заданы требования к АЧХ ПФ (табл. 4 и пример 3 в [6]). Значения amax = 0,4455 дБ и amin = 40 дБ (rp и rs) были вычислены в [6] в примере 1. Синтезировать оптимальные БИХ–фильтры Баттерворта, Чебышева I и II рода и Золотарева – Кауэра методом билинейного Z–преобразования.

  Параметры WDp и WDs представляют собой векторы из двух элементов:

  Выведем рассчитанные значения порядков R1, R2, R3, и R4 соответственно оптимальных ПФ Баттерворта, Чебышева I и II рода и Золотарева – Кауэра:

  Наименьший порядок имеет ФНЧ Золотарева – Кауэра.

  Графики АЧХ БИХ–фильтров (рис. 3) ПФ Баттерворта, Чебышева I и II рода и Золотарева – Кауэра строятся так же, как для ФНЧ (пример 2).

  Анализ БИХ–фильтра

  В состав MATLAB входит программа GUI FVTool (Filter Visualization Tool — средства визуализации фильтра), предназначенная для анализа характеристик синтезированных ЦФ в окне Figure : Filter Visualization Tool, обращение к которой производится с помощью функции fvtool:

  📎📎📎📎📎📎📎📎📎📎