开发者

computing fft and ifft with fftw.h in C

开发者 https://www.devze.com 2023-03-01 12:32 出处:网络
Hi all I am using the fftw C libraries to compute the frequency spectrum for some signal processing applications on embedded systems. However, in my project I have run into a slight hinderence.

Hi all I am using the fftw C libraries to compute the frequency spectrum for some signal processing applications on embedded systems. However, in my project I have run into a slight hinderence.

Below is a simple program I wrote to ensure I am implementing the fftw functions correctly. Basically I want to calculate the fft of a sequence of 12 numbers, then do the ifft and obtain the same sequence of numbers again. If you have fftw3 and gcc installed this program should work if you compile with:

gcc -g -lfftw3 -lm fftw_test.c -o fftw_test

Currently my fft length is the same size as the input array.

#include <stdio.h>
#include <stdlib.h>
#include <sndfile.h>
#include <stdint.h>
#include <math.h>
#include <fftw3.h>

int main(void)
{
double array[] = {0.1, 0.6, 0.1, 0.4, 0.5, 0, 0.8, 0.7, 0.8, 0.6, 0.1,0};
//double array2[] = {1, 6, 1, 4, 5, 0, 8, 7, 8, 6, 1,0};
double *out;
double *err;
int i,size = 12;

fftw_complex *out_cpx;

fftw_plan fft;
fftw_plan ifft;
out_cpx = (fftw_complex*) fftw_malloc(sizeof(fftw_complex)*size);
out = (double *) malloc(size*sizeof(double));
err = (double *) malloc(size*sizeof(double));

fft = fftw_plan_dft_r2c_1d(size, array, out_cpx, FFTW_ESTIMATE);  //Setup fftw plan for fft
ifft = fftw_plan_dft_c2r_1d(size, out_cpx, out, FFTW_ESTIMATE);   //Setup fftw plan for ifft

fftw_execute(fft);
fftw_execute(ifft);

//printf("Input:    \tOutput:    \tError:\n");
printf("Input:    \tOutput:\n");
for(i=0;i<size;i++)
{
err[i] = abs(array[i] - out[i]);    
printf("%f\t%f\n",(array[i]),out[i]);
//printf("%f\t%f\t%f\n",(array[i]),out[i],err[i]);
}

fftw_destroy_plan(fft);
fftw_destroy_plan(ifft);
fftw_free(out_cpx);
free(err);
free(out);
return 0;
}

Which Produces the following output:

Input:      Output:
0.100000    1.200000
0.600000    7.200000
0.100000    1.200000
0.400000    4.800000
0.500000    6.000000
0.000000    0.000000
0.800000    9.600000
0.700000    8.400000
0.800000    9.600000
0.600000    7.200000
0.100000    1.200000
0.000000    0.000000

So obviously the ifft is producing some scaled up result. In the fftw docs found here: fftw docs about 开发者_StackOverflow中文版scaling. It mentions about some scaling, however I am using the "r2c" and "c2r" transforms rather than the FFT_FORWARD and FFT_BACKWARD. Any insight would be appreciated.


Looking at the great documentation for the functions you use, you will see you are using FFT_FORWARD and FFT_BACKWARD, and exactly where it is intended. Therefore, the scaling information you found previously also applies here.


Sorry to be pedantic, but your size for out_cpx is incorrect. instead of being size long, it should be size/2 + 1. This is because FFT of a real signal is Hermitian. You can verify what I say by initializing out_cpx to some random number (all 3.14159). Run both the forward and backward and then print out out_cpx from size/2 + 1 to size. It will not have changed.

http://www.fftw.org/fftw3_doc/Real_002ddata-DFT-Array-Format.html#Real_002ddata-DFT-Array-Format


r2c and c2r do essentially the same as the regular Fourier transform. The only difference is that both the input and output array need to hold half of the numbers. Please take a look at the last paragraph of the manual of FFTW r2c and c2r. So the normalization factor is precisely the number of elements of the real array, or the variable size (== 12) in your case.

0

精彩评论

暂无评论...
验证码 换一张
取 消

关注公众号