How can I port code that uses numpy.fft.rfft from python to C++?

I have code written in python. It computes positive part of FFT of real input using numpy. I need to port this code to C++.

import numpy as np
interp=[131.107, 133.089, 132.199, 129.905, 132.977]
res=np.fft.rfft(interp)
print res

Result of rfft is [ 659.27700000+0.j, 1.27932533-1.4548977j, -3.15032533+2.1158917j]

I tried to use OpenCV for 1D DFT:

std::vector<double> fft;
std::vector<double> interpolated = {131.107, 133.089, 132.199, 129.905, 132.977};
cv::dft( interpolated, fft );
for( auto it = fft.begin(); it != fft.end(); ++it ) {
    std::cout << *it << ' ';
}
std::cout << std::endl;

Result of cv::dft is {1.42109e-14, -127.718, -94.705, 6.26856, 23.0231}. It is much different from numpy.fft.rfft. It looks strange that DC value (zero element) is near zero on all inputs after OpenCV's dft computed.

Usage of FFTW3 library gave me the same results as OpenCV's results:

std::vector<double> interpolated = {131.107, 133.089, 132.199, 129.905, 132.977};
fftw_complex* out = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * 3 );
fftw_plan plan = fftw_plan_dft_r2c_1d( interpolated.size( ), interpolated.data( ), out, FFTW_ESTIMATE );
fftw_execute(plan);
fftw_destroy_plan(plan);
for( size_t i = 0; i < interpolated.size( ); ++i ) {
    std::cout << " (" << out[ i ][ 0 ] << ", " << out[ i ][ 1 ] << ")";
}
fftw_free(out);

This code gives me the same results as OpenCV. It prints: (1.42109e-14, 0) (-127.718, -94.705) (6.26856, 23.0231).

Why do I get different results of dft in C++ and in python? What am I doing wrong?

Thanks!

Solution

I'm using gcc 4.6 at the moment, which doesn't have C++11, so I tried this version of your code, using OpenCV 2.4.8:

#include <iostream>
#include "opencv2/core/core.hpp"

int main(int argc, char *argv[])
{
    double data[] = {131.107, 133.089, 132.199, 129.905, 132.977};
    std::vector<double> interpolated (data, data + sizeof(data) / sizeof(double));
    std::vector<double> fft;

    cv::dft(interpolated, fft);

    for (std::vector<double>::const_iterator it = fft.begin(); it != fft.end(); ++it) {
        std::cout << *it << ' ';
    }
    std::cout << std::endl;
}

The output

659.277 1.27933 -1.4549 -3.15033 2.11589

agrees with numpy and with the cv2 python module:

In [55]: np.set_printoptions(precision=3)

In [56]: x
Out[56]: array([ 131.107,  133.089,  132.199,  129.905,  132.977])

In [57]: np.fft.rfft(x)
Out[57]: array([ 659.277+0.j   ,    1.279-1.455j,   -3.150+2.116j])

In [58]: cv2.dft(x)
Out[58]: 
array([[ 659.277],
       [   1.279],
       [  -1.455],
       [  -3.15 ],
       [   2.116]])

I don't know why your code is not working, so I guess this is more of a long comment than an answer.