+ FFT1::FFT1()\r
+ {\r
+ construct_default();\r
+ }\r
+ FFT1::FFT1(int wid)\r
+ {\r
+ construct_default();\r
+ set(wid);\r
+ }\r
+ FFT1::FFT1(const Matrix &source)\r
+ {\r
+ construct_default();\r
+ set(source);\r
+ }\r
+ FFT1::~FFT1()\r
+ {\r
+ release();\r
+ }\r
+ void FFT1::construct_default()\r
+ {\r
+ width_=0;\r
+ img_spc=0; frq_spc=0;\r
+ }\r
+ void FFT1::release()\r
+ {\r
+ if(img_spc!=0) { free(img_spc); img_spc = 0; }\r
+ if(frq_spc!=0) { free(frq_spc); frq_spc = 0; }\r
+ width_=0;\r
+ }\r
+ void FFT1::set(int wid)\r
+ {\r
+ if(img_spc==0 || frq_spc==0 || (width_ != wid) ) {\r
+ release();\r
+ img_spc = (fftw_complex*)malloc(sizeof(fftw_complex) * wid);\r
+ frq_spc = (fftw_complex*)malloc(sizeof(fftw_complex) * wid);\r
+ }\r
+ width_ = wid;\r
+ x_zero = wid + 1;\r
+ left = (wid)/2;\r
+ }\r
+ void FFT1::set(const Matrix &source)\r
+ {\r
+ if(source.getRows()==1) {\r
+ set(source.getCols());\r
+ for(int x=0; x<width_; x++) {\r
+ img_spc[x][REAL] = source(1, x+1);\r
+ img_spc[x][IMAGINARY] = 0.0;\r
+ }\r
+ }\r
+ }\r
+\r
+ double FFT1::getDC()\r
+ {\r
+ return sqrt(frq_spc[0][REAL]*frq_spc[0][REAL]+frq_spc[0][IMAGINARY]*frq_spc[0][IMAGINARY]);\r
+ }\r
+ double FFT1::setDC(double l)\r
+ {\r
+ frq_spc[0][REAL] = l;\r
+ frq_spc[0][IMAGINARY] = 0;\r
+ return l;\r
+ }\r
+\r
+\r
+ void FFT1::fft()\r
+ {\r
+ fftw_plan plan = fftw_plan_dft_1d(width_, img_spc, frq_spc, FFTW_FORWARD, FFTW_ESTIMATE);\r
+ fftw_execute(plan);\r
+ normalizeFFT();\r
+ }\r
+ void FFT1::ifft()\r
+ {\r
+ fftw_plan plan = fftw_plan_dft_1d(width_, frq_spc, img_spc, FFTW_BACKWARD, FFTW_ESTIMATE);\r
+ fftw_execute(plan);\r
+ }\r
+ void FFT1::normalizeFFT()\r
+ {\r
+ double num = width_;\r
+ for (int i=0; i<num; i++) {\r
+ frq_spc[i][REAL] /= num;\r
+ frq_spc[i][IMAGINARY] /= num;\r
+ }\r
+ }\r
+ void FFT1::getImage(Matrix &target)\r
+ {\r
+ double rl, im;\r
+ int h, v;\r
+ if( (target.getCols()!=width_)) {\r
+ target.release();\r
+ target.set(1, width_);\r
+ }\r
+ for (int x=0; x<width_; x++) {\r
+ rl = img_spc[x][REAL];\r
+ im = img_spc[x][IMAGINARY];\r
+ target(1, x+1) = sqrt(rl*rl+im*im);\r
+ }\r
+ }\r
+ void FFT1::getSpectrum(Matrix &target, double gamma)\r
+ {\r
+ const int right = width_%2==0 ? left : left+1;\r
+ double rl, im;\r
+ Color col;\r
+ if( (target.getCols()!=width_) ) {\r
+ target.release();\r
+ target.set(1, width_);\r
+ }\r
+ int h, v = 1;\r
+ for (int x=0; x<right; x++) {\r
+ h = x+left + 1;\r
+ rl = frq_spc[x][REAL];\r
+ im = frq_spc[x][IMAGINARY];\r
+ target(v, h) = pow( sqrt(rl*rl+im*im), gamma );\r
+ }\r
+ for (int x=right; x<width_; x++) {\r
+ h = x-right + 1;\r
+ rl = frq_spc[x][REAL];\r
+ im = frq_spc[x][IMAGINARY];\r
+ target(v, h) = pow( sqrt(rl*rl+im*im), gamma );\r
+ }\r
+ }\r
+\r
+\r