This is the mail archive of the gsl-discuss@sources.redhat.com mailing list for the GSL project.


Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]
Other format: [Raw text]

Having fft-problems


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


Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]