Практичне застосування перетворення фур'є обробки сигналів. Цифрова обробка сигналів на ARM7-мікроконтролерах Перетворення фур'є avr

Теорема Фур'є говорить, що будь-який сигнал можна розкласти в ряд за ортонормованим набором періодичних функцій (наприклад, синусами і косинусами) з частотами кратними частоті періодичного сигналу. Таким чином, в основі спектрального аналізу сигналу лежить пошук вагових коефіцієнтів (загалом комплексних), модуль яких відповідає частки потужності коливань відповідної гармоніки, що вноситься в загальну суперпозицію всіх гармонік.

Швидке перетворення Фур'є

Швидке перетворення Фур'є - це алгоритм обчислення, який успішно використовує властивості періодичності тригонометричних функцій для того, щоб уникнути непотрібних обчислень у дискретному перетворенні Фур'є (ДПФ), за рахунок чого здійснюється пошук коефіцієнтів у розкладанні Фур'є. Основна відмінність від дискретного перетворення полягає лише у методі обчислення числових значень (алгоритмі), а чи не у обробці сигналу. І у разі БПФ, і у разі ДПФ результат обчислень один і той же. Єдиною вимогою алгоритму БПФ є розмір вибірки, кратний N = 2L, де L - будь-яке позитивне ціле число. Найбільш поширеними алгоритмами БПФ на основі 2 є: з проріджуванням за часом і з проріджуванням за частотою.

У цій роботі реалізований алгоритм БПФ на основі 2 з проріджуванням за часом (алгоритм Кулі - Тьюкі). Його легко одержати, досліджуючи деякі закономірності ДПФ. Введемо так званий поворотний коефіцієнт:

В цьому випадку в ДПФ коефіцієнти Фур'є для ряду значень сигналу (f0, f1, ..., fN-1) виражаються співвідношенням:

Розглянемо ряд сигналу із 4 значень: (f0, f1, f2, f3). Представимо перетворення Фур'є в матричному вигляді (нормувальний коефіцієнт 1/N внесений у вектор-стовпець Сij у правій частині виразу):

Розписавши поворотні коефіцієнти за формулою Ейлера і визначивши їх значення при k = 0, 1, 2,.. 9, можна побудувати діаграму (Рис. 2), з якої видно закономірність коефіцієнтів, що повторюються.

Малюнок 2. Ступіньний ряд w для N=4

Підставляючи чисельні значення (4) отримаємо:

Тобто значення w, починаючи з w4, дорівнюють відповідному значенню від w0 до w3. Далі перепишемо матричне рівняння (4) у нестандартному вигляді (подібні позначення введені для наочності подальших операцій):

Поміняємо місцями стовпці матриці, розділивши її на дві групи: за парними f0, f2 і непарними f1, f3 індексами:

Врахуємо, що wk+1 = wkw1, тоді вираз (6) перепишеться у вигляді:

Використовуючи співвідношення:

Отримуємо шукані коефіцієнти розкладання у вигляді вектора-стовпця зі значеннями осередків:

Графічне зображення алгоритму (Рис. 3) схоже на метелика з розкритими крилами, тому цей метод обчислення називають «метеликом».

Малюнок 3. Граф «Метелик» для ряду з 4 членів

Отже, першому кроці алгоритму йде розбиття по парним і непарним індексам членів низки значень сигналу. Потім виконується граф «метелик», він складається з двох ступенів, їх кількість дорівнює ступеню двійки обсягу вибірки (N = 4 = 22). На кожному щаблі виконується по дві «метелики» та їх загальна кількість незмінна. Кожній операції "метелик" відповідає одна операція множення. Для порівняння: у ДПФ з вибіркою (f0, f1, f2, f3) операцію множення необхідно було б виконати 4Ч4 = 16 разів, а у разі БПФ лише 4 рази.

Цей ряд може бути записаний у вигляді:

(2),
де, k-я комплексна амплітуда.

Зв'язок між коефіцієнтами (1) та (3) виражається такими формулами:

Зазначимо, всі ці три уявлення низки Фур'є цілком рівнозначні. Іноді під час роботи з рядами Фур'є буває зручніше використовувати замість синусів і косінусів експоненти уявного аргументу, тобто використовувати перетворення Фур'є у комплексній формі. Але нам зручно використовувати формулу (1), де ряд Фур'є представлений у вигляді суми косінусоїд з відповідними амплітудами та фазами. У будь-якому разі неправильно говорити, що результатом перетворення Фур'є дійсного сигналу будуть комплексні амплітуди гармонік. Як правильно говориться в Вікі «Перетворення Фур'є (?) - операція, яка зіставляє однієї функції речової змінної іншу функцію, також речової змінної.»

Разом:
Математичною основою спектрального аналізу сигналів є перетворення Фур'є.

Перетворення Фур'є дозволяє представити безперервну функцію f(x) (сигнал), визначену на відрізку (0, T) у вигляді суми нескінченного числа (нескінченного ряду) тригонометричних функцій (синусоїд та косинусоід) ​​з певними амплітудами і фазами, також розглядаються на відрізку (0, T). Такий ряд називається поряд Фур'є.

Зазначимо ще деякі моменти, розуміння яких потрібно правильного застосування перетворення Фур'є до аналізу сигналів. Якщо розглянути ряд Фур'є (суму синусоїд) на всій осі Х, то можна побачити, що поза відрізком (0, T) функція представлена ​​поруч Фур'є буде періодично повторювати нашу функцію.

Наприклад, на графіці рис.7 вихідна функція визначена на відрізку (-T2, + T2), а ряд Фур'є представляє періодичну функцію, визначену на всій осі х.

Це тому, що синусоїди самі є періодичними функціями, відповідно і їх сума буде періодичною функцією.


рис.7 Подання неперіодичної вихідної функції поруч Фур'є

Таким чином:

Наша вихідна функція – безперервна, неперіодична, визначена на деякому відрізку довжиною T.
Спектр цієї функції – дискретний, тобто представлений у вигляді нескінченного ряду гармонійних складових – низки Фур'є.
За фактом, поряд Фур'є визначається деяка періодична функція, що збігається з нашою на відрізку (0, T), але для нас ця періодичність не суттєва.

Періоди гармонійних складових кратні величині відрізка (0, T), на якому визначено вихідну функцію f(x). Інакше кажучи, періоди гармонік кратні тривалості вимірювання сигналу. Наприклад, період першої гармоніки низки Фур'є дорівнює інтервалу Т, у якому визначено функцію f(x). Період другої гармоніки ряду Фур'є дорівнює інтервалу Т/2. І так далі (див. мал. 8).


рис.8 Періоди (частоти) гармонійних складових ряду Фур'є (тут Т=2?)

Відповідно, частоти гармонійних складових кратні величині 1/Т. Тобто частоти гармонійних складових Fk дорівнюють Fk= к\Т, де до пробігає значення від 0 до?, наприклад до=0 F0=0; к=1 F1=1\T; к=2 F2=2\T; к = 3 F3 = 3 \ T; ... Fk = к \ Т (при нульовій частоті - постійна складова).

Нехай наша вихідна функція є сигнал, записаний протягом Т=1 сек. Тоді період першої гармоніки дорівнюватиме тривалості нашого сигналу Т1=Т=1 сек і частота гармоніки дорівнює 1 Гц. Період другої гармоніки дорівнюватиме тривалості сигналу, поділеної на 2 (Т2=Т/2=0,5 сек) і частота дорівнює 2 Гц. Для третьої гармоніки Т3=Т/3 с і частота дорівнює 3 Гц. І так далі.

Крок між гармоніками у разі дорівнює 1 Гц.

Таким чином, сигнал тривалістю 1 сек можна розкласти на гармонійні складові (отримати спектр) з роздільною здатністю по частоті 1 Гц.
Щоб збільшити розподільну здатність в 2 рази до 0,5 Гц - треба збільшити тривалість вимірювання в 2 рази - до 2 сек. Сигнал тривалістю 10 с можна розкласти на гармонійні складові (отримати спектр) з роздільною здатністю по частоті 0,1 Гц. Інших способів збільшити роздільну здатність за частотою немає.

Існує спосіб штучного збільшення тривалості сигналу шляхом додавання нулів до масиву відліків. Але реальну роздільну здатність за частотою він не збільшує.

3. Дискретні сигнали та дискретне перетворення Фур'є

З розвитком цифрової техніки змінилися способи зберігання даних вимірювань (сигналів). Якщо раніше сигнал міг записуватися на магнітофон і зберігатися на стрічці у аналоговому вигляді, нині сигнали оцифровуються і зберігаються у файлах у пам'яті комп'ютера як набору чисел (отсчетов).

Звичайна схема вимірювання та оцифрування сигналу виглядає наступним чином.


рис.9 Схема вимірювального каналу

Сигнал з вимірювального перетворювачанадходить на АЦП протягом періоду часу Т. Отримані за час Т відліки сигналу (вибірка) передаються на комп'ютер і зберігаються у пам'яті.


рис.10 Оцифрований сигнал - N відліків отриманих за час Т

Які вимоги висуваються до параметрів оцифрування? Пристрій, що перетворює вхідний аналоговий сигнал на дискретний код ( цифровий сигнал) називається аналого-цифровий перетворювач (АЦП, англ. Analog-to-digital converter, ADC) (Wiki).

Однією з основних параметрів АЦП є максимальна частота дискретизації (чи частота семплювання, англ. sample rate) - частота взяття відліків безперервного у часі сигналу за його дискретизації. Вимірюється у герцах. ((Wiki))

Згідно з теоремою Котельникова, якщо безперервний сигнал має спектр, обмежений частотою Fмакс, то він може бути повністю і однозначно відновлений за його дискретними відліками, взятими через інтервали часу , тобто. з частотою Fd? 2*Fмакс, де Fd – частота дискретизації; Fмакс – максимальна частота спектра сигналу. Іншими словами, частота оцифрування сигналу (частота дискретизації АЦП) повинна як мінімум у 2 рази перевищувати максимальну частоту сигналу, який ми хочемо виміряти.

А що буде, якщо ми братимемо відліки з меншою частотою, ніж потрібно за теоремою Котельникова?

У цьому випадку виникає ефект "аліасингу" (він же стробоскопічний ефект, муаровий ефект), при якому сигнал високої частоти після оцифрування перетворюється на сигнал низької частоти, якого насправді не існує. На рис. 5 червона синусоїда високої частоти – це реальний сигнал. Синя синусоїда нижчої частоти - фіктивний сигнал, що виникає внаслідок того, за час взяття відліку встигає пройти більше, ніж півперіоду високочастотного сигналу.


Мал. 11. Поява помилкового сигналу низької частоти за недостатньо високої частоти дискретизації

Щоб уникнути ефекту аліасингу перед АЦП ставлять спеціальний антиаліасинговий фільтр - ФНЧ (фільтр нижніх частот), який пропускає частоти нижче половини частоти дискретизації АЦП, а вищі частоти зарізає.

Для того щоб обчислити спектр сигналу за його дискретними відліками використовується дискретне перетворення Фур'є (ДПФ). Зазначимо ще раз, що спектр дискретного сигналу «за визначенням» обмежений частотою Fмакс меншою половиною частоти дискретизації Fd. Тому спектр дискретного сигналу може бути представлений сумою кінцевого числа гармонік, на відміну від нескінченної суми для ряду безперервного Фур'є сигналу, спектр якого може бути необмежений. Відповідно до теореми Котельникова максимальна частота гармоніки повинна бути такою, щоб на неї припадало щонайменше два відліки, тому число гармонік дорівнює половині числа відліків дискретного сигналу. Тобто якщо у вибірці є N відліків, то число гармонік у діапазоні дорівнюватиме N/2.

Розглянемо тепер дискретне перетворення Фур'є (ДПФ).

Порівнюючи з рядом Фур'є

Бачимо, що вони збігаються, за винятком того, що час у ДПФ має дискретний характер і кількість гармонік обмежена величиною N/2 – половиною числа відліків.

Формули ДПФ записуються у безрозмірних цілих змінних k, s, де k – номери відліків сигналу, s – номери спектральних складових.
Величина s показує кількість повних коливань гармоніки на періоді Т (тривалості вимірювання сигналу). Дискретне перетворення Фур'є використовується знаходження амплітуд і фаз гармонік чисельним методом, тобто. "на комп'ютері"

Повертаючись до результатів, отриманих на початку. Як було зазначено вище, при розкладанні у ряд Фур'є неперіодичної функції (нашого сигналу), отриманий ряд Фур'є фактично відповідає періодичної функції з періодом Т. (рис.12).


рис.12 Періодична функція f(x) з періодом Т0, з періодом вимірювання Т>T0

Як видно на рис.12 функція f(x) періодична з періодом Т0. Однак через те, що тривалість вимірювальної вибірки Т не збігається з періодом функції Т0, функція, що отримується як ряд Фур'є, має розрив у точці Т. В результаті спектр цієї функції міститиме велику кількість високочастотних гармонік. Якби тривалість вимірювальної вибірки Т збігалася з періодом функції Т0, то в отриманому після перетворення Фур'є спектрі була б присутня тільки перша гармоніка (синусоїда з періодом рівним тривалості вибірки), оскільки функція f(x) являє собою синусоїду.

Іншими словами, програма ДПФ «не знає», що наш сигнал є «шматком синусоїди», а намагається представити у вигляді ряду періодичну функцію, яка має розрив через нестикування окремих шматків синусоїди.

У результаті спектрі з'являються гармоніки, які мають у сумі зобразити форму функції, включаючи цей розрив.

Таким чином, щоб отримати «правильний» спектр сигналу, що є сумою кількох синусоїд з різними періодами, необхідно щоб на період вимірювання сигналу вкладалося ціле число періодів кожної синусоїди. На практиці цю умову можна виконати за досить великої тривалості вимірювання сигналу.


Рис.13 Приклад функції та спектра сигналу кінематичної похибки редуктора

За меншої тривалості картина виглядатиме «гірше»:


Рис.14 Приклад функції та спектру сигналу вібрації ротора

Насправді буває складно зрозуміти, де «реальні складові», а де «артефакти», викликані некратностью періодів складових і тривалості вибірки сигналу чи «стрибками і розривами» форми сигналу. Звичайно слова «реальні складові» та «артефакти» не дарма взяті у лапки. Наявність на графіку спектра безлічі гармонік не означає, що наш сигнал насправді з них «складається». Це все одно що вважати, що число 7 «складається» з чисел 3 і 4. Число 7 можна представити у вигляді суми чисел 3 і 4 - це правильно.

Так і наш сигнал… а точніше навіть не «наш сигнал», а періодичну функцію, складену шляхом повторення нашого сигналу (вибірки) можна у вигляді суми гармонік (синусоїд) з певними амплітудами та фазами. Але в багатьох важливих для практики випадках (див. малюнки вище) дійсно можна пов'язати отримані в спектрі гармоніки і з реальними процесами, що мають циклічний характер і значний внесок у форму сигналу.

Деякі підсумки

1. Реальний виміряний сигнал, тривалістю T сек, оцифрований АЦП, тобто представлений набором дискретних відліків (N штук), має неперіодичний дискретний спектр, представлений набором гармонік (N/2 штук).

2. Сигнал представлений набором дійсних значень та його спектр представлений набором дійсних значень. Частоти гармонік позитивні. Те, що математикам буває зручніше уявити спектр у комплексній формі з використанням негативних частот не означає, що так правильно і так завжди треба робити.

3. Сигнал, виміряний на відрізку часу Т визначено тільки на відрізку часу Т. Що було до того, як ми почали вимірювати сигнал, і що буде після того – науці це невідомо. І в нашому випадку – нецікаво. ДПФ обмеженого у часі сигналу дає його «справжній» спектр, у тому сенсі, що за певних умов дозволяє обчислити амплітуду та частоту його складових.

Використані матеріали та інші корисні матеріали.

Як пристрій відображення використовується дворядковий символьний РК індикатор. p align="justify"> Основним моментом при реалізації даного проекту є не апаратна частина, а програмна, точніше реалізація дискретного перетворення Фур'є (ДПФ) на 8-розрядному мікроконтролері. Відразу слід зазначити, що автор не є експертом у цій галузі і тому почав з основ – з простого дискретного перетворення Фур'є. Алгоритм швидкого перетворення Фур'є не лише швидким, а й досить складним.

Дискретне перетворення Фур'є (в англомовній літературі DFT, Discrete Fourier Transform) - це одне з перетворень Фур'є, які широко застосовуються в алгоритмах цифрової обробки сигналів (його модифікації застосовуються в стисканні звуку в MP3, стисненні зображень в JPEG та ін.), а також в інших областях, пов'язаних з аналізом частот у дискретному (наприклад, оцифрованому аналоговому) сигналі. Дискретне перетворення Фур'є вимагає як вход дискретну функцію. Такі функції часто створюються шляхом дискретизації (вибірки значень безперервних функцій).

Принципова схема аналізатора спектра звукового сигналудуже проста та умовно її можна розділити на цифрову частину та аналогову.

Цифрова частина утворена мікроконтролером та підключеним до нього РК індикатором. Мікроконтролер тактується від кварцового резонатора 16 МГц, як опорну напругу для АЦП мікроконтролера використовується напруга живлення +5 В.
Шина даних РК індикатора підключена до порту C мікроконтролера (лінії введення/виведення PC0-PC3), шина управління підключена до порту D(PD5, PD6) мікроконтролера. Індикатор працює у 4-бітному режимі. Змінний резистор номіналом 4.7 ком використовується для регулювання контрастності. Для роботи з індикатором були створені символи для відображення 8 горизонтальних стовпчиків аналізатора, ці символи займають всі 64 Байта ОЗУ ЖК індикатора.

Мікроконтролер працює від зовнішнього кварцового резонатора 16 МГц.

Аналогова частина пристрою - це найважливіша частина і є попереднім підсилювачем сигналу електретного мікрофона, вихід якого підключається до каналу ADC0 вбудованого в мікроконтролер АЦП. Рівень нуля на вході АЦП необхідно встановити рівним точно половині опорного напруги, тобто. 2.5 У. ​​У разі ми зможемо використовувати позитивну і негативну напівхвилю сигналу, та її амплітуда має перевищувати встановлений межа, тобто. коефіцієнт посилення має бути точно налаштований для запобігання навантаженню. Всім вищезгаданим умовам відповідає поширена мікросхема низькоспоживаючого операційного підсилювача.

Алгоритм ДПФ дещо повільніший у порівнянні зі швидким перетворенням Фур'є. Але наш аналізатор спектру не вимагає високої швидкості, і якщо він здатний забезпечити швидкість оновлення близько 30 кадрів на секунду, цього буде більш ніж достатньо візуалізації спектру звукового сигналу. У будь-якому випадку, у нашому варіанті можна досягти швидкості 100 кадрів в секунду, але це вже занадто високе значення параметра для дворядкового символьного РК індикатора і воно не рекомендується. Частота дискретизації дорівнює 20 кГц для 32 точкового дискретного перетворення Фур'є і оскільки результат перетворення симетричний, потрібно використовувати лише першу половину, тобто. перші 16 результатів. Отже, ми можемо відображати частотний спектр у діапазоні до 10 кГц і роздільна здатність аналізатора становить 10 кГц/16 = 625 Гц.

Автором конструкції було здійснено спроби збільшення швидкості обчислення ДПФ. Якщо це перетворення має N точок, ми повинні знайти N2/2 значень синуса і косинуса. Для нашого 32 точкового перетворення необхідно знайти 512 значень синуса та косинуса. Але, перш ніж знайти їх нам необхідно обчислити кут (градуси), що займе деякий процесорний час, тому було вирішено використовувати ці обчислення таблиці значень. При розрахунках у програмі мікроконтролера не використовуються обчислення з плаваючою точкою та числа подвійної точності (double), оскільки це займе більше часу на обробку на 8-розрядному мікроконтролері. Замість цього значення в таблицях пошуку використовуються 16-розрядні дані цілого чисельного типу (integer), помножені на 10000. Потім після виконання перетворення результати діляться на 10000. При такому підході є можливість виконувати 120 32-точкових перетворень за секунду, що більш ніж достатньо для нашого пристрої.

Демонстрація роботи аналізатора спектру на мікроконтролері ATmega32

Завантаження

Вихідний код (програма мікроконтролера, таблиці даних синуса, косинуса та кута) -

  • Зрозуміло, що на АВР-ці далі світломузики складно виїхати, не ті параметри. Але 120 32-точкових перетворень на секунду для більшості завдань може бути достатньо. А вибірку 625Гц, можна й іншу взяти, по точніше втративши частоту оновлення. Варто зазначити, що МК погано почуватиметься, у плані продуктивності мало що на нього ще навішаєш. Але тут можна організувати видачу результату з апаратних методів передачі. Тоді це буде допоміжний МК, а основний тільки прийматиме з нього даний і оброблятиме сумісно з іншими процесами. За великим рахунком все ж таки в частоту праці впирається. Колись виходило розганяти мегу вище 20 МГц, але для цих завдань напевно отримаємо глюки на високих частотах. Ідея хороша, аби більше мат частини розписано було б... саме її реалізація на МК
  • я і цікавіше аналізатори робив: You Tube або варіант на кольоровому РКІ: You Tube в основі знаменита бібліотека Чена:)
  • "Нам необхідно обчислити кут (градуси)" А може хтось докладніше розповісти як розраховуються значення для цих таблиць?
  • З таблицею синусів та косинусів все зрозуміло. Не зрозуміло, як розраховуються значення в таблиці degree_lookup?

Всі сигнали, незалежно від того, ви їх придумали або спостерігали у Всесвіті, насправді просто сума простих синусоїд різних частот.

Я зробив невеликий аудіо аналізатор спектру (0 - 10 кГц) з РК-дисплея 16x2 та мікроконтролера ATmega32. Я почав із простих ДПФ (Дискретне Перетворення Фур'є). БПФ (Швидке Перетворення Фур'є) відрізняється від ДПФ тільки більшою швидкістю і трохи складнішим алгоритмом, я не став його використовувати, можливо, я додам його пізніше.

ДПФ повільний проти БПФ. Мій РК аналізатор спектру не вимагає великої швидкості, яку може забезпечити БПФ, і якщо зображення на екрані змінюватиметься з частотою близько 30 кадрів/сек, то це більш ніж достатньо для візуалізації звукового спектру. Але я можу досягти частоти близько 100 кадрів/сек, проте для РК-дисплея не рекомендується занадто висока частота оновлення. Звук із частотою дискретизації 20 кГц дає 32 точки ДПФ. Оскільки результат перетворення є симетричним, мені потрібно використовувати лише перші 16 результатів. Відповідно, максимальна частота 10 кГц. Отже, 10кГц/16 = 625Гц.

Я намагався збільшити швидкість обчислення ДПФ. Якщо є точка N ДПФ, то необхідно знайти синус та косинус (N^2)/2. Для 32-точкового ДПФ, необхідно знайти синус та косинус 512. Перш ніж шукати синус та косинус, нам потрібно знайти кут (градуси), який займає деякий час процесора. Для цього я зробив таблиці для синуса та косинуса. Я зробив синус і косинус 16-бітними змінними, помноживши значення синуса і косинуса на 10000. Після перетворення я повинен розділити кожен результат на 10000. Тепер я можу розрахувати 120 32-точкових ДПФ на секунду, що більш ніж достатньо для спектра аналізатора.

Дисплей

Я використовував символи для РК-дисплея завантажені в 64 Байт вбудованої пам'яті РК-дисплея. В інтернеті я побачив відео, де РК-дисплей 16х2 використовується як дисплей аналізатора спектру і використовував цю ідею.

Аудіо вхід

Однією з найважливіших частин аналізатора спектра є отримання сигналу з електретного мікрофона. Особлива увага має бути приділена розробці попереднього підсилювачадля мікрофону. Нам необхідно встановити нульовий рівень на вході АЦП і максимальний рівень дорівнює половині напруги живлення, тобто. 2,5В. На нього може подаватися напруга від -2,5 до +2,5В. Підсилювач має бути налаштований так, щоб не перевищувати ці межі. Я використовував операційний підсилювач LM324 як попередній підсилювач для мікрофона.

Список радіоелементів

Позначення Тип Номінал Кількість ПриміткаМагазинМій блокнот
Дисплей
МК AVR 8-біт

ATmega32

1 До блокноту
Конденсатор22 пФ2 До блокноту
Конденсатор0.1 мкФ1 До блокноту
Електролітичний конденсатор100 мкФ1 До блокноту
Резистор

100 Ом

1 До блокноту
Підстроювальний резистор4.7 ком1 До блокноту
Кварцовий резонатор16 МГц1 До блокноту
LCD-дисплей16х21 До блокноту
Блок живлення5 В1 До блокноту
Аудіо вхід
U1 Операційний посилювач

LM324

1 До блокноту
З 1 Конденсатор1 мкФ1 До блокноту
С8 Конденсатор0.01 мкф1 До блокноту
R1 Резистор

220 ком

1 До блокноту
R2, R3 Резистор

10 ком

2 До блокноту
R4, R9 Резистор

1 ком

2 До блокноту
R5 Резистор

Для цифрової обробки сигналів (DSP) існує безліч спеціалізованих процесорів, це DSP від ​​Texas Instruments серії TMS320, яка включає в себе і прості цілочисленні ядра, і такі монстри як сімейство сімейства C6000, що обробляють дані з плаваючою точкою. Є ціла серія ADSP від ​​Analog Devices (до якої входить більш-менш універсальний BlackFin), є і більше прості рішеннявід MicroChip – dsPIC.

Однак, спеціалізований DSP - це добре, але чи завжди він так необхідний? Так, при величезному потоці інформації він просто незамінний, проте є і простіші завдання обробки. Конкретно мене цікавило завдання подвійного перетворення - робиться згортка аудіосигналу, тим самим виходить спектр, потім над спектром можна провести будь-які операції та виконати зворотне перетворення, отримавши тим самим оброблений сигнал. Все це потрібно провести в реальному часі і отримати якість не нижче за телефонний.

Зараз не 2000 рік, існують і суттєво впали в ціні однокристальні рішення на продуктивних ядрах ARM7/Cortex-M3, вони є 32 бітними, мають апаратну реалізацію операції 32-бітного множення (мало того, майже DSP-операцію множення з накопиченням і 64 бітним результатом), а Cortex-M3 включає ще й апаратний поділ.

Відразу хочу попередити, що обробка сигналів не є моїм профілем, майже всі знання (а точніше розуміння принципів) збереглися ще з інституту, зараз же захотілося просто перевірити та реалізувати їх. Це я до того, що можуть зустрічатися неточності в описі, підміни понять та інше. Власне, академічна точність мене не дуже хвилювала.

Майже для будь-якого DSP, озвучене вище завдання є простим і підйомним. Але як поведеться на ній RISC-ядро загального застосування? Якщо розглядати AVR чи PIC, їх навряд чи вистачить. Дається взнаки 8-бітність і низька тактова частота. Хоча, у Elm-Chan є конструкції, в яких він на AVR проводить БПФ і малює спектр сигналу. Однак, у цьому випадку в реальному часі робиться просто візуалізація (з мінімальною точністю обробки), а не повна обробка сигналу з прийнятною аудіо-якістю.

Як піддослідний чіп був обраний LPC2146, заснований на ядрі ARM7TDMI-S і має максимальну тактову частоту 60МГц (на практиці, працює на 72 і навіть на 84МГц). Чому? По-перше, я маю налагоджувальну плату йому, по-друге, на борту є АЦП і ЦАП, тобто. потрібне мінімальне зовнішнє обважування.

Трохи теорії

Насамперед, цікаво було оцінити продуктивність БПФ (Швидкого Перетворення Фур'є) на ARM-мікроконтролерах. За цією оцінкою можна зробити висновок: чи вистачить у нього швидкості обробки потоку аудіоданих, і сигнал з якою частотою дискретизації і на скільки каналів можна буде обробити на такому мікроконтролері.

На основі перетворення Фур'є можна споруджувати хитрі фільтри (при цьому з дуже привабливими характеристиками). Мене ж насамперед цікавили завдання зміни тональності сигналу (підняття та опускання спектру) та "відображення" спектру. Останнє потрібно в SDR-радіоприймачі для прослуховування радіопередачі в нижній бічній смузі LSB.

Не буду завантажувати теорією і пояснювати, що таке Перетворення Фур'є, матеріалів на цю тему досить багато, з того, що використав я: вікі та глава з дуже гарної та інформативної книги.

Програмна реалізація

Програмних реалізацій БПФ існує дуже багато, проте я писав власну. Основна мета, яку я мала: оптимізація коду на конкретну архітектуру. По-перше, я відразу орієнтувався на 32-бітність, по-друге були потрібні тільки цілочисленні обчислення і бажано було уникнути операції поділу. Під ці вимоги підібрати щось готове проблемно.

Усі константи, які можна було розрахувати заздалегідь, було розраховано і в таблиці (переважно це значення тригонометричних функцій). Це основна оптимізація алгоритму, в іншому він майже повністю повторює описаний алгоритм.

Найважливішим же є вимога до цілочисленних обчислень. У процесі реалізації була навіть допущена помилка, через яку відбувалося переповнення в одній з 32 бітних змінних циклу. Мало того, вона виникала не на всіх тестових даних, тому завдавала пристойного головного болю, поки не була знайдена.

Усі вихідні тексти зібрав у одному архіві . Реалізація не остаточна, є дублюючі обчислення (не враховується симетрія спектру та фаз), і потрібна оптимізація використання буферів, тому що зараз для розрахунків використовується занадто багато оперативної пам'яті(Майже 12к для перетворення з 1024 пікселів).

Перші тести

Барабанний дріб: тестую швидкість роботи перетворення для вибірки з 1024 пікселів, частота процесорного ядра 60МГц. Тестування проводилося в емуляторі, так що не претендує на 100% точність, проте цей результат можна використовувати як показник (за моїм попереднім досвідом емулятор хоч і брехав, але не сильно). Тест першої версії коду, компілятор RealView MDK, опція оптимізації O3, ARM-режим генерації коду.

Отже, що ми бачимо:

По 6мс на кожне перетворення, всього трохи більше 12мс на перетворення "туди і назад". Виходить, що при частоті дискретизації 44100Гц (стандартна для аудіо) та вибірках роздільною здатністю до 16 біт, чисті обчислення займатимуть ~44*12мс = 528мс. І це на проміжній версії мікропрограми, коли ще не закінчено деякі оптимізації коду (за прикидками алгоритм може бути прискорений ще чи не в 2 рази)! По-моєму, просто відмінний показник.

Загалом, завантаження ядра передбачається в районі 50%, ще 50% залишається на перетворення над спектром і накладні витрати при роботі з АЦП, ЦАП та інші передачі даних. Якщо ж знизити частоту вибірок до "телефонного" рівня (близько 4800-9600Гц), то завантаження ядра буде ще нижче (близько 15-30%).

Отже, з математичною частиною більш-менш зрозуміло. Можна розпочати конкретну реалізацію.

Залізо

Для тестової платформи взято налагоджувальну плату Keil MCB2140 з динаміком. Напаяний шнурок Mini-Jack для підключення до лінійного виходу пристрою і зібраний найпростіший вхідний ланцюжок. Як було зазначено, на платі вже встановлено динамік, підключений до аналогового виходу мікроконтролера і є елементи управління (кнопка і потенціометр).

Ось малюнок вхідного ланцюга:


Налагодження софту відбувалося поетапно:

  1. Налагодження всієї необхідної периферії: АЦП, ЦАП, таймери, світлодіодна індикація.
  2. Тест з оцифруванням сигналу: оцифровую дані з необхідною швидкістю і кладу в буфер, потім витягую дані з буфера і відтворюю сигнал. Тобто. простий зсув сигналу у часі, без будь-яких перетворень. На цьому етапі тестується механізм роботи з двома буферами, необхідний для подальшої роботи.
  3. До попередньої версіїдодається пряме та зворотне перетворення Фур'є. Цей тест остаточно перевіряє коректність роботи коду БПФ, а також перевірка, що код вкладається в доступну продуктивність.
  4. Після цього основний кістяк програми зроблено, можна перейти до практичних застосувань.

Проблема виникла після додавання до коду БПФ: у сигналі з'явилися сторонні шуми та свисти. Взагалі, ця поведінка видалася мені досить дивною, т.к. без перетворення сигнал проходячи цифровий тракт залишався досить чистим. Перша причина цього: після аналогового ланцюга, амплітуда сигналу на АЦП була не повною (0-3.3В), а тільки в межах 0.5-2В на максимальній гучності плеєра, друга: досить сильний шум через цілих обчислень (+-1 одиниця, що виявилося достатнім для появи чутних перешкод).

Для боротьби з першою проблемою було вирішено зайнятися підстроюванням аналогової частини. А для вирішення питання з шумами – спробувати використати lowpass-фільтр.

Прикладне застосування 1: змінюємо тональність сигналу

На платі передбачений потенціометр ( змінні резистор), його можна використовувати для керування. В даному випадку він задає зсув спектру сигналу вгору-вниз, цілком достатньо щоб "перетворити" улюблені композиції.

Ось що відбувається у частотній області:


При цьому результат перетворення міститься у 2 буферах. Один буфер - справжня частина, а інший - уявна. Фізичний зміст отриманих у них чисел: у дійсній частині містяться значення гармонік, уявної - зрушення фаз даних гармонік. Мало того, очевидно, початковий сигнал описується N-значеннями, а після перетворення виходить 2N-значень. Кількість інформації не змінюється, а збільшення у 2 рази кількості інформації відбувається через те, що дані буфера мають надмірність у вигляді дублювання значень.




Top