FFT C++ issue

[pajczur]

Hello,
I'm tryung to implement Fast Fourier Transformation. And I'm not sure if I understand it correctly.

As I understand when I transform (for example) simple sine wave by FFT, and than make inversion I should get exactly the same sine wave. Am I right?

But in some way it doesn't work for me. I found that code with the friend google:

#include <complex>
#include <iostream>
#include <valarray>
#include <QVector>
using namespace std;
 
const double PI = 3.141592653589793238460;
 
typedef std::complex<double> Complex;
typedef std::valarray<Complex> CArray;
 
// Cooley–Tukey FFT (in-place, divide-and-conquer)
// Higher memory requirements and redundancy although more intuitive
void fft(CArray& x)
{
    const size_t N = x.size();
    if (N <= 1) return;
 
    // divide
    CArray even = x[std::slice(0, N/2, 2)];
    CArray  odd = x[std::slice(1, N/2, 2)];
 
    // conquer
    fft(even);
    fft(odd);
 
    // combine
    for (size_t k = 0; k < N/2; ++k)
    {
        Complex t = std::polar(1.0, -2 * PI * k / N) * odd[k];
        x[k    ] = even[k] + t;
        x[k+N/2] = even[k] - t;
    }
}
 
// inverse fft (in-place)
void ifft(CArray& x)
{
    // conjugate the complex numbers
    x = x.apply(std::conj);
 
    // forward fft
    fft( x );
 
    // conjugate the complex numbers again
    x = x.apply(std::conj);
 
    // scale the numbers
    x /= x.size();
}
 
int main()
{
    const Complex test[] = { 1.0, 1.0, 1.0, 1.0, 0.0, 0.0, 0.0, 0.0 };
    CArray data(test, 8);
 
    // forward fft
    fft(data);
 
    std::cout << "fft" << std::endl;
    for (int i = 0; i < 10; ++i)
    {
        std::cout << data[i] << std::endl;
    }
 
    // inverse fft
    ifft(data);
 
    std::cout << std::endl << "ifft" << std::endl;
    for (int i = 0; i < 10; ++i)
    {
        std::cout << data[i] << std::endl;
    }
    return 0;
}

#include <complex>
#include <iostream>
#include <valarray>
#include <QVector>
using namespace std;

const double PI = 3.141592653589793238460;

typedef std::complex<double> Complex;
typedef std::valarray<Complex> CArray;

// Cooley–Tukey FFT (in-place, divide-and-conquer)
// Higher memory requirements and redundancy although more intuitive
void fft(CArray& x)
{
    const size_t N = x.size();
    if (N <= 1) return;

    // divide
    CArray even = x[std::slice(0, N/2, 2)];
    CArray  odd = x[std::slice(1, N/2, 2)];

    // conquer
    fft(even);
    fft(odd);

    // combine
    for (size_t k = 0; k < N/2; ++k)
    {
        Complex t = std::polar(1.0, -2 * PI * k / N) * odd[k];
        x[k    ] = even[k] + t;
        x[k+N/2] = even[k] - t;
    }
}

// inverse fft (in-place)
void ifft(CArray& x)
{
    // conjugate the complex numbers
    x = x.apply(std::conj);

    // forward fft
    fft( x );

    // conjugate the complex numbers again
    x = x.apply(std::conj);

    // scale the numbers
    x /= x.size();
}

int main()
{
    const Complex test[] = { 1.0, 1.0, 1.0, 1.0, 0.0, 0.0, 0.0, 0.0 };
    CArray data(test, 8);

    // forward fft
    fft(data);

    std::cout << "fft" << std::endl;
    for (int i = 0; i < 10; ++i)
    {
        std::cout << data[i] << std::endl;
    }

    // inverse fft
    ifft(data);

    std::cout << std::endl << "ifft" << std::endl;
    for (int i = 0; i < 10; ++i)
    {
        std::cout << data[i] << std::endl;
    }
    return 0;
}

To copy to clipboard, switch view to plain text mode 

So as I understand by ifft(data) I should receive the same as it was on the beggining, but the results are quite different. Can anyone help me?
For any help thanks in advance.
Best Regards

[ars]

Hello,

running your program I get the following output (reducing the number of lines to the size of the input data, i.e. 8 instead of 10):

fft
(4,0)
(1,-2.41421)
(0,0)
(1,-0.414214)
(0,0)
(1,0.414214)
(0,0)
(1,2.41421)
 
ifft
(1,-0)
(1,-5.55112e-17)
(1,2.4895e-17)
(1,-5.55112e-17)
(5.55112e-17,0)
(5.55112e-17,5.55112e-17)
(0,-2.4895e-17)
(5.55112e-17,5.55112e-17)

fft
(4,0)
(1,-2.41421)
(0,0)
(1,-0.414214)
(0,0)
(1,0.414214)
(0,0)
(1,2.41421)

ifft
(1,-0)
(1,-5.55112e-17)
(1,2.4895e-17)
(1,-5.55112e-17)
(5.55112e-17,0)
(5.55112e-17,5.55112e-17)
(0,-2.4895e-17)
(5.55112e-17,5.55112e-17)

To copy to clipboard, switch view to plain text mode 

Your input was a vector of 4 ones and 4 zeros (real numbers). And that's (up to rounding errors) "exactly" what you get as output of the ifft call: You get a sequence of complex numbers and you should read the x.yze-17 numbers as 0.

Compiled with gcc6.

Best regards
ars

[pajczur]

Hello,
thanks for your reply.
Yes you are right. When I try it with the sine wave, with much more sampling points, the iFFT gives me really close results to the input. I think I understand it now.