[tocha 0]

Делая FFT для реального сигнала длиной 2*N, мы получаем комплексный результат, причём выход имеет complex-conjugate симметрию. И реализовать его можно через комплексное FFT длиной N + небольшие манипуляции и на этом сэкономить вычислительные ресурсы.

Нужно сделать обратное преобразовние IFFT для комплексного входа длиной 2*N, который имеет свойство complex-conjugate симметрии. Соответствено результат будет реальным.

Правильно то, шо если мы сделали RFFT длиной 2*N упрощённым методом(CFFT длиной N), то и complex conjugate IFFT длиной 2*N мы можем вычислить через FFT длиной N?

 

 

 

[ivanpa 0]

похоже на правду...

а почему бы просто не проверить в матлабе?=)

[Alexey Lukin 0]

Да, это верно. Соответствующие процедуры имеются во всех распространённых библиотеках FFT. Называются real-valued FFT/IFFT.

[Xenia 35]

Делая FFT для реального сигнала длиной 2*N, мы получаем комплексный результат, причём выход имеет complex-conjugate симметрию. И реализовать его можно через комплексное FFT длиной N + небольшие манипуляции и на этом сэкономить вычислительные ресурсы.

Нужно сделать обратное преобразовние IFFT для комплексного входа длиной 2*N, который имеет свойство complex-conjugate симметрии. Соответствено результат будет реальным.

Правильно то, шо если мы сделали RFFT длиной 2*N упрощённым методом(CFFT длиной N), то и complex conjugate IFFT длиной 2*N мы можем вычислить через FFT длиной N?

Не совсем поняла суть вашего вопроса в части употребления обозначений N и 2*N. Например, если реального сигнал у вас длиной 2*N, то что вы принимаете за N? Это осталось непонятным. Тем не менее, я постараюсь кратко описать эту ситуацию.

 

Если мы пользуемся обычным классическим FFT, которое является по своей сути комплексным, то мы заполняем своим реальным сигналом действительную часть комплексного массива длиной N, а мнимую заполняем нулями. После преобразования получим на этом же месте две зеркальные половинки, из которых для дела достаточно первых половин. Зеркальная же часть здесь содержит избыточную информацию из-за того, что число коэффициентов (или частот) в результате FFT вдвое меньше числа точек. Т.е. на середине отрезка (N/2) находится частота Найквиста, выше которой частоты определить теоретически невозможно.

 

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

 

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

 

Частота Найквиста только не влезает, и ее приходится запихивать в пустующую мнимую часть при постоянной составляющей (частоты 0). Это всегда возможно сделать, поскольку частота Найквиста всегда действительна, т.к. мнимую/синусную части она никогда не содержит. Всю эту операцию обычно называют реальным FFT или RFFT. Но в книжках стараются не афишировать, как она вычисляется :), ограничиваясь примерами комплексного FFT.

 

RFFT экономнее, чем CFFT, как по занимаемой памяти (не получаем зеркального мусора и не выделяем для него места в памяти), так и по времени, но зато код становится очень противным :).

 

Обратное преобразование IRFFT устроено подобным же образом. Там тоже складывают пополам и прокручивают вместе на мясорубке за один проход :). Самое простое - не вникать во все эти премудрости, а взять готовый код. А то мне вникать пришлось, когда была задача написать RFFT и IRFFT на ассемблере x87, чтобы все расчеты были на регистрах сопроцессора и никакого другого обращения к памяти, кроме массива, не было. Как вспомню, так вздрогну :).

[tocha 0]

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

Именно так и сделал.

 

Обратное преобразование IRFFT устроено подобным же образом. Там тоже складывают пополам и прокручивают вместе на мясорубке за один проход :). Самое простое - не вникать во все эти премудрости, а взять готовый код. А то мне вникать пришлось, когда была задача написать RFFT и IRFFT на ассемблере x87, чтобы все расчеты были на регистрах сопроцессора и никакого другого обращения к памяти, кроме массива, не было. Как вспомню, так вздрогну :).

А теперь хочу сделать это (IRFFT), то есть irfft симметричного комплексного сигнала у которого будет реальный выход. Макимально эффективно с точки зрения вычислительных затрат.

Если можете ткнуть в готовый код, буду благодарен (в пятницу навскидку не нашёл - может плохо искал...).

То есть, если я правильно понял, вы такое (IRFFT) делали: обратное преобразование симметричного сигнала длиной 2*N с помощью N-point ICFFT плюс немного вычислений и этим самым экономили ресурсы? Если скажете да, буду тоже благодарен.

 

 

[Xenia 35]

А теперь хочу сделать это (IRFFT), то есть irfft симметричного комплексного сигнала у которого будет реальный выход. Максимально эффективно с точки зрения вычислительных затрат.

Если можете ткнуть в готовый код, буду благодарен (в пятницу навскидку не нашёл - может плохо искал...).

То есть, если я правильно понял, вы такое (IRFFT) делали: обратное преобразование симметричного сигнала длиной 2*N с помощью N-point ICFFT плюс немного вычислений и этим самым экономили ресурсы? Если скажете да, буду тоже благодарен.

IRFFT как бы прилагается к RFTP в качестве обратной функции, а потому желательно, чтобы они были "от одного производителя". Лучше поищите вы в том же месте, где нашли прямую RFFT.

 

Тут есть один нюанс: если действительный массив, длиной 2*N, обработать RFTP, а к результату этой обработки сразу применить обратное IRFTP, то может получиться два варианта, в зависимости от реализации алгоритма. В первом варианте вы получите назад свой исходный массив с вычислительной точностью. А во втором варианте вы получите его увеличенным по амплитуде в N раз. Вот из-за этого множителя и возникает "несовместимость" прямого и обратного RFFTP у разных реализаторов. Причем это касается не только RFFT, но и простого CFFT - увеличение в N раз у обоих.

 

А сама проблема сводится к тому, за чей счет подавление этого множителя N относить, если уж решено сделать так, чтобы после обратного преобразования возвращался точно исходный результат. И тут кто во что горазд. Например, в MatLab'е прямое FFT делают, не вводя поправки на множитель, а внутри процедуры обратного IFFT делают деление всех элементов массива на N в начале или конце процедуры (эффект от места делёжки не зависит). Такое решение кажется несправедливым по отношении к IFFT и сильно похоже на подгон результата. Можно было бы разделить коэффициент N и по "справедливости", если делить после каждого преобразования, (как прямого, так и обратного) на квадратный корень из N. Но это противно, т.к. надо вычислять корень.

 

Я же задалась целью выяснить, как должно быть на самом деле :). Для чего создала на всю длину массива синусоиду на полный период с известной мне амплитудой 1000 единиц, и посмотрела, какие коэффициенты дает прямая FFT. очевидно, что из такой синусоиды по всем канонам должен получиться единственный мнимый коэффициент, равной 1000, во втором по порядку элементе массива. А оно не получилось... Тогда пришлось модифицировать процедуры так: после FFT умножаю на множитель 2/N, а после обратной IFFT делю на 2. Такой прием дает как правильную амплитуду гармоник, так и поглощает множитель при обратном преобразовании, т.к. (2/N)/2 = 1/N.

 

Очень возможно, что в этом деле я что-то недопонимаю и изобретаю велосипед, но мне удивительно, что люди почему-то над этими вещами не задумываются, а вполне механически применяют процедуры из книг или библиотек функций. Может быть я тут где-то и пронеслась, но вам я советую ОБЯЗАТЕЛЬНО проверить прямое и обратное FFT на тестовом примере из синусоиды (или косинусоиды) с известной вам амплитудой, чтобы проверить FFT-результат на предмет множителя.

 

Подробнее обсуждать алгоритмы мне с вами сложно, т.к. вы не пишите ни процессор, на котором собираетесь гонять эти процедуры, ни язык, на котором собираетесь их программировать. А мои наработки для Pentium могут оказаться для ваш излишней информацией, т.к. там есть своя специфика. Например, получать значения "опорных" синусов и косинусов на Пентиуме быстрее через тригонометрическую функцию sincos, входящую в набор инструкций сопроцессора (!), которая за один проход дает оба результата - синус и косинус. Это оказалось быстрее, чем получать синусы и косинусы половинных углов по тригонометрической формуле для половинного угла. А, скажем, для DSP-процессора, у которого тригонометрических функций в наборе инструкций нет, это было бы уже плохим решением. А вот в тех случаях, когда N фиксировано, то набор синусов для N/4 (от 0 до 90 градусов) лучше вообще зашить в ПЗУ и не заниматься их вычислением. Иными словами, здесь есть очень много нюансов, которые зависят от сопутствующих обстоятельств.

[tocha 0]

IRFFT как бы прилагается к RFTP в качестве обратной функции, а потому желательно, чтобы они были "от одного производителя". Лучше поищите вы в том же месте, где нашли прямую RFFT.

 

Тут есть один нюанс:...

Надо для 32-bit fixed point, коеффициенты - в табличке. Вопросы масштабирования пока не волнуют. Будем искать....Вам спасибо за ответ, Ксения.