[task.solver 0]

Есть функция делающая БПФ (FFT) для комплексного float входа (для длин степени двойки).

 

У меня есть два вопроса. Искал ответы в интернете, но не нашёл. Нашёл только те же упоминания возможности. Особенно интересует второй вопрос.

 

1) Как можно сделать за один её вызов два БПФ для действительного входа той же длины?

Знаю что один вход надо записать в действительную часть, а второй в мнимую. Но что делать после вызова? Как разделить результаты?

(Кое что было написано тут, но непонятно: http://electronix.ru/forum/index.php?showtopic=71731)

 

2) Как можно сделать за один её вызов один БПФ для действительного входа двойной длины?

Этот вопрос важнее, чем первый.

 

Если можно напишите формулы.

 

Собственно второй вопрос навеян найденными высказываниями в интернете о подобной возможности. И вопросом о том, что действительное БПФ длины N можно реализовать эффективней чем комплексное БПФ с нулевыми мнимыми частями так же длины N.

 

Изменено 17 января, 2014 пользователем Task Solver

[thermit 1]

Есть функция делающая БПФ (FFT) для комплексного float входа (для длин степени двойки).

 

У меня есть два вопроса. Искал ответы в интернете, но не нашёл. Нашёл только те же упоминания возможности. Особенно интересует второй вопрос.

 

1) Как можно сделать за один её вызов два БПФ для действительного входа той же длины?

Знаю что один вход надо записать в действительную часть, а второй в мнимую. Но что делать после вызова? Как разделить результаты?

(Кое что было написано тут, но непонятно: http://electronix.ru/forum/index.php?showtopic=71731)

 

Что именно здесь непонятно?

2) Как можно сделать за один её вызов один БПФ для действительного входа двойной длины?

Этот вопрос важнее, чем первый.

 

Если можно напишите формулы.

 

Собственно второй вопрос навеян найденными высказываниями в интернете о подобной возможности. И вопросом о том, что действительное БПФ длины N можно реализовать эффективней чем комплексное БПФ с нулевыми мнимыми частями так же длины N.

 

N=1024;
x=randn(1,N);
ref=fft(x);
%Эти 2 бпф выполняются за Одно как в пункте 1
s1=fft(x(1:2:end));
s2=fft(x(2:2:end));

s=[s1+s2.*exp(-1j*2*pi/N*(0:N/2-1)) s1-s2.*exp(-1j*2*pi/N*(0:N/2-1))];

plot(abs(ref-s));

 

 

 

 

[ViKo 1]

В книжке Р. Лайонс "Цифровая обработка сигналов" про это написано.

 

 

[task.solver 0]

Что именно здесь непонятно?

Это не на Си, и не математические формулы. Непонятны записи. Что например вот это: x(2:2:end) ?

[thermit 1]

Task Solver:

Это не на Си, и не математические формулы. Непонятны записи. Что например вот это: x(2:2:end) ?

 

Это матлаб.

x - вектор

x(2:2:end) - вектор из x(2) x(4) ... x(end)

x(1:2:end) - вектор из x(1) x(3) ... x(end-1)

[task.solver 0]

Спасибо. Буду разбираться, пока до конца всё равно не понятно. Книжку я скачал только что. Тоже спасибо.

[V_G 11]

До применения математики важно понять принцип.

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

[Xenia 45]

1) Как можно сделать за один её вызов два БПФ для действительного входа той же длины?

Знаю что один вход надо записать в действительную часть, а второй в мнимую. Но что делать после вызова? Как разделить результаты?

(Кое что было написано тут, но непонятно: http://electronix.ru/forum/index.php?showtopic=71731)

2) Как можно сделать за один её вызов один БПФ для действительного входа двойной длины?

Этот вопрос важнее, чем первый.

 

FFT в своем алгоритмическом воплощении оставляет "мусор" из-за того, что в нем используются преобразования типа "бабочка". Поэтому в старшей части массива мы имеем второе крыло, зеркально отраженное от первого. Вот оно и сослужит нам добрую службу.

 

Рассмотрим самый примитивный пример (N=8), когда мы преобразуем комплексную функцию, разложенную на два массива:

double Re[8]; // это действительная часть исходной функции

double Im[8]; // это мнимая часть исходной функции

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

 

После FFT-преобразования результат получится в той же паре массивов, а исходные значения испортятся. Получится вот что:

Re[0] - постоянная составляющая, Im[0] - обычно остается нулем.

Re[1] и Im[1] - амплитуды сos и sin 1-ой гармоники.

Re[2] и Im[2] - амплитуды сos и sin 2-ой гармоники.

Re[3] и Im[3] - амплитуды сos и sin 3-ей гармоники.

Re[4] - амплитуда сos 4-ой гармоники. Это частота Найквиста, т.к. более высокие таблично не представимы. Im[4] - равна нулю, т.к. у частоты Найквиста не может быть sin-составляющей.

Re[5] и Im[5] - мусор, копия Re[3] и -Im[3], соответственно.

Re[6] и Im[6] - мусор, копия Re[2] и -Im[2], соответственно.

Re[7] и Im[7] - мусор, копия Re[1] и -Im[1], соответственно.

 

Т.е. мусор обладает той особенностью, что амплитуда cos-гармоник отражается зеркально, сохраняя свою амплитуду. А амплитуды sin-гармоник при отражении изменяют свои знаки. Это происходит потому, что функция cos четная (симметричная относительно вертикальной оси), а функция sin нечетная (симметричная относительно точки).

 

Если у нас есть две действительные функции, то мы можем преобразовать их вдвоем сразу, за один прогон алгоритма FFT. Для этого надо первую действительную функцию запихнуть в Re[], как это и было раньше, а вторую записать в Im[], который раньше был обнулен.

 

После преобразования получим результат преобразования, соответствующий входной функции F1+i*F2 (множитель i из-за того, что F2 в мнимую часть была записана).

Однако, используя свойство зеркальности, можем выделить из суммы результаты, относящиеся в F1 и F2 по отдельности. Делается это так:

 

Re[0] – постоянная составляющая F1.

Im[0] – постоянная составляющая F2.

 

(Re[1]+Re[7])/2 - cos-амплитуда F1 для 1-ой гармоники, т.е. средняя арифметическая результата и зеркального мусора.

Почему так? А потому что cos-амплитуда F1 симметрична и от этой операции и останется прежней. Этой операцией она очистится от вклада sin-амплитуды F2, которая имеет противоположные по знаку амплитуды в Re[1] и Re[7].

(Re[1]-Re[7])/2 - sin-амплитуда F1 для 1-ой гармоники, тут уже симметричные косинусы взаимно уничтожились, а синусы остались живы.

 

Точно так же обрабатываем зеркальные пары:

(Re[2]+Re[6])/2 - cos-амплитуда F1 для 2-ой гармоники.

(Re[2]-Re[6])/2 - sin-амплитуда F1 для 2-ой гармоники.

 

(Re[3]+Re[5])/2 - cos-амплитуда F1 для 3-ой гармоники.

(Re[3]-Re[5])/2 - sin-амплитуда F1 для 4-ой гармоники.

 

Частоты Найквиста зеркального отражения не имеет, с ней проще:

Re[4] - cos-амплитуда F1 для частоты Найквиста.

Im[4] - cos-амплитуда F2 для частоты Найквиста.

 

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

[task.solver 0]

Про зеркальный спектр от действительных данных я знал. Благодаря этой особенности и линейности FFT от суммы двух сигналов я понял как вывести формулы для разделения спекторов. В общем с первым вопросом я разобался. Осталось разобраться со вторым. Всем, кто помогает - спасибо.

 

Появлился новый вопрос. Что такое частота Найквиста? Что она характеризует? И что будет, если в предыдущем примере Re[4] обнулить и вычислить обратное преобразование? Тоже получится чисто действительный сигнал? (С нулевой мнимой частью.)

Изменено 18 января, 2014 пользователем Task Solver

[Xenia 45]

Появлился новый вопрос. Что такое частота Найквиста? Что она характеризует? И что будет, если в предыдущем примере Re[4] обнулить и вычислить обратное преобразование? Тоже получится чисто действительный сигнал? (С нулевой мнимой частью.)

 

Частота Найквиста выглядит как знакопеременный шум:

..., +1, -1, +1, -1, +1, -1, +1, -1, +1, -1, +1, -1, +1, -1, ...

Оттого и формально считается косинусом, имеющим в первой половине периода +1, а во второй -1. Соответствующей ей синус был бы обязан начинаться с нуля, как и все синусы. И таковым он был бы и во второй части периода.

 

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

 

Частоту Найквиста можно вырезать и без FFT, если усреднять по парам соседних точек, т.к. именно этот вклад она и вносит.

 

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

 

Поэтому частота Найквиста вполне реальна и пренебрегать ею не следует.

[thermit 1]

 

 

 

Изменено 18 января, 2014 пользователем thermit

[task.solver 0]

Какие хорошие формулы. Теперь всё понятно. :rolleyes:

[ViKo 1]

Странно, что люди изобрели БПФ, который делает вдвое больше работы, чем нужно. :rolleyes:

[task.solver 0]

Наверно с комплексными числами формулы были проще. В оригинале там вообще корни из единицы любого кольца.

[utherVV 0]

 

Numerical Recipes in C

http://apps.nrbook.com/c/index.html

 

страница 511

"go to page: 511"

 

всё, что нужно