Re: Having fft-problems

[Anders Wallander <baw at biosci dot ki dot se>]

Sorry!
Somehow one always tend to find the errors when telling others how stupid one are ;-)
Adjusted some more parameters to FOURIER_WORKSPACE/2 and suddenly it returned the correct values, or a least close enough.

Cheers!
/Anders

On Sat, 5 Jul 2003 17:01:58 +0200
Anders Wallander <baw@biosci.ki.se> wrote:

> 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;
> }


--