Matlab数字图像处理实验指导
实验目的：
通过实验，深入理解和掌握图像处理的基本技术，提高动手实践能力。
实验环境：
Matlab变成
实验一图像的几何变换
实验内容：设计一个程序，能够实现图像的各种几何变换。
实验要求：读入图像，打开图像，实现图像的平移变换、比例缩放、转宜变换、镜像变换、旋转变换等操作。
实验原理：
图像儿何变换乂称为图像空间变换，它将一幅图像中的朋标位置映射到另一?幅图像中的新处标位査。学习儿何变换的关键就是要确定这种空间映射关系，以及映射过程小的变化参数。
几何变换不改变图像的像素值，只是在图像平面上进行像素的重新安排。一个几何变换需要两部分运算：首先是空间变换所盂的运算，如平移、镜像和旋转等，需要川它來表示输出图像与输入图像Z间的(像素)映射关系；此外，述需耍使用灰度插值算法，因为按照这种变换关系进行计算，输出图像的像索可能被映射到输入图像的非整数坐标上。
设原图像f(xO,yO)经过几何变换产生的口标图像为g(xlYL),则该空间变换(映射)关系可表示为：
xl=s(xO,yO)
yl=t(xO,yO)
其中,s(xO,yO)和t(xO,yO)为由f(xO,yO)到g(xl,yl)的坐标换变换函数。
一、图像平移
图像平移就是将图像屮所有-的点按照指定的平移最水平或者垂肓移动。
二、 图像镜像
镜像变换又分为水平镜像和垂肓镜像。水平镜像即将图像左半部分和右半部分以图像竖 直中轴线为中心轴进行对换；而竖直镜像则是将图像上半部分和下半部分以图像水平中轴线 为中心轴进行对换。
三、 图像转置
图像转置是将图像像索的x 坐标和y 坐标呼唤。图像的大小会随Z 改变——高度和宽度 将呼唤。
原图像 100
150
200
250
50
100 150 200 250 50
100 150 200 250 图像平移
50 100 150 200 250
原图像 水平镜像
原图像图像转蛊
四、 图像的缩放
图像缩放是指将图像大小按照指定的比率放大或者缩小。
图像缩放函数imresizef);调用格式如卜：
B=imresize(A /Scale,method);
参数A 为要进行缩放的原始图像。
Scale 为统一的缩放比例。如果希望在x 和y 方向上以不同比例进行缩放，可用如下调 用形式。B= imresize(A,[mrows ncols],method);向量参数[mrows ncols]指明了变换后目标图像 B 的具休行数和列数。
可选参数method 允许imresize()函数指定的插值方法。其合法值同imtransform()函数, 但默认为最近邻插值。
五、 图像旋转
图像旋转一般是指将图像围绕某一指定点旋转一定的角度。旋转通常也会改变图像的人 小。
1)
以原点为中心的图像旋转 2) 以任意点为中心的图像旋转
图像缩放
图像旋转变换的效果受具体插值方法的影响较为明显，木实验给出的旋转采川最近邻插 值。
1) 围绕图像中心的旋转变换函数imrotateO,调用方式如卜?：
B=imrotate(A,a ngle, method,'crop');
A 是要旋转的图像。
Angle 为旋转和度，单位为度，如为其指定一个正肓，则按逆时针方向旋转图像。 可选参数method 允许imrototef)函数指定的插值方法。
'crop ，选项会裁剪旋转后增大的图像，使得到的图像和原图像人小一致。
2) 以任意点为中心的图像旋转(自学)
通过调用imtransformO 函数来实现以任意点为屮心的图像旋转。
提不：
读取图像文件:l=imread ('filename.fmf)
写入图像:imwrite (l/filename.fmt 7)
显示图像:imshow (I)
显示子图：subplot(m,m,p)
二维空间变换：B=imtransform(A,TFORM,method);
参数A 为要进行几何变换的图像。
空间变换结构TFORM 指定了具体的变换类型。可以通过两种方法来创建TFORM 结构，即使 用maketform()W 数和卬2tform()函数。这里给出使用maketform()函数获得TFORM 结构的方 法。T=maketform(transformtype,Matrix);参数 transformtype 指定了变换的类型，如常见 的Qffine ，为二维或多维仿射变换，包括平移、旋转、比例、拉伸和错切等。Matrix 为相应的 仿射变换矩阵。
可选参数method 允许imtransform()函数选择的插值方法。其合法值'bicubic'-三次插
值，原團像 逆时针旋转30度
'bilinear7-双线性插值，’nearest，■最近邻插值。默认吋为'bilinear'-双线性插值。函数输出
B为经imtransform()变换后的口标图像。
插值算法
实现几何运算时，有两种方法。第一种为向前映射法，其原理是将输入图像的灰度一个像素一个像素地转移到输岀图像屮，即从原图像坐标计算出目标图像坐标：
g(xl,y l)=f(s(xO,yO),t(xO,yO)) o
第二种为向后映射法，它是向前映射变换的逆，即输出像素一个一个地映射回输入图像中。如果一个输出像素映射到的不是输入图像的采样栅格的整数坐标处的像素点，则其灰度值就需要基于整数坐标的灰度值进行推断，这就是插值。
内容：设计一个或多个程序，实现二维傅里叶变换和逆变换，并对图像进行傅里叶变换平移。要求：读入图像，打开图像，实现二维傅里叶变换和逆变换，等操作。
实验原理：
图像变换在图像处理和图像分析屮片有重要的地位。将图像从空间域变换到其它域(如频域)的数学变换，便于进行图像处理和图像分析。常用的变换冇：傅里叶变换、离散余弦变换、K?L变换等。
一、傅立叶变换
若f(x)为一维连续实函数，则它的傅里叶变换可定义为
= r f 傅立叶逆变换定义如下: 函数f(x)和F(u)被称为傅立叶变换对.即对于任一函数f(x),其傅立叶变换F(u)是惟一的; 反之，对于任一函数F(u),其傅立叶逆变换f(x)也是惟一的. F(u)可以表示为如下形式： =彳(")+ jl (iz) ［斤2(") + l\uW Flu) E(u) =\尸(〃)|2 E (u)称为函数f(x)的能最谱或功率谱。 二、离散余弦变换(DCT变换) 傅立叶变换的一个最人的问题是：它的参数都是复数，在数据的描述上相当于实数的两倍。为此，我们希望有一种能够达到相同功能但数据量又不人的变换。在此期望下，产生了DCT变换。余弦变换实际上是傅立叶变换的实数部分.余弦变换主要用于图像的压缩，如H 前的国际压缩标准的JPEG格式中就用到了DCT变换.具体的做法与DFT相似.给高频系数大间隔量化，低频部分小间隔量化。 三、K?L变换 K?L变换(Karhunen-Loeve Transform)是数字图像处理屮具有广泛应用的一类重要变换; 乂称为特征向虽变换、主分：S变换或霍特林变换;K?L变换既冇连续形式的变换也冇离散形式的变换;它是完全从图像的统计性质出发实现的变换。 数字图像中主要应用离散K-L变换，优点是去相关性好，该变换在数据压缩、图像旋转、遥感多光谱图像的特征选择和统计识别等方面具有重耍意义。 K?L变换形式：y = A(X - /zz v) 操作提示： 1, 傅里叶变换 读入图像:imread(); 傅里叶变换要求输入的矩阵为双精度浮点类型,需进行数据类型转换:im2doubel(); 快速傅里叶变换；fft2()；傅里叶逆变换；ifft2(); 把傅里叶变换的零频率部分移到频谱的中间，使用fftshift函数，调用格式为：Y= fftshift(X);把fft函数、fft2函数输出的结果的零频率部分移到数组的中间，对于观察傅里叶变换频谱屮间零频率部分十分有效。 在同一窗口显示3幅图像;subplot();imshow(); 2, 离散余弦变换 二维离散余弦变换；dct2(); 二维离散余弦逆变换；idct2(); 实验三灰度变换和直方图均衡处理 内容：设计一个程序，达到空域图像增强和频域图像增强的ri的。 要求：读入图像，打开图像，实现空域图像增强和频域图像增强等操作。 操作提示: Matlab傅立叶变换、余弦变换和小波变换 1. 离散傅立叶变换的Matlab实现 Matlab函数fft、fft2和fftn分别可以实现一维、二维和N维DFT算法；而函数ifft、 ifft2和ifftn则用来计算反DFT。这些函数的调用格式如下: A=fft(X,N,DIM) 其中，X表示输入图像；N表示采样间隔点，如果X小于该数值，那么Matlab将会对X 进行零填充，否则将进行截取，使之长度为N ；DIM表示要进行离散傅立叶变换。 A=fft2(X,MR0WS, NCOLS) 其中，MROWS和NCOLS指定对X进行零填充后的X大小。别可以实现一维、二维和N 维DFT A=fftn(X,SIZE) 其中，SIZE是一个向量，它们每一个元素都将指定X相应维进行零填充后的长度。 函数ifft、ifft2和ifftn的调用格式于对应的离散傅立叶变换函数一致。 别可以实现一维、二维和N维DFT 例子：图像的二维傅立叶频谱 1.离散傅立叶变换的Matlab实现％读入原始图像 I =imread('le na.bmp'); imshow(l) %求离散傅立叶频谱 J=fftshift(fft2(l)); figure;别可以实现一维、二维和N维DFF imshow( log(abs(J) ),[8,10]) 2离散余弦变换的Matlab实现 2.1. dct2 函数 功能：二维DCT变换Matlab 格式：B=dct2(A) B=dct2(代m,n) B=dct2(A,[m,n])函数fft、fft2 和fftn 分 说明：B=dct2(A)计算A的DCT变换B , A与B的人小相同；B=dct2(代m,n)和 B=dct2(A,[m,n])通过对A补0或剪裁，使B的大小为mxn。 2.2. dict2 函数 功能：DCT反变换 格式:B=idct2(A) B=idct2(A,m,n)?J可以实现一维、二维和N维DFT B=idct2(A,[m,n]) 说明:B=idct2(A)计算A的DCT反变换B , A与B的大小相同；B=idct2(A,m,n)和 B=idct2(A,[m,n])通过对A补0或剪裁，使B的大小为mm。 Mat lab 2.3. det mt x 函数 功能：计算DCT变换矩阵 格式:D=dctmtx(n) 说明：D=dctmtx(n)返冋一个nx n的DCT变换矩阵，输出矩阵D为double类型。 3. 图像小波变换的Matlab实现函数 3.1 —维小波变换的Mat lab实现 (1) dwt 函数Mat lab 功能：一维离散小波变换 格式：[cA,cD] = dwt(X,wname,) [cA,cD] = dwt(X,Lo_D,Hi_D)别可以实现一维、二维和N 维DFT 说明：[cA,cD] = dwt(X,,wname,)使用指定的小波基函数’wname'对信号X进行分解，cA、cD 分别为近似分量和细节分量；[cA,cD] = dwt(X,Lo_D,Hi_D)使用指定的滤波器组Lo_D、Hi_D 对信号进行分解。 (2) idwt 函数 功能：一维离散小波反变换 格式：X= idwt(cA,cD,'wname,) X=idwt(cA,cD,Lo_R,Hi_R) X= idwt(cA,cD,'wname,,L)函数 X= idwt(cA,cD,Lo_R,Hi_R,L) 说明：X=idwt(cA,cD,wname')由近似分量cA和细节分量cD经小波反变换重构原始信号 X O 'wname'为所选的小波函数 X=idwt(cA,cD,Lo_R,Hi_R)用指定的重构滤波器Lo_R和Hi_R经小波反变换重构原始信号X O X= idwt(cA,cD,'wname',L)和X= idwt(cA,cD,Lo_R,Hi_R,L)指定返回信号X 中心附近的L 个点。 1 ?离散傅立叶变换的Matlab实现 3.2二维小波变换的Matlab实现 二维小波变换的函数别可以实现一维、二维和N维DFT dwt2二维离散小波变换 wavedec2二维信号的多层小波分解 idwt2二维离散小波反变换Matlab waverec2二维信号的多层小波重构 wrcoef2由多层小波分解重构某一层的分解信号upcoef2由多层小波分解重构近似分虽或细节分:fi 1.离散傅立叶变换的Matlab实现detcoef2提収二维信号小波分解的细节分量appcoef2提取二维信号小波分解的近似分量upwlev2二维小波分解的单层重构1 ?离散傅立叶变换的Matlab实现dwtpet2二维周期小波变换 idwtper2二维周期小波反变换 (1) wcodemat 函数 功能：对数据矩阵进行伪彩色编码 格式：Y= wcodemat(X, NB,OPT,ABSOL) Y= wcod e m at (X, N B, OPT) Y=wcodemat(X,NB) Y= wcodemat (X) 说明：Y= wcodemat(X,NB,OPT,ABSOL)返回数据矩阵X的编码矩阵丫；NB伪编码的最人值，即编码范围为0?NB,缺省值NB=16； OPT指定了编码的方式(缺省值为’mat'),即：別可以实现一维、二维和N维DFT OPT='row* ,按彳亍编码 OPT=，cor ,按列编码 OFT^'mat1 ,按整个矩阵编码 ABSOL是函数的控制参数(缺省值为T),即： ABSOL-0时，返回编码矩阵 ABSOL-1时，返回数据矩阵的绝对值ABS(X)1.离散傅立叶变换的Matlab实现 ⑵dwt2函数 功能：二维离散小波变换 格式：[cA,cH,c=cD] = dwt2(X,'wname') [cA, cH,cM cD] = dwt2( X, Lo_D, Hi_D) 说明：[cA,cH,cUcD] = dwt2(X,'wname‘)使用指定的小波基函数’wname'对二维信号X进行二维离散小波变幻；cA, cH.cVcD分别为近似分量、水平细节分量、垂直细节分量和对角细节分[cA,cH,cVcD] = dwt2(X,Lo_D5Hi_D)使用指定的分解低通和高通滤波器Lo_D 和Hi_D分解信号X o 1.离散傅立叶变换的Matlab实现 ⑶wavedec2函数 功能：二维信号的多层小波分解1?离散傅立叶变换的Matlab实现 格式：[C,S] = wavedec2(X,N,'wname') [C,S] = wavedec2( X, N,Lo_D, Hi_D) 说Rhl：[C,S] = wavedec2(X,N,,wname,)使用小波基函数，wname，对二维信号X进行N层分解；[C,S] = wavedec2(X,N,Lo_D,Hi_D)使用指定的分解低通和高通滤波器Lo_D和Hi_D 分解信号X o别可以实现一维、二维和N维DFT (4) idwt2 函数 功能：二维离散小波反变换 格式：X= idwt2(cA,cH,cVcD,,wname') X= idwt2(cA,cH,cVcD,Lo_R,Hi_R) X= idwt2(cA,cH,cVcD；wname,,S)别可以实现一?维、二维和N 维DFT X= idwt2(cA,cH,cVcD,Lo_R,Hi_R,S) 说明：X=idwt2(cAcH,cVcD,wname')由信号小波分解的近似信号cA和细节信号cH、cH、eV、cD经小波反变换重构原信号X ； X= idwt2(cA,cH,cMcD,Lo_R,Hi_R)使用指定的重构低通和高通滤波器Lo_R和Hi_R重构原信号X ；X=idwt2(cA,cH,cV,cD,'wname',S)和X= idwt2(cA,cH,cV,cD,Lo_R,Hi_R,S)返回中心附近的S 个数据点。 ⑸waverec2函数 说明：二维信号的多层小波重构 格式：X=waverec2(C,S,,wname,) X= waverec2(C, S, Lo_R, Hi_R) 说明：X=waverec2(C,S；wname*)由多层二维小波分解的结果C、S重构原始信号X , ?wname'为使川的小波基函数；X= waverec2(C,Lo_R,Hi_R)使用重构低通和高通滤波器Lo_R和Hi_R 重构原信号。 Allnodes计算树结点 appcoef提取一维小波变换低频系数 appcoef2提取二维小波分解低频系数bestlevt计算完整敲佳小波包树别可以实现一维、二维和N维DFF besttree计算最佳(优)树 * biorfilt双止交样条小波滤波器组 biorwavf双正交样条小波滤波器Matlab * centfrq求小波中心频率cgauwavf Complex Gaussian 小波 cmorwavf coiflets小波滤波器1.离散傅立叶变换的Matlab实现cwt 一维连续小波变换dbaux Daubechies小波滤波器计算 dbwavf Daubechies 小波滤波器dbwavf(W) W='dbN' N=1,2,3,...,50 别可以实现一维、 二维和N维DFT ddencmp获取默认值阈值(软或硬)爛标准depo2ind将深度■位置结点形式转化成索引结点形式 detcoef提収一维小波变换高频系数Matlab detcoef2提取二维小波分解高频系数 disp显示文木或矩阵 drawtree画小波包分解树(GUI)别可以实现一维、二维和N维DFT dtree 构造DTREE 类dwt单尺度一维离散小波变换 dwt2单尺度二维离散小波变换别可以实现一维、二维和N维DFT dwtmode离散小波变换拓展模式 * dyaddown 二元取样 * dyadup二元插值1.离散傅立叶变换的Mat lab实现entrupd更新小波包的嫡值fbspwavf B样条小波gauswavf Gaussian 小波Mat lab get获取对象属性值idwt单尺度一维离散小波逆变换idwt2单尺度二维离散小波逆变换ind2depo将索引结点形式转化成深度一位置结点形式 * i nt wave积分小波数 isnode判断结点是否存在Matlab istnode判断结点是否是终结点并返回排列值 iswt —维逆SWT(Stationary Wavelet Transform)变换iswt2 二维逆SWT 变换Matlab leaves Determine terminal no des mexihat墨西哥帽小波 meyer Meyer小波别町以实现一维、二维和N维DFT meyeraux Meyer小波辅助函数 morlet Morlet 小波 nodease计算上溯结点nodedesc计算下溯结点（子结点）nodejoin重组结点 nodepar寻找父结点别可以实现一维、二维和N维DFT nodespit分割份解）结点 no leaves Determi ne no ntermi nal no des ntnode Number of terminal nodes 函数fft> fft2 和fftn 分 ntree Constructor for the class NTREE orthfilt正交小波滤波器组 plot绘制向量或矩阵的图形qmf镜像二次滤波器 rbiowavf Fteverse biorthogonal spline wavelet filters read读取二进制数据 readtree读取小波包分解树 scal2frq Scale to frequency set Mat lab shanwavf Shannon wavelets swt 一维SWT（Stationary Wavelet Transform）变换swt2二维SWT变换 symaux Symlet wavelet filter computation, symwavf ^mlets小波滤波器 thselect信号消噪的阈值选择 thodes Ftefere nces treedpth求树的深度treeord求树结构的叉数upcoef 一维小波分解系数的直接重构 upcoef2二维小波分解系数的直接重构upwlev单尺度一维小波分解的重构upwlev2单尺度二维小波分解的重构wavedec单尺度一维小波分解 wavedec2多尺度二维小波分解Matlab wavedemo小波工具箱函数demo wavefun小波函数和尺度函数 wavefun2二维小波函数和尺度函数别可以实现一维、二维和N维DFT wavemenu小波工具箱函数menu图形界面调用函数 wavemngr小波管理函数 waverec多尺度?维小波重构1.离散傅立叶变换的Matlab实现waverec2多尺度二维小波重构 wbmpe n Pen al i zed threshold for wavelet 1-D or 2-D de?noising wcodemat对矩阵进行量化编码1.离散傅立叶变换的Matlab实现wdcbm Thresholds for wavelet 1?D using Birge-Massart strategy wdcbm2 Thresholds for wavelet 2-D using Birge-Massart strategy wden用小波进行一维信号的消噪或压缩 wdencmp De-noising or compression using wavelets wentropy计算小波包的爛 wextend Extend a vector or a matrix * wfilters小波滤波器 wkeep提取向量或矩阵屮的一部分 * wmaxlev计算小波分解的最大尺度1.离散傅立叶变换的Matlab实现wnoise产生含噪声的测试函数数据 wnoisest估计一维小波的系数的标准偏差 wp2wtree从小波包树中提取小波树 1.离散傅立叶变换的Matlab实现 wpcoef计算小波包系数 wpcutree剪切小波包分解树 wpdec 一维小波包的分解wpdec2二维小波包的分解wpdencmp用小波包进行信号的消 噪或压缩 wpfun小波包两数 wpjoin 重组小波包 wprcoef小波包分解系数的重构 wprec ?维小波包分解的重构1 .离散傅立叶变换的Matlab实现wprec2二维小波包分解的重构 wpsplt分割(分解)小波包 wpthcoef进行小波包分解系数的阈值处理wptree 显示小波包树结构 wpviewcf Plot the colored wavelet packet coefficients. wrcoef对一维小波系数进行单支垂构别口J以实现一维、二维和N维DFT wrcoef2对二维小波系数进行单支重构 wrev向录逆序 write向缓冲区内存写进数据 wtbo Constructor for the class WTBO wthcoef 一维信号的小波系数阈值处理wthcoef2二维信号的小波系数阈值处理1 .离散傅立叶变换的Matlab实现wthresh进行软阈值或硬阈值处理 wthrmngr阈值设置管理wtreemgr 管理树结构Matlab> = 50 1. cwt函数 函数功能：实现一维连续小波变换的函数。 语法格式：OOEFS=cwt(S, SCALES, 'wname') COEFS=cwt(S, SCALES, 'wname', 'plot') COEFS=cwt(S, SCALES, ,wname,, 'PLOTMODE) COEFS=cwt(S, SCALES, 'wname1, 1 PLOTMODE, XLIM) 使用说明：cwt为一维小波变换的两数。 格式COEFS=cwt(S, SCALES, 'wname')釆用'wname1小波，在正、实尺度SCALES K计 算向量一维小波系数。 格式COEFS=cwt(S, SCALES,'wnamejplot')除了计算小波系数外，还加以图形显示。格式COEFS=cwt(S, SCALES, 'wname*, 'PLOTMODE)计算并画出连续小波变换的系 数，并使用PLOTMODE 对图形着色。 格式 COEFS=cwt(S, SCALES, ,wname ,, 'plot')和当于 格式 COEFS=cwt(S, SCALES, 'wname*. ,PLOTMODE)屮的语法 OOEFS=cwt(S, SCALES, *wname', bbsglb') 格式 COEFS=cwt(S, SCALES ,wname ,J PLOTMODE, XLIM)能够计算并画出连续小波 变换的系数。系数使用PLOTMODE 和XLIM 进行着色。其中：XLIM=[x1 ,x2],并且冇如下 关系：1 v = x1 v = x2v = length(S)。 MODE 值 含义 scale-by-scale 着色模式 考虑所有尺度的着色模式 使用系数绝对值的scale-by-scale 着色模式 使用系数绝对值并考虑所冇尺度的着色模式 例子：cwt 函数的应用举例 t=linspace(-1:1:512); s=1-abs(t); c=cwt(s,1:32,'cgau4,); c=cwt(s,[64 32 16:-2:2]); c=cwt(s,[3 18 12.9 7 1.5],'db2); c= cwt(s, 1:64/sym4','abslvr,[100:400]); 2. pat2cwav 函数 函数功能：由模式构造小波。 语法格式：[PSI,XVAL,NC] = pat2cwav(YPAT METHOD, POLDEGREE, REGULARITY) 使川说明：该函数计算由XVAL 和PSI 给定并川于连续小波变换的小波函数，该小波向量 YPAT 定义的模式构造，方羌为1。其屮模式隐含的x 值是xpat=linspace(051,length(YPZ^T))o 常数NC 的选取应保证通过以下方式的最小二乘拟合，NCPSI 在区间[0,1]±近似于 YPAT ： 当 METHD 等于'polynomial'时，为 POLDEGREE 多项式； 当METHD 等于Sthconst'时，为正交两数空间的投影； 参数REGULARITY 定义了在0点和1点的边界约束，可以是'continuous', ?differentiable* 或'none'。 当 METHOD 为'polynomiar 时： 当 REGUUVRITY 等于'continuous',则必须 POLDEGREE> = 3； 当 REGULARITY 等于?differentiable',则必须 POLDEGREE 函数名称：im2dowbel 把图像数据类型转换为双精度浮点类型。I2 = im2double(⑴ 把灰度图像11的数据类型转换成转换成双精度浮点类型，返回。(但11并不发生改变) RGB2 = im2double(RGB1) 把真彩色图像的数据类型转换成双粹度浮点类型的，并返回。 Ivl' ?gib* 'abslvl'或'Ivlabs' 'absglb 或 glbabs' I = im2double(BW) 把二值图像BW的数据格式转换成双精度浮点类型的。 X2 = im2double(X1,'indexed') 函数名称：im2uint8 函数功能：把图像数据类型转换为无符号八位整型。如果输入图像是无符号八位整型的，返回的图像和源图像相同。如果源图像不是无符号八位整烈的，该函数将返四和源图像相同但数据类型为uint8的图像(必要时对图像进行调整)。 调用格式： I2 = im2uint8(l1) 把灰度图像11的数据类型转换成无符号八位整型的后，返回。(但11并不发牛改变) RGB2 = im2uint8(RGB1) 把具彩色图像的数据类型转换成无符号八位整型的,并返冋。注意这电并不是将真彩色图像转换成为黑口图像，这里仅仅改变了图像数据的类型。 I = im2uint8(BW) 把「?值图像BW的数据格式转换成uint8格式的，由于二值图像的图像数据仅仅取0和1, 所以这里的转换效果仅仅是把1映射成255o 函数名称：im2bw matlab中DIP工具箱函数im2bw使用阈值(threshold)变换法把灰度图像(grayscale image) 转换成一?值图像。所谓二值图像，一?般意义上是指只有纯黑(0)、纯H (255)两种颜色的图像。当然，也可以是其他任意两种颜色的组合。 语法 BW = im2bw(l, level) BW = im2bw(X, map, level) BW = im2bw(RGB, level) 其中level就是设置阈值的。level取值范围［0,1］。 小波的基函数有哪些？ 常见的小波基有：Haar小波基、db系列小波基、Biorthogonal(biorNr.Nd)小波系＞ CoiHet(coifN) 小波系、SymletsA(symN)小波系、Molet(morl)小波、Mexican Hat (mexh)小波、Meyer 小波 function x=MyIFFT_FB(y) %MyIFFT_TB:My Inverse Fast Fourier Transform Time Based %按频率抽取基2-傅里叶逆变换算法 %input: % y -- 傅里叶正变换结果,1*N的向量 %output: % x -- 逆变换结果，1*N的向量 %参考文献： % https://www.doczj.com/doc/d92176962.html,/view/fea1e985b9d528ea81c779ee.html N=length(y); x=conj(y); %求共轭 x=MyFFT_FB(x);%求FFT x=conj(x);%求共轭 x=x./N;%除以N end %% 内嵌函数====================================================== function y=MyFFT_FB(x,n) %MYFFT_TB:My Fast Fourier Transform Frequency Based %按频率抽取基2-fft算法 %input: % x -- 输入的一维样本 % n -- 变换长度，缺省时n=length(x) 当n小于x数据长度时，x数据被截断到第n个数据% 当n大于时，x数据在尾部补0直到x 含n个数据 %output: % y -- 1*n的向量，快速傅里叶变换结果 %variable define: % N -- 一维数据x的长度 % xtem -- 临时储存x数据用 % m,M -- 对N进行分解N=2^m*M,M为不能被2整除的整数 % two_m -- 2^m % adr -- 变址，1*N的向量 % l -- 当前蝶形运算的级数 % W -- 长为N/2的向量，记录W(0,N),W(1,N),...W(N/2-1,N) % d -- 蝶形运算两点间距离 % t -- 第l级蝶形运算含有的奇偶数组的个数 % mul -- 标量，乘数 % ind1,ind2 -- 标量，下标 % tem -- 标量，用于临时储存 %参考文献： % https://www.doczj.com/doc/d92176962.html,/view/fea1e985b9d528ea81c779ee.html %% 输入参数个数检查 function [Spec,Freq]=STFT(Sig,nLevel,WinLen,SampFreq) %计算离散信号的短时傅里叶变换； % Sig 待分析信号； % nLevel 频率轴长度划分（默认值512）； % WinLen 汉宁窗长度（默认值64）； % SampFreq 信号的采样频率（默认值1）； if (nargin <1), error('At least one parameter required!'); end; Sig=real(Sig); SigLen=length(Sig); if (nargin <4), SampFreq=1; end if (nargin <3), WinLen=64; end if (nargin <2), nLevel=513; end nLevel=ceil(nLevel/2)*2+1; WinLen=ceil(WinLen/2)*2+1; WinFun=exp(-6*linspace(-1,1,WinLen).^2); WinFun=WinFun/norm(WinFun); Lh=(WinLen-1)/2; Ln=(nLevel-1)/2; Spec=zeros(nLevel,SigLen); wait=waitbar(0,'Under calculation,please wait...'); for iLoop=1:SigLen, waitbar(iLoop/SigLen,wait); iLeft=min([iLoop-1,Lh,Ln]); iRight=min([SigLen-iLoop,Lh,Ln]); iIndex=-iLeft:iRight; iIndex1=iIndex+iLoop; iIndex2=iIndex+Lh+1; Index=iIndex+Ln+1; Spec(Index,iLoop)=Sig(iIndex1).*conj(WinFun(iIndex2)); end; close(wait); Spec=fft(Spec); Spec=abs(Spec(1:(end-1)/2,:)); Matlab傅里叶变换傅里叶逆变换 %% 信号经过傅里叶变换然后进行傅里叶逆变换后信号的变化 clear all;clc; %------Author&Date------ %Author: %Date: 2013/07/31 %========================================================================== Fs=8e3; %采样率 t=0:1/Fs:1; %采样点 len=length(t); %采样长度 f1=10; %频率1 f2=100; %频率2 f3=1000; %频率3 A1=1; %幅度1 A2=0.8; %幅度2 A3=0.3; %幅度3 MaxS=A1+A2+A3; %信号幅度的最大值 signal=A1*sin(2*pi*f1*t)+A2*sin(2*pi*f2*t)+A3*sin(2*pi*f3*t); X=fft(signal,len); %傅里叶变换 magX=abs(X); %信号的幅度 angX=angle(X); %信号的相位 Y=magX.*exp(1i*angX); %信号的频域表示 y=ifft(Y,len); %信号进行傅里叶逆变换 y=real(y); er=signal-y; %原始信号和还原信号的误差 subplot(311);plot(t,signal);axis([0 1 -MaxS MaxS]);xlabel('时间');ylabel('振幅');title('原始信号'); subplot(312);plot(t,y);axis([0 1 -MaxS MaxS]);xlabel('时间');ylabel('振幅');title('还原信号'); subplot(313);plot(t,er);xlabel('时间');ylabel('振幅');title('误差'); % End Script %傅里叶变换 clc;clear all;close all; tic Fs=128;%采样频率，频谱图的最大频率 T=1/Fs;%采样时间，原始信号的时间间隔 L=256;%原始信号的长度，即原始离散信号的点数 t=(0:L-1)*T;%原始信号的时间取值范围 x=7*cos(2*pi*15*t-pi)+3*cos(2*pi*40*t-90*pi/180)+3*cos(2*pi*30*t-90*pi/ 180); z=7*cos(2*pi*15*t-pi)+3*cos(2*pi*40*t-90*pi/180); z1=6*cos(2*pi*30*t-90*pi/180); z1(1:L/2)=0; z=z+z1; y=x;%+randn(size(t)); figure; plot(t,y) title('含噪信号') xlabel('时间（s）') hold on plot(t,z,'r--') N=2^nextpow2(L);%N为使2^N>=L的最小幂 Y=fft(y,N)/N*2; Z=fft(z,N)/N*2;%快速傅里叶变换之后每个点的幅值是直流信号以外的原始信号幅值的N/2倍（是直流信号的N倍） f=Fs/N*(0:N-1);%频谱图的频率取值范围 A=abs(Y);%幅值 A1=abs(Z); B=A; %让很小的数置零. B1=A1; A(A<10^-10)=0; % A1(A1<10^-10)=0; P=angle(Y).*A./B; P1=angle(Z).*A1./B1; P=unwrap(P,pi);%初相位值，以除去了振幅为零时的相位值 P1=unwrap(P1,pi); figure subplot(211) plot(f(1:N/2),A(1:N/2))%函数ffs返回值的数据结构具有对称性，因此只取前一半 hold on plot(f(1:N/2),A1(1:N/2),'r--') title('幅值频谱') 1.请用MATLAB编写程序，实现任意两个有限长度序列的卷积和。要求用图 形显示两个序列及卷积结果。 解：y(n)=∑x(i)h(n-i) 假设x(n)={1,2,3,4,5}; h(n)={3,6,7,2,1,6}; y(n)=x(n)*h(n) 验证：y[n]=[1,12,28,46,65,72,58,32,29,30] 【程序】 N=5 M=6 L=N+M-1 x=[1,2,3,4,5] h=[3,6,7,2,1,6] y=conv(x,h) nx=0:N-1 nh=0:M-1 ny=0:L-1 subplot(131);stem(nx,x,'*b');xlabel('n');ylabel('x(n)');grid on subplot(132);stem(nh,h,'*b');xlabel('n');ylabel('h(h)');grid on subplot(133);stem(ny,y,'*r');xlabel('n');ylabel('y(h)');grid on 【运行结果】 2.已知两个序列x[n]=cos(n*pi/2), y[n]=e j*pi*n/4x[n],请编写程序绘制 X(e jw)和Y(e jw)和幅度和相角,说明它们的频移关系。 –提示：用abs函数求幅度，用angle求相角。 【程序】 n=0:15; x=cos(n*pi/2); y=exp(j*pi*n/4).*x; X=fft(x); Y=fft(y); magX=abs(X); angX=angle(X); magY=abs(Y); angY=angle(Y); subplot(221);stem(n,magX,'*r');xlabel('频率');ylabel('幅度');grid on; subplot(222);stem(n,angX,'*b');xlabel('频率');ylabel('相位');grid on; subplot(223);stem(n,magY,'*r');xlabel('频率');ylabel('幅度');grid on; subplot(224);stem(n,angY,'*b');xlabel('频率');ylabel('相位');grid on; Matlab数字图像处理实验指导 实验目的： 通过实验，深入理解和掌握图像处理的基本技术，提高动手实践能力。 实验环境： Matlab变成 实验一图像的几何变换 实验内容：设计一个程序，能够实现图像的各种几何变换。 实验要求：读入图像，打开图像，实现图像的平移变换、比例缩放、转置变换、镜像变换、旋转变换等操作。 实验原理： 图像几何变换又称为图像空间变换，它将一幅图像中的坐标位置映射到另一幅图像中的新坐标位置。学习几何变换的关键就是要确定这种空间映射关系，以及映射过程中的变化参数。 几何变换不改变图像的像素值，只是在图像平面上进行像素的重新安排。一个几何变换需要两部分运算：首先是空间变换所需的运算，如平移、镜像和旋转等，需要用它来表示输出图像与输入图像之间的（像素）映射关系；此外，还需要使用灰度插值算法，因为按照这种变换关系进行计算，输出图像的像素可能被映射到输入图像的非整数坐标上。 设原图像f(x0,y0)经过几何变换产生的目标图像为g(x1,y1),则该空间变换（映射）关系可表示为： x1=s(x0,y0) y1=t(x0,y0) 其中，s(x0,y0)和t(x0,y0)为由f(x0,y0)到g(x1,y1)的坐标换变换函数。 一、图像平移 图像平移就是将图像中所有的点按照指定的平移量水平或者垂直移动。 二、图像镜像 镜像变换又分为水平镜像和垂直镜像。水平镜像即将图像左半部分和右半部分以图像竖直中轴线为中心轴进行对换；而竖直镜像则是将图像上半部分和下半部分以图像水平中轴线为中心轴进行对换。 三、图像转置 图像转置是将图像像素的x坐标和y坐标呼唤。图像的大小会随之改变——高度和宽度将呼唤。 一、傅立叶变化的原理； （1）原理 正交级数的展开是其理论基础！将一个在时域收敛的函数展开成一系列不同频率谐波的叠加，从而达到解决周期函数问题的目的。在此基础上进行推广，从而可以对一个非周期函数进行时频变换。 从分析的角度看，他是用简单的函数去逼近（或代替）复杂函数，从几何的角度看，它是以一族正交函数为基向量，将函数空间进行正交分解，相应的系数即为坐标。从变幻的角度的看，他建立了周期函数与序列之间的对应关系；而从物理意义上看，他将信号分解为一些列的简谐波的复合，从而建立了频谱理论。 当然Fourier积分建立在傅氏积分基础上，一个函数除了要满足狄氏条件外， 一般来说还要在积分域上绝对可积，才有古典意义下的傅氏变换。引入衰减因子e^(-st)，从而有了Laplace变换。（好像走远了）。 （2）计算方法 连续傅里叶变换将平方可积的函数f（t）表示成复指数函数的积分或级数形式。 这是将频率域的函数F(ω)表示为时间域的函数f（t）的积分形式。 为 连续傅里叶变换的逆变换 (inverse Fourier transform) 即将时间域的函数f（t）表示为频率域的函数F(ω)的积分。 一般可称函数f（t）为原函数，而称函数F(ω)为傅里叶变换的像函数，原函数和像函数构成一个傅里叶变换对（transform pair）。 二、傅立叶变换的应用； DFT在诸多多领域中有着重要应用，下面仅是颉取的几个例子。需要指出 的是，所有DFT的实际应用都依赖于计算离散傅里叶变换及其逆变换的快速算 法，即快速傅里叶变换（快速傅里叶变换（即FFT ）是计算离散傅里叶变换及其逆变换的快速算法。）。（1）、频谱分析DFT 是连续傅里叶变换的近似。因此可以对连续信号x(t)均匀采样并截断以得到有限长的离散序列，对这一序列作离散傅里叶变换，可以分析连续信号x(t)频谱的性质。前面还提到DFT 应用于频谱分析需要注意的两个问题：即采样可能导致信号混叠和截断信号引起的频谱泄漏。可以通过选择适当的采样频率（见奈奎斯特频率）消减混叠。选择适当的序列长度并加窗可以抑制频谱泄漏。（2）、数据压缩由于人类感官的分辨能力存在极限，因此很多有损压缩算法利用这一点将语音、音频、图像、视频等信号的高频部分除去。高频信号对应于信号的细节，滤除高频信号可以在人类感官可以接受的范围内获得很高的压缩比。这一去除高频分量的处理就是通过离散傅里叶变换完成的。将时域或空域的信号转换到频域，仅储存或传输较低频率上的系数，在解压缩端采用逆变换即可重建信号。（3）、OFDM OFDM （正交频分复用）在宽带无线通信中有重要的应用。这种技术将带宽为N 个等间隔的子载波，可以证明这些子载波相互正交。尤其重要的是，OFDM 调制可以由IDFT 实现，而解调可以由DFT 实现。OFDM 还利用DFT 的移位性质，在每个帧头部加上循环前缀（Cyclic Prefix ），使得只要信道延时小于循环前缀的长度，就能消除信道延时对传输的影响。三、傅里叶变换的本质； 傅里叶变换的公式为dt e t f F t j ?+∞∞--=ωω)()(可以把傅里叶变换也成另外一种形式： t j e t f F ωπ ω),(21)(=可以看出，傅里叶变换的本质是内积，三角函数是完备的正交函数集，不同频率的三 角函数的之间的内积为0，只有频率相等的三角函数做内积时，才不为0。)(2,21)(2121Ω-Ω==?Ω-ΩΩΩπδdt e e e t j t j t j 目录 用Matlab 对信号进行傅里叶变换 (2) Matlab 的傅里叶变换实例 (5) Matlab 方波傅立叶变换画出频谱图 (7) 用 Matlab 对信号进行傅里叶变换 1. 离散序列的傅里叶变换 DTFT(Discrete Time Fourier Transform) 代码： %原离散信号有 8 点 %原信号是 1行 8列的矩阵 %构建原始信号，为指数信号 %频域共-800 +800 的长度（本应是无穷， 高 %求 dtft 变换，采用原始定义的方法，对复指 7 subplot(311) 8 stem(n,xn); 9 title('原始信号(指数信号 )'); 10 subplot(312); 11 plot(w/pi,abs(X)); 12 title('DTFT 变换 ') 结果： 分析：可见，离散序列的 dtft 变换是周期的，这也符合 Nyquist 采样 定理的描述， 连续时间信号经周期采样之后， 所得的离散信号的频谱 是原连续信号频谱的周期延拓。 2. 离散傅里叶变换 1 N=8; 2 n=[0:1:N-1] 3 xn=0.5.^n; 4 5 w=[-800:1:800]*4*pi/800; 频分量很少，故省去) 6 X=xn*exp(-j*(n'*w)); 数分 量求和而得 与 1 中 DTFT 不一样的是， DTFT 的求和区间是整个频域，这对 N=8; % 原离散信号有 8 点 n=[0:1:N-1] %原信号是 1行 8列的矩阵 xn=0.5.^n; %构建原始信号，为指数信号 w=[-8:1:8]*4*pi/8; %频域共 -800 +800 的长度(本应是无穷， 高频分量很少， 故省去) X=xn*exp(-j*(n'*w)); %求 dtft 变换，采用原始定义的方法，对复指数分量求和而得 subplot(311) stem(n,xn); w1=[-4:1:4]*4*pi/4; X1=xn*exp(-j*(n'*w1)); title(' 原始信号 (指数信号 )'); subplot(312); stem(w/pi,abs(X)); title(' 原信号的 16 点 DFT 变换 ') subplot(313) stem(w1/pi,abs(X1)); title(' 原信号的 8 点 DFT 变换 ') 计算机的计算来说是不可以实现的， DFT 就是序列的有限傅里叶变换。 实际上， 1 中代码也只是对频域的 -800 +800 中间的 1601 结果图： 分析： DFT 只是 DTFT 的现实版本，因为 DTFT 要求求和区间无穷， 而 DFT 只在有限点内求和。 3. 快速傅里叶变换 FFT ( Fast Fourier Transform ) 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 1 利用FFT 计算连续时间信号的傅里叶变换 设()x t 是连续时间信号，并假设0t <时()0x t =，则其傅里叶变换由下式给出 0()()i t X x t e dt ωω∞ -=? 令Γ是一个固定的正实数，N 是一个固定的正整数。当,0,1,2,,1k k N ω=Γ=-L 时，利用FFT 算法可计算()X ω。 已知一个固定的时间间隔T ，选择T 足够小，使得每一个T 秒的间隔(1)nT t n T ≤<+内，()x t 的变化很小，则式中积分可近似为 (1)0 ()()()n T iwt nT n X e dt x nT ω∞+-==∑? (1)01[ ]()i t t n T t nT n e x nT i ωω ∞-=+==-=∑ 0 1()i T i nT n e e x nT i ωωω-∞-=-=∑ （27） 假设N 足够大，对于所有n N ≥的整数，幅值()x nT 很小，则式（27）变为 1 01()()i T N i nT n e X e x nT i ωωωω---=-=∑ （28） 当2/k NT ωπ=时，式（28）两边的值为 2/2/12/0211()()[]2/2/i k N i k N N i nk N n k e e X e x nT X k NT i k NT i k NT ππππππ----=--==∑ （29） 其中[]X k 代表抽样信号[]()x n x nT =的N 点DFT 。最后令2/NT πΓ=，则上式变为 2/1()[]0,1,2,,12/i k N e X k X k k N i k NT ππ--Γ==-L （30） 首先用FFT 算法求出[]X k ，然后可用上式求出0,1,2,,1k N =-L 时的()X k Γ。 应该强调的是，式（28）只是一个近似表示，计算得到的()X ω只是一个近似值。通过取更小的抽样间隔T ，或者增加点数N ，可以得到更精确的值。如果B ω>时，幅度谱()X ω很小，对应于奈奎斯特抽样频率2s B ω=，抽样间隔T 选择/B π比较合适。如果已知信号只在时间区间10t t ≤≤内存在，可以通过对1nT t >时的抽样信号[]()x n x nT =补零，使N 足够大。 例1 利用FFT 计算傅里叶变换 一、实验目的 1.利用MATLAB 编写FFT 快速傅里叶变换。 2.比较编写的myfft 程序运算结果与MATLAB 中的FFT 的有无误差。 二、实验条件 PC 机，MATLAB7.0 三、实验原理 1. FFT （快速傅里叶变换）原理： 将一个N 点的计算分解为两个N/2点的计算，每个N/2点的计算再进一步分解为N/4点的计算，以此类推。根据DFT 的定义式，将信号x[n]根据采样号n 分解为偶采样点和奇采样点。设偶采样序列为y[n]=x[2n]，奇采样序列为z[n]=x[2n+1]。 上式中的k N W -为旋转因子N k j e /2π-。下式则为y[n]与z[n]的表达式： 2. 蝶形变换的原理： 下图给出了蝶形变换的运算流图，可由两个N/2点的FFT （Y[k]和Z[k]得出N 点FFT X[k]）。同理，每个N/2点的FFT 可以由两个N/4点的FFT 求得。按这种方法，该过程可延迟后推到2点的FFT 。 下图为N=8的分解过程。图中最右边的为8个时域采样点的8点FFTX[k]，由偶编号采样点的4点FFT 和奇编号采样点的4点得到。这4点偶编号又由偶编号的偶采 样点的2点FFT 和奇编号的偶采样点的2点FFT 产生。相同的4点奇编号也是如此。依次往左都可以用相同的方法算出，最后由偶编号的奇采样点和奇编号的偶采样点的2点FFT 算出。图中没2点FFT 成为蝶形，第一级需要每组一个蝶形的4组，第二级有每组两个蝶形的两组，最后一级需要一组4个蝶形。 四、实验内容 1.定义函数disbutterfly ，程序根据FFT 的定义：]2[][][N n x n x n y + +=、n N W N n x n x n z -+-=])2 [][(][，将序列x 分解为偶采样点y 和奇采样点z 。 function [y,z]=disbutterfly(x) N=length(x); n=0:N/2-1; w=exp(-2*1i*pi/N).^n; x1=x(n+1); x2=x(n+1+N/2); y=x1+x2; z=(x1-x2).*w; 2.定义函数rader ，纠正输出序列的输出顺序。 function y=rader(x,N) n=[0:N-1]; bn=dec2bin(n); rbn=fliplr(bn); rn=bin2dec(rbn); y=x(rn+1); 3.定义函数myfft ，程序中套了两个循环。 function X=myfft(x) N=length(x); h=log2(N); %h=3 for i=1:h %第一次i=1;第二次i=2 s=[]; for j=1:2^(i-1);%i=1时，j=1;i=2时,j=1:2 M=2^(h-i+1);%M:M=8；M=4 xj=x([1:M]+(j-1)*M);%xj=x([1:8]+(1-1)*8)=x(1)+x(2)...+x(8); %j=1:xj=x([1:4]);j=2:xj=x([1:4]+4) [y,z]=disbutterfly(xj); s=[s,y,z]; end x=s; 陕西科技大学实验报告 班级信工142 学号22 姓名何岩实验组别实验日期室温报告日期成绩报告内容：(目的和要求，原理，步骤，数据，计算，小结等) 1.求信号的离散时间傅立叶变换并分析其周期性和对称性； 给定正弦信号x(t)=2*cos(2*pi*10*t),fs=100HZ,求其DTFT。 (a)代码： f=10;T=1/f;w=-10:0.2:10; t1=0:0.0001:1;t2=0:0.01:1; n1=-2;n2=8;n0=0;n=n1:0.01:n2; x5=[n>=0.01]; x1=2*cos(2*f*pi*t1); x2=2*cos(2*f*pi*t2); x3=(exp(-j).^(t2'*w)); x4=x2*x3; subplot(2,2,1);plot(t1,x1); axis([0 1 1.1*min(x2) 1.1*max(x2)]); xlabel('x(n)');ylabel('x(n)'); title('原信号x1'); xlabel('t');ylabel('x1'); subplot(2,2,3);stem(t2,x2); axis([0 1 1.1*min(x2) 1.1*max(x2)]); title('原信号采样结果x2'); xlabel('t');ylabel('x2'); subplot(2,2,2);stem(n,x5); axis([0 1 1.1*min(x5) 1.1*max(x5)]); xlabel('n');ylabel('x2'); title('采样函数x2'); subplot(2,2,4);stem(t2,x4); axis([0 1 -0.2+1.1*min(x4) 1.1*max(x4)]); xlabel('t');ylabel('x4'); title('DTFT结果x4'); （b）结果： 实验四 离散信号的频域分析 一、 实验目的 1. 掌握序列的傅里叶变换、离散傅里叶级数、离散傅里叶变换、快速傅里叶变换的Matlab 实现； 2. 学习用FFT 对连续信号和离散信号进行谱分析的方法，了解可能出现的分析误差及其原因，以便在实际中正确应用FFT 。 二、 实验内容及步骤 1. 计算序列的DTFT 和DFT ，观察栅栏效应 设)()(4n R n x =，要求用MATLAB 实现： （1）计算)(n x 的傅里叶变换)(ωj e X ，并绘出其幅度谱； （2）分别计算)(n x 的4点DFT 和8点DFT ，绘出其幅度谱。并说明它们和)(ωj e X 的关系。 （提示：DFT 变换可用MA TLAB 提供的函数fft 实现，也可以自己用C 语言或matlab 编写） 2.计算序列的FFT ，观察频谱泄漏 已知周期为16的信号)1612cos()1610cos()(n n n x π π +=。 （1） 截取一个周期长度M=16点，计算其16点FFT ，并绘出其幅度谱； （2） 截取序列长度M=10点，计算其16点FFT ，绘出其幅度谱，并与（1）的结果进行比 较，观察频谱泄漏现象，说明产生频谱泄漏的原因。 三、 实验报告要求 1. 结合实验中所得给定典型序列幅频特性曲线，与理论结果比较，并分析说明误差产生的原因以及用FFT 作谱分析时有关参数的选择方法。 2. 总结实验所得主要结论。 1. 计算序列的DTFT 和DFT ，观察栅栏效应 设)()(4n R n x =，要求用MATLAB 实现： （1）计算)(n x 的傅里叶变换)(ωj e X ，并绘出其幅度谱； （2）分别计算)(n x 的4点DFT 和8点DFT ，绘出其幅度谱。并说明它们和)(ωj e X 的关系。 （1）代码： n=0:3; M=10; fs=51.2; N=1024; n=0:N-1; t=n/fs; x=0.5-0.5*sign(t-1); Y=fft(x,N); mag=abs(Y); Y1=fftshift(Y); mag1=abs(Y1); fn2=(-N/10.24:N/10.24)*fs/N; subplot(2,1,1) plot(fn2,mag1((N/2-N/10.24+1):(N/2+N/10.24+1))); set(gca,'XTick',(-5:0.5:5)); set(gca,'YTick',(0:10:60)); xlabel('频率/Hz');ylabel('振幅'); title('图1：矩形函数的FFT结果，N=512，fs=51.2Hz'); grid on; f=linspace(-5,5,1000); y=sqrt(2-2*cos(2*pi*f))./abs((2*pi*f)); subplot(2,1,2) plot(f,y); set(gca,'XTick',(-5:0.5:5)); set(gca,'YTick',(0:0.2:1)); xlabel('频率/Hz');ylabel('振幅'); title('图2：矩形函数傅里叶变换的理论结果'); grid on -5-4.5-4-3.5-3-2.5-2-1.5-1-0.500.51 1.52 2.53 3.54 4.55010203040 50 60 频率/Hz 振幅 图1：矩形函数的FFT 结果，N=512，fs=51.2Hz -5-4.5-4-3.5-3-2.5-2-1.5-1-0.500.51 1.52 2.53 3.54 4.55 00.2 0.4 0.60.8 1频率/Hz 振幅图2：矩形函数傅里叶变换的理论结果 function y=MyFFT_TB(x,n) %MYFFT_TB:My Fast Fourier Transform Time Based %按时间抽取基2-fft算法 %input: % x -- 输入的一维样本 % n -- 变换长度，缺省时n=length(x) 当n小于x数据长度时，x数据被截断到第n个数据% 当n大于时，x数据在尾部补0直到x 含n个数据 %output: % y -- 1*n的向量，快速傅里叶变换结果 %variable define: % N -- 一维数据x的长度 % xtem -- 临时储存x数据用 % m,M -- 对N进行分解N=2^m*M,M为不能被2整除的整数 % two_m -- 2^m % adr -- 变址，1*N的向量 % l -- 当前蝶形运算的级数 % W -- 长为N/2的向量，记录W(0,N),W(1,N),...W(N/2-1,N) % d -- 蝶形运算两点间距离 % t -- 第l级蝶形运算含有的奇偶数组的个数 % mul -- 标量，乘数 % ind1,ind2 -- 标量，下标 % tem -- 标量，用于临时储存 %参考文献： % https://www.doczj.com/doc/d92176962.html,/view/fea1e985b9d528ea81c779ee.html %% 输入参数个数检查 msg=nargchk(1,2,nargin); error(msg); %% 输入数据截断或加0 N=length(x); if nargin==2 if N M A t l a b傅里叶变换实 验报告 文档编制序号：[KKIDT-LLE0828-LLETD298-POI08] 班级信工142 学号 22 姓名何岩实验组别 实验日期室温报告日期成绩 报告内容：(目的和要求，原理，步骤，数据，计算，小结等) 1.求信号的离散时间傅立叶变换并分析其周期性和对称性； 给定正弦信号x(t)=2*cos(2*pi*10*t),fs=100HZ,求其DTFT。(a)代码： f=10;T=1/f;w=-10::10; t1=0::1;t2=0::1; n1=-2;n2=8;n0=0;n=n1::n2; x5=[n>=]; x1=2*cos(2*f*pi*t1); x2=2*cos(2*f*pi*t2); x3=(exp(-j).^(t2'*w)); x4=x2*x3; subplot(2,2,1);plot(t1,x1); axis([0 1 *min(x2) *max(x2)]); xlabel('x(n)');ylabel('x(n)'); title('原信号x1'); xlabel('t');ylabel('x1'); subplot(2,2,3);stem(t2,x2); axis([0 1 *min(x2) *max(x2)]); title('原信号采样结果x2'); xlabel('t');ylabel('x2'); subplot(2,2,2);stem(n,x5); axis([0 1 *min(x5) *max(x5)]); xlabel('n');ylabel('x2'); title('采样函数x2'); subplot(2,2,4);stem(t2,x4); axis([0 1 +*min(x4) *max(x4)]); xlabel('t');ylabel('x4'); title('DTFT结果x4'); （b）结果： 2.用以下两个有限长序列来验证DTFT的线性、卷积和共轭特性； 傅里叶变换m a t l a b代 码 Document serial number【UU89WT-UU98YT-UU8CB-UUUT-UUT108】 %傅里叶变换 clc;clear all;close all; tic Fs=128;%采样频率，频谱图的最大频率 T=1/Fs;%采样时间，原始信号的时间间隔 L=256;%原始信号的长度，即原始离散信号的点数 t=(0:L-1)*T;%原始信号的时间取值范围 x=7*cos(2*pi*15*t-pi)+3*cos(2*pi*40*t-90*pi/180)+3*cos(2*pi*30*t-90*pi/180); z=7*cos(2*pi*15*t-pi)+3*cos(2*pi*40*t-90*pi/180); z1=6*cos(2*pi*30*t-90*pi/180); z1(1:L/2)=0; z=z+z1; y=x;%+randn(size(t)); figure; plot(t,y) title('含噪信号') xlabel('时间（s）') hold on plot(t,z,'r--') N=2^nextpow2(L);%N为使2^N>=L的最小幂 Y=fft(y,N)/N*2; Z=fft(z,N)/N*2;%快速傅里叶变换之后每个点的幅值是直流信号以外的原始信号幅值的N/2倍（是直流信号的N倍） f=Fs/N*(0:N-1);%频谱图的频率取值范围 A=abs(Y);%幅值 A1=abs(Z); B=A; %让很小的数置零. B1=A1; A(A<10^-10)=0; % A1(A1<10^-10)=0; P=angle(Y).*A./B; P1=angle(Z).*A1./B1; P=unwrap(P,pi);%初相位值，以除去了振幅为零时的相位值 P1=unwrap(P1,pi); figure subplot(211) plot(f(1:N/2),A(1:N/2))%函数ffs返回值的数据结构具有对称性，因此只取前一半 hold on Fs = 1000; % Sampling frequency采样频率 T = 1/Fs; % Sample time采样周期 L = 1000; % Length of signal信号长度(点的个数) t = (0:L-1)*T; % Time vector时间向量(序列)(用来画图) % Sum of a 50 Hz sinusoid and a 120 Hz sinusoid一个50赫兹正弦加上120赫兹正弦x = 0.7*sin(2*pi*50*t) + sin(2*pi*120*t); y = x + 2*randn(size(t)); % Sinusoids plus noise正弦之和加上正态噪声 figure(1); plot(Fs*t(1:50),y(1:50)) title('Signal Corrupted with Zero-Mean Random Noise')%被零均值噪声祸害的信号xlabel('time (milliseconds)') %It is difficult to identify the frequency components by looking at the original signal. Converting to the frequency domain, the discrete Fourier transform of the noisy signal y is found by taking the fast Fourier transform (FFT): %从时域上直接看原信号难以确定各频率分量. 通过使用快速傅里叶变换FFT, 实现了含噪声信号Y的离散傅里叶变换,从而把信号转换到频域上(确定频率分量) NFFT = 2^nextpow2(L); % Next power of 2 from length of y Y = fft(y,NFFT)/L; f = Fs/2*linspace(0,1,NFFT/2+1); % Plot single-sided amplitude spectrum. figure(2); plot(f,2*abs(Y(1:NFFT/2+1))) title('Single-Sided Amplitude Spectrum of y(t)') xlabel('Frequency (Hz)') ylabel('|Y(f)|') A=xlsread('SHUJU2.xlsx') ; %这里的SHUJU1.xlsx是数据 FS=1000; % FS是采样率 T=1/FS; L=length(A) t=(0:L-1)*T; figure(3); 应用MATLAB对信号进行频谱分析及滤波 设计目的 要求学生会用MATLAB语言进行编程，绘出所求波形，并且运用FFT求对连续信号进行分析。 一、设计要求 1、用Matlab产生正弦波，矩形波，并显示各自的时域波形图； 2、进行FFT变换，显示各自频谱图，其中采样率、频率、数据长度自选，要求注明； 3、绘制三种信号的均方根图谱； 4、用IFFT回复信号，并显示恢复的正弦信号时域波形图。 二、系统原理 用FFT对信号作频谱分析是学习数字信号处理的重要内容。经常需要进行频谱分析的信号是模拟信号和时域离散信号。频谱分辨率直接和FFT的变换区间N有关，因为FFT能够实现频率分辨率是2π/N。 x(n)是一个长度为M的有限长序列，则x(n)的N点离散傅立叶变换为： N?1?2?kn)(nx j?W W NN e?0?n N X(k)=DFT[x(n)]=,k=0,1,...,N-1N?11?kn?)(WXk N N0?n x(n) =IDFT[X(k)]= 逆变换：,k=0,1,...,N-1 但FFT是一种比DFT更加快速的一种算法，提高了DFT的运算速率,为数字信号处理技术应用于各种信号处理创造了条件，大大提高了数字信号处理技术的发展。本实验就是采用FFT，IFFT对信号进行谱分析。 三、程序设计 fs=input('please input the fs:');%设定采样频率 N=input('please input the N:');%设定数据长度 t=0:0.001:1; f=100;%设定正弦信号频率 %生成正弦信号 x=sin(2*pi*f*t); figure(1); subplot(211); plot(t,x);%作正弦信号的时域波形 axis([0,0.1,-1,1]); title('正弦信号时域波形'); z=square(50*t); subplot(212) plot(t,z) axis([0,1,-2,2]); title('方波信号时域波形');grid; 电子1004 黄佳卿201081163 图像傅里叶变换代码如下： I = imread('C:\Users\jiaqing\Desktop\matlab homework\zhen.jpg'); % 读入图像 subplot(2,2,1),imshow(I); title('zhen'); FI = abs((fft2(I))); NFI = 255*mat2gray(FI); %归一化 SFI = fftshift(NFI); imgray = rgb2gray(SFI); %灰阶 subplot(2,2,2),imshow(imgray); title('fft of zhen'); I = imread('C:\Users\jiaqing\Desktop\matlab homework\cat.jpg'); % 读入图像 subplot(2,2,3),imshow(I); title('cat'); FI = abs((fft2(I))); NFI = 255*mat2gray(FI); %归一化 SFI = fftshift(NFI); imgray = rgb2gray(SFI); %灰阶 subplot(2,2,4),imshow(imgray); title('fft of cat'); 运行结果如下： 图像添加高斯噪声与去除代码如下： I = imread('C:\Users\jiaqing\Desktop\matlab homework\zhen.jpg'); % 读入图像 %给图像添加噪声 grayI = rgb2gray(I); nI=imnoise(grayI,'gaussian',0.05); subplot(2,2,1);imshow(nI); title('加入高斯噪声后的图像'); %加入噪声后图像的傅立叶变换 FI = abs((fft2(nI))); NFI = 255*mat2gray(FI); %归一化 SFI = fftshift(NFI); subplot(2,2,2),imshow(SFI); title('加入噪声后的傅里叶变换'); xd=wiener2(nI,[3 3]); subplot(2,2,3);imshow(xd); %去除噪声后图像的傅立叶变换按频率抽取基2-快速傅里叶逆变换算法_MATLAB代码
短时傅里叶变换matlab程序
Matlab傅里叶变换傅里叶逆变换-FFT-IFFT
傅里叶变换matlab代码
matlab-离散信号傅里叶变换
MATLAB数字图像处理几何变换傅里叶变换
【免费下载】matlab实现傅里叶变换
用Matlab对信号进行傅里叶变换实例
傅里叶变换的应用,matlab程序,C语言程序
利用MATLAB编写FFT快速傅里叶变换
MAtlab-傅里叶变换-实验报告
离散信号变换的matlab实现
傅里叶变换MATLAB程序
快速傅里叶变换_蝶形运算_按时间抽取基2-fft算法_MATLAB代码
MAtlab傅里叶变换实验报告
傅里叶变换matlab代码
傅里叶变换matlab程序
MATLAB的离散傅里叶变换的仿真
matlab图像傅里叶变换图像
相关主题
文本预览