实验三图像的傅立叶变换

1、实验三、图像的傅立叶变换一、实验目的1了解图像变换的意义和手段；2熟悉傅里叶变换的基本性质；3熟练掌握FFT的方法及应用；4通过实验了解二维频谱的分布特点；5通过本实验掌握编程实现数字图像的傅立叶变换。二、实验原理1应用傅立叶变换进行图像处理傅里叶变换是线性系统分析的一个有力工具，它能够定量地分析诸如数字化系统、采样点、电子放大器、卷积滤波器、噪音和显示点等的作用。通过实验培养这项技能，将有助于解决大多数图像处理问题。对任何想在工作中有效应用数字图像处理技术的人来说，把时间用在学习和掌握博里叶变换上是很有必要的。2傅立叶（Fourier）变换的定义对于二维信号，二维Fourier变换定义为：二

2、维离散傅立叶变换为：图像的傅立叶变换与一维信号的傅立叶变换变换一样，有快速算法，具体参见参考书目，有关傅立叶变换的快速算法的程序不难找到。实际上，现在有实现傅立叶变换的芯片，可以实时实现傅立叶变换。三、实验内容1.根据二维离散Fourier变换的定义编写程序2.实现图象的变换3.画出图象的频谱图。2提交实验的原始图像和结果图像。3报告内容主要是将原理与实验现象相结合，利用原理简要解释观察到的现象。六、思考题1傅里叶变换有哪些重要的性质?2图像的二维频谱在显示和处理时应注意什么?1. 源程序：/ MyDIPView.cpp : implementation of the CMyDIPView c

3、lass/#include "stdafx.h"#include "MyDIP.h"#include "MyDIPDoc.h"#include "MyDIPView.h"#include "math.h"#include<complex>using namespace std;#define PI 3.14159#ifdef _DEBUG#define new DEBUG_NEW#undef THIS_FILEstatic char THIS_FILE = _FILE_;#endif

4、/* * * 函数名称： * FFT() * * 参数: * complex<double> * TD- 指向时域数组的指针 * complex<double> * FD- 指向频域数组的指针 * r2的幂数，即迭代次数 * * 返回值: * 无。 * * 说明: * 该函数用来实现快速付立叶变换。 * */VOID WINAPI FFT(complex<double> * TD, complex<double> * FD, int r)/ 付立叶变换点数LONGcount;/ 循环变量inti,j,k;/ 中间变量intbfsize,p;/ 角

5、度doubleangle;complex<double> *W,*X1,*X2,*X;/ 计算付立叶变换点数count = 1 << r;/ 分配运算所需存储器W = new complex<double>count / 2;X1 = new complex<double>count;X2 = new complex<double>count;/ 计算加权系数for(i = 0; i < count / 2; i+)angle = -i * PI * 2 / count;Wi = complex<double> (c

6、os(angle), sin(angle);/ 将时域点写入X1memcpy(X1, TD, sizeof(complex<double>) * count);/ 采用蝶形算法进行快速付立叶变换for(k = 0; k < r; k+)for(j = 0; j < 1 << k; j+)bfsize = 1 << (r-k);for(i = 0; i < bfsize / 2; i+)p = j * bfsize;X2i + p = X1i + p + X1i + p + bfsize / 2;X2i + p + bfsize / 2 =

7、(X1i + p - X1i + p + bfsize / 2) * Wi * (1<<k);X = X1;X1 = X2;X2 = X;/ 重新排序for(j = 0; j < count; j+)p = 0;for(i = 0; i < r; i+)if (j&(1<<i)p+=1<<(r-i-1);FDj=X1p;/ 释放内存delete W;delete X1;delete X2;/* * * 函数名称： * Fourier() * * 参数: * LPSTR lpDIBBits - 指向源DIB图像指针 * LONG lWidt

8、h - 源图像宽度（象素数） * LONG lHeight - 源图像高度（象素数） * * 返回值: * BOOL - 成功返回TRUE，否则返回FALSE。 * * 说明: * 该函数用来对图像进行付立叶变换。 * BOOL WINAPI Fourier(LPSTR lpDIBBits, LONG lWidth, LONG lHeight)/ 指向源图像的指针unsigned char*lpSrc;/ 中间变量doubledTemp;/ 循环变量LONGi;LONGj;/ 进行付立叶变换的宽度和高度（2的整数次方）LONGw;LONGh;intwp;inthp;/ 图像每行的字节数LONG

9、lLineBytes;/ 计算图像每行的字节数lLineBytes = WIDTHBYTES(lWidth * 8);/ 赋初值w = 1;h = 1;wp = 0;hp = 0;/ 计算进行付立叶变换的宽度和高度（2的整数次方）while(w * 2 <= lWidth)w *= 2;wp+;while(h * 2 <= lHeight)h *= 2;hp+;/ 分配内存complex<double> *TD = new complex<double>w * h;complex<double> *FD = new complex<doub

10、le>w * h;/ 行for(i = 0; i < h; i+)/ 列for(j = 0; j < w; j+)/ 指向DIB第i行，第j个象素的指针lpSrc = (unsigned char*)lpDIBBits + lLineBytes * (lHeight - 1 - i) + j;/ 给时域赋值TDj + w * i = complex<double>(*(lpSrc), 0);for(i = 0; i < h; i+)/ 对y方向进行快速付立叶变换FFT(&TDw * i, &FDw * i, wp);/ 保存变换结果for(i

11、 = 0; i < h; i+)for(j = 0; j < w; j+)TDi + h * j = FDj + w * i;for(i = 0; i < w; i+)/ 对x方向进行快速付立叶变换FFT(&TDi * h, &FDi * h, hp);/ 行for(i = 0; i < h; i+)/ 列for(j = 0; j < w; j+)/ 计算频谱dTemp = sqrt(FDj * h + i.real() * FDj * h + i.real() + FDj * h + i.imag() * FDj * h + i.imag() /

12、 100;/ 判断是否超过255if (dTemp > 255)/ 对于超过的，直接设置为255dTemp = 255;/ 指向DIB第(i<h/2 ? i+h/2 : i-h/2)行，第(j<w/2 ? j+w/2 : j-w/2)个象素的指针/ 此处不直接取i和j，是为了将变换后的原点移到中心/lpSrc = (unsigned char*)lpDIBBits + lLineBytes * (lHeight - 1 - i) + j;lpSrc = (unsigned char*)lpDIBBits + lLineBytes * (lHeight - 1 - (i<

13、h/2 ? i+h/2 : i-h/2) + (j<w/2 ? j+w/2 : j-w/2);/ 更新源图像* (lpSrc) = (BYTE)(dTemp);/ 删除临时变量delete TD;delete FD;/ 返回return TRUE;/ CMyDIPViewIMPLEMENT_DYNCREATE(CMyDIPView, CScrollView)BEGIN_MESSAGE_MAP(CMyDIPView, CScrollView)/AFX_MSG_MAP(CMyDIPView)ON_COMMAND(ID_MENUITEM32778, OnMenuitem32778)/AFX_M

14、SG_MAP/ Standard printing commandsON_COMMAND(ID_FILE_PRINT, CScrollView:OnFilePrint)ON_COMMAND(ID_FILE_PRINT_DIRECT, CScrollView:OnFilePrint)ON_COMMAND(ID_FILE_PRINT_PREVIEW, CScrollView:OnFilePrintPreview)END_MESSAGE_MAP()/ CMyDIPView construction/destructionCMyDIPView:CMyDIPView()/ TODO: add const

15、ruction code hereCMyDIPView:CMyDIPView()BOOL CMyDIPView:PreCreateWindow(CREATESTRUCT& cs)/ TODO: Modify the Window class or styles here by modifying/ the CREATESTRUCT csreturn CScrollView:PreCreateWindow(cs);/ CMyDIPView drawingvoid CMyDIPView:OnDraw(CDC* pDC) CMyDIPDoc* pDoc = GetDocument();ASS

16、ERT_VALID(pDoc);if(pDoc->m_hDIB = NULL)return ;/ TODO: add draw code for native data hereint i,j; unsigned char *lpSrc;LPSTR lpDIB = (LPSTR) :GlobalLock(HGLOBAL) pDoc->m_hDIB);int cxDIB = (int) :DIBWidth(lpDIB); / Size of DIB - xint cyDIB = (int) :DIBHeight(lpDIB); / Size of DIB - yLPSTR lpDIB

17、Bits=:FindDIBBits (lpDIB);/ 计算图像每行的字节数long lLineBytes = WIDTHBYTES(cxDIB * 8);/ 每行for(i = 0; i < cyDIB; i+)/ 每列for(j = 0; j < cxDIB; j+)/ 指向DIB第i行，第j个象素的指针lpSrc = (unsigned char*)lpDIBBits + lLineBytes * (cyDIB - 1 - i) + j;/ 计算新的灰度值/*(lpSrc) = BYTE(255-*lpSrc);:GlobalUnlock(HGLOBAL) pDoc->

18、m_hDIB);CRect rect(0,0,cxDIB,cyDIB), rcDIB(0,0,cxDIB,cyDIB);:PaintDIB(pDC->m_hDC, &rect, pDoc->m_hDIB, &rcDIB, pDoc->m_palDIB);/ CMyDIPView printingBOOL CMyDIPView:OnPreparePrinting(CPrintInfo* pInfo)/ default preparationreturn DoPreparePrinting(pInfo);void CMyDIPView:OnBeginPrintin

19、g(CDC* /*pDC*/, CPrintInfo* /*pInfo*/)/ TODO: add extra initialization before printingvoid CMyDIPView:OnEndPrinting(CDC* /*pDC*/, CPrintInfo* /*pInfo*/)/ TODO: add cleanup after printing/ CMyDIPView diagnostics#ifdef _DEBUGvoid CMyDIPView:AssertValid() constCScrollView:AssertValid();void CMyDIPView:

20、Dump(CDumpContext& dc) constCScrollView:Dump(dc);CMyDIPDoc* CMyDIPView:GetDocument() / non-debug version is inlineASSERT(m_pDocument->IsKindOf(RUNTIME_CLASS(CMyDIPDoc);return (CMyDIPDoc*)m_pDocument;#endif /_DEBUGvoid CMyDIPView:OnMenuitem32778() / 图像付立叶变换/ 获取文档CMyDIPDoc* pDoc = GetDocument()

21、;/ 指向DIB的指针LPSTRlpDIB;/ 指向DIB象素指针LPSTR lpDIBBits;/ 锁定DIBlpDIB = (LPSTR) :GlobalLock(HGLOBAL) pDoc->GetHDIB();/ 找到DIB图像象素起始位置lpDIBBits = :FindDIBBits(lpDIB);/ 判断是否是8-bpp位图（这里为了方便，只处理8-bpp位图的付立叶变换，其它的可以类推）if (:DIBNumColors(lpDIB) != 256)/ 提示用户MessageBox("目前只支持256色位图的付立叶变换！", "系统提示&qu

22、ot; ,MB_ICONINFORMATION | MB_OK);/ 解除锁定:GlobalUnlock(HGLOBAL) pDoc->GetHDIB();/ 返回return;/ 更改光标形状BeginWaitCursor();/ 调用Fourier()函数进行付立叶变换if (:Fourier(lpDIBBits, :DIBWidth(lpDIB), :DIBHeight(lpDIB)/ 设置脏标记pDoc->SetModifiedFlag(TRUE);/ 更新视图pDoc->UpdateAllViews(NULL);else/ 提示用户MessageBox("分配内存失败！", "系统提示" , MB_ICONINFORMATION | MB_OK);/ 解除锁定:GlobalUnlock(HGLOBAL) pDoc->GetHDIB();/ 恢复光标EndWaitCursor();2.程序运行结果截图实验原始图像：一次傅里叶变换第二次傅里叶变换第三次傅里叶变换3.傅立叶变换的性质及意义傅立叶变换是数字信号处理领域一种很重要的算法。要知道傅立叶变换算法的意义，首先要了解傅立叶原