fft.rar_C Builder_fft

#include "stdafx.h" #include "fft.h" #include "Complex.h" #include<math.h> #define M_PI 3.141592654 // The mathematical constant for PI void FFT1D(Complex<double> *fx, Complex<double> *Fu, int twoK, int stride) { // This function performs the 1D fast fourier transform from input data fx, // it outputs the frequency coefficients in the array Fu. // Input: fx - pointer to an array of data in complex number format // twoK - the total number of data in fx // stride - the number of data to be skipped when reading the next data in fx // e.g. fx = [0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15], twoK = 4 // if stride = 1, the data to be processed is [0 1 2 3] // if stride = 2, the data to be processed is [0 2 4 8] // if stride = 4, the data to be processed is [0 4 8 12] // Output: Fu - pointer to an array buffer for transformed cofficients in complex number format // Code to be added here ... Complex<double> *F_even,*F_odd, e; int u,K; if(twoK==2) { Fu[0]=(fx[0]+fx[1*stride])*0.5; Fu[1]=(fx[0]-fx[1*stride])*0.5; return; } K=twoK/2; F_even=new Complex<double>[K]; F_odd=new Complex<double>[K]; FFT1D(fx,F_even,K,2*stride); FFT1D(fx+stride, F_odd, K, 2*stride); for(u=0;u<K;u++) { e.Exp(-M_PI*((double)u/(double)K)); Fu[u]=(F_even[u]+e*F_odd[u])*0.5; Fu[u+K]=(F_even[u]-e*F_odd[u])*0.5; } delete []F_even; delete []F_odd; } void FFT2D(Complex<double> *f, Complex<double> *F, int width, int height) { // This function performs the 2D Fast Fourier Transform for a given image f // You can safely assume that the width = 2^m, height = 2^n for some integers m and n, // and the input image are already properly padded // Input: f - An 2D Image represented in Complex Numbers // Output: F - The transformed coefficients, also represented in Complex Numbers // Code to be added here ... int i,j; Complex<double> **s = new Complex<double>*[height]; for(int j=0; j<height; j++){ s[j] =new Complex<double>[width]; } for(j=0;j<height;j++) for(i=0;i<width;i++) s[j][i]=f[(j*width)+i]; Complex<double> *row; row=new Complex<double>[width]; for(j=0;j<height;j++) { for(i=0;i<width;i++) row[i]=s[j][i]; FFT1D(row,row,width,1); for(i=0;i<width;i++) s[j][i]=row[i]; } Complex<double> *col; col=new Complex<double>[height]; for(i=0;i<width;i++) { for(j=0;j<height;j++) col[j]=s[j][i]; FFT1D(col,col,height,1); for(j=0;j<height;j++) s[j][i]=col[j]; } for(j=0;j<height;j++) for(i=0;i<width;i++) F[(j*width)+i]=s[j][i]; delete []row; delete []col; delete s; } void IFFT2D(Complex<double> *F, Complex<double> *f, int width, int height) { // This function performs the Inverse 2D Fast Fourier Transform for a given image f (in one color channel R, G or B) // You can safely assume that the width = 2^m, height = 2^n for some integers m and n // Input: F - The transformed coefficients, also represented in Complex Numbers // Output: f - An 2D Image represented in Complex Numbers // Code to be added here ... int i; for(i=0;i<width*height;i++) { F[i]=Complex<double>(F[i].Real(),-F[i].Imag()); } FFT2D(F,f,width,height); for(i=0;i<width*height;i++) { f[i]=Complex<double>(f[i].Real(),-f[i].Imag()); f[i]=f[i]*width*height; } } void ScaleFFT(Complex<double> *F, Complex<double> *PS, int width, int height) { // This function takes a 2D coefficients in F, and compute the power spectrum of F. // This function also scales the power spectrum, such that its values lie within the range [0, 255] // for proper display. // Input: F - 2D coefficients of the FFT image // width - width of the F // height - height of the F // Output: PS - scaled power spectrum (with value normalized to the range of [0, 255]) // Code to be added here ... int i; double min,max,diff; for(i=0;i<width*height;i++) { PS[i] = Complex <double>(pow((F[i].Mod()*F[i].Mod()),0.2)); } min=max=PS[0].Real(); for(i=1;i<width*height;i++) { if(min > PS[i].Real()) min = PS[i].Real(); if(max < PS[i].Real()) max =PS[i].Real(); } diff = max-min; for(i=0;i<width*height;i++) { PS[i] = Complex <double>(((PS[i].Real()-min)/diff)*255.0); } } void BLPF(Complex<double> *F, Complex<double> *LF, int width, int height) { // This function takes a 2D coefficients in F, and perform a second other ButterWorth // low pass filtering on the coefficents. The filtered coefficients are put back into LF // Input: F - 2D coefficients of the FFT image // width - width of the F // height - height of the F // Output: LF - ButterWorth filtered coefficients // The coded below is a low pass filter implementation, however, it cannot filter the image // properly. Modify the code below such that the filter becomes a butterworth filter of order 2. /*int i, j ; int Df ; int centerX, centerY ; centerX = width / 2 ; centerY = height / 2 ; Df = 32 ; for(j = 0; j < height; j++) { for(i = 0; i < width; i++) { if((i - centerX) * (i - centerX) + (j - centerY) * (j - centerY) > Df * Df) { LF[i + j * width] = 0.0 ; } else { LF[i + j * width] = F[i + j * width] ; } } }*/ int i, j,D; double T,H; int Df ; int centerX, centerY ; centerX = width / 2 ; centerY = height / 2 ; Df = 32 ; for(j = 0; j < height; j++) { for(i = 0; i < width; i++) { D=(i - centerX) * (i - centerX) + (j - centerY) * (j - centerY); T=(D/(Df*Df))*(D/(Df*Df)); H=1/(1+T); LF[i+j*width]=F[i+j*width]*H; } } }