This is the mail archive of the
gsl-discuss@sources.redhat.com
mailing list for the GSL project.
Having fft-problems
- From: Anders Wallander <baw at biosci dot ki dot se>
- To: gsl-discuss at sources dot redhat dot com
- Date: Sat, 5 Jul 2003 17:01:58 +0200
- Subject: Having fft-problems
- Organization: KI
Hi!
In the manuals it's written:
"...a call to gsl_fft_complex_forward followed by a call to gsl_fft_complex_inverse should return the original data (within numerical errors)."
That is only true for the real-values in my case. Someone have a clue what I'm doing wrong here?
I've gone down to this level now in pure frustration. At the moment I'm trying to write proper 2D_fft-routines taking gsl_matrixes as arguments. This should be as simple as doing fft_forward on all rows and then all columns. Same thing should work for fft_inverse on the matrixes, fft_inverse on all rows and columns. The reason for doing this is that I'm trying to implement a cross-correlation routine taking 2 matrixes as arguments.
Regards,
Anders
-----------
#include <stdio.h>
#include <stdlib.h>
#include <gsl/gsl_math.h>
#include <gsl/gsl_fft_complex.h>
#define FOURIER_WORKSPACE 32
int
main()
{
size_t i;
double data[FOURIER_WORKSPACE];
gsl_fft_complex_wavetable *wave_table;
gsl_fft_complex_workspace *work_space;
wave_table = gsl_fft_complex_wavetable_alloc(FOURIER_WORKSPACE);
work_space = gsl_fft_complex_workspace_alloc(FOURIER_WORKSPACE);
for (i = 0; i < (FOURIER_WORKSPACE/2); i++){
data[2*i] = sin(i/(2*M_PI));
data[2*i+1] = 0.0;
}
printf("\n Data:\n ");
for (i = 0; i < (FOURIER_WORKSPACE/2); i++)
printf("%f+i%f\n", data[2*i], data[2*i+1]);
gsl_fft_complex_forward (data, 1, FOURIER_WORKSPACE, wave_table, work_space);
printf("\n\n transformerad data: \n ");
for (i = 0; i < (FOURIER_WORKSPACE/2) ; i++)
printf("%f+i%f\n", data[2*i], data[2*i+1]);
gsl_fft_complex_inverse (data, 1, FOURIER_WORKSPACE, wave_table, work_space);
printf("\n\n inverse trasnformerad data: \n ");
for (i = 0; i < (FOURIER_WORKSPACE/2) ; i++)
printf("%f+i%f\n", data[2*i], data[2*i+1]);
gsl_fft_complex_workspace_free(work_space);
gsl_fft_complex_wavetable_free(wave_table);
return 0;
}