基于Matlab的频谱分析

1、数字信号处理课程设计题目： 基于Matlab的频谱分析 学院： 专业： 班级： 学号： 姓名： 指导教师： 目 录摘要1第1章 绪论41.1 课题背景41.2 研究意义41.3 本文研究内容4第2章 频谱分析背景知识52.1 时域抽样定理52.2 离散傅立叶变换（DFT）52.3 快速傅立叶变换（FFT）6第3章 设计内容与步骤83.1 考察序列x(n)有限个样本的频谱83.2 用DFT计算序列的幅度频谱样值103.3 采样含三个频谱分量的连续信号113.4 处理被噪声污染的信号15第4章 结论16第5章 收获与体会17参考文献18摘要MATLAB是一个数据分析和处理功能十分强大的工程实用软件

2、，运用它来进行信号的分析和处理相当便捷。 本文介绍了在Matlab环境中如何采集语音信号和语音信号采集后的频谱分析处理，并通过实例析了语音信号处理的Matlab。 并以理想正弦分布的纹理图像的 Matlab仿真和正弦纹理图像的频谱特征分析为例，叙述了基于傅立叶能量谱的纹理图像分析的基本原理和基本过程。编写了程序，获得了具有理想正弦分布的空域纹理图像，并对其频谱特征进行了分析。 并且根据离散傅里叶变换DFT的定义，推导出一种用 DFT计算离散信号幅值谱的方法，通过 MATLAB 仿真验证了该方法的有效性。 关键词：MATLAB；频谱分析；傅里叶变换；信号分析 Abstract MATLAB is

3、 very powerful and practical engineering software which is used in data analysis and processing fraction and use it to analysis and process voice signal which is very convenient. This paper describes the Matlab environment how to collect and process speech signal and through example analysis of the

4、speech signal processing in Matlab. Then take the Matlab simulation of ideal sinusoid texture image and analysis of its spectrum characteristics for example; basic principles and process of analyzing the texture image based on Fourier energy spectrum were described. The program was written for obtai

5、ning spatial domain texture image with ideal sinusoid distribution. And its spectrum characteristics were analyzed. Moreover, according to the definition of DFT, a calculation method of amplitude spectrum for periodic signal is derived in this paper. Through simulation of amplitude spectrum calculat

6、ion in MATLAB, the efficiency of the method is validated to satisfy the need for project. In the end, the problems of spectrum leakage and picket fence effect are explained and the corresponding solution is proposed. Keywords: MATLAB; Spectrum Analysis FFT; signal analysis第1章 绪论1.1 课题背景随着信息时代和数字世界的到

7、来，数字信号处理己成为当今一门极其重要的学科和技术领域，数字信号处理在通信、语音、图像、自动控制、医疗和家用电器等众多领域得到了广泛的应用。任意一个信号都具有时域与频域特性，信号的频谱完全代表了信号，因而研究信号的频谱就等于研究信号本身。通常从频域角度对信号进行分析与处理，容易对信号的特性获得深入的了解。因此，信号的频谱分析是数字信号处理技术中的一种较为重要的工具。 1.2 研究意义 信号处理几乎涉及到所有的工程技术领域，而频谱分析又是信号处理中一个非常重要的分析手段。一般的频谱分析都依靠传统频谱分析仪来完成，价格昂贵，体积庞大，不便于工程技术人员的携带。虚拟频谱分析仪改变了原有频谱分析仪的整

8、体设计思路，用软件代替了硬件，使工程技术人员可以用一部笔记本电脑到现场就可轻松完成信号的采集、处理及频谱分析。 1.3 本文研究内容 信号的频谱分析就是利用傅里叶分析的方法，求出与时域描述相对应的频域描述,从中找出信号频谱的变化规律，以达到特征提取的目的。不同信号的傅里叶分析理论与方法，在有关专业书中都有介绍。 但实际的待分析信号一般没有解析式，直接利用公式进行傅里叶分析非常困难. DFT是一种时域和频域均离散化的傅里叶变换，适合数值计算且有快速算法，是分析信号的有力工具。 DFT及FFT是数字信号处理的重要内容。DFT是FFT的基础，FFT是DFT的快速算法，在 MAT LAB中可以利用函数

9、 FFT来计算序列的离散傅里叶变换DFT。基于此首先介绍了Matlab软件；其次给出了基于Matlab软件的DFT和 FFT频谱分析的方法，利用 Matlab软件方法，使得设计方便、快捷，大大减轻了工作量；最后结合实例给出了仿真结果。 本文将重点介绍基于MATLAB的频谱分析设计。第2章 频谱分析背景知识2.1 时域抽样定理时域抽样定理给出了连续信号抽样过程中信号不失真的约束条件：对于基带信号，信号抽样频率 大于等于2倍的信号最高频率 ，即 。时域抽样是把连续信号 变成适于数字系统处理的离散信号 。对连续信号 以间隔T抽样，则可得到的离散序列为 。 图1 连续信号抽样的离散序列若 ， 则信号

10、与 的频谱之间存在：其中： 的频谱为 ， 的频谱为 。由可见，信号时域抽样导致信号频谱的周期化。 (rad/s)为抽样角频率， 为抽样频率。数字角频率与模拟角频率的关系为：=T。2.2 离散傅立叶变换（DFT） 离散傅里叶变换（DFT），是连续傅里叶变换在时域和频域上都离散的形式，将时域信号的采样变换为在离散时间傅里叶变换（DTFT）频域的采样。在形式上，变换两端（时域和频域上）的序列是有限长的，而实际上这两组序列都应当被认为是离散周期信号的主值序列。即使对有限长的离散信号作DFT，也应当将其看作经过周期延拓成为周期信号再作变换有限长序列的离散傅立叶变换（DFT）为逆变换为2.3 快速傅立叶变

11、换（FFT）在各种信号序列中，有限长序列占重要地位。对有限长序列可以利用离散傅立叶变换(DFT)进行分析。DFT不但可以很好的反映序列的频谱特性，而且易于用快速算法(FFT)在计算机上进行分析。有限长序列的DFT是其z变换在单位圆上的等距离采样，或者说是序列傅立叶的等距离采样，因此可以用于序列的谱分析。FFT是DFT的一种快速算法，它是对变换式进行一次次分解，使其成为若干小数据点的组合，从而减少运算量。MATLAB为计算数据的离散快速傅立叶变换，提供了一系列丰富的数学函数，主要有Fft、Ifft、Fft2 、Ifft2, Fftn、ifftn和Fftshift、Ifftshift等。当所处理的

12、数据的长度为2的幂次时，采用基-2算法进行计算，计算速度会显著增加。所以，要尽可能使所要处理的数据长度为2的幂次或者用添零的方式来添补数据使之成为2的幂次。Fft函数调用方式：、Y=fft(X)；、Yfft(X,N)；、Yfft(X,dim)或Yfft(X,N,dim)。 函数Ifft的参数应用与函数Fft完全相同。时间抽取算法令信号序列的长度为N=2,其中M是正整数，可以将时域信号序列x(n)分解成两部分，一是偶数部分x（2n），另一是奇数部分x（2n+1），其中。于是信号序列x(n）的离散傅里叶变换可以用两个N/2抽样点的离散傅里叶变换来表示和计算。因此，一个抽样点数为N 的信号序列x(n

13、）的离散傅里叶变换，可以由两个 N/2抽样点序列的离散傅里叶变换求出。依此类推，这种按时间抽取算法是将输入信号序列分成越来越小的子序列进行离散傅里叶变换计算，最后合成为N点的离散傅里叶变换。 频率抽取算法按频率抽取的 FFT算法是将频域信号序列X(k）分解为奇偶两部分，但算法仍是由时域信号序列开始逐级运算，同样是把N点分成N/2点计算FFT，可以把直接计算离散傅里叶变换所需的N次乘法缩减到次。 实际上，频率抽取算法与时间抽取算法的信号流图之间存在着转置关系，如将流图适当变形，可以得出多种几何形状。 实际上，频率抽取算法与时间抽取算法的信号流图之间存在着转置关系，如将流图适当变形，可以得出多种几

14、何形状。 除了基2的FFT算法之外，还有基4、基8等高基数的FFT算法以及任意数为基数的FFT算法。 第3章 设计内容与步骤 3.1 考察序列x(n)有限个样本的频谱 为了说明高密度频谱和高分辨率频谱之间的区别，考察序列x(n)有限个样本的频谱。当时，求x(n)的DFT.当时，求x(n)的DFT.程序:function DFTX1 %定义DFTX 1N=input('请输入离散信号的长度：');n=0:N-1; k =0:N-1; %序列范围xn=cos(0.48*n*pi)+cos(0.52*n*pi); %定义序列 WN=exp(-j*2*pi/N); %设定wn因子 nk

15、=n'*k; %赋值nk转置为n行1列矩阵，与k相乘 WNnk=WN.nk; %矩阵点幂 Xk=xn*WNnk; %DFT结果 subplot(2,1,1) stem(n,xn); %绘制离散图 title(' x(n):'); grid on; subplot(2,1,2) stem(k,abs(Xk); %绘制离散图 title(' x(n)的DFT:'); grid on;>> DFTX1请输入离散信号的长度：10>> DFTX1请输入离散信号的长度：100频谱图：图2 0n10时的DFT的频谱图图3 0n100时的DFT的

16、频谱图3.2 用DFT计算序列的幅度频谱样值已知一模拟信号，现在采样率进行采样。用DFT计算当序列长度分别为L=100,L=20时，N=200点的幅度频谱样值并通过作图与理论上准确的频谱样值进行比较。 程序：L1=100;L2=20; %对L1 L2赋值t1=0:L1-1;x1=exp(-t1/20);t2=0:L2-1;x2=exp(-t2/20);N=0:199; %定义N的范围WN=exp(-j*2*pi/length(N); %设定wn因子nk1=t1'*N;nk2=t2'*N; %赋值nk1 nk2转置为t1 t2行1列矩阵，与N相乘WNnk1=WN.nk1;WNnk

17、2=WN.nk2; %矩阵点幂Xk1=x1*WNnk1; %DFT结果Xk2=x2*WNnk2; %DFT结果subplot(4,1,1)stem(t1,x1);title('序列长度L=100时x(n)图:'); grid on;subplot(4,1,3)stem(t2,x2);title('序列长度L=20时x(n)图:'); grid on;subplot(4,1,2)stem(N,abs(Xk1);title('序列长度L=100时x(n)的DFT图:'); grid on;subplot(4,1,4)stem(N,abs(Xk2);t

18、itle('序列长度L=20时x(n)的DFT图:'); grid on;实验图：图4 序列长度L=100时的DFT的频谱图图5 序列长度L=20时的DFT的频谱图3.3 采样含三个频谱分量的连续信号 一个连续信号含三个频谱分量，经采样得以下序列：（1）N64，df分别为、1/64，观察其频谱；程序：function dftx3(df) %定义DFTX 3N=64; k=0:127; %定义k的范围n=0:N-1; %定义n的范围 xn=sin(2*pi*0.15*n)+cos(2*pi*(0.15+df)*n)+cos(2*pi*(0.15+2*df)*n); %定义序列WN

19、=exp(-j*2*pi/length(k); %设定wn因子 nk=n'*k; %赋值nk转置为n行1列矩阵，与k相乘WNnk=WN.nk; %矩阵点幂Xk=xn*WNnk; %DFT结果subplot(2,1,1)plot(n,xn);title('df=1/16,长度为N=64的有限序列x(n)图:'); grid on;subplot(2,1,2)plot(abs(Xk(1:64);title('取上图x(n)做128点的DFT图:'); grid on;>> dftx3(1/16)>> dftx3(1/64)>&g

20、t;实验图：图6 df=1/16,N=64，当x(n)做128点的DFT的图图7 df=1/64,N=64，当x(n)做128点的DFT的图（2）N64、128，df为1/64，做128点得FFT，其结果有何不同？程序：n1=0:63; n2=0:127; df=1/64;x1=sin(2*pi*0.15*n1)+cos(2*pi*(0.15+df)*n1)+cos(2*pi*(0.15+2*df)*n1); %定义x1序列x2=sin(2*pi*0.15*n2)+cos(2*pi*(0.15+df)*n2)+cos(2*pi*(0.15+2*df)*n2); %定义x2序列y1=fft(x1

21、,128); %对x1序列进行128点ffty2=fft(x2,128); %对x2序列进行128点fftsubplot(3,1,1)plot(n2,x2); title('有限序列x(n):'); grid on;subplot(3,1,2)plot(abs(y1(1:64);title(' 取有限序列x(n)长度N=64的128点FFT图:'); grid on;subplot(3,1,3)plot(abs(y2(1:64);title(' 取有限序列x(n)长度N=128的128点FFT图:'); grid on;实验图：图8 df=1/1

22、6,当N=64和N=128的128点的DFT的图3.4 处理被噪声污染的信号被噪声污染的信号，较难看出所包含得频率分量，如一个由50Hz和120Hz正弦信号构成的信号，受零均值随机噪声的干扰，数据采样率为1000Hz，试用FFT函数来分析其信号频率成分，要求：（1）画出时域波形；（2）分析信号功率谱密度。注：在MATLAB中，可用函数rand（1，N）产生均值为0，方差为1，长度为N的高斯随机序列。t=0:0.001:1; %定义t的范围以及采样间隔 x=sin(2*pi*50*t)+sin(2*pi*120*t)+rand(size(t); %定义x序列subplot(2,1,1)plot(x(1:50);title('由50Hz和120Hz正弦信号构成且受噪声污染的信号图:');grid on; Y=fft(x,51