基于matlab的功率谱分析方法研究毕业论文.doc

基于matlab的功率谱分析方法研究摘 要数字信号处理(DSP)重要的应用领域之一，是建立在周期信号和随机信号基础上的功率谱估计。在实际应用中往往不能获得具体信号的表达式,需要根据有限的数据样本来获得较好的谱估计效果，因而谱估计被广泛的应用于各种信号处理中1。本论文研究了功率谱估计的几种常用的方法，包括经典谱估计和现代谱估计的各种方法，且对每种方法的估计质量做了数学推导，并给出仿真程序及仿真图。经典法主要包括周期图法、自相关法，但这两种方法都存在缺陷，即认为观测数据之外的数据都为零，所以对经典法中的周期图法进行了加窗、平均等修正，因此提出了周期图法的改进方法；现代谱估计的方法分类比较多，AR模型法,MA模型法和ARMA模型法是现代功率谱估计中最主要的参数模型，本论文着重讨论了AR模型参数法2。同时论文将通过对经典谱估计和现代谱估计的实现方法及仿真图的比较，得出经典功率谱估计方法的方差性较差，分辨率较低，而现代谱估计的目标正是在于努力改善谱估计的分辨率，因此能得到较好的谱估计效果，为此应用更为广泛3。关键字：数字信号处理，功率谱估计，周期图法，自相关法，AR模型法 ABSTRACTDigital signal processing (DSP) important application of one of the field. Actually, we cant get the expression of a specific signal, so we need to estimate the power spectral of a signal according to some sample data sequences.So spectrum estimation which is widely used in various signal processing. In this thesis, some common methods of Power Spectral Estimation, such as classical spectral estimation and modern spectral estimation, are studied. The quality of each estimation method is derived, simulation program and simulation figure is given. Classical methods of Power Spectral Estimation mainly include the Periodogram and the BT method. But both of them have a common drawback: the data sequences, beyond the area of the observed sequences, are all presumed to zero. So the Windows and the average method are introduced to improve the quality of the Periodogram. Therefore the improvement of The Periodogram estimation method is proposed. The classification of modern spectral estimation methods are more , AR，MA, and ARMA is the most important parameters of modern spectral estimation. This thesis will focus on discussion of AR model parameters method. At the same time , It can be seen from the comparison and realization of classical spectral estimation and modern spectral estimation, classical power spectrum estimation variance is poor, low resolution .The goal of modern spectral estimation is working to improve the resolution of spectral estimation, better results of the estimation of the power spectrum can be obtained, so it is applied more widely. Keywords: digital signal processing, Power Spectrum Estimation, The Periodogram, the BT methods，AR model目 录摘 要IABSTRACTII1 绪论1.1功率谱估计的发展11.2功率谱估计的方法11.3功率谱估计的提出21.4经典谱估计21.5现代谱估计31.6功率谱估计应用及用途32 谱估计中的变量2.1随机信号简介52.2平稳随机信号72.3估计质量的评价标准103 经典功率谱估计3.1谱估计与相关函数123.2 周期图法153.3 自相关法183.4 直接法和间接法的关系183.5谱估计仿真与比较183.6 本章小结254 现代谱估计4.1平稳随机信号的参数模型274.2 AR模型的正则方程与参数计算284.3 AR模型谱估计的实现及性质314.4 MA模型谱估计334.5 ARMA模型谱估计354.6 小结364 论文总结37参考文献38致 谢401 绪论1.1功率谱估计的发展功率谱估计技术渊源流长，在过去的几十年获得了飞速的发展。功率谱估计涉及信号与系统、随机信号分析、概率统计、随机过程、矩阵代数等一系列的基础科学，广泛应用于雷达、声纳、通信、地址勘探、天文、生物医学工程等众多领域，其内容、方法不断更新，是一个具有强大生命力的研究领域4。功率谱估计(PSD)是用有限长的数据来估计信号的功率谱, 它对于认识一个随机信号或其它应用方面来讲都是极其重要的, 是数字信号处理的重要研究内容之一，在军事、生物医学、通信等领域得到了较为广泛的应用5。“谱”最早是由英国科学家牛顿提出来的，后来法国工程师傅里叶提出了著名的傅里叶谐波分析理论，该理论至今仍然是我们进行信号分析和处理的理论基础。傅里叶级数首先在观察自然界中的周期现象得到应用，但傅里叶的计算比较复杂，促使人们研制相应的机器来计算傅里叶级数。在19世纪末，Schuster提出傅里叶系数的平方，并命名为周期图，这是经典谱估计的最早提出法，至今仍被人们沿用6。后来，鉴于周期图的起伏剧烈，提出了平均周期图的概念，并提出了在对有限长数据计算傅里叶系数时所存在的边瓣问题，这就是后来我们所熟悉的窗函数的影响。周期图较差的方差性能促使人们研究另外的分析方法。Yule在1927年提出了用线性回归方程来模拟一个时间序列，从而发现隐含在该时间序列中的周期,从而发现了现代谱估计中最重要的方法参数模型法。Walker利用Yule的分析方法研究了衰减正弦时间序列，并得出了在对最小二乘分析中经常应用的Yule-Walker方程。7Yule的工作使人们重新想起了早在1795年Prony提出的指数拟合法，从而Prony方法形成了现代谱估计的又一重要内容。之后又陆续提出了Wiener-khintchine定理、谱估计自相关法BT法等。所有这些都为现代谱估计的发展打下了良好的基础7。1.2功率谱估计的方法功率谱估计可以分为经典谱估计(非参数估计) 和现代谱估计(参数估计)。经典谱估计的方法主要方法有自相关估计法和周期图法以及对周期图的改进方法; 现代谱估计的内容极其丰富，涉及的学科及应用领域也相当广泛，方法大致可分为参数模型谱估计和非参数模型谱估计，前者有AR 模型法(最大熵谱分析法)、MA模型，ARMA模型、Prony 指数模型等；后者有最小方差法，多分量的MUSIC方法等8。其中周期图法和AR 模型法是用得较多且最具代表性的方法。从信号的来源分，又可分为一维谱估计、二维谱估计及多维谱估计。从使用的统计量来分，目前大部分工作是建立在二阶矩基础上的，但由于功率谱密度是频率的实函数，缺少相位信息，因此，建立在高阶矩基础上的谱估计方法正引起人们的注意。从信号的特征来分，在这之前所说的方法都是对平稳随机信号而言，其谱分量不随时间变化，对非平稳随机信号，其谱是时变的，近20年来，以wigner分析为代表的时域分析引起了人们的广泛兴趣，形成了现代谱估计的一个新的研究领域。1.3功率谱估计的提出在通信系统中，往往需要研究具有统计特性的随机信号。由于随机信号是一持续时间无限长，具有无限大能量的功率信号，它不满足傅里叶变换条件，而且也不存在解析表达式，因此就不能够应用确定信号的频谱计算方法去分析随机信号的频谱9。然而，虽然随机信号的频谱不存在，但其相关函数是可以确定的。如果随机信号是平稳的，那么其相关函数的傅里叶变换就是它的功率谱密度函数，简称功率谱。功率谱反应了单位频带内随机信号的一个样本信号来对该随机过程的功率谱密度函数做出估计。1.4经典谱估计直接法：又称周期图法，它是把随机序列的N个观测数据视为一能量有限的序列，直接计算的离散傅立叶变换，得，然后再取其幅值的平方，并除以N，作为序列真实功率谱的估计。周期图这一概念早在1989年就提出了，但由于点数N一般比较大，该方法的计算量过大在当时无法使用，在1965年FFT出现后，此法才变成谱估计的一个常用方法。周期图法包含了下列二条假设：认为随机序列是广义平稳且各态遍历的，可以用其一个样本中的一段来估计该随机序列的功率谱,这当然必带来误差。由于对有限序列采用DFT，就默认此有限序列时域是周期的，以及该有限序列在频域是周期的。这种方法把随机序列样本看成是截得一段的有限序列的周期延拓，这也就是周期图法这个名字的来历。间接法：也叫相关法。间接法先由序列估计出自相关函数，然后对进行傅立叶变换，便得到的功率谱估计。周期图法与相关法相比，相关法在求相关函数时将有限长序列以外的数据看做是零，因此相关法认为除有限长序列外是全零序列，这种处理方法显然和周期法不一样。但是，当相关法被引入基于FFT的快速变换后，相关法和周期图法开始融合。改进的周期图法：对于直接法的功率谱估计，当数据长度N太大时，谱曲线起伏加剧，若N太小，谱的分辨率又不好，因此需要改进。其又有以下几种方法：1.Bartlett法Bartlett平均周期图的方法是将N点的有限长序列分段求周期图再平均。2. Welch法Welch法对Bartlett法进行了两方面的修正:一是选择适当的窗函数，并再周期图计算前直接加进去，加窗的优点是无论什么样的窗函数均可使谱估计非负。二是在分段时，可使各段之间有重叠，这样会使方差减小。不同窗函数的welch谱估计在选择窗函数时，一般有如下要求：1) 窗口宽度M要远小于样本序列长度N，以排除不可靠的自相关值；2) 当平稳信号为实过程时，为保证平滑周期图和真实功率谱也是实偶函数，平滑窗函必须是实偶对称的；3) 平滑窗函数应当在是峰值，并且m随绝对值增加而单调下降，使可靠的自相关值有较大的权值；4) 功率谱是频率的非负函数且周期图是非负的，因而要求窗函数的fourier变换是非负的。在经典谱估计中，无论是周期图法还是其改进方法，都存在着频率分辨率低、方差性能不好的问题，原因是谱估计时需要对数据加窗截断，用有限个数据或其自相关函数来估计无限个数据的功率谱，这其实是假设了窗以外的数据或自相关函数全为零，这种假设是不符合实际的，正是由于这些不符合实际的假设造成了经典谱估计分辨率较差。另外，经典谱估计的功率谱定义中既无求均值运算又无求极限运算，因而使得谱估计的方差性能较差，当数据很短时，这个问题更为突出，如何选取最佳窗函数、提高频率分辨率，如何在数据情况下提高信号谱估计质量，还需要进一步研究10。1.5现代谱估计现代谱估计与经典谱估计的主要区别就在于，现代谱估计一般采用信号模型法，信号模型法将原始信号视为白噪声通过一系统的输出信号，通过对输出信号的观测，按照一定的准则，求出相应的系统函数，这样再由输入白噪声和以求得的系统函数就很容易得到输出信号的功率谱。由已知白噪声和系统函数求得的输出序列，实际上是对原始观测到的输出信号的两端进行了估计或延拓。数据长度加宽以后，频谱分辨率会得到改善！因此现代谱估计优于经典谱估计。1.6功率谱估计应用及用途功率谱估计有着极其广泛的应用，不仅在认识一个随机信号时，需要估计它的功率谱。它还被广泛的应用于各种信号处理中。在信号处理的许多场所，要求预先知道信号的功率谱密度。例如，在最佳线性过滤问题中，要设计一个维纳滤波器就首先要求知道信号与噪声的功率谱密度，根据信号与噪声的功率谱才能设计出能够尽量不失真的重现信号，而把噪声最大限度抑制的维纳滤波器常常利用功率谱估计来得到线性系统的参数估计。例如，当我们要了解某一系统的幅频特性时，可用一白色噪声通过该系统，再从该系统的输出样本估计功率谱密度，故通过估计输出信号的PSD，可以估计出系统的频率特性。从宽带噪声中检测窄带信号。这是功率谱估计在信号处理中的一个重要用途。但是这要求功率谱估计有足够好的频率的分辨率，否则就不一定能够清楚地检测出来。所谓谱估计的分辨率可以粗略的定义为能够分辨出的二个分立的谱分量间的最小频率间隙，提高谱估计的分辨率已成为目前谱估计研究中的一个重要方向。功率谱估计就是通过信号的相关性估计出接受到信号的功率随频率的变化关系，实际用途有滤波、信号识别、信号分离、系统辨识等。谱估计技术是现代信号处理的一个重要部分，还包括空间谱估计、高阶谱估计等11。2 谱估计中的变量2.1随机信号简介2.1.1 随机变量随机变量（random variable）表示随机现象（在一定条件下，并不总是出现相同结果的现象称为随机现象）各种结果的变量（一切可能的样本点）。例如某一时间内公共汽车站等车乘客人数，电话交换台在一定时间内收到的呼叫次数等等，都是随机变量的实例。随机变量在不同的条件下由于偶然因素影响，其可能取各种不同的值，具有不确定性和随机性，但这些取值落在某个范围的概率是一定的，此种变量称为随机变量。随机变量可以是离散型的，也可以是连续型的。如分析测试中的测定值就是一个以概率取值的随机变量，被测定量的取值可能在某一范围内随机变化，具体取什么值在测定之前是无法确定的，但测定的结果是确定的，多次重复测定所得到的测定值具有统计规律性。随机变量与模糊变量的不确定性的本质差别在于，后者的测定结果仍具有不确定性，即模糊性。按照随机变量可能取得的值，可以把它们分为两种基本类型：离散型随机变量，即在一定区间内变量取值为有限个，或数值可以一一列举出来。例如某地区某年人口的出生数、死亡数，某药治疗某病病人的有效数、无效数等。连续型随机变量，即在一定区间内变量取值有无限个,或数值无法一一列举出来。例如某地区男性健康成人的身长值、体重值，一批传染性肝炎患者的血清转氨酶测定值等。1.随机变量的分布函数设X是随机变量，对任意实数，事件X<x的概率称为随机变量X的分布函数。记为，即，易知，对任意实数a,b,。分布函数的性质(1) 单调不减性：若, 则；(2) 归一性：对任意实数，且 ，(3) 左连续性：对任意实数x，2.数学期望、方差、标准差定义： ,为的数学期望值，或简称为均值。；以上分别称为X的标准差和方差。若为离散型随机变量，则上述的求均值运算将有积分改为求和。例如，式中的是取值为时的概率12。3.随机向量在某些实际问题中，往往需要同时用两个或两个以上的随机变量来描述试验的结果。设E是一个随机试验, 样本空间是，设和是定义在上的随机变量， 由它们构成的一个向量叫做二维随机向量或二维随机变量。（注: 二维随机向量 的性质不仅与和有关,而且还依赖于这两个随机变量的相互关系。）4.概率密度函数概率密度函数是为了表示瞬时数据落在指定幅值范围的概率。其定义为：。瞬时值小于或等于某值x的概率定义为概率分布函数或累计概率分布函数5.相关函数表征了一个随机过程自身在不同时刻的状态间，或者两个随机过程在某个时刻状态间线性依从关系的数字特征。相关函数是两随机变量之积的数学期望，称为相关性。统计学中用相关系数xy来描述变量x，y之间的相关性，函数的相关系数，简称相关函数：2.1.2随机信号的特征随机信号具有不重复性、不确定性，通常用概率与统计方法研究其中是否存在某些重复、确定的成分。随机过程在某一时刻的均值（一阶矩）可将总体中各样本函数在的瞬时值相加，然后除以样本函数的个数而得到。自相关函数即为随机过程两不同时刻之值的相关性，又称二阶矩。用和两时刻瞬时值乘积的总体平均值得到。自相关函数的性质：(1) 自相关函数是 的偶函数 ；(2) 当 时，自相关函数具有最大值，；(3) 周期信号的自相关函数仍然是同频率的周期信号，但不保留原信号的相位信息。(4) 当随机信号中含有周期信号时，中也必定有周期性分量，且周期相同。(5) 对变化迅速的信号（宽带随机过程），相关的程度在很小时就完全丧失。2.2平稳随机信号2.2.1平稳随机信号的定义平稳信号分严平稳和宽平稳，严平稳的条件在信号处理中太严格，不实用，一般所说的平稳是指宽平稳，满足三个条件：(1) 均值为与时间无关的常数；(2) 均方有界；(3) 自相关函数与信号时间的起始点无关，只和时间差有关（宽平稳信号的方差和均方也是与时间无关的）。(1) 平稳随机过程的定义:如果对于任意n和以及有则称为严平稳随机过程，或称狭义平稳随机过程。(2) 平稳随机过程的数字特征：1) ,平稳随机过程的数学期望与时间无关；2) ，平稳随机过程的方差与时间无关；3) 其中：；4) 。平稳随机过程的数学期望及方差与无关，它的自相关函数和协方差函数只与时间间隔有关；随机过程的这种“平稳”数字特征，有时就直接用来判断随机过程是否平稳。即若一个随机过程的数学期望及方差与时间无关，而其相关函数仅与有关，即我们就称这个随机过程是广义平稳的。(3) 宽平稳随机过程（广义平稳）:若的数学期望为常数，且自相关函数只与有关，则称为宽平稳随机过程，或称广义平稳随机过程。不难看出，严平稳过程一定是宽平稳过程，反之不一定。但对于正态随机过程两者是等价的。本论文若不加特别说明，平稳过程均指宽平稳过程13。2.2.2 平稳随机信号的自相关函数实随机信号的自相关函数定义：，由于平稳随机信号的统计特性与时间的起点无关，设, 则有。所以，平稳随机信号的自相关函数是时间间隔的函数，记为。平稳随机信号自相关函数的性质：设为平稳随机过程，其自相关函数为，自协方差函数，则它们有如下性质：（1）时的自相关函数等于均方差，自协方差函数等于方差， 即 (2) 当平稳随机信号是实函数时，其相关函数是偶函数，即： (3) =0时的自相关函数、自协方差函数取最大值，即 (4) 若X(t)=X(t+T)，则其自相关函数也是周期为T的周期函数，即 (5) 若均值，当时，与相互独立，有 ，即对于零均值的平稳随机信号，当时间间隔很大时，与相互独立，互不相关。2.2.3 平稳随机信号的功率谱1. 功率谱密度定义：设,<<是均方连续的随机过程，称为的平均功率。称为的功率谱密度，简称谱密度。2. 功率谱密度的性质（1）若 ，则 是 的傅里叶变换；(2) SX()是的非负实函数；(3) 实平稳过程的谱密度是偶函数；当 是的有理函数时，其形式必为，其中,为常数，且, ，分母无实根。3.随机序列的功率谱随机序列，它的相关函数满足其功率谱密度 具有如下式子：2.3估计质量的评价标准1.无偏性对于待估参数，不同的样本值就会得到不同的估计值，这样，要确定一个估计量的好坏，就不能仅仅依据某次抽样的结果来衡量，而必须由大量抽样的结果来衡量对此，一个自然而基本的衡量标准是要求估计量无系统偏差。也就是说，尽管在一次抽样中得到的估计值不一定恰好等于待估参数的真值，但在大量重复抽样时，所得到的估计值平均起来应与待估参数的真值相同换句话说，我们希望估计量的均值（数学期望）应等于未知参数的真值，这就是所谓无偏性(Unbiasedness)的要求。定义:设来自总体X的一个样本，是总体参数 的一个估计量，若 ,则称是 的无偏估计量(Unbiased Estimator)。一个估计量如果不是无偏的就称它是有偏估计量。称为估计量的偏差。无偏估计的实际意义就是无系统偏差，估计量是否无偏是评价估计量好坏的一个重要标准，若 ，但有，则称是的渐近无偏估计。 2.有效性比较两个无偏估计量优劣的一个重要标准就是观察它们哪一个取值更集中于待估参数的真值附近，即哪一个估计量的方差更小，这就是下面给出的有效性(Effectiveness)概念。 定义:设与都是总体参数的无偏估计，若 ,则称比 更有效。在的所有无偏估计量中，如果存在一个估计量，它的方差最小，则此估计量应当最好，并称此估计量为 的最小方差无偏估计，也称其为最有效的3.相合性估计量的无偏性和有效性都是在样本容量n固定的情况下讨论的。由于估计量和样本容量n有关，我们自然希望当很大时，一次抽样得出的的值能以很大的概率充分接近被估参数，这就提出了相合性(Consistency)（一致性）的要求。定义：设 是总体参数的估计量，如果对任意都有, 则称是的相合估计量（或一致估计量）。是的相合估计就意味着依概率收敛于.根据大数定律，无论总体X服从什么分布，只要其阶原点矩 存在，则对任意 都有，所以样本的阶原点矩始终是总体阶原点矩的相合估计。 进一步地, 可以证明：只要相应的总体矩存在，矩估计必定是相合估计。特别地， 总是 的相合估计, 样本方差 和样本的二阶中心矩都是总体方差 的相合估计和又都是的相合估计。由相合性定义可以看出，若是的相合估计，当样本容量很大时，一次抽样得到的值便可作为的较好近似值14。3 经典功率谱估计经典谱估计方法是以傅里叶变换为基础的方法，主要有两类：周期图法和自相关法（布莱克曼图基法，简称BT法）。这两类方法都与相关函数有着密切的联系，由维纳辛钦定理可知，功率谱和相关函数之间的关系是一对傅里叶变换，因而可以从观测数据直接估计相关函数，根据估计出来的相关函数，求它的傅立叶变换，就可以得到功率谱的估计值。3.1谱估计与相关函数 3.1.1 相关函数和功率谱若 常数，即，则称为广义平稳序列。若和均为广义平稳序列，且即：，则称和为广义联合平稳序列。广义平稳随机序列的相关函数和它的功率谱密度之间是傅立叶变换对的关系，即这一关系式常称为维纳辛钦定理。由自相关函数和功率谱密度的定义，不难得出它们的一些基本性质，主要有：1、当为复序列时，；若为实序列，则相关函数为偶函数，即。2、相关函数的极大值出现在处，即。3、若含有周期性分量，则也含有同一周期的周期性分量，否则，当时，。4、当为实序列时，为非负实对称函数，即和。5、平稳随机序列的自相关函数是实的且为正，而且对任一序列和任一，自相关函数（ACF）满足：，这个函数称为半正定的。自相关函数（ACF）和互相关函数（CCF）的变换定义为：；，若令为归一化频率，频率区间为基本周期。则上述两式功率谱密度又可分别表示为：其中，是实的，且非负。当一平稳随机序列通过一个脉冲响应为的线性非时变系统时，其输出序列也是一平稳随机序列。它的自相关函数为：若为实系统，则。令，得到相应的功率谱表达：或，上述关系对以后讨论谱估计问题是很有用的。为输出过程的平均功率。经常会遇到的一种过程是离散白噪声，它的自相关函数（ACF）定义为： ，其中是离散冲激函数。这就是说，各样本之间彼此是不相关的。所以 ， 这表明它在各频率上是完全平坦的。换句话说，白噪声的所有频率分量均具有相同的功率15。3.1.2 相关函数的估计1、自相关函数的各态历经性一般说来，严格各态历经过程允许我们用时间平均来代替系综平均（集合平均或统计平均），用时间平均作为广义平稳随机过程均值的估计。2.相关函数的估计实际所能得到的随机序列的样本数总是有限的，由有限个样本通过某种运算求出的序列的均值和自相关函数统计特征值叫做它们的估计值。下面讨论随机序列有限个样本的相关函数的估计问题。设为实随机序列的一批样本，共有N个值。有时简称之为长度为N的随机序列。方法一：根据假定的自相关函数的各态历经性（或遍历性），可用下式估计它的自相关函数，即当时，因此是相关函数的无偏估计且是渐近一致的，即当为有限值时，是的一致估计。方法二：有限长度序列的相关函数的另一种估计方法可表示为（3-1）可见，它是相关函数的有偏估计。但是，当，估计值是渐近无偏的。当时，即（3-1）式的也是的一致估计。（3-1）式所定义的相关函数取傅立叶变换求功率谱估计时，在计算上有某些方便之处，以后的讨论中，如不作特别申明，将采用这种有偏估计表示式求相关函数的估计式。功率谱密度的另一个定义:可以证明，功率谱密度（PSD）的一个近似等效的定义是（3-2）上式定义的PSD与维纳一辛钦定理（3-3）是等效的。由（3-2）式和由（3-3）式维纳一辛钦定理给出的PSD的等效定义将作为经典谱估计方法的基础16。3.2 周期图法3.2.1周期图法的定义周期图法，它是把随机序列的N个观测数据视为一能量有限的序列，直接计算的离散傅立叶变换，得，然后再取其幅值的平方，并除以N，作为序列真实功率谱的估计。周期图谱估计定义为： 可以证明，周期图等于估计出的自相关序列的傅里叶变换，或其中是有偏自相关函数的估计值，定义为:3.2.2 周期图的性能周期图的期望值是:（3-4）式中是Bartlett窗（三角形窗）的傅里叶变换。由式（3-4）可知周期图的均值是真实PSD和Bartlett窗傅里叶变换的卷积，在平均意义上得到真实功率谱密度（PSD）的平滑形式。因此对有限记录数据，周期图一般有偏的，但是当时，它是无偏的。 这是由于收敛到狄拉克函数。对于高斯白噪声的特殊情况，结果为:对于白噪声情况，即使有限记录数据，周期图也是无偏的。对于白高斯过程，方差对任何非白过程，只要记录数据足够长,对于不靠近0或的频率，且时，上式近似地退化为：可以看出，方差与数据长度N无关，即方差不随N的增大而减小至零。由此可得出一个重要的看法：周期图估计器是不可靠的，因为标准差和均值一样大，因而周期图不是一致估计而其均值近似地等于要估计的量值。上述论证表明，我们不能寄希望于直接用周期图方法获得良好的谱估计，必须采用适当的修正措施减小估计方差，才能使之成为一种实用的方法17。3.2.3周期图法改进措施1、加窗周期图周期图法只用了N个样本，这可以看作是用一长度为N的矩形窗函数与原来无限长的序列相乘的结果，我们知道，时域中两函数相乘对应于频域中它们的傅立叶变换的卷积。由此可以想到，当用周期图方法作谱估计时，它的谱分辨率约与N成反比，且和信号本身的特征，例如信噪比等无关。此外，如果序列是由多个正弦波信号组成的，而各分量强度不等，则弱信号分量可能淹没在强信号谱的旁瓣中而无法发现。这种所谓信号能量（向旁瓣）泄漏现象如果不设法消除，也将妨碍周期图谱估计法的应用。因此提出了周期图的改进方法：周期图改进的方法之一是将长度为N的序列乘以同一长度的数据窗。数据加窗后的周期图谱估计值的数学期望值等于谱的真实值与窗谱函数的平方的卷积。显然它不等于功率谱的真实值，因而是有偏估计。若序列为单频信号，则为函数，这样，数据加窗后的谱估计值的均值与窗谱函数的平方形状相同，因此选用低旁瓣的数据窗可使得杂散响应减少。但旁瓣的降低必然使主瓣加宽，而且降低了分辨率。数据加窗后，周期图谱估计值的方差大于或近似等于谱估计值的平方，且不随数据长度的增大而减小到0。从以上的分析可知，数据加窗用于周期图谱估计可以降低谱估计值的旁瓣，但要降低谱估计的分辩率，而用数据加窗的办法不能减小估计方差，因而无法降低分辨率。2、平均周期图平均周期图的思想：对一个随机变量进行观测，得到L组独立记录数据，用每一组数据求其均值，然后将L个均值加起来求平均。这样得到的均值，其方差将是用一组数据得到的均值的方差的1/L。为了改进周期图的统计特性，可以近似地用对一组周期图进行平均的方法完成期望运算。假定在区间上有组独立记录数据，并且都是同一随机过程的现实。平均周期图估计器定义为： 其中是第个数据组的周期图。因为，所以方差将减少K倍。又因为 ，所以谱估计器具有样本的分辨率。显然，为了获得最大分辨率，L将尽可能选得大一些或就选，此时，平均周期图就成为标准周期图估计。但有时为了减小方差，应该把选得大一些，或等效地把L取得小一些。因为这两个目的不可能同时实现，所以必须调整L在方差和偏差之间进行折衷18。3.3 自相关法自相关法是先估计自相关函数，再进行傅里叶变换得到功率谱。下面介绍自相关函数的估计。利用观测到的实随机序列，估计自相关函数的两种方法是：无偏自相关函数估计和有偏自相关函数估计。BT法功率谱估计采用有偏自相关函数估计法。自相关法的理论基础是维纳-辛钦定理，即先由估计出自相关函数，然后对求傅里叶变换得到的功率谱，因为由这种方法求出的功率谱是通过自相关函数间接得到的，所以称为间接法。3.4 直接法和间接法的关系由以下式子可以得出：m=k-n，则k=n+m，得 可见周期图法就是BT法中M=N-1时的功率谱估计19。3.5谱估计仿真与比较1.周期图法以下是利用matlab对周期图法进行仿真，并得到仿真图（3-1）。程序如下：clear;Fs=1000;%采样频率n=0:1/Fs:1;%产生含有噪声的序列xn=cos(2*pi*40*n)+3*cos(2*pi*100*n)+randn(size(n);nfft=1024;window=boxcar(length(n);%加矩形窗Pxx,f=periodogram(xn,window,nfft,Fs);%直接法 plot(f,10*log10(Pxx);图3-1 周期图法仿真图2.平均周期图以下是利用matlab对平均周期图法进行仿真，并得到仿真图（3-2）和（3-3）。程序如下：clear;Fs=1000;n=0:1/Fs:1;xn=cos(2*pi*40*n)+3*cos(2*pi*100*n)+randn(size(n);nfft=1024;window=boxcar(length(n); %加矩形窗noverlap=0; %数据无重叠 p=0.9; %置信概率Pxx,Pxxc=psd(xn,nfft,Fs,window,noverlap,p);index=0:round(nfft/2-1);k=index*Fs/nfft;plot_Pxx=10*log10(Pxx(index+1);plot_Pxxc=10*log10(Pxxc(index+1);figure(1) plot(k,plot_Pxx); pause; figure(2)plot(k,plot_Pxx plot_Pxx-plot_Pxxc plot_Pxx+plot_Pxxc);图3-2 周期图法仿真图图3-3平均周期图法仿真图3.加窗周期图以下是利用matlab对加窗周期图法进行仿真，并得到仿真图（3-4），（3-5）和（3-6）。程序如下：clear;Fs=1000;n=0:1/Fs:1;xn=cos(2*pi*40*n)+3*cos(2*pi*100*n)+randn(size(n);nfft=1024;window=boxcar(100);%矩形窗 window1=hamming(100);%汉明窗window2=blackman(100); %blackman窗noverlap=20; %数据无重叠range='half'%频率间隔为0 Fs/2，只计算一半的频率Pxx,f=pwelch(xn,window,noverlap,nfft,Fs,range);Pxx1,f=pwelch(xn,window1,noverlap,nfft,Fs,range);Pxx2,f=pwelch(xn,window2,noverlap,nfft,Fs,range);plot_Pxx=10*log10(Pxx);plot_Pxx1=10*log10(Pxx1);plot_Pxx2=10*log10(Pxx2);figure(1)plot(f,plot_Pxx);pause;figure(2)plot(f,plot_Pxx1);pause;figure(3)plot(f,plot_Pxx2);图3-4 加矩形窗仿真图图3-5 加汉明窗仿真图图3-6 加blackman窗仿真图从图中可以看出，采用周期图法估计得到的功率谱很不平滑，相应的估计协方差比较大，而且采用增加采样点的办法也不能使周期图变得更加平滑，这是周期图法的缺点。周期图法得出的估计谱方差特性不好：当数据长度N太大时，谱线的起伏加剧20；N太小时，谱的分辨率又不好。对其改进的主要方法有两种，即平均和平滑，平均就是将截取的数据段再分成L个小段，分别计算功率谱取功率谱的平均，这种方法使估计的方差减少，但偏差加大，分辨率下降。平滑是用一个适当的窗函数与算出的功率谱进行卷积，使谱线平滑，这种方法得出的谱线是无偏的，方差也小，但分辨率下降。4.自相关法以下是利用matlab对自相关法进行仿真，并得到仿真图（3-7)。程序如:clear;Fs=1000; %采样频率n=0:1/Fs:1;%产生含有噪声的序列xn=cos(2*pi*40*n)+3*cos(2*pi*100