Matlab方波频谱分析课程设计报告

ZhenRihouInstituteof^eronautkalIndustryMana^eineiit通信系统建模与仿真课程设计2010级通信工程专业1013072班级题目基于Matlab/Simulink的信号频谱的估计姓名学号指导教师 胡娟，王丹，王娜，闫利超 1任务书用Matlab编程方式产生一个100Hz的方波，画出其波形。并用fft指令计算其频谱，做出幅度谱和相位谱，与理论结果进行对比。用Simulink方式重做上题，并通过统计模块在时域和频域同时计算信号的功率，看两者计算结果是否一致，验证帕萨瓦尔定理。2理论分析方波的一个周期可用1,0<t<；T 八-1,-—<t<0I2依据周期信号傅里叶级数系数的定义，有F=T\TI2f(t)e-jn^tdt兀兀 兀因此，方波信号的fT(t)的傅里叶级数展开式为f((t)=£Fjt登▲(1-cos…n=-s n=一8根据周期信号傅里叶级数同傅里叶变换之间的关系：f(t)= FejnQt £2时(.-n。)n=-sn=-sn=-s可知，方波信号的傅里叶变换是f(t)sH(加)=Uj^~(cosnn-1)8(①一n。)T nn=-s显然，当n为偶数时， ，因此方波信号中只存在奇次谐波其功率谱c号4A2, 〜=乙 (cosn兀一1)28(①一n。)n2n=—s化为以频率为自变量表示的功率密度谱，得到H(f)|2=£-A^(cosn冗一1)28(f-n)n侦2 Tn=—s可见，方波在几次谐波处存在冲激谱线，其功率谱谱线冲激强度为数列C，n取奇数，C为常数。n2离散时间信号的帕斯瓦尔定理：对于N点的离散序列及其离散傅里叶变换f(n)〜F(k)，其时域能量等于频域能量，即E=T习f(n)n=E=T习f(n)n=0k=0时域和频域的平均功率关系为P=E=E=—切If(n)2=—习If(k)|2其中，T为米avLNTN' N21 1n=0 k=0样时间间隔；N为离散时间序列的点数；L=NT为离散时间序列的时间长度。Matlab代码详述clear;clc;fs=1e6;

%计算时间范围ft二square(2*pi*100*t1,50);subplot(3,1,1);plot(t1,ft);axis([00.1-1.21.2]);T1=0.01;w1=2*pi/T1;n=-59:2:59;W=w1.*n;F_w=-4*j./n;subplot(3,1,2);stem(W,abs(F_w));holdon;%时域波形%信号周期%信号角频率%时域波形%信号周期%信号角频率%奇次谐波数%数字角频率%频谱理论结果%频域幅度谱T=pi/w_m; %采样间隔L=5.9;t=0:T:L; %时域截断x_t=square(2*pi*100*t,50); %信号序列N=length(x_t); %序列长度(点数)%FFT计算%FFT计算w0=2*pi/(N*T); %离散频率间隔kw=2*pi/(N*T).*[0:N-1]; %离散频率样点X_kw=T.*X_k; %乘以T得到连续傅里叶变换频谱的样值plot(kw-w_m,abs(fftshift(X_kw)),'.','MarkerSize',10);%做出数值计算的幅度谱点subplot(3,1,3);stem(W,angle(F_w)); %频域相位谱holdon;plot(kw-w_m,angle(fftshift(X_kw))); %做出数值计算的相位谱点SIMULINK各模块说明由于Simulink中FFT模块只接受2的整数幕次点数数据，故设计变换数据采样率为2048样值/秒。FFT变换数据长度到2048,对应时间长度为1秒。因此，频率分辨率为1Hz。由式TOC\o"1-5"\h\zF(f)|412=T•DFT[f(n)]|2=T2•\F(k)|2，=NTs S S确定功率谱估计值。其中Ts=1/2048秒，N=2048。依据教材式(3.34)确定频域、时域平均功率，即 P=1如f(k)2=—云IF(k)|2avN' N21 1n=0 k=0图1为测试模型。其中，仿真步长为固定的1/2048秒。各个参数设置如下：Zero-OrderHold采样得到离散时间信号以便进行离散傅里叶变换，设置采样时间间隔也为1/2048秒。Buffer模块设置缓存长度为2048，刚好能缓存1秒的数据。采用FFT模块进行快速傅里叶变换，变换点数取决于数据帧长度。以Abs模块、乘法器模块、Mean平均模块以及增益模块等实现对频域、时域平均功率的计算。Display模块显示时域功率计算结果，Display1模块显示时域功率计算结果，显然两者应当相同，为1W（幅度为1V的方波理论计算功率为1W）。用VectorScope观察功率谱，同时用示波器观测时域波形。模型图如下所示：5仿真结果分析编程输出结果如图1所示，时域仿真时间0.1s，但是该信号时域是无限长的，因此其频谱也是无限宽的。经理论公式计算画图显示后发现，当频率大于30000Hz后，幅度谱值接近于零。因此做FFT变换时选择频域截断区为30000Hz。图1编程仿真结果Scope-昌首I I费通圈旧■点图2Simulink仿真输出的方波时域波形-1-O.S-O.e-0.4 -0.2 0 0.2 0.4O.BFrame:2. Frequency(kH^)-1-O.S-O.e-0.4 -0.2 0 0.2 0.4O.BFrame:2. Frequency(kH^)0.8 166432JIcJcJQO.O.0.号三dE<nch3exampIe3/VectorScopeli■=■iieir^~FileViewAxesChannels-WindowHelp -声■GGM图3Simulink仿真输出的幅度谱0.06-1 -0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8Frame-2 Frequency(kHz)0.06-1 -0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8Frame-2 Frequency(kHz)14535251

O.3O.2O..1O

ooogpm=dE<rch3example3/V&ctorScop&FileViewAxesChannelsWindowHelp~3i图4Simulink仿真输出的功率谱6遇到的问题及解决的方法在产生方波公式时，一直无从下手。数学表达式很简单，可是无法在matlab脚本文件编程准确实现。总有不尽人意的地方。后来改用matlab程序库函数square函数输出方波，以前没用过这个函数，又查了查这个函数的用法，以及它的内容与数学公式的差别，收获挺大的。在计算连续周期函数方波的频谱时，遇到了一系列的问题。由于方波是由函数输出的，所以在理论计算时显然用不到。所以改用方波的数学表达式来计算频谱。其中应用到了《信号与系统（上）》相关章节，由于这门课是去年修完的，记忆不太完整，所以又翻了翻这本书，查找了相关知识，从而正确计算出了计算方波频谱的数学表达式。在写程序时发现数学表达式和程序的表达还是有差别的，又认真学习了程序的表达，终于写出相对简单容易理解的理论编程程序。对于理论与实践收获还是相当大的。7结束语结本文介绍了基于MATLAB/Simulink的信号频谱的估计和使用MATLAB/.M