Хабр Курсы для всех
РЕКЛАМА
Практикум, Хекслет, SkyPro, авторские курсы — собрали всех и попросили скидки. Осталось выбрать!
Хабр доступен 24/7 благодаря поддержке друзей
Комментарии 23
Для увеличения разрешающей способности можно так же применить super-resolution FFT методы, например WOLA FFT (по-сути filter banks). Интересно, вы сравнивали эти подходы?
Да, WOLA FFT (или полифазные схемы БПФ) я применял в других проектах. «Супер-разрешение» достаточно опасное понятие в этом плане, поскольку отражается только в визуальном улучшении качества спектра, а не в настоящей разрешающей способности по частоте.
На результирующих спектрах после WOLA в некоторых случаях хорошо видны усредненные гармоники с сравнении с обычным БПФ, оно и понятно — увеличение времени наблюдения дает больше информации. Но разрещающая способность по частоте fs/N остается прежней! А для очень близких гармоник может произойти их наложение из-за применения оконной функции (расширяется главный лепесток спектра). Здесь подробно обсуждались нюансы WOLA в задачах Super-Resolution, а тут — теория.
Все эти свойства я описывал в своей последней лекции по ЦОС. Вот графики:
Эффекты:
1) Потеря первой гармоники из-за некратности фазы при разбиении сигнала N на M пачек длины K (FFT lenght = K).
2) С фазой все хорошо -график приличный.
3) Переход к наложению оконной функции, если сигналы далеко — все хорошо.
4) Неверный выбор оконной функции. Слипание из-за наложения окна.
5) Без оконной функции, но с фазами снова все хорошо.
6) Сдвигаем еще ближе, выполняя условие по фазам. Ура. Super-resolution.
Кроме того, WOLA — не честный реал-тайм и выдает спектр с периодом M = N / K.
Применять этот метод следует очень аккуратно и поставленную задачу он не решает.
На результирующих спектрах после WOLA в некоторых случаях хорошо видны усредненные гармоники с сравнении с обычным БПФ, оно и понятно — увеличение времени наблюдения дает больше информации. Но разрещающая способность по частоте fs/N остается прежней! А для очень близких гармоник может произойти их наложение из-за применения оконной функции (расширяется главный лепесток спектра). Здесь подробно обсуждались нюансы WOLA в задачах Super-Resolution, а тут — теория.
Все эти свойства я описывал в своей последней лекции по ЦОС. Вот графики:
Эффекты:
1) Потеря первой гармоники из-за некратности фазы при разбиении сигнала N на M пачек длины K (FFT lenght = K).
2) С фазой все хорошо -график приличный.
3) Переход к наложению оконной функции, если сигналы далеко — все хорошо.
4) Неверный выбор оконной функции. Слипание из-за наложения окна.
5) Без оконной функции, но с фазами снова все хорошо.
6) Сдвигаем еще ближе, выполняя условие по фазам. Ура. Super-resolution.
Кроме того, WOLA — не честный реал-тайм и выдает спектр с периодом M = N / K.
Применять этот метод следует очень аккуратно и поставленную задачу он не решает.
1) Да, расстояние между бинами остается прежним, тут super-resolution никакого нет.
2) Мне не совсем понятно, зачем в WOLA использовать оконные функции, ведь принцип действия заключается в том, что бы каждый бин отфильтровать «хорошим» фильтром, т.е. уйти от фильтра с прямоугольной импульсной характеристикой.
3) Я не знаком с ограничениями WOLA на реал-тайм, можно хоть каждый входной отсчет учитывать. Может быть я не правильно понял Вашу формулировку?
4) Требование разрешения до единиц Гц — это продуктовая задача или академическая? Ниже прочитал Ваш ответ про астрономию. Почитаю.
2) Мне не совсем понятно, зачем в WOLA использовать оконные функции, ведь принцип действия заключается в том, что бы каждый бин отфильтровать «хорошим» фильтром, т.е. уйти от фильтра с прямоугольной импульсной характеристикой.
3) Я не знаком с ограничениями WOLA на реал-тайм, можно хоть каждый входной отсчет учитывать. Может быть я не правильно понял Вашу формулировку?
4) Требование разрешения до единиц Гц — это продуктовая задача или академическая? Ниже прочитал Ваш ответ про астрономию. Почитаю.
2) Если не задать окно или выбрать его неправильно — то возможна ситуация потери гармоник.
3) Возможно, я некорректно выразился. WOLA может работать в реал-тайме, то есть может получать новые отсчеты на каждом такте. Но выдавать результат пользователю будет в соответствии с разбиением исходной последовательности с периодом M = N / K. В алгоритме сверхдлинного БПФ — результат приходит и уходит на каждом такте.
4) Продуктовая. Описанный в статье метод применяется в реальных задачах. К сожалению, не могу раскрывать детали из-за NDA.
3) Возможно, я некорректно выразился. WOLA может работать в реал-тайме, то есть может получать новые отсчеты на каждом такте. Но выдавать результат пользователю будет в соответствии с разбиением исходной последовательности с периодом M = N / K. В алгоритме сверхдлинного БПФ — результат приходит и уходит на каждом такте.
4) Продуктовая. Описанный в статье метод применяется в реальных задачах. К сожалению, не могу раскрывать детали из-за NDA.
2) Мне все равно пока не очевидно применение окон. Ладно, попробую на досуге поэкспериментировать.
3) WOLA может работать в потоковом режиме — результат так же может быть получен на каждом такте. Тут пример реализации WOLA, D может быть равен 1. Мы обсуждаем одно и тоже?
3) WOLA может работать в потоковом режиме — результат так же может быть получен на каждом такте. Тут пример реализации WOLA, D может быть равен 1. Мы обсуждаем одно и тоже?
3) Да, схема правильная.
Смотрите: на вход поступает D = K * M отсчетов. Они разбиваются на последовательности M по K сэмплов. Затем сумируются и вычисляется БПФ (поворот фазы опустим, он необязателен). На выходе получаете пачку длины K. Чтобы получить вторую пачку данных длины K нужно снова набрать и обработать D = K * M.
Смотрите: на вход поступает D = K * M отсчетов. Они разбиваются на последовательности M по K сэмплов. Затем сумируются и вычисляется БПФ (поворот фазы опустим, он необязателен). На выходе получаете пачку длины K. Чтобы получить вторую пачку данных длины K нужно снова набрать и обработать D = K * M.
При D=1 схема перестанет работать?
Да, верно. Алгоритму нечего взвешивать и суммировать.
А как же измененный x[]? Тот факт, что вектор обновился ни на что не повлияет? Это схема является разверткой фильтрации на каждой поднесущей. Один новый отсчет на входе, как в FIR фильтре, обновляет результат на выходе.
Прошу прощения, я неправильно понял назначение D на этой картинке. Тут D отвечает за степень параллелизма подачи входных данных, если D = 1, то все будет работать, просто будет поступать 1 отсчет за такт, а не D.
Но это дополнительная полифазная обработка (например, когда у вас частота дискретизации 4 ГГц, а в ПЛИС вы можете работать только на 400 МГц, тогда D = 10). То есть это просто аспект реализации на скорость обработки в FPGA.
Таким образом, на add стадию после overlap этот D не влияет: вам все равно придется ждать всю последовательность, чтобы произвести суммирование.
Но это дополнительная полифазная обработка (например, когда у вас частота дискретизации 4 ГГц, а в ПЛИС вы можете работать только на 400 МГц, тогда D = 10). То есть это просто аспект реализации на скорость обработки в FPGA.
Таким образом, на add стадию после overlap этот D не влияет: вам все равно придется ждать всю последовательность, чтобы произвести суммирование.
Хорошая работа.
Если пойти дальше и собрать длиное Фурье не в двухмерное, а трёх-, четырёх-, N- (нет N мало, возьмём K, пока размер не станет 2), то можно дойти до того же "оптимального" соотношения логики(умножителей) к памяти сколько есть в наличии у ПЛИС или максимальная частота при дальнейшем разложении тоже начнёт проседать?
И гораздо интереснее, что это за задачи такие, где надо иметь спектр сигнал с полосой в сотню МГц с разрешением меньше Гц причём во всей полосе сразу?
По сути дальнейшее разложение приведет к «разворачиванию» схемы до элементарной бабочки (Radix-2, например). То есть изначальная конвейерная схема соединения «бабочка — узел перестановки» превратится в многомерную (параллельную?) форму, где бабочки будут только Radix-2, а узлы кросс-коммутации разложатся по памяти, но это выглядит нереализуемым по ресурсам ПЛИС.
Для примера взять БПФ N = 8. В одномерном случае это 3 стадии бабочек. На каждой стадии выполнится четыре Radix-2 операции. В сумме — 12 операций. В двумерном случае N = K1 * K2 = 4 * 2. Для K1 надо 8 операций, для K2 = 4. Итого 12 операций. Для трехмерного БПФ получаем элементарный куб K1 * K2 * K3 = 2 * 2 * 2, который тоже выполняет 12 операций. Но в многомерном случае приходится где-то хранить всю выборку N после прохода по одному измерению (и столько раз по числу измерений).
Например, в астрономии в задачах спектроскопии.
Для примера взять БПФ N = 8. В одномерном случае это 3 стадии бабочек. На каждой стадии выполнится четыре Radix-2 операции. В сумме — 12 операций. В двумерном случае N = K1 * K2 = 4 * 2. Для K1 надо 8 операций, для K2 = 4. Итого 12 операций. Для трехмерного БПФ получаем элементарный куб K1 * K2 * K3 = 2 * 2 * 2, который тоже выполняет 12 операций. Но в многомерном случае приходится где-то хранить всю выборку N после прохода по одному измерению (и столько раз по числу измерений).
Например, в астрономии в задачах спектроскопии.
Я бы сначал сделал filter bankом.
Еще для очень длинных FFT просто изумительно подходят GPU вычисления на Nvidia.
Дешево и сердито.
Я сам разработчик на ПЛИС — а похвалил видеокарту :).
Еще для очень длинных FFT просто изумительно подходят GPU вычисления на Nvidia.
Дешево и сердито.
Я сам разработчик на ПЛИС — а похвалил видеокарту :).
Некоторе время назад проектировал коррелятор для приёма широкополосного сигнала с кодовым разделением и большой базой (DSSS), полоса сигнала около 5 МГц, пришлось разбивать на небольшие кусочки и искать свёртку в каждым канале. Вычислить БПФ от всего сигнала было очень проблематично было. Теперь по-новому взглянул на построение коррелятора. Спасибо за статью
Удивительная статья! Но как вы вычисляете сами множители для БПФ, чтобы не страдала точность?
Возьмём экстремальный случай — 256 млн (точнее, 256*2^20 = 2^28 = 268435456) отсчётов. Для 1-й гармоники множитель сдвига фазы на 1 отсчёт равен:
— действительная часть буквально находится на границе точности формата double: в этом формате она равна 0x3FEFFFFFFFFFFFFE — вся мантисса, кроме самого младшего бита, состоит из единиц. Сколько бит нужно данному алгоритму, чтобы оперировать такими множителями?
Возьмём экстремальный случай — 256 млн (точнее, 256*2^20 = 2^28 = 268435456) отсчётов. Для 1-й гармоники множитель сдвига фазы на 1 отсчёт равен:
exp(2*pi*i/2^28) =
0.999999999999999726063448749220403437928237007258173656276+
2.340668926827455275950549341903484403788620722377879033845e-8*i— действительная часть буквально находится на границе точности формата double: в этом формате она равна 0x3FEFFFFFFFFFFFFE — вся мантисса, кроме самого младшего бита, состоит из единиц. Сколько бит нужно данному алгоритму, чтобы оперировать такими множителями?
Вы правы — для таких длин БПФ нужна очень высокая точность промежуточных вычислений, в частности для поворачивающих множителей. С одной стороны, для этого можно использовать нестандартные самописные floating-point вычислители, у которых мантисса шире, чем у double. Или все вычисления проводить в fixed-point с линейным ростом разрядности от бабочки к бабочке (а в FPGA это сделать проще, чем в floating-point). С другой стороны — я не могу придумать реальную задачу, где действительно была бы нужна такая точность, поэтому промежуточные результаты вычислений и хранимые поворачивающие множители можно округлять.
Александр, просветите: для каких задач необходимо сверхдлинное БПФ? Ведь частотное разрешение считается по формуле Fs/N, где N = 2^FftOrder — число точек. Если у Вас N = 2^28, а Fs, к примеру, 128 МГц = 2^27 Гц, то Вы получите две точки на герц! Зачем нужна столь подробная детализация?
Найти прикладное применение для таких длинных БПФ непросто, но оно есть: астрономия (в задачах спектрометрии частоты дискретизации начинаются от гигагерц), медицина (ЭЭГ, например). В классической радиотехнике тоже можно использовать. Я делал сверхдлинный БПФ для задачи радиоприема (к сожалению, ответить подробнее не могу).
Мне тут подсказывают
А не БПФ-алгоритм ли это Кули-Тьюки ?
Зарегистрируйтесь на Хабре, чтобы оставить комментарий
Сверхдлинное преобразование Фурье на FPGA