На сайте http://www.zetms.ru/catalog/analyzers/a19.php есть описание анализатора спектра.

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

1) Частота входного сигнала до 100кГц.
2) Частота дискретизации до 250кГц.
3) Количество разрядов АЦП 16.
4) Для отображения сигнала в частотную область использую БПФ (на ADSP 21375).
5) Динамический диапазон то же - 85дБ.

Возникли вопросы по алгоритму:

1) Если сигнал случайный (музыка, речь и.т.д), то как определить частоту гармоники с наибольшей амплитудой в спектре?
2) Как наиболее достоверно (погрешность от 2 до 5%) отобразить спектр?
3) Какую фильтрацию использовать?
4) Можно ли обойтись без октавных фильтров, а использовать скажем оконную фильтрацию (при условии, что накапливаю несколько периодов сигнала основной частоты)?

Сообщение отредактировал Nukers - Mar 7 2010, 11:16

...мало того, можно его отобразить с достоверностью 0.01%
... меньше чем разрешение в пикселах на дисплее будет... в чем проблема то?


...ну... можно сказать и так... любое отображение одного множества на другое - своего рода фильтрация...



А откуда октавные фильтры? И опять эта "основная частота"...


P.S.: а чем Вы раньше занимались, и с чего вдруг интерес к такого рода заданию?

P.P.S.: Я думаю будет нелишне прочесть, что такое спектр вообще, в чем суть спектра мощностей и об алгоритмах DFT и FFT(к Вашему процу есть готовые примеры) ...

P.P.P.S.: что-бы не выглядеть только критикующим - Вам так или иначе нужно делать FFT(алгоритм заведомо блочный), потом выделять магнитуду как sqrt(I*I+Q*Q). Для каждого частотного "бина"... Может оказаться(в зависимости от того сколько точек взято для FFT), что спектр будет меняться быстрее чем успевает это воспринять оператор - тогда для каждого "бина" по оси времени можно применить скользящее среднее или блочное среднее. Если разрешение выходного БПФ больше чем отведенный массив пикселей на экранчиГе - то блочно масштабировать разрешение. Децимировать ни в коем случае нельзя - это безвозвратные потери полосы по Котельникову-Найквисту - возможны различные значительные артефакты.

Знаю, что можно. Что на это влияет?

Интересуют оконные КИХ фильтры. Точнее правильное расположение окна - чтобы снизить погрешность.



Вот код:


Думаю, понятно что он делает.
Так вот там есть так называемая "основная частота" - т.е. та частота, что имеет наибольшую амплитуду (возможно я неправильно подобрал определение просто).

Вопрос не стоит в отображении как таковом. Устройство лишь должно произвести БПФ, получить спектр, затем отправить пакетик через USB или Ethernet в компьютер. Далее работа программы верхнего уровня, которую нет смысла затрагивать. На данном этапе важно лишь получить вектор, содержащий информацию об амплитуде (фазу пока тоже не затрагиваем) (см. вышеописанный код).

Сэмпл, к вышеприведенному коду:  speech_dft.wav ( 214.95 килобайт ) Кол-во скачиваний: 118


Есть под ADSP 21375 EZ-KIT Lite примеры выполнения FFT Radix2 и Radix4 ? Буду благодарен.

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

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

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


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


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


Нет... ее там нет и быть не может... Есть понятие "основной тон речи"... Но это вообще не из этой степи, и Вам нужно на время забыть это определение...





...вряд-ли он чем-то мне поможет.


...они на сайте AD лежат...


Фирмы в открытую не афишируют свои блок-схемы... Читайте и думайте значит еще - все это выплывает из природы вещей которые мы обсуждаем... То, что студент - похвально! Я в Ваши годы вообще ничего из этой лабуды почти не знал ...

Фильтра там нет пока. Но будет. Выбираю.

"Фильтр" или "оконная функция"? Вот Вы приходите в магазин, а на полке книга - "нейронные сети" - не удивляет что она называется не "фильтрация"... Хотя где-то в глубине души - это было-бы верно...

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

Возьмите кольцевой буфер на 4096 точек. Разбейте его на 4 равных интервала - 1024 точек. Примените трапецевидное окно. 0...1024 - прямая проходящая через начало координат[0,0] и точку [1024,1] у Вас же float, так? Далее два интервала по 1024 точек просто прямая параллельная оси Ox - т.е. с этими данными ничего не делаем. Последний интервал опять прямая - [3072,1]..[4096,0].
Делаем так - с приходом каждой новой порции 1024-х отсчетов производим наложение окна и делаем FFT всего буфера на 4096 бинов. Т.е. за обход всего буфера у нас получается 4-е FFT на 4096 точек - при этом у окон по временной линии частичное перекрытие - на четверть длинны интервала рассмотрения. Далее эти 4-е результата FFT по-точечно суммируем и получаем одно FFT длинной 4096 бинов. Что логично ... Далее - все как-обычно... Для 4096 точек увеличение ДД по частоте будет существенное - я немного отошел от этих дел и сейчас не скажу - но для float ДД ожидается в районе 90-100 dB.
Если бы это была AVR-ка и расход памяти был бы весом - как функцию окна можно было бы применить трапецию в четвертьсинусоидальными скатами (гармоническая функция по-любому будет храниться в памяти для вычисления FFT).

P.S.: кстати - Вам нужно применять RealFFT. Вход то у Вас, по-идее, один...

P.P.S.: для диплома(если речь о нем) это все очень знаетели жирно... Напомните руководителю о истории возникновения лицензии BSD, о том что в развитых странах для того что-бы полюбоваться студенческими работами не нужна санкция прокурора, и о том что порядочные люди зарабатывают работая, а не эксплуатируя...

Тут я постарался графически отобразить мое понятие временного наложения кольцевого буфера (Оно конечно абстрактно)

Себе я выстроил такой алгоритм (поправьте, если заблуждаюсь):
1) Принимаем порцию из 1024 отсчетов.
2) Накладываем окно (т.е. копируем эту порцию 4 раза).
3) Делаем 4096 точечное БПФ.
4) Удаляем порцию из 1024 отсчетов №4.
5) Принимаем новую порцию из 1024 отсчетов.
6) Идем на шаг 3.
(Выходит пока буфер "протолкнет" последний отсчет, принятый в самом начале пройдет 4 цикла БПФ - усвоил).
и так пока солнце не сгорит.... Верно?

В MatLab можно с помощью утилиты Filter Design & Analysis Tool можно рассчитать функцию окна. Там же можно узнать АЧХ этой функции и кучу разных параметров, вообщем вы знаете. Если эту оконную функцию экспортировать в Simulink, то будет построен стандартный прямой КИХ фильтр (вот почему я упоминал слово фильтр).

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

Не пытаюсь спорить с вами, просто хотелось бы знать: Реализованный КИХ фильтр на основе весовой функции в MatLab и весовая функция на основе кольцевых буферов дадут одинаковые результаты?

Да, использовать буду RealFFT.

P.S. Спасибо за переживания. Никто меня не эксплуатирует. Тема интересная, потому и взял. Времени мало, так что, скорее всего, прийдется ограничиться микропроцессорной системой, кусками кода и алгоритмом. Все же знают, что у студентов на бумаге все работает. Так что после университета будем развивать тему дальше

Сообщение отредактировал Nukers - Mar 7 2010, 21:03

Нет... не так...


...да...
...нет-нет-нет! Какой ужас ... Мы ничего никуда не копируем - его принцип неприкосновенен... кольцевой буфер живет своей жизнью - самый новый сэмпл затирает самый старый. Данные все поступают и поступают... И так без конца... После каждых новых 1024-х сэмплов мы его полностью переписываем в отдельный массив(4096 слов) - накладываем функцию окна(уже на этот массив) и делаем FFT(прям на месте этого массива обычно и получим результат - так работает большинство алгоритмов).
И храним до получения еще 3-х подобных...

над массивом...
...да...

ну они сами "затрутся"... на то буфер и кольцевой...


ну да... даже правильнее говорить "дождемся поступления еще 1024-х отсчетов"... Они могут как блочно приходить так и по-одному... зависит от ситуации "на местах"... ... Вообще за этим должна следить переменная а-ля семафор(unsigned int). Тем или иным способом в памяти появилось 1024 новых сэмпла(DMA, через процессор, другой процессор, как угодно...) - она меняет свое состояние заставляя систему работать. Пусть она считает по-блочно. Тогда это должно выглядеть как 3-2-1-0-3-2-1-0-... С каждым ее изменением мы производим перерасчет FFT, а с каждым "обнулением" происходит "дополнительное действие" - сложение 4-х предыдущих результатов FFT в одно...


Нет... на шаг "1"...


Верно-то оно верно - но... так как нарисовано у меня ...



Ммм...

Ээээ... Вы опять где-то что-то недочитали... Нет там никакого сдвига... Это не FIR алгоритм... Вот подумайте откуда берется прямоугольная функция окна по-умолчанию для FFT?... Да это же просто интервал рассмотрения - блок... Можно сказать, что все сэмплы на входе FFT умножены на "1"-цу. А вот при определенной оконной функции - не все на "1"-у ...


Каждый бин FFT по природе соей фильтр - поскольку это значение корреляции с синусоидой на определенном промежутке. По-умолчанию фйильтр этот с sinc(w)-подобной АЧХ. Наложив гладкую функцию окна на интервал рассмотрения Вы домножите АЧХ функции рассмотрения на АЧХ каждого бина в отдельности. Вот так они связаны. Не на прямую... И они не то что результаты разные дадут - они вообще разные вещи в данном контексте...


Это гут... Но их на сайте нет ...

Не могу не поучаствовать.
При наложении окна на выборку окно ничего не копирует 4 раза. Просто окно как бы сглаживает выборку покраям, уменьшая разрывы функции при FFT.
В матлаб для реализации оконной функции есть Window Design & Analysis Tool. И сглаживание, о котором я тебе говорил, не имеет отношение к КИХ-фильтрации. Просто окно может применяться для построения КИХ-фильтров.

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

P.S. Не путай окно как функцию спадающую покраям (Хэмминг, Ханн, Блэкмэн...) и временное окно (время за которое берётся определённое количество выборок).
Стучи в ICQ!)

----------------------------------------------------------------------------------------------------------------------------------------------

То, насколько я знаю проблему, позволяет мне сформулировать вопрос.
Сигнал поступающий на осциллограф может иметь полосу до 100кГц.
Остаётся открытым вопрос о длине выборки. Так как сигнал заведомо неизвестен.
Допустим сигнал содержит составляющие 100Гц, 105Гц и 302 Гц.
если N = 4096 то частотное разрешение в спектре составит 200000/4096 = 48,82 Гц.
И мы просто не увидим с пектре гармоники 105, 302 Гц. Они просто переползут на ближайшие бины БПФ. и сильно могут искозить то, что должно приходиться на эти бины.
Можно сделать вывод, что нужно сначала задаться наименьшим возможным шагом гармоник в сигнале?
Как это решается в современных анализаторах спектра?

Могу предположить, что в выборки сигнала можно вставить нули....

Сообщение отредактировал TigerSHARC - Mar 8 2010, 07:29

Кто-нибудь может дать ссылку, где можно найти информацию по времени выполнения реального БПФ для ADSP-21375 при количестве точек 256, 512, 1024, 2048 и 4096?

Заранее огромное спасибо

Надуманная проблема. А вот если изменения звукового давления - не нарушая Найквиста-Котельникова - будут приходиться не на целые дискреты времени от начала выборки(уже смешно ) а на дробные? Перезапустить прибор? И так пока не угадаем?

То что Вы описываете на самом деле не проблема. Временной сэмпл или частотный бин - однозначно-обратимо-разделимо размазывается по соседним выборкам.
Но для отдельных случаев, например извечная тема - измерение Кг 50(60) Hz-овой синусоиды, это не очень удобно. Потому там и мухлюют с этим периодом измерений и числом точек FFT. Для любых "чуть более стохастических" сигналов с этим уже не заморачиваются...
Единственное что - иногда при частоте дискретизации не кратной степени двойки (например кратной 44100) хотят получить кратный десяткам(удобно для восприятия человеком) шаг частотных бинов. Тогда нужно писать специфическое FFT - не для степеней 2-ки...


Никак... Такой проблемы нет.


Куда? Интересно бы мне было узнать...


В листинге самих исходников приведены такты... ...Кстати "реального" FFT, насколько помню там нет(а оно таки есть)... Но для 1024 точек RealFFT - заведомо быстрее ComplexFFT для 1024 точек. Прикиньте - не ошибетесь...

Возьмите RealFFT из этой библиотеки. Далее есть 2-а пути - так и оставить - сомневаюсь что ассемблерный результат будет колоссально выше. Второй путь - убрать внутреннее комплексное БПФ (а оно там есть по-любому) и вставить ассемблерное с сайта AD...

Спаасибо огромное. Но я нашел лишь для 21364. А для 375-ого нет. Не ткнете пальцем?

ADSP-2136x Application Code Examples...

Кстати - Single-Channel Real Rad2 FFT - там есть...

Ну думаю теперь всё решено.

Уважаемый DRUID3. Есть к вам вопрос.
Уолт Кестер в Главе 5, стр. 22 Быстрое Преобразование Фурье http://www.analog.com.ru/Public/5.pdf пишет, что максимальное время выполнения БПФ при реализации в реальном масштабе времени рассчитывается как:

Tmax=N/fs

Таким образом при N=1024 и fs=250кГц получаем Tmax=4096мкс.

Я понял так:

1) Включили прибор
2) Подождали неск. тактов пока закончатся переходные процессы (в осн. для аналоговой части).
3) Копим первые 1024 отсчета с fs=250кГц - это и занимает 4096мкс. (В это время подпрограмма выполнения БПФ просто не вызывается=процессор простаивает).
4) Обрабатываем 1024 полученных отсчета по алгоритму БПФ. (Для ADSP-21375 при 266МГц и 1024 точках это занимает ~27мкс).
5) На пункт 3.
...и.т.д.

Скажите, я не ошибся? Верно ли, что пока копятся первых 1024 отсчета у меня есть 4096мкс, чтобы выполнить БПФ на 1024 точки, а так же на сопутствующие операции?

Сообщение отредактировал Nukers - Mar 15 2010, 13:10