MATLAB для DSP
В. Анохін, А. Ланне
MATLAB для DSP. Застосування багатошвидкісних фільтрів у задачах вузькосмугової фільтрації
Багатошвидкісні фільтри та банки фільтрів (фільтри та банки фільтрів з багаточастотною дискретизацією) визначили самостійний напрямок у теорії та практиці цифрової обробки сигналів (ЦГЗ). У зв'язку з цим досить згадати квадратурно-дзеркальні фільтри — новий клас багатосмугових фільтрів, банки багатошвидкісних фільтрів для реалізації вейвлет-перетворень, поліфазні фільтри. Багатошвидкісна фільтрація знайшла широке практичне застосування у завданнях компресії мови, звуку та зображень, побудови ефективних систем фільтрації, очищення сигналу від перешкод, оптимізації обчислювальних ресурсів під час реалізації алгоритмів ЦОС.
Для вирішення завдань аналізу (моделювання) та синтезу (проектування) систем із багаточастотною дискретизацією пакет MATLAB надає широкі можливості. В рамках проекту "MATLAB для DSP" розглянемо приватне, але практично важливе завдання проектування та моделювання вузькосмугового фільтра. Для цих цілей, як це і обумовлено в проекті, використовуються два найбільш дружні для користувача інструменти MATLAB - Simulink і GUI. При цьому передбачається, що читач знайомий із предметною областю. У наступних розділах на прикладі вирішення конкретної задачі - синтез та моделювання вузькосмугового ФНЧ - викладаються правила та практичні рекомендації щодо використання інструментів MATLAB.
Основні співвідношення теорії багаточастотної дискретизації.
Основи теорії багатошвидкісної фільтрації (багатошвидкісної обробки сигналів) та опис додатків можна знайти в [1-4].
В якостібазових при багаточастотній дискретизації використовуються:
- M-кратний дециматор (рис. 1), де YD(n) = x(Mn);
- L-кратний інтерполятор (рис. 2), де
Таким чином, частота дискретизації f s сигналу x(n) пов'язана з частотою дискретизації f's співвідношенням Mf's = f's або співвідношенням періодів дискретизації T' = MT. Частоти дискретизації сигналів x(n) (частота f s ) та y l (n) (частота f s ) при інтерполяції пов'язані співвідношеннями f's = Lf s або T = LT'.
Приклад децимації сигналу x(n) при M = 2: відліки сигналу до (а) та після (б) децимації.
Приклад інтерполяції сигналу x(n) при L = 2: відліки сигналу до (а) та після (б) інтерполяції
При децимації та інтерполяції сигналу відбувається деформація спектрів. Для децимації
де нормована частота w = wT', а T' – період дискретизації після децимації. Він буде в M разів більшим, так що в шкалі ненормованих частот отримаємо
Враховуючи, що T' = MT можна записати остаточне співвідношення між спектрами децимованого сигналу Y D (e jw ) і вихідного X(e jw ):
Для інтерполяції Y I (z) = X (z L) або Y I (e jw) = X (e jwL) або Y I (e jw) = X (e jwL). Таким чином, спектр децимованого сигналу є зваженою сумою вихідного спектра X(e jw ) та її (M–1) зрушених по частоті копій (відбитків) з кроком 2Пw s /M. Спектр інтерполованого сигналу є спектром вихідного сигналу зі зміненим періодом за частотою. Період збільшується в L раз. Враховуючи згадані вище властивості спектру, необхідно перед децимацією ставити фільтр децимації, щоб виключити накладання спектру, а для інтерполяції - фільтр інтерполяції для усунення відбитків, тобто додаткових компонент спектру, які потрапили в робочу смугу [ 0, F N] за рахунокзбільшення періоду спектра. Чудові тотожності. При побудові систем з багаточастотною дискретизацією дуже корисні перетворення, зображені малюнки. Вони корисні в багатьох випадках при реалізації фільтрів, що буде показано в наступному розділі.
Поліфазне розбиття та поліфазні фільтри. Передатна функція нерекурсивного (КІХ - кінцевої імпульсної характеристики) фільтра
може бути представлено сумою
H(z) = h(0) + h(2)z-2+. + h (1) z-1 + h (3) z-3 +. або H(z) = h(0) + h(2)z-2 + . + z-1 (H (1) + h (3) z-2 + .) = E0 (z2) + E1 (z2).
Сенс многочленів E 0 і E 1 зрозумілий з контекста.Если E 0 і E 1 розглядати як передавальні функції КІХ-фільтрів, то неважко помітити, що базовим елементом затримки таких фільтрів є z -2, що забезпечує затримку на два такти. Отже, фільтри E 0 і E 1 можуть працювати на частоті дискретизації, удвічі меншою за вихідну. Якщо використовувати розкладання за ступенями z-3 або z-4 і так далі, то можна отримати блоки фільтрів, що працюють на більш низьких частотах дискретизації. Підсумкова схема фільтра, коли H(z) = E0(z2) + z-1E1(z2) показана на малюнку 7.
Розглянуте розбиття називається поліфазним, а схеми, що його реалізують, - поліфазними фільтрами. Як приклад, що ілюструє побудову поліфазного фільтра і використання чудових тотожностей, покажемо, як можна ефективно реалізувати КІХ-фільтр дециматор.
У схемі рис. 8 для обчислення кожного відліку необхідно виконати N+1 множень та N додавань. У той же час, за рахунок децимації (M = 2) половину результатів відрахунків ми відкидаємо і, отже, використовуємо ресурси обчислювача неефективно.чудовою тотожністю (рис. 5). Послідовність перетворень у своїй показано на рис. 9.
В результаті фільтри E 0 (z 2 ) і E 1 (z 2 ), що працюють на частоті f s замінюються фільтрами E 0 (z) і E 1 (z), що працюють на частоті f s /2. Таким чином ми отримаємо дворазову економію у швидкості обчислень. При коефіцієнті децимації M застосування поліфазних фільтрів дозволяє отримати економію M разів. Аналогічні результати виходять при поліфазному побудові фільтрів інтерполяторів [1,2]. Таким чином, поліфазні реалізації разом з раціональними перетвореннями схем фільтрів дозволяють будувати ефективні обчислювачі в завданнях ЦГЗ.
Розрахунок вузькосмугового низькочастотного фільтра.
У процесі цифрової обробки сигналів нерідко виникає завдання фільтрації сигналу дуже вузькому частотному діапазоні. Прикладом такої задачі може бути ситуація, коли сигнал містить кілька складових на близьких частотах, і потрібно виділити одну з них. Для цього необхідно спроектувати цифровий фільтр із дуже вузькими відносними смугами пропускання та переходу. Залежно від розташування складових сигналу, це може бути фільтр нижніх, верхніх частот, смуговий або фільтр, що загороджує. Пряме проектування таких пристроїв часто призводить до фільтрів дуже високих порядків, практична реалізація яких або нераціональна, або неможлива. Розглянемо наступний приклад. Нехай є дискретний сигнал
u(n) = sin(2 π f 1 nT) + sin(2 π f 2 nT) + e(n) ,
де f 1 = 90 Гц і f 2 = 105 Гц – частоти складових; T = 1/f s = 1/8000 с – період дискретизації; f s = 8000 Гц – частота дискретизації; e(n) - гаусівський шум з нульовим середнім та дисперсією σ 2 = 0,1. Потрібно виділити синусоїдальну складову з частотою f1. Дляцього необхідно побудувати фільтр нижніх частот з граничними значеннями частот смуги пропускання f pb і смуги затримування fsb, що задовольняють умову f 1 pb sb 2 . Оскільки частота Найквіста f N = f s /2 = 4000 Гц, для нормалізованих значень f pb = f pb / f N і f sb = f sb / f N ця умова буде виглядати так:
Отже, перехідна смуга f амплітудно-частотної характеристики фільтра (АЧХ) не повинна перевищувати величину f s1 = 2 · 90 = 180 Гц. Для виділення складової з частотою f s виконаємо фільтрацію сигналу, використовуючи два поліфазні фільтри-дециматори і два поліфазні фільтри-інтерполятори. Нехай перший фільтр має смугу пропускання
[0, 1/20-0,025] · F N = [0, 100] Гц
та смугу затримування
[1/20, 1] · F N = [200, 4000] Гц.
АЧХ у смузі пропускання цього фільтра має знаходитися в межах [0,995, 1,005], а в смузі затримування - не перевищувати величину 10-4. Виходячи з цих вимог, коефіцієнт децимації М покладемо рівним 20. При цьому частота дискретизації на виході першого фільтра дециматора дорівнюватиме 8000/20 = 400 Гц. Розрахунок коефіцієнтів першого фільтра можна виконати знову з використанням GUI sptool. Введені дані для розрахунку:
Результати розрахунку дають приблизну оцінку порядку фільтра n 1 = 270, однак отримані значення Actual Rp = 0,1436 і Actual Rs = 75,58 (вони відображаються у правій частині вікна Filter Designer) вказують на необхідність збільшити порядок фільтра та повторити обчислення. У нашому прикладі ми вибрали значення n 1 = 289, що відповідає Actual Rp = 0,08281 і Actual Rs = 80,16. Другий фільтр дециматор характеризується такими параметрами:
Результати обчислень дають попередню оцінку порядку фільтра n 2 = 270 та фактичні значення Rp та Rs для цьогопорядку, які знову ж таки не вкладаються у задані вимоги (пульсації неприпустимо великі). Вибір n 2 = = 275 забезпечує виконання вимог. Щоб коефіцієнти спроектованих фільтрів стали доступні для подальшого використання, їх треба експортувати в робочий простір MATLAB. Для цього необхідно відкрити меню File головної панелі sptool та вибрати розділ Export… Якщо експортовані структури мають імена filt1 та filt2, вибрати відповідні коефіцієнти можна за допомогою команд
Застосування двох фільтрів дециматорів, що реалізуються у вигляді поліфазних структур, вимагає виконання
290/2 MACs/с · 8000 відліків/с · 1/20 + 276/2 MACs/с · 400 відліків/с · = 85600 MACs/с
і зберігання 145 + 138 = 283 коефіцієнтів. Для відновлення вихідної частоти дискретизації 8000 Гц слід до виходу другого фільтра послідовно підключити два фільтри-інтерполятори, причому коефіцієнт інтерполяції першого фільтра дорівнює 2, а другого - 20. фільтри виконують стільки ж операцій на секунду, як і попередні фільтри- дециматори, загальна кількість операцій множення з накопиченням на секунду дорівнюватиме 85600·2 = 171200, а число коефіцієнтів фільтрів - 283.
Моделювання вузькосмугового низькочастотного фільтра
Описана процедура може бути легко реалізована у вигляді моделі Simulink, для чого нам знадобляться такі блоки:
Правила побудови моделей докладно викладені в довідковій системі MATLAB, з ними також можна ознайомитися, звернувшись до відповідних джерел [6-8]. Побудована модель, де блоки (F+D)1 і (F+D)2 - це фільтри -дециматори, а блоки (I+F)1 та (I+F)2 - фільтри-інтерполятори. Перед тим, як її запустити, необхідно встановити параметри блоків і режиму роботи моделі. У нашому прикладіпараметри блоків мають такі значення:
Параметри режиму роботи моделі встановлюються у вікні Simulation Parameters.
Розглянутий приклад наочно демонструє можливості багаточастотної дискретизації на вирішення важливих практичних завдань. Справді, щодо вузькосмугові фільтри є найважливішим елементом радіосистем і у сенсі мають велике самостійне значення. Крім цього, багаточастотна дискретизація дуже ефективна для побудови банків фільтрів та для вирішення завдань обробки сигналів за допомогою техніки вейвлет-розкладів.