[PriBoris 0]

Какой эффективный алгоритм для длин такого типа (2^n+2^(n-1)) выбрать ?

Вопрос возник при попытке реализации фильтрации с помощью быстрой линейной свёртки. Длина фильтра почти 256 точек, данные поступают в реальном времени блоками по 512.

Изменено 27 июля, 2010 пользователем PriBoris

[Xenia 45]

Вопрос возник при попытке реализации фильтрации с помощью быстрой линейной свёртки.

Быть может здесь было бы проще дополнить массив 768 нулями до массива 1024 и действовать обычным образом?

[Связист 0]

Какой эффективный алгоритм для длин такого типа (2^n+2^(n-1)) выбрать ?

Вопрос возник при попытке реализации фильтрации с помощью быстрой линейной свёртки. Длина фильтра почти 256 точек, данные поступают в реальном времени блоками по 512.

 

Не совсем понятно что Вы делаете?

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

 

Про эффективные алгоритмы лучше всего посмотреть книжку Блейхута "Быстрые алгоритмы цифровой обработки сигналов".

Алгоритм Кули-Тьюки наверное подойдёт для Вашей задачи. Получится табличка 3*256 точек. 3-точечные БПФ делаем по плгоритму Винограда, а 256-точечное БПФ делаем бабочками (Кули-Тьюки по основанию 2).

Или алгоритм Гуда-Томаса тоже подойдёт, но его понять сложнее, поэтому разбираться надо дольше. Там тоже получится табличка 3*256 точек.

[PriBoris 0]

Быть может здесь было бы проще дополнить массив 768 нулями до массива 1024 и действовать обычным образом?

Да,я так и хочу сделать. Просто на всякий случай спросил, вдруг есть возможность сэкономить.

 

[DRUID3 0]

Да,я так и хочу сделать. Просто на всякий случай спросил, вдруг есть возможность сэкономить.

...это ерунда, причем полная. Вон товарищ связист дело говорит...

[thermit 1]

768 = 3*2^8

последняя стадия по основанию 3, остальные - по основанию 2. или 16.

C перестановками, правда будет возня. Может, действительно проще будет взять на 1к точек?

[Xenia 45]

Да,я так и хочу сделать. Просто на всякий случай спросил, вдруг есть возможность сэкономить.

Если вы достигните своего желания получить FFT на 768-точечном массиве, то у вас получится циклическая свертка. А вам оно надо? В большинстве практических случаев нужна простая свертка, а не циклическая, поскольку образ, находящийся в массиве 768, вряд ли периодический.

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

Более того, даже если бы ваш массив был длиной 512, что позволило бы легко преобразовать его в FFT-образ, то и в том случае стоило бы подумать о целесообразности дополнения его нулями до 1024, чтобы избавиться от циклической свертки. А у вас длина массива просто идеальна для получения нормальной свертки, т.к. 768+256=1024.

[Связист 0]

Если вы достигните своего желания получить FFT на 768-точечном массиве, то у вас получится циклическая свертка. А вам оно надо? В большинстве практических случаев нужна простая свертка, а не циклическая, поскольку образ, находящийся в массиве 768, вряд ли периодический.

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

Более того, даже если бы ваш массив был длиной 512, что позволило бы легко преобразовать его в FFT-образ, то и в том случае стоило бы подумать о целесообразности дополнения его нулями до 1024, чтобы избавиться от циклической свертки. А у вас длина массива просто идеальна для получения нормальной свертки, т.к. 768+256=1024.

 

Кажется Вы пытаетесь сказать новое слово в цифровой обработке сигналов...

Насколько мне известно, раньше пользовались окнами для устранения краевых эффектов.

А вот дополнение нулями это конечно интересно, но требует количественной оценки побочных эффектов и сравнение его с окнами.

 

[Xenia 45]

Насколько мне известно, раньше пользовались окнами для устранения краевых эффектов.

Цикличность свертки это и есть причина краевых эффектов :)

[bahurin 0]

Какой эффективный алгоритм для длин такого типа (2^n+2^(n-1)) выбрать ?

Длина фильтра почти 256 точек, данные поступают в реальном времени блоками по 512.

 

если надо поставить фильтр, например ФНЧ, то быстрее всего это БИХ - фильтр рассчитать и использовать. В этом случае выши 256 точек сведутся к БИХ фильтру порядка не больше 10.

 

Если же очень хочется ких фильтр, то тогда 512 точек сигнала переводите в спектр и потом уже спектр умножаете на частотную характеристику фильтра, и переводите через ifft. Как-то я не нашел в вашем вопросе 738 точек, поскольку:

данные поступают в реальном времени блоками по 512.

 

Изменено 27 июля, 2010 пользователем bahurin

[Xenia 45]

Кажется Вы пытаетесь сказать новое слово в цифровой обработке сигналов...

Попытка не пытка! Тогда уж я еще одно новое слово скажу :)

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

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

Однако FFT-алгоритм нам на зло станет рассматривать этот массив, как свернутый в кольцо, а потому усмотрит в нем чудовичный разрыв непрерывности в том месте, где происходит склейка в кольцо. Из этого "уступа" в частотном образе сигнала появятся большие вклады не только высокочастностных составляющих, но и низкочастотных. И все лишь для того, чтобы этот уступ воспроизвести. Можно, конечно, окнами с ними бороться - это традиционный путь. А вот я "изобрела" совсем простой путь, но исключительно эффективнный. Он состоит в том, что надо дополнить массив (примерно в два раза до следующего размера, кратного степени числа 2), но только не нулями! Вместо нулей следует в этом месте провести НАКЛОННУЮ ПРЯМУЮ ЛИНИЮ, которая одним свои концом исходит из последней точки исходных данных, а другим концом должна дойти до конца рассширенного массива именно на тот уровнь, соотвествующий первой точке исходных данных. Тогда при склейке в кольцо окажется, что переход получился максимально плавным.

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

[PriBoris 0]

Как-то я не нашел в вашем вопросе 738 точек

Для вычисления фильтрованного значения первой точки свежего блока используются 255 (256-1) точек предыдущего блока.

[bahurin 0]

А вот я "изобрела" совсем простой путь, но исключительно эффективенный. Он состоит в том, что надо дополнить масссив (примерно в два раза до следующего размера, кратного степени числа 2), но только не нулями! Вместо нулей следует в этом месте провести НАКЛОННУЮ ПРЯМУЮ ЛИНИЮ, которая одним свои концом исходит из последней точки исходных данных, а другим концом должна дойти до конца рассширенного массива именно на тот уровнь, соотвествующий первой точке исходных данных. Тогда при склейке в кольцо окажется, что переход получился максимально плавный.

 

1. вы нелинейно исказили ваш сигнал. Если вы например возьмете синусоиду и целое количество отсчетов на период этой синусоиды, то не увидите эффектов растекания спектра и получите одну палку на нужной частоте. Применив свой чудо метод вы получите вместо синусоиды синус с разрывами и спектр будет отличаться от истинного из-за нелинейного искажения исходного сигнала, другими словами вы получите кучу ложных гармоник и никак не сможете их контроллировать.

 

2. Если уж вы хотите сказать свое слово, то это слово должно быть обоснованным. А так получается что вы усилием воли решили а давайте сделаем линейную интерполяцию. А почему собственно линейную? А может лучше будет если кубическим полиномом соединить, тогда можно еще и производные выровнять и устранить точки перелома? А если полином 5 степени или 7-ой? И почему вы удваиваете длину выборки, а не увеличиваете ее в 4 раза например или в 8? Короче не убедительно это все. Кроме того, я бы еще проанализировал, что будет если просто удвоить количество отсчетов и наложить сверху весовое окно? Мне кажется спектр чище будет, чем с вашей интерполяцией.

 

3. В отличии от вашего метода умножение на весовое окно если линейная операция, поскольку спектр в этом случае сворачивается со спектром окна.

 

Для вычисления фильтрованного значения первой точки свежего блока используются 255 (256-1) точек предыдущего блока.

если у вас предыдущие 512 отсчетов и следующие 512 во времени не разрываются, тогда бих фильтр вам предпочтительнее. Если они имеют дырку между собой, то брать предыдущие отсчеты нельзя. Если же все таки хотите через fft, то вам необходимо применить обработку с перекрытием. тогда надо делать fft на 512 отсчетов, умножение на чх фильтра, затем ifft, потом сдвиг на 256 отсчетов и снова fft на 512. Т.о. вы получите аналог ких фильтра, но суть в том, что надо делать fft на удвоенную длительность всегда, чтобы было перекрытие.

Изменено 27 июля, 2010 пользователем bahurin

[DRUID3 0]

Попытка не пытка! Тогда уж я еще одно новое слово скажу :)

В чем, собственно, недостаток FFT в практическом смысле? В прежде всего в его циклическом характере базисных функций и самого способа построения ортонормированного базиса. FFT преобразует массив так, как будто он закольцован, хотя на самом деле этого нет. - не так FFT, как то, что мы делаем всю математику подразумевая sin(t)/cos(t) - а они бесконечные циклические функции.

Т.е. проблему решит не какое-то дополнение нулями - ну же, Xenia, у меня блок где замерян постоянный ток. 11111...111 (V) - как его не дополняй - спектр станет заметно отличным от простой палки на 0-ой частоте... Или вон пример который привел товарищ bahurin - тоже очевиден для мыслеэксперимента.

Проблему решит:

1) Выбор других базисных функций и метода построения базиса.

2) Выбор "иного" математического мира - и так между прочим тоже делают...

 

Остальное прочел... Здорово. Но это частные случаи.

 

P.S.: был задан конкретный вопрос и товарищ связист дал на него самый компетентный ответ - смысле на самый первый вопрос(уточнение аФФтАр топика дал просто чумовое ). Дальше понеслось - и дополнение 0-ми, и IIR... Сейчас пойдут теоретико-числовые алгоритмы и кольца Галуа... :laughing:

 

P.P.S.: Xenia, перепишите DCT(FCT) для своего mp3 плеера или JPEG просмотрщика - так "грамотно". С "не кольцевой" сверткой. Послушайте(посмотрите) его минут 10. А выводы опубликуйте на форуме...

 

P.P.P.S.:

Для вычисления фильтрованного значения первой точки свежего блока используются 255 (256-1) точек предыдущего блока.

Мда... Это диагноз... :( PriBoris, без обид, но Вы жжоте...

Вопрос возник при попытке реализации фильтрации с помощью быстрой линейной свёртки. Длина фильтра почти 256 точек, данные поступают в реальном времени блоками по 512.

Почти это как? 255.5?

 

делайте по 2-а FFT по 256 точек...

или одно по 512 на блок...

 

Физический смысл быстрой свертки в том, что у Вашего фильтра есть определенные АЧХ/ФЧХ(которое однозначно связано с его "импульсной"). Мы определяем его и накладываем операцией умножения(свертка-фильтрация над временным рядом соответствует умножению в частотном отображении) на частотный спектр сигнала. Размер блока может быть любым - впринципе, даже меньше длинны фильтра - но тогда его нужно будет пересемплировать с другой частотой Ничего из уже обработанного блока хранить не надо. Это не FIR. А все ускорение за счет того, что переход в частотную область с помощью FFT, посемпловое умножение, и обратный переход уже где-то начиная с 128 точек требует меньше операций чем свертка "в лоб".

 

Кстати блоки по 256 точек, 32-х битная система и тип данных int(если у Вас, PriBoris, так) - это и в самом деле лучше на кольцах Галуа. Зачем нам комплексность и округленные значения sin(t)/cos(t) !? :laughing:

 

В чем вопрос? :laughing:

[Xenia 45]

1. вы нелинейно исказили ваш сигнал. Если вы например возьмете синусоиду и целое количество отсчетов на период этой синусоиды, то не увидите эффектов растекания спектра и получите одну палку на нужной частоте. Применив свой чудо метод вы получите вместо синусоиды синус с разрывами и спектр будет отличаться от истинного из-за нелинейного искажения исходного сигнала, другими словами вы получите кучу ложных гармоник и никак не сможете их контроллировать.

Где ж это в практических случаях бывает, чтобы изменяемая гармоника да уложилась бы что в целое число периодов в массив данных? Вы привели выгодный для себя пример :). А давайте я вам приведу пример, когда гармоника (предположим, косинусоида) заканчивается "неудачно". Скажем начала косинусоида значение +1, а занончилась -1 (минус один), пол периода ей не хватило. Знаете, как гадость тогда получается? Одной палкой даже не пахнет :). А вот мой метод с такими случаями отлично справляется. И пусть там не чистая палка, но боковых липествов посчи нет.

И что из того, что нелинейно? Откуда у вас такая приверженность к линейности и отвращение к нелинейности? Линейным преобразованием из сигнала даже большую гадость можно сотворить, чем линейным :). Опять же сейчас речь идет не об усилительных каскадах, где нелинейные искажения почитаются грехом, а о совершенно иной задаче.

 

2. Если уж вы хотите сказать свое слово, то это слово должно быть обоснованным. А так получается что вы усилием воли решили а давайте сделаем линейную интерполяцию. А почему собственно линейную? А может лучше будет если кубическим полиномом соединить, тогда можно еще и производные выровнять и устранить точки перелома? А если полином 5 степени или 7-ой? И почему вы удваиваете длину выборки, а не увеличиваете ее в 4 раза например или в 8? Короче не убедительно это все. Кроме того, я бы еще проанализировал, что будет если просто удвоить количество отсчетов и наложить сверху весовое окно? Мне кажется спектр чище будет, чем с вашей интерполяцией.

Основания моих слов имеются, причем они примерно такие же, как причины выбора формы окна при фильтровании в частотной области. Ведь никто же не режет отсечкой с некоторой частоты, поскольку знает, что за этим поледует. А последует "зашумленность" той самой граничной частой (при обратном FFT), по которй обрезали. Поэтому если и давят частоты, то осторожно - так, чтобы коэффициент подавления сходил к нулю плавно. И тут кто во что горазд: и треугольное окно, и косинусоидное, и гауссовый профиль, и т.п.

Дополнительный треугольник, достоенный в исходной области данных - наименьшее зло, по сравнению со всеми остальными случаями, когда краевые точки массива сильно отличаются по величине. Причем сам треугольник ведет себя на редкость смирно. Можно построить "домик" (линейное возрастание до середины массива с убыванием его во второй половине массива к начальному значению). Тут "мусор", порождаемый таким треугольником слаб, а амплитуды частот изменяются непрерывно, без каких либо максимумов. Анализировать спектральные данные (наложенные на этот треугольник) это не мешает. Опять же естественные процессы линейного дрейфа сигнала есть ни что иное, как добавление аналогичного трегольника. Поэтому такую форму искажения тоже можно считать естественной, присущей реальным задачам обсчета временных трендов.