high/low-pass filter z użyciem FFTW

Witam.

Chcę zastosować filtr górno i dolnoprzepustowy podpierając się biblioteką FFTW.
Dane wczytuję z pliku tekstowego.
Doszedłem do miejsca kiedy stosuję FFT (szybka transformata fouriera) jednak nie wiem co dalej, proszę o pomoc.
Jak wygląda sprawa obcięcia? Czytałem że wykorzystywane jest także IFFT jednak nie wiem do czego.
Proszę o odrobinę ludzko brzmiącej teorii z zastosowaniem praktycznym w tym przypadku ;)

 
#include <iostream>
#include <fstream>
#include <sstream>
#include <string>
#include <vector>
#include <stdio.h>
#include <stdlib.h>
#include "fftw3.h"

using namespace std;

int main( int argc, char** argv )
{
	bool cut_up = true;
	bool cut_down = false;
	int cut_freq = 20; 


	  vector<string> text_file;

	  ifstream ifs( "data.txt" );
	  string temp;
	  int nol = 0;
	  while( getline( ifs, temp ) )
	  {
		  nol++;
	  }

	  vector <vector <string> > file_data;
	  ifstream infile( "data.txt" );

	  while (infile)
	  {
	    string s;
	    if (!getline( infile, s )) break;

	    istringstream ss( s );
	    vector <string> record;

	    while (ss)
	    {
	      string s;
	      if (!getline( ss, s, ',' )) break;
	      record.push_back( s );
	    }

	    file_data.push_back( record );
	  }
	  if (!infile.eof())
	  {
	    cerr << "Fooey!\n";
	  }


	  const int N = nol;
	  //cout << N <<" " << nol <<endl;

		      fftw_complex    *data, *fft_result, *ifft_result;
		      fftw_plan       plan_forward, plan_backward;
		      int             i;

		      data        = ( fftw_complex* ) fftw_malloc( sizeof( fftw_complex ) * SIZE );
		      fft_result  = ( fftw_complex* ) fftw_malloc( sizeof( fftw_complex ) * SIZE );
		      ifft_result = ( fftw_complex* ) fftw_malloc( sizeof( fftw_complex ) * SIZE );


		      /* populate input data */
		      double val, val1;
		      for( i = 0 ; i < N ; i++ ) {

		    	  val = atof(file_data[i][0].c_str());
		          data[i][0] = val;

		          val1 = atof(file_data[i][1].c_str());
		          data[i][1] = val1;
		      }


		      plan_forward  = fftw_plan_dft_1d( N, data, fft_result, FFTW_FORWARD, FFTW_ESTIMATE );
		      plan_backward = fftw_plan_dft_1d( N, fft_result, ifft_result, FFTW_BACKWARD, FFTW_ESTIMATE );



		      /* print initial data */
		      for( i = 0 ; i < N ; i++ ) {
		          fprintf( stdout, "data[%d] = { %2.2f, %2.2f }\n",
		                      i, data[i][0], data[i][1] );
		      }

		      fftw_execute( plan_forward );


		      /* print fft result */
		      for( i = 0 ; i < N ; i++ ) {
		          fprintf( stdout, "fft_result[%d] = { %2.2f, %2.2f }\n",
		                      i, fft_result[i][0], fft_result[i][1] );

		      }
		      
		      
		      if (cut_up) {
		    	  for( i = 0 ; i < N ; i++ ){
		    		  if  (fft_result[i][0] < 20)
		    		  {
		    			  fft_result[i][0] = 0 ;
		    		  }
		    	  }
		      }
		      
		      if (cut_down) {
		    	  for( i = 0 ; i < N ; i++ ){
		    		  if  (fft_result[i][0] > 20)
		    		  {
		    			  fft_result[i][0] = 0 ;
		    		  }
		    	  }
		      }

		      /*
		      // write to file
		      ofstream results;
		      results.open("results.txt");
		      results << fft_result[i][0]<<","<<fft_result[i][1]<<"\n";
		      results.close();
		      */



		      fftw_destroy_plan( plan_forward );
		      fftw_destroy_plan( plan_backward );

		      fftw_free( data );
		      fftw_free( fft_result );
		      fftw_free( ifft_result );



return 0;
}

Przykładowe dane:
| 1 | 0,2 | |
| 1 | 20 | |
0 | 0 | 0 | SUMA
0,02 | 0,019998667 | 0,077883668 | 0,097882335
0,04 | 0,039989334 | 0,143471218 | 0,183460552
0,06 | 0,059964006 | 0,186407817 | 0,246371824
0,08 | 0,079914694 | 0,199914721 | 0,279829415

Nie wiem czy dobrze to rozumiem jednak w pliku z danymi (data.txt) mam je zapisane jako
0.02,0.097882335
0.04,0.183460552
...

Wynikiem FFT jest sygnał w dziedzinie częstotliwości, a nie czasu. Żeby odfiltrować niskie i wysokie częstotliwości, wyzeruj odpowiednie wartości w sygnale przedstawionym w dziedzinie częstotliwości (po FFT), a następnie na tym zastosuj IFFT (Inverse Fast Fourier Transform), aby ponownie otrzymać sygnał w dziedzinie czasu, z odfiltrowanymi odpowiednimi częstotliwościami.

Czy mógłbym prosić o odniesienie się bardziej do mojego kodu- co jest w nim nie tak?
Z tego co rozumiem poprawnie dokonałem transformaty i zerowania jednak należy to teraz jeszcze przepuścić przez IFFT?

fftw_execute( plan_forward );

for(int i=cos; i<n; ++i) {
    fft_result[i][0] = 0;
    fft_result[i][1] = 0;
}

plan_backward = fftw_plan_dft_1d( N, fft_result, ifft_result, FFTW_BACKWARD, FFTW_ESTIMATE );
fftw_execute( plan_backward );

Żeby odfiltrować niskie i wysokie częstotliwości, wyzeruj odpowiednie wartości w sygnale przedstawionym w dziedzinie częstotliwości (po FFT), a następnie na tym zastosuj IFFT (Inverse Fast Fourier Transform), aby ponownie otrzymać sygnał w dziedzinie czasu, z odfiltrowanymi odpowiednimi częstotliwościami.
To jest rozwiązanie mało wydajne (mnóstwo obliczeń) i o marnej dokładności.

Do takich zastosowań jak filtrowanie częstotliwości używa się nie FFT tylko FIR.
Po wyznaczeniu współczynników samo przetwarzanie jest obliczeniowo banalne.

MarekR22 napisał(a):

fftw_execute( plan_forward );

for(int i=cos; i<n; ++i) {
    fft_result[i][0] = 0;
    fft_result[i][1] = 0;
}

plan_backward = fftw_plan_dft_1d( N, fft_result, ifft_result, FFTW_BACKWARD, FFTW_ESTIMATE );
fftw_execute( plan_backward );

czym są w tym przypadku zmienne "cos" oraz "n"?