28 #include <Base/Debug/Error.hh>
52 if (_Size!=0 || _InputSize!=0 || _OutNum!=0 || _out || _in || _p){
53 BIASERR (
"FFT: call Release() before calling Init() a second time");
54 BIASEXIT(EXIT_FAILURE,
"call Release() before calling Init() a second time");
58 _InputSize=_Size*
sizeof(double);
60 _out =(fftw_complex*)fftw_malloc(
sizeof(fftw_complex) * _OutNum);
61 _in = (
double*)fftw_malloc(
sizeof(
double) * _Size);
63 _p = fftw_plan_dft_r2c_1d(_Size, _in, _out, FFTW_MEASURE);
68 if (_p) fftw_destroy_plan(_p);
69 if (_in) fftw_free(_in);
70 if (_out) fftw_free(_out);
74 _Size=_InputSize=_OutNum=0;
79 memcpy((
void *)_in, (
void *)in, _InputSize);
86 if ((
int)in.
Size()!=_Size){
87 BIASERR (
"called FFT::Compute() with vector of wrong size "<<in.
Size());
88 BIASEXIT(EXIT_FAILURE,
"called FFT::Compute() with vector of wrong size "<<in.
Size());
96 complex<double> *p=
reinterpret_cast<complex<double> *
>(_out), *e=p+_OutNum;
107 complex<double> *p=
reinterpret_cast<complex<double> *
>(_out), *e=p+_OutNum;
118 if ((
int)phase.
Size()!=_OutNum){
121 complex<double> *p=
reinterpret_cast<complex<double> *
>(_out);
122 for (
int i=0; i<_OutNum; i++){
129 if ((
int)mag.
Size()!=_OutNum){
132 complex<double> *p=
reinterpret_cast<complex<double> *
>(_out);
133 for (
int i=0; i<_OutNum; i++){
143 os <<
"#result of fft: "<<endl;
144 os <<setw(4)<<
"#num"<<setw(15)<<
"freq"<<setw(15)<<
"mag"<<setw(15)
146 for (
int i=0; i<_OutNum; i++){
148 <<setw(15)<<(double)i/(
double)_Size
150 <<setw(15)<<phase[i]<<
"\n";
void GetMagnitude(double *mag)
returns the resulting magnitude, mag must be of length size/2+1 (rounded downwards) ...
void Release()
cal this before a second call to Init()
unsigned int Size() const
length of the vector
Vector< T > & newsize(Subscript N)
void Init(int size)
initializes the internal variables.
T * GetData() const
get the pointer to the data array of the vector (for faster direct memory access) ...
void DumpResult(std::ostream &os=std::cout)
void Compute()
input must be filled before with GetInput
void GetPhase(double *phase)
returns the resulting phase, phase must be of length size/2+1 (rounded downwards) ...