傅里叶级数连续时间信号分析MATLAB课程

1、课程设计任务书学生姓名： 专业班级： 指导教师： 工作单位： 题 目: 连续时间信号的傅利叶变换及MATLAB实现 初始条件：MATLAB软件，微机要求完成的主要任务: 利用MATLAB强大的图形处理功能，符号运算功能和数值计算功能，实现连续时间非周期信号频域分析的仿真波形；1、用MATLAB实现典型非周期信号的频域分析；2、用MATLAB实现信号的幅度调制；3、用MATLAB实现信号傅立叶变换性质的仿真波形；4、写出课程设计报告。时间安排：学习MATLAB语言的概况 第1天学习MATLAB语言的基本知识 第2天学习MATLAB语言的应用环境，调试命令，绘图能力 第3、4天课程设计 第5-9天

2、答辩 第10天指导教师签名： 年 月 日系主任（或责任教师）签名： 年 月 日目 录摘要IABSTRACTI绪论I1傅里叶变换原理概述 11.1 傅里叶变换及逆变换的MATLAB实现 22 用MATLAB实现典型非周期信号的频域分析 32.1 单边指数信号时域波形图、频域图 32.2 偶双边指数信号时域波形图、频域图 42.3 奇双边指数信号时域波形图、频域图 42.4 直流信号时域波形图、频域图 52.5 符号函数信号时域波形图、频域图 52.6 单位阶跃信号时域波形图、频域图 62.7 单位冲激信号时域波形图、频域图 62.8 门函数信号时域波形图、频域图 73 用MATLAB实现信号的幅

3、度调制 83.1 实例1 83.2 实例2 104 实现傅里叶变换性质的波形仿真114.1 尺度变换特性 114.2 时移特性 144.3 频移特性 164.4 时域卷积定理 18 4.5 对称性质 204.6 微分特性 22心得体会25参考文献26附录27 摘要MATLAB和Mathematica、Maple并称为三大数学软件。MATLAB在数学类科技应用软件中在数值计算方面首屈一指。Simulink是MATLAB软件的扩展，它是实现动态系统建模和仿真的一个软件包。MATLAB具有强大的图形处理功能、符号运算功能和数值计算功能。其中系统的仿真（Simulink）工具箱是从底层开发的一个完整的

4、仿真环境和图形界面。在这个环境中，用户可以完成面向框图系统仿真的全部过程，并且更加直观和准确地达到仿真的目标1。本文主要介绍基于MATLAB的一阶动态电路特性分析。关键字：MATLAB；仿真；图形处理；一阶动态电路。AbstractMATLAB, and Mathematica, Maple, and known as the three major mathematical software. It is the application of technology in mathematics classes in numerical computing software, second t

5、o none. Simulink is an extension of MATLAB software, which is the realization of dynamic system modeling and simulation of a package. MATLAB has a powerful graphics processing capabilities, symbolic computing and numerical computing functions. One system simulation (Simulink) toolbox from the bottom

6、 of the development of a complete simulation environment and the graphical interface. In this environment, the user can complete system simulation block diagram for the entire process and achieve a more intuitive and accurate simulation of goal1.In this paper, MATLAB-based first-order characteristic

7、s of dynamic circuits.Keywords: MATLAB；Simulation；Graphics；First Order Circuit。 绪论 在科学技术飞速发展的今天，计算机正扮演着愈来愈重要的角色。在进行科学研究与工程应用的过程中，科技人员往往会遇到大量繁重的数学运算和数值分析，传统的高级语言Basic、Fortran 及C 语言等虽然能在一定程度上减轻计算量，但它们均用人员具有较强的编程能力和对算法有深入的研究。MATLAB 正是在这一应用要求背景下产生的数学类科技应用软件。MATLAB 是matrix 和laboratory 前三个字母的缩写，意思是“矩阵实验室”

8、，是Math Works 公司推出的数学类科技应用软件2。MATLAB 具有以下基本功能：（1）数值计算功能；（2）符号计算功能；（3）图形处理及可视化功能；（3）可视化建模及动态仿真功能。本文介绍了如何利用MATLAB强大的图形处理功能、符号运算功能以及数值计算功能，实现连续时间系统频域分析。本次课程设计介绍了用MATLAB实现典型非周期信号的频谱分析，用MATLAB实现信号的幅度调制以及用MATLAB实现信号傅里叶变换性质的仿真波形。291傅里叶变换原理概述 设有连续时间周期信号，它的周期为T，角频率，且满足狄里赫利条件，则该周期信号可以展开成傅里叶级数，即可表示为一系列不同频率的正弦或复

9、指数信号之和。傅里叶级数有三角形式和指数形式两种3。 1. 三角形式的傅里叶级数2： 式中系数，称为傅里叶系数，可由下式求得： 2. 指数形式的傅里叶级数2： 式中系数称为傅里叶复系数，可由下式求得： 周期信号频谱具有三个特点1：（1） 离散性，即谱线是离散的；（2） 谐波性，即谱线只出现在基波频率的整数倍上；（3） 收敛性，即谐波的幅度随谐波次数的增高而减小。 周期信号的傅里叶分解用Matlab进行计算时，本质上是对信号进行数值积分运算。在Matlab中有多种进行数值积分运算的方法，我们采用quadl函数，它有两种其调用形式。 (1) yquadl(func, a, b)。 其中func是一

10、个字符串，表示被积函数的.m文件名（函数名）；a、b分别表示定积分的下限和上限。 (2) yquadl(myfun, a, b)。其中“”符号表示取函数的句柄，myfun表示所定义函数的文件名。1.1 傅里叶变换及逆变换的MATLAB实现MATLAB 的Symbolic Math Toolbox 提供了能直接求解傅里叶变换及逆变换的函数Fourier()及Fourier()4。 1.1 fourier 变换1. (1) F=fourier(f); (2) F=fourier(v); (3) F=fourier(f,u,v);说明：(1) F=fourier(f)是符号函数f 的Fourier

11、变换，缺省返回是关于的函数。如果 f=f(),则fourier 函数返回关于t 的函数。 (2)F=fourier(f,v)返回函数F 是关于符号对象v 的函数，而不是缺省的(3)F=fourier(f,u,v)对关于u 的函数f 进行变换，返回函数F 是关于v 的函数。 1.2 fourier 逆变换 1. (1) f=ifourier(F); (2) f=ifourier(F,u); (3) f=ifourier(F,v,u);说明：(1) f=ifourier(F)中输入参量F是傅里叶变换的符号表达式，缺省为符号变量w的函数，输出参量f是F的傅里叶逆变换的符号表达式，缺省为符号变量x的函

12、数。 (2)f=ifourier(F,u)中输入参量F是傅里叶变换的符号表达式，缺省为符号变量w的函数，输出参量f是F的傅里叶逆变换的符号表达式，为指定符号变量u的函数(3)f=ifourier(F,v,u)中输入参量F是傅里叶变换的符号表达式，为指定符号变量v的函数，输出参量f是F的傅里叶逆变换的符号表达式，缺省为符号变量u的函数。2 用MATLAB实现典型非周期信号的频域分析2.1单边指数信号时域波形图、频域图 的时域波形图和频谱图如图2.1.1 ：图2.1.1单边指数信号2.2偶双边指数信号时域波形图、频域图偶双边指数信号时域波形图、频域图如下图图2.2.1：图2.2.1偶双边指数信号2

13、.3奇双边指数信号时域波形图、频域图奇双边指数信号时域波形图、频域图如下图图2.3.1：图2.3.1奇双边指数信号2.4 直流信号时域波形图、频域图直流信号f（t）=A,不满足绝对可积条件，但傅里叶变换却存在。可以把单位直流信号看做双边指数信号当a趋于0时的极限。直流信号时域波形图、频域图如下图2.4.1：图2.4.1直流信号2.5 符号函数信号时域波形图、频域图符号函数信号时域波形图、频域图如下图2.5.1：图2.5.1符号函数信号图5符号函数信号波形图2.6 单位阶跃信号时域波形图、频域图单位阶跃函数信号时域波形图、频域图如下图2.6.1：图2.6.1单位阶跃函数信号2.7 单位冲激信号时

14、域波形图、频域图单位冲激函数信号时域波形图、频域图如下图2.7.1：图2.7.1单位冲激函数信号2.8 门函数信号时域波形图、频域图门函数信号时域波形图、频域图如下图2.8.1：图2.8.1门函数信号3用MATLAB实现信号的幅度调制设信号f (t) 的频谱为F( jw) ,现将f (t) 乘以载波信号cos (w0t) ，得到高频的已调信号y(t ) ，即：y(t ) = f (t) cos (w0t)从频域上看，已调制信号y(t ) 的频谱为原调制信号f (t) 的频谱搬移到0 ±w 处，幅度降为原F( jw) 的1/2，即上式即为调制定理，也是傅里叶变换性质中“频移特性”的一种

15、特别情形。MATLAB 提供了专门的函数modulate（）用于实现信号的调制。调用格式为：y=modulate(x,Fc,Fs,'method')y,t=modulate(x,Fc,Fs)其中，x 为被调信号，Fc 为载波频率，Fs 为信号x 的采样频率，method 为所采用的调制方式，若采用幅度调制、双边带调制、抑制载波调制，则method为am或amdsd-sc。其执行算法为y=x*cos(2*pi*Fc*t)其中y 为已调制信号，t 为函数计算时间间隔向量。涉及到一个函数，暂时不容易理解，因此查阅工具书，特在此说明：MATLAB 的“信号处理工具箱函数”中的估计信号的

16、功率谱密度函数psd(),其格式是：Pxx，f=psd(x，Nfft，Fs，window，noverlap，dflag)其中，x 是被调制信号（即本例中的f (t) ），Nfft 指定快速付氏变换FFT 的长度，Fs为对信号x 的采样频率。后面三个参数的意义涉及到信号处理的更深的知识，在此暂不介绍。3.1 信号调制实例1例1： f(t)=sin(100t) f=400Hz，绘出原信号f(t)以及调制信号y(t)=f(t)coswt的实域波形图、频铺图以及功率谱。程序如下：Fm=50;Fc=400;%载波频率Fs=1000;%信号x的抽样频率N=1000;k=0:N-2;%采样点t=k/Fs;%

17、采样时间x=sin(2*pi*Fm*t);%被调制信号subplot(221);plot(t,x);%画出被调制信号的波形xlabel('t(s)');ylabel('x');title('被调制信号的波形');axis(0 0.1 -1 1);%坐标系范围 t取值范围不能大，因为采样频率很高，不便于观察Nfft=1024;window=hamming(512);noverlap=256;dflag='none'Pxx,f=psd(x,Nfft,Fs,window,noverlap,dflag);%估算被调信号的功率谱密度Nfft

18、是快速傅里叶变换的长度subplot(222);plot(f,Pxx);%画出被调信号的功率谱密度的波形ylabel('功率谱(x)');xlabel('f(hz)');%axis(0 600 0 100);%坐标系的范围title('被调信号的功率谱密');gridy=modulate(x,Fc,Fs,'am');%得到调制信号subplot(223);plot(t,y);%会出调制信号的波形xlabel('t(s)');ylabel('y');title('已调信号');axis

19、(0 0.1 -1 1);%坐标系的范围t取值范围不能大，因为采样频率很高，不便于观察Pxx,f=psd(y,Nfft,Fs,window,noverlap,dflag);%估算被调信号的功率谱密度Nfft是快速傅里叶变换的长度subplot(224);plot(f,Pxx);%画出被调信号的功率谱密度的波形ylabel('功率谱(y)');xlabel('f(hz)');%axis(0 600 0 100);%坐标系的范围title('已调信号功率谱');grid图3.1.1调制信号与被调信号3.2 信号调制实例2例2：设 ,绘出原信号f(t)

20、以及调制信号y(t)=f(t)coswt的实域波形图、频谱图以及功率谱。解：n=0.005;t=-1.5:n:1.5;f=Heaviside(t+1)-Heaviside(t-1);ft=f.*cos(10*pi*t);%FT为已调信号，要满足矩阵相乘规则，点乘，.wsubplot(221);plot(t,f);%画出被调制信号波形xlabel('t');ylabel('f(t)');title('被调制信号波形');subplot(222);plot(t,ft);%画出已调制信号波形xlabel('t');ylabel('

21、;ft(t)');title('已调制信号波形');w1=40;N=1000;k=-N:N;w=w1*k/N;Fw=f*exp(-j*t'*w)*n;%得到被调制信号频谱Ftw=ft*exp(-j*t'*w)*n;%得到已调制信号频谱Fwr=real(Fw);%热挡?Ftwr=real(Ftw);subplot(223);plot(w,Fwr);%画出被调制信号频谱xlabel('w');ylabel('F(jw)');title('被调制信号频谱');subplot(224);plot(w,Ftwr);%

22、画出已调制信号频谱xlabel('w');ylabel('Ft(jw)');title('已调信号频谱');图3.2.1 原信号f (t) 、调制信号ft( t) 的波形及其频谱F( jw) 、Ft( jw)4用MATLAB实现信号傅立叶变换性质的仿真波形4.1 傅里叶变换的尺度变换特性若f (t) « F( jw) ，则傅里叶变换的尺度变换特性为5：例1: 设 ，即门宽为=2 的门信号，用MATLAB 求 的频谱Y ( jw) ，并与f (t) 的频谱F( jw) 进行比较。 %尺度变换n=0.02;%采样间隔t=-2:n:2;%采样

23、范围f=Heaviside(t+1)-Heaviside(t-1);%脉宽为2的门信号h=Heaviside(2*t+1)-Heaviside(2*t-1);%脉宽为1的门信号w1=5*2*pi;N=500;k=-N:N;w=k*w1/N;F=f*exp(-j*t'*w)*n;%求出FwH=h*exp(-j*t'*w)*n;%求出Hwsubplot(221);plot(t,f);%画出脉宽为2的门信号的时域波形xlabel('t');ylabel('f(t)');title('脉宽为2的门信号的时域波形');axis(-2.5 2

24、.5 0 1.1);subplot(222);plot(t,h);%画出脉宽为1的门信号的时域波xlabel('t');ylabel('h(t)');title('脉宽为1的门信号的时域波形');axis(-2.5 2.5 0 1.1);subplot(223);plot(w,F);%画出脉宽为2的门信号的频域波xlabel('w');ylabel('F(w)');title('脉宽为2的门信号的频域波形');axis(-5*pi 5*pi -0.5 2.1);subplot(224);plot(w

25、,H);%画出脉宽为1的门信号的频域波xlabel('w');ylabel('H(w)');title('脉宽为1的门信号的频域波形');axis(-5*pi 5*pi -0.5 2.1);图4.1.1 傅里叶变换的尺度变换特性由图4.1.1，y(t ) 信号相当于原信号f (t) 在时域上压缩一倍，即y(t ) = f (2t) ，a = 2 ，按式，Y ( jw) 的频域宽度应是F( jw) 的两倍，而幅度下降为F( jw) 的一半。4.2 傅里叶变换的时移变换特性若f (t) « F( jw) ，则傅里叶变换的时移特性为：例2：

26、设,试用MATLAB 绘出f (t-t0),f (t+t0) 及其频谱（幅度谱及相位谱）。t0=0.2;%时移大小n=0.02;%采样间隔t=-5:n:5;%采样范围f1=1/2*exp(-2*t).*Heaviside(t);%定义函数f1f2=1/2*exp(-2*(t-t0).*Heaviside(t-t0);%定义函数f2，时域右移t0f3=1/2*exp(-2*(t+t0).*Heaviside(t+t0);%定义函数f1，时域左移t0subplot(311);plot(t,f1);%画出f1,f2,f3的时域波形xlabel('t');ylabel('f(t

27、)');hold onplot(t,f2,'-.');plot(t,f3,':');axis(-6 6 0 0.7);legend('f1(t)','f2(t)','f3(t)');标注f1f2f3N=300;w1=5*pi*2;%频谱宽度k=-N:N;w=k*w1/N;F1=n*f1*exp(-j*t'*w);%f1的傅里叶变换F1f=abs(F1);%f1的幅度频谱F1a=angle(F1);%f1的相位频谱F2=n*f2*exp(-j*t'*w);%f2的傅里叶变换F2f=abs(F2

28、);%f2的幅度频谱F2a=angle(F2);%f2的相位频谱F3=n*f3*exp(-j*t'*w);%f3的傅里叶变换F3f=abs(F3);%f3的幅度频谱F3a=angle(F3);%f3的相位频谱subplot(312);plot(w,F1f);xlabel('w');ylabel('F(w)');hold onplot(w,F2f,'-.');plot(w,F3f,':');%画出f1f2f3的幅度谱axis(-6 6 0 0.7);legend('F1f(w)','F2f(w)

29、9;,'F3f(w)');subplot(313);plot(w,F1a*180/pi);xlabel('w');ylabel('P(度)');hold onplot(w,F2a*180/pi,'-.');plot(w,F3a*180/pi,':');%画出f1f2f3的相位axis(-6 6 -200 200); 图4.2 .1f (t), f (t-t0),f (t+t0)及其幅频特性与相频特性从图4.2.1可以看出，信号时移后其幅度频谱并没有改变，只是相位频谱发生了改变，增加0.2w或减少0.2w。4.3 傅

30、里叶变换的频移变换特性若f (t)的傅里叶变换为 F( jw) ，则傅里叶变换的频移特性为：例4:设f (t) = e (t +1) - e (t - 1) ，试用MATLAB 绘出的频谱F1( jw) 及F 2 ( jw) ， 并与f (t) 的频谱F( jw) 进行比较。w0=30;n=0.02;t=-1:n:1;g=Heaviside(t+1)-Heaviside(t-1);gt=1/2*g.*(exp(j*w0*t)+exp(-j*w0*t);subplot(221);plot(t,g);xlabel('t');ylabel('g');title(

31、9;脉宽为2的门函数的时域波形','FontSize',15);axis(-1.5 1.5 0 1.2);subplot(222);plot(t,gt);xlabel('t');ylabel('gt');title('调制信号的时域波形','FontSize',15);axis(-1.5 1.5 0 1.2);N=500;w1=2*w0;k=-N:N;w=k*w1/N;gw=n*g*exp(j*t'*w);gtw=n*gt*exp(j*t'*w);gwf=abs(gw);subplot(22

32、3);plot(w,gwf);xlabel('w');ylabel('gwf');title('门函数的幅度谱','FontSize',15);axis(-60 60 -0.5 2.5);gtwf=abs(gtw);subplot(224);plot(w,gtwf);xlabel('w');ylabel('gtwf');title('调制信号的幅度谱','FontSize',15);axis(-60 60 -0.5 2.5); 图4.4.1 傅里叶变换的频移特性由图4

33、.4.1可见，对比的结果可知 g(jw) 及) gt ( jw) 是将F( jw) 分别搬移到w 0= -20及w0 = 20 处的频谱。4.4 傅里叶变换的时域卷积定理变换的时域卷积定理如下：若信号f1( t) , f2( t) 的傅里叶变换分别为, F1( jw) F2( jw) ，则： f1( t) * f2( t) «F1( jw) × F2 (jw) 例5：设f (t) =u(t +1) - u (t -1), y(t ) = f (t) * f (t )试用MATLAB 给出f(t)、y(t)、F(j)、F(j)?F(j)及Y（j)的图形，验证式（9-13）的时

34、域卷积定理。n=0.005;t=-2:n:2;f=Heaviside(t+1)-Heaviside(t-1);%脉宽为2的门函数subplot(221);plot(t,f);xlabel('t');ylabel('f(t)');axis(-2 2 0 1.1);title('脉宽为2的门函数','FontSize',15);y=conv(f,f)*n;%对f f求积dt=-4:n:4;%y的取值范围为f取值范围的2倍subplot(222);plot(dt,y);xlabel('t');ylabel('f(

35、t)');axis(-4 4 0 2.3);title('y=f(t)*f(t)','FontSize',15);N=500;w1=2*pi*5;k=-N:N;w=k*w1/N;F=n*f*exp(j*t'*w);F=abs(F);F1=F.*F;subplot(223);plot(w,F);xlabel('w');ylabel('F(w)的幅度频谱');title('F(w)的幅度频谱','FontSize',15);subplot(224);plot(w,F1);%画出F1(w)

36、的幅度频xlabel('w');ylabel('F1(w)的幅度频谱');title('F1(w)的幅度频谱','FontSize',15);图4.5.1 时域卷积定理由图4.5.1可见，F（j）与F(j)´ F(j)的图形一致，而y(t)的波形正是我们熟知的t)*f(t)的波形，F（j）也是熟知的y(t)的付氏变换，从而验证时域卷积定理。4.5 傅里叶变换的对称性例5 设f (t ) = Sa(t ), 已知信号f (t) 的傅里叶变换为： 用MATLAB 求f2( t) = pg（t）的傅里叶变换F1 (jw) ,

37、并验证对称性。解：MATLAB程序为：n=0.1;t=-20:n:20;f1=sinc(t/pi);%抽样信号f2=pi*(Heaviside(t+1)-Heaviside(t-1);%门函数subplot(221);plot(t,f1);%画出抽样信号的时域波形xlabel('t');ylabel('f1');title('抽样信号的时域波形','FontSize',15);axis(-21 21 -1 1.5);subplot(222);plot(t,f2);%画出门函数的时域波形xlabel('t');yla

38、bel('f2');title('门信号的时域波形','FontSize',15);axis(-2 2 -3 5);N=500;%采样数k=-N:N;w1=20*pi;w=k*w1/N;F1=f1*exp(-j*t'*w)*n;%抽样信号的傅里叶变换F2=f2*exp(-j*t'*w)*n;%门函数的傅里叶变换subplot(223);plot(w,F1);%画出抽样信号的频谱xlabel('w');ylabel('F1(w)');title('抽样信号的频域','FontS

39、ize',15);axis(-2 2 -3 7);subplot(224);plot(w,F2);%画出门函数的频谱xlabel('w');ylabel('F2(w)');title('门函数信号的频域','FontSize',15);axis(-20 20 -3 7);图4.5.1 傅里叶变换对称性实例由图4.5.1可见，f (t ) = Sa(t ) 的傅里叶变换为的傅里叶变换为考虑到Sa(w) 是w 的偶函数，因此我们有： F(t)=2*pi*f (jw) ，即验证了傅里叶变换的对称性。4.6 傅里叶变换的时域微分性

40、质傅里叶变换的时域微分特性为：若f (t)的傅里叶变换为F(jw)，则: 例6: 已知f (t) 的波形如图9.13 所示，试用MATLAB求f (t) 及df (t)/ dt的傅里叶变换，F(jw) 及F'(jw)，并验证时域微分特性。 图4.7 f(t)的波形解：在MATLAB中，有专门的三角波形生成函数sawtooth()，其格式为：f = sawtooth(t, width)其中width(0<width1的标量)用于确定最大值的位置，即当t从0到2´ width变化时，f从-1上升到+1，然后当t从2´ width至4时f(t)又线性地从+1下降到-

41、1，周而复始。当width=0.5时，可产生一对称的标准三角波。利用此三角波与一门信号g2(t)相乘，再进行必要的幅度调整（乘系数2/ ），并时移（左移）可得到f(t)：又设f 1 ( t) = df( t)/dt，其波形为：f1( t) 可用阶跃函数Heaviside()生成：即验证：n=0.02;t=-6:n:6;h1=Heaviside(t+pi)-Heaviside(t-pi);f1=Heaviside(t+pi)-2*Heaviside(t)+Heaviside(t-pi);%三角波信号的一介倒数f=pi/2*(sawtooth(t+pi,0.5)+1).*h1;%三角波函数表达式N

42、=500;%采样点数k=-N:N;w1=10*pi;%采样频谱宽度w=k*w1/N;F=n*f*exp(-j*t'*w);%三角波函数的傅里叶变换F1=n*f1*exp(-j*t'*w);%三角波一介函数的傅里叶变换F2=F1./(j*w);%三角波一介函数的傅里叶变换除以jwsubplot(221);plot(t,f);%画出三角波函数的时域波形set(gca,'box','off');xlabel('t');ylabel('f(t)');title('三角波时域波形','FontSize',15);subplot(222);plot(t,f1);%画出三角波一介导数的时域波形set(gca,'box','off');xlabel('t