中科大信号与系统实验报告3.docx
《中科大信号与系统实验报告3.docx》由会员分享,可在线阅读,更多相关《中科大信号与系统实验报告3.docx(18页珍藏版)》请在冰豆网上搜索。
中科大信号与系统实验报告3
信号与系统
实验报告
学号:
姓名:
信息科学技术学院
电子科学与技术系
一、实验目的
1.熟悉傅里叶变换的性质
2.熟悉常见信号的傅里叶变换
3.了解傅里叶变换的MATLAB实现方法
二、实验原理
傅里叶变换是信号分析的最重要的内容之一。
从已知信号
求出相应的频谱函数
的数学表示为:
的傅里叶变换存在的充分条件是
在无限区间内绝对可积,即
满足下式:
但上式并非傅里叶变换存在的必要条件。
在引入广义函数概念之后,使一些不满足绝对可积条件的函数也能进行傅里叶变换。
傅里叶反变换的定义为:
。
在这一部分的学习中,大家都体会到了这种数学运算的麻烦。
在MATLAB语言中有专门对信号进行正反傅里叶变换的语句,使得傅里叶变换很容易在MATLAB中实现。
在MATLAB中实现傅里叶变换的方法有两种,一种是利用MATLAB中的SymbolicMathToolbox提供的专用函数直接求解函数的傅里叶变换和傅里叶反变换,另一种是傅里叶变换的数值计算实现法。
下面分别介绍这两种实现方法的原理。
1.直接调用专用函数法
在MATLAB中实现傅里叶变换的函数为:
F=fourier(f) 对f(t)进行傅里叶变换,其结果为F(w)
F=fourier(f,v) 对f(t)进行傅里叶变换,其结果为F(v)
F=fourier(f,u,v) 对f(u)进行傅里叶变换,其结果为F(v)
傅里叶反变换
f=ifourier(F) 对F(w)进行傅里叶反变换,其结果为f(x)
f=ifourier(F,U) 对F(w)进行傅里叶反变换,其结果为f(u)
f=ifourier(F,v,u) 对F(v)进行傅里叶反变换,其结果为f(u)
由于MATLAB中函数类型非常丰富,要想了解函数的意义和用法,可以用mhelp命令。
如在命令窗口键入:
mhelpfourier回车,则会得到fourier的意义和用法。
注意事项
(1)在调用函数fourier()及ifourier()之前,要用syms命令对所有需要用到的变量(如t,u,v,w)等进行说明,即要将这些变量说明成符号变量。
对fourier()中的f及ifourier()中的F也要用符号定义符sym将其说明为符号表达式。
(2)采用fourier()及fourier()得到的返回函数,仍然为符号表达式。
在对其作图时要用ezplot()函数,而不能用plot()函数。
(3)fourier()及fourier()函数的应用有很多局限性,如果在返回函数中含有δ(ω)等函数,则ezplot()函数也无法作出图来。
另外,在用fourier()函数对某些信号进行变换时,其返回函数如果包含一些不能直接表达的式子,则此时当然也就无法作图了。
这是fourier()函数的一个局限。
另一个局限是在很多场合,尽管原时间信号f(t)是连续的,但却不能表示成符号表达式,此时只能应用下面介绍的数值计算法来进行傅氏变换了,当然,大多数情况下,用数值计算法所求的频谱函数只是一种近似值。
例①求门函数
的傅里叶变换,并画出幅度频谱图
MATLAB程序如下:
symstw %定义两个符号变量t,w
Gt=sym('Heaviside(t+1)-Heaviside(t-1)'); %产生门宽为2的门函数
Fw=fourier(Gt,t,w); %对门函数作傅氏变换求F(jw)
FFw=maple('convert',Fw,'piecewise'); %数据类型转换,转为分段函数,此处可以去掉
FFP=abs(FFw); %求振幅频谱|F(jw)|
ezplot(FFP,[-10*pi10*pi]);grid; %绘制函数图形,并加网格
axis([-10*pi10*pi02.2]) %限定坐标轴范围
运行结果:
Fw=exp(i*w)*(pi*Dirac(w)-i/w)-exp(-i*w)*(pi*Dirac(w)-i/w)
%Dirac(w)为δ(ω),即傅立叶变换结果中含有奇异函数,故绘图前需作函数类型转换
FFw=-i*exp(i*w)/w+i*exp(-i*w)/w%FFw为复数
FFP=abs(-i*exp(i*w)/w+i*exp(-i*w)/w)%求FFw的模值
例②求函数
的傅里叶反变换f(t)
MATLAB程序如下:
symstw %定义两个符号变量t,w
Fw=sym('1/(1+w^2)'); %定义频谱函数F(jw)
ft=ifourier(Fw,w,t); %对频谱函数F(jw)进行傅氏反变换
运行结果:
ft=
1/2*exp(-t)*Heaviside(t)+1/2*exp(t)*Heaviside(-t)
2、傅里叶变换的数值计算实现法
严格说来,如果不使用symbolic工具箱,是不能分析连续时间信号的。
采用数值计算方法实现连续时间信号的傅里叶变换,实质上只是借助于MATLAB的强大数值计算功能,特别是其强大的矩阵运算能力而进行的一种近似计算。
傅里叶变换的数值计算实现法的原理如下:
对于连续时间信号f(t),其傅里叶变换为:
其中τ为取样间隔,如果f(t)是时限信号,或者当|t|大于某个给定值时,f(t)的值已经衰减得很厉害,可以近似地看成是时限信号,则上式中的n取值就是有限的,假定为N,有:
若对频率变量ω进行取样,得:
通常取:
,其中
是要取的频率范围,或信号的频带宽度。
采用MATLAB实现上式时,其要点是要生成f(t)的N个样本值
的向量,以及向量
,两向量的内积(即两矩阵的乘积),结果即完成上式的傅里叶变换的数值计算。
注意事项
时间取样间隔τ的确定,其依据是τ必须小于奈奎斯特(Nyquist)取样间隔。
如果f(t)不是严格的带限信号,则可以根据实际计算的精度要求来确定一个适当的频率
为信号的带宽。
例③用数值计算法实现上面门函数
的傅里叶变换,并画出幅度频谱图.
分析:
该信号的频谱为
,其第一个过零点频率为π,一般将此频率认为是信号的带宽。
但考虑到
的形状(为抽样函数),假如将精度提高到该值的50倍,即取
,则据此确定的Nyquist取样间隔为:
。
MATLAB程序如下:
R=0.02; %取样间隔τ=0.02
t=-2:
R:
2; %t为从-2到2,间隔为0.02的行向量,有201个样本点
ft=[zeros(1,50),ones(1,101),zeros(1,50)]; %产生f(t)的样值矩阵(即f(t)的样本值组成的行向量)
W1=10*pi; %取要计算的频率范围
M=500;k=0:
M;w=k*W1/M; %频域采样数为M,w为频率正半轴的采样点
Fw=ft*exp(-j*t'*w)*R; %求傅氏变换F(jw)
FRw=abs(Fw); %取振幅
W=[-fliplr(w),w(2:
501)]; %由信号双边频谱的偶对称性,利用fliplr(w)形成负半轴的点,%w(2:
501)为正半轴的点,函数fliplr(w)对矩阵w行向量作180度反转
FW=[fliplr(FRw),FRw(2:
501)]; %形成对应于2M+1个频率点的值
Subplot(2,1,1);plot(t,ft);grid; %画出原时间函数f(t)的波形,并加网格
xlabel('t');ylabel('f(t)'); %坐标轴标注
title('f(t)=u(t+1)-u(t-1)'); %文本标注
subplot(2,1,2);plot(W,FW);gridon; %画出振幅频谱的波形,并加网格
xlabel('W');ylabel('F(W)'); %坐标轴标注
title('f(t)的振幅频谱图'); %文本标注
运行结果如下:
三、实验内容
1.编程实现求下列信号的幅度频谱
【实验步骤】:
1.将所需表达式通过符号表达式表达出来
2.用matlab实验输出
【实验原理】:
1.几种绘图函数总结:
plot:
绘制连续函数,可以表示成显函数的函数
fplot:
绘制连续的隐函数,符号表达式的函数
stem:
绘制离散函数
(1) 求出
的频谱函数F1(jω),请将它与上面门宽为2的门函数
的频谱进行比较,观察两者的特点,说明两者的关系。
【理论值的计算】:
Fw1=(cos(w/2)*1i+sin(w/2))/w-(cos(w/2)*1i-sin(w/2))/w
Fw2=-(cos(w)*1i-sin(w))/w+(sin(w)+cos(w)*1i)/w
【程序如下】:
图片如下:
【绘图如下】:
【结论】:
通过观察可以知道,f2(t)是对f1(t)的一个x轴的2倍拉伸,这种变换使得Fw2相对于Fw1是一种1/2倍的x轴压缩效应。
(2)三角脉冲
【理论值的计算】:
Fw=2/w^2-(exp(-w*1i)*1i)/w+(exp(w*1i)*1i)/w-exp(-w*1i)/w^2-exp(w*1i)/w^2+(sin(w)+cos(w)*1i)/w-(cos(w)*1i-sin(w))/w
【程序如下】:
【作图如下】:
(3)单边指数信号
【理论值的计算】:
Fw=1/(1+w*1i)
【程序如下】:
【程序注意点】:
1.在使用fourier进行傅立叶变换时,自变量一定要是符号函数
2.matlab可以通过subplot把一个窗口分成N*M个窗口,编号从左到右,从上到下
3.fplot可以通过在函数后面加[a,b]实现图像生成的范围
【作图如下】:
(4)高斯信号
【理论值的计算】:
Fw=pi^(1/2)*exp(-w^2/4)
两种实现方法:
方法一:
【程序的数值法实现】
方法二:
【程序的库函数实现】
【作图如下】:
2.利用ifourier()函数求下列频谱函数的傅氏反变换
【实验步骤】:
1.将所需表达式通过符号表达式表达出来
2.将表达式进行反傅立叶变换
3.用matlab实验输出
(1)
【理论值的计算】:
ft=exp(-4*abs(t))*sign(t)
【程序如下】:
【作图如下】:
(2)
【理论值的计算】:
ft=(2*pi*dirac(t)-3*pi*exp(-t)*(sign(t)+1)+2*pi*exp(-5*t)*(sign(t)+1)-(pi*exp(-t)*dirac(t))/2+(pi*exp(-5*t)*dirac(t))/2)/(2*pi)
【程序如下】:
【作图如下】:
四、实验收获
1.学会了matlab中傅立叶变换的两种方法:
直接调用库函数和使用定义计算
2.熟悉了matlab中求反傅立叶函数的方法
3.了解了原函数拉伸变换后傅立叶变换后的函数的变换
4.熟悉了各种基本函数的傅立叶变换
5.了解了fplot,ezplot,plot间的差别和使用
【二维图像matalb实现总结】:
1.plot函数
plot函数有很多重载方法,这里只做简单的介绍
1.1plot(Y)
1.1若Y是向量,绘制向量Y对其索引值的曲线。
1.2若Y是实数矩阵,绘制矩阵的每列对应于行数的曲线集合。
1.3若Y是复数矩阵,等价于plot(real(Y),imag(Y));
1.2plot(X,Y)
1.2.1若X,Y,均为向量,绘制向量Y对应向量X的曲线此时X的长度跟Y的长度必须相等。
1.2.2若X为向量Y为矩阵,则X的长度与矩阵Y的行数或列数必须相等:
1.2.2.1X的长度与Y的列数相等,或X的长度与矩阵Y的行数和列数均相等(Y为方阵),绘制矩阵Y的每列对应X的曲线集合
1.2.2.2X的长度与Y的行数相等,绘制矩阵Y的每行对应X的曲线集合
1.2.2若X为矩阵Y为向量,则Y的长度与矩阵X的行数或列数必须相等,绘制方法与1.2.1类似,不在过多介绍
1.2.3若X,Y均为矩阵,则X,Y的大小必须相等,绘制矩阵Y的每列对应X的每列的曲线。
P.S如果矩阵是复数矩阵,会自动忽略掉复数的虚部
1.3plot(X,Y,LineSpec)
绘制Y对应于X的曲线集合,并指定曲线的LineSpec,比如线型,标记符号,和颜色或其任意组合,API上关于LineSpec的讲解很详细。
1.4plot(X,Y,‘属性’,属性值)
跟1.3相比,只是属性值得不同。
API上有很多
P.S1.有一种调用的方法h=plot(x....)返回和曲线的句柄,其实感觉就像是指针
2.semilogxsemilogy和loglog等函数的使用方法和plot基本类似,只是在曲线的外观上有所不同
2.plotyy函数
这个函数又叫双Y轴函数(上面的例子都是单Y轴函数,只有左侧一个Y轴)。
plotyy是为了满足:
对函数值变化范围较大的两组数据同事绘图(这个时候用holdon,会很难从图形中辨识函数值变化范围较小的那组数据变化趋势的细节信息的,因为坐标系太大,所以才有这样的双坐标轴绘图),
3.polar函数
前面的函数都是在直角坐标系的绘图的,但是有时候需要在极坐标或者柱坐标中绘图,于是就有了poltar(极坐标绘图函数)
>>polar(theta,rho,LineSpec)其中theta是极脚,rho是极径,其他的方法跟其他绘图函数基本相同。
4.fplot函数
如果不太了解某个函数随自变量的变化趋势,随便的就取定自变量的范围,很可能就用为自变量取值范围不好而是的绘制出的图片失真,为解决此问题可以使用fplot。
这个函数,据说可以通过其内部的自适应算法动态的决定自变量的间隔。
比如,函数值得变化较为剧烈,那么自变量的取值间隔就小。
从而保证绘制出的图形的质量和效率。
绘制出的图形的质量和效率。
基本的调用函数形式为:
fplot(fun,limits)其中,fun为要绘制的函数,fun可以为可执行字符串,M文件,inline或者匿名函数。
fun为要绘制的函数,limits为指定的范围,可以使二维向量平[xmin,xmax]或四维向量[xmin,xmax,ymin,ymax]
5.ezplot函数
ezplot函数可以直接绘制一元函数如y=f(x)参数方程
构成的函数y=f(x),以及隐函数f(x,y)=0的图形。
它的调用跟fplot类似,基本的调用方法是:
ezplot(fun,limits)参见4.
6.子图
为了突出的对比几个相似图形,一般可以使用子图(subplot),这样就可以在一个figure中做若干个图。
使用方法:
subplot(m,n,p)
其中,m表示是图排成m行,n表示图排成n列,也就是整个figure中有n个图是排成一列的,一共m行,如果m=2就是表示2行图。
p表示图所在的位置,p=1表示从左到右从上到下的第一个位置。
网上找了个例子:
7.交互式绘图
交互式绘图的常用方法是:
ginput,gtext
在绘图前调用,然后在调用绘图函数,即可。
比如用ginput可以方便的通过鼠标来读取二维平面图形的任意一个点的坐标值。
当调用时,如ginput,会是当前的图形从后台调到前台,然后咱们可以移动鼠标点击想要选取的点,完成后(到可N值,或按下space),会在命令窗口看见点的坐标