快速傅立叶变换程序(更新版)

1、1加图标ID_FFT，类向导在视类添加相应函数2BMP文件中增加：BITMAPINFO* lpDIB = NULL; /原时域图像BITMAPINFO* lpDIB_FFT = NULL; /频谱图像static DWORD LineBytes = 0; /图像每行像素所占字节数static LONG size = 0; /图像实际数据所占字节数注意：在LoadBmpFile函数中/删除LONG size= sizeof(BITMAPINFOHEADER) + NumColors * sizeof(RGBQUAD) + ImgSize;中的LONG/删除 DWORD LineBytes = (

2、bi.biWidth * bi.biBitCount + 31)/32 * 4;中的DWORD/在“if (NULL = (lpBitsInfo = (LPBITMAPINFO) malloc(size)return FALSE;”后面给lpDIB分配内存“if (NULL = (lpDIB= (LPBITMAPINFO) malloc(size)return FALSE; “/在“fread(char*)lpBitsInfo, size, 1, fp);”后面添加“lpDIB = lpBitsInfo”在视类OnDraw前中增加：extern BITMAPINFO *lpDIB_FFT;ex

3、tern BITMAPINFO *lpDIB;OnDraw中增加：if (lpDIB_FFT)lpBits = (LPVOID)&lpDIB_FFT-bmiColorslpDIB_FFT-bmiHeader.biClrUsed;StretchDIBits( pDC-GetSafeHdc(),600,0,lpDIB_FFT-bmiHeader.biWidth,lpDIB_FFT-bmiHeader.biHeight,0,0,lpDIB_FFT-bmiHeader.biWidth,lpDIB_FFT-bmiHeader.biHeight,lpBits,lpDIB_FFT, / bitmap dat

4、a DIB_RGB_COLORS,SRCCOPY);3在视类添加void Fourier();void CFFTView:OnFft() if (lpDIB_FFT)free(lpDIB_FFT);Fourier();Invalidate();void CFFTView:OnUpdateFft(CCmdUI* pCmdUI) pCmdUI-Enable(lpDIB != NULL);4BMP文件前面添加#include using namespace std;（PPT中Copy）#define PI 3.4BMP文件中添加函数void Fourier()5BMP文件中添加函数一维傅里叶变换函数

5、FFT()(红色加粗是需要自己写的，其他的为PPT中Copy)void FFT(complex * TD, complex * FD, int r)/ 计算付立叶变换点数LONG count = 1 r;/ 计算加权系数int i;double angle;complex* W = new complexcount / 2;for(i = 0; i count / 2; i+)angle = -i * PI * 2 / count;Wi = complex (cos(angle), sin(angle);/ 将时域点写入X1complex* X1 = new complexcount;memc

6、py(X1, TD, sizeof(complex) * count);/ 采用蝶形算法进行快速付立叶变换，输出为频域值X2complex* X2 = new complexcount; int k,j,p,size;complex* temp;for (k = 0; k r; k+)for (j = 0; j 1 k; j+)size = 1 (r-k);for (i = 0; i size/2; i+)p = j * size;X2i + p = X1i + p + X1i + p + size/2;X2i + p + size/2 = (X1i + p - X1i + p + size/

7、2) * Wi * (1k);temp = X1;X1 = X2;X2 = temp;/ 重新排序（码位倒序排列）for (j = 0; j count; j+)p = 0;for (i = 0; i r; i+)if (j & (1i)p += 1bmiHeader.biWidth;int height = lpDIB-bmiHeader.biHeight;/指向图像数据指针BYTE* lpBits = (BYTE*)&lpDIB-bmiColorslpDIB-bmiHeader.biClrUsed;/ FFT宽度（必须为2的整数次方）int FFT_w = 1;/ FFT宽度的幂数，即迭代

8、次数int wp = 0;while(FFT_w * 2 = width)FFT_w *= 2;wp +;/ FFT高度（必须为2的整数次方）int FFT_h = 1;/ FFT高度的幂数，即迭代次数int hp = 0;while(FFT_h * 2 = height)FFT_h *= 2;hp +;/ 分配内存complex* TD = new complexFFT_w * FFT_h;complex* FD = new complexFFT_w * FFT_h;int i, j;BYTE* pixel;for(i = 0; i FFT_h; i+) / 行for(j = 0; j FF

9、T_w; j+) / 列/ 指向DIB第i行，第j个象素的指针pixel = lpBits + LineBytes * (height - 1 - i) + j;/ 给时域赋值TDj + FFT_w * i = complex(*pixel, 0);for(i = 0; i FFT_h; i+)/ 对y方向进行快速付立叶变换FFT(&TDFFT_w * i, &FDFFT_w * i, wp);/ 保存中间变换结果for(i = 0; i FFT_h; i+)for(j = 0; j FFT_w; j+)TDi + FFT_h * j = FDj + FFT_w * i;for(i = 0; i bmiColorslpDIB_FFT-bmiHeader.biClrUsed;double temp;for(i = 0; i FFT_h; i+