если результаты отличаются от исходного сигнала

Question

если результаты отличаются от исходного сигнала

БПФ работает нормально, но когда я хочу взять IFFT, я всегда вижу один и тот же график из его результатов. Результаты сложны, и график всегда одинаков независимо от исходного сигнала.

в реальной части графа -син с периодом = размер кадра

в мнимой части это -кос с тем же периодом

Где может быть проблема?

исходный сигнал:
оригинальный сигнал

IFFT реальная стоимость (на фото только половина кадра):

Алгоритм БПФ, который я использую.

double** FFT(double** f, int s, bool inverse) {
if (s == 1) return f;
int sH = s / 2;

double** fOdd = new double*[sH];
double** fEven = new double*[sH];
for (int i = 0; i < sH; i++) {
int j = 2 * i;
fOdd[i] = f[j];
fEven[i] = f[j + 1];
}

double** sOdd = FFT(fOdd, sH, inverse);
double** sEven = FFT(fEven, sH, inverse);

double**spectr = new double*[s];

double arg = inverse ? DoublePI / s : -DoublePI / s;
double*oBase = new double[2]{ cos(arg),sin(arg) };
double*o = new double[2]{ 1,0 };

for (int i = 0; i < sH; i++) {
double* sO1 = Mul(o, sOdd[i]);

spectr[i] = Sum(sEven[i], sO1);
spectr[i + sH] = Dif(sEven[i], sO1);

o = Mul(o, oBase);
}

return spectr;
}

1

c++fft ifft

Решение

Другие решения

Что касается памяти, вы можете использовать std::vector инкапсулировать динамически распределяемые массивы и гарантировать их освобождение, когда выполнение выходит из области видимости. Вы могли бы использовать unique_ptr<double[]> но прирост производительности не стоит ИМО, и вы теряете безопасность at() метод.

(На основе ответа @ Робба)

Несколько других советов:

Избегайте загадочных идентификаторов — программы должны быть читаемыми, а имена похожи на «f» а также «s«сделать вашу программу труднее читать и поддерживать.
Венгерская нотация на основе типов не одобряется, поскольку современные редакторы автоматически отображают информацию о типах, что добавляет ненужные сложности именам идентификаторов.
использование size_t для индексов, не int
STL ваш друг, используйте его!
Превентивно предотвращать ошибки с помощью const для предотвращения случайного изменения данных только для чтения.

Вот так:

#include <vector>

using namespace std;

vector<double> fastFourierTransform(const vector<double> signal, const bool inverse) {

if( signal.size() < 2 ) return signal;
const size_t half = signal.size() / 2;

vector<double> lower; lower.reserve( half );
vector<double> upper; upper.reserve( half );

bool isEven = true;
for( size_t i = 0; i < signal.size(); i++ ) {
if( isEven ) lower.push_back( signal.at( i ) );
else         upper.push_back( signal.at( i ) );

isEven = !isEven;
}

vector<double> lowerFft = fastFourierTransform( lower, inverse );
vector<double> upperFft = fastFourierTransform( upper, inverse );

vector<double> result;
result.reserve( signal.size() );

double arg = ( inverse ? 1 : -1 ) * ( DoublePI / signal.size() );

// Ideally these should be local `double` values passed directly into `Mul`.
unique_ptr<double[]> oBase = make_unique<double[]>( 2 );
oBase[0] = cos(arg);
oBase[1] = sin(arg);

unique_ptr<double[]> o = make_unique<double[]>( 2 );
o[0] = 0;
o[1] = 0;

for( size_t i = 0; i < half; i++ ) {

double* lower1 = lower.at( i );
double* upper1 = Mul( o, upper.at( i ) );

result.at( i        ) = Sum( lower1, upper1 );
result.at( i + half ) = Dif( lower1, upper1 );

o = Mul( o, oBase );
}

// My knowledge of move-semantics of STL containers is a bit rusty - so there's probably a better way to return the output 'result' vector.
return result;
}

1

Источник

Accepted Answer

Часть «бабочка» неправильно применяет коэффициенты:

for (int i = 0; i < sH; i++) {
double* sO1 = sOdd[i];
double* sE1 = Mul(o, sEven[i]);

spectr[i] = Sum(sO1, sE1);
spectr[i + sH] = Dif(sO1, sE1);

o = Mul(o, oBase);
}

Примечание:

Я сохранил ваши записи, но это сбивает с толку:

fOdd имеет индексы 0, 2, 4, 6, … так и должно быть fEven

fEven имеет индексы 1, 3, 5, 7, … так и должно быть fOdd

действительно sOdd должно быть sLower а также sEven должно быть sUpper так как они соответствуют 0:s/2 а также s/2:s-1 элементы спектра соответственно:

sLower = FFT(fEven, sH, inverse); // fEven is 0, 2, 4, ...
sUpper = FFT(fOdd, sH, inverse); // fOdd is 1, 3, 5, ...

Тогда бабочка становится:

for (int i = 0; i < sH; i++) {
double* sL1 = sLower[i];
double* sU1 = Mul(o, sUpper[i]);

spectr[i] = Sum(sL1, sU1);
spectr[i + sH] = Dif(sL1, sU1);

o = Mul(o, oBase);
}

Когда написано так, легче сравнивать с этим пример псевдокода в википедии.

И @Dai это правильно, вы собираетесь потерять много памяти

2