利用经典谱估计法估计信号的功率谱(随机信号)

随机信号功率谱估计作者：杨太张荣龙袁晓华来源：《中国新通信》2012年第14期1引言功率谱估计（PSD）是利用给定的一组样本数据估计一个平稳随机信号的功率谱密度，它能给出被分析对象的能量随频率分布的情况。

在雷达信号处理中可以根据回波信号的功率谱密度、谱峰的宽度、高度和位置可以确定运动目标的位置、辐射强度和运动速度。

功率谱估计是数字信号处理的重要研究内容之一。

功率谱估计可以分为经典谱估计（非参数估计）和现代谱估计（参数估计）。

经典的功率谱估计有2种：一种是直接法；另一种是间接法。

直接法就是先计算N个数据的傅里叶变换，然后取频域和其共轭的乘积得到功率谱；间接法则是先计算N个样本数据的估计自相关函数，然后再计算自相关数据的傅里叶变换得到功率谱。

间接法的主要方法有最大熵谱分析法（AR模型法）、Pisarenko谐波分解法、Prony提取极点法、Prony 谱线分解法以及Capon最大似然法。

其中周期图法和最大熵谱分析法是信号处理中最具代表性的方法。

2周期图法2.1周期图谱估计周期图法是直接将信号的采样数据x（n）进行Fourier变换求取功率谱密度估计的方法。

它的特点是：先取信号序列的离散傅里叶变换，然后取其幅频特性的平方并除以序列长度，由于序列的离散傅里叶变换具有周期性，因而这种功率谱也具有周期性。

假定有限长随机信号序列为x（n）。

它的Fourier变（2）其中，FFT[x（n）]为对序列x（n）的Fourier变换，由于FFT[x（n）]的周期为N，求得的功率谱估计以N为周期，因此这种方法称为周期图法。

2.2周期图法的性质差的分析，以此来说明周期图法来作为功率谱估计的不可靠性。

在信号处理中通常用偏差和方差这两种基本的统计量来度量估计子的性能，设a为一个带估计量的变量，a赞是它的一个估计。

估计的均方差（MSE）为下面是对一随机信号分别取长度为256和1024的功率谱估计：从图1中可以看出，采用周期图法估计得到的功率谱很不平滑，也就是估计的协方差比较大，而且采用增加采样点的方法也不能使周期图变得更加平滑，这是周期图法的缺点。

随机信号利用经典谱估计法估计信号的功率谱作业综述：给出一段信号“asd.wav”,利用经典谱估计法的原理，通过不同的谱估计方法，求出信号的功率谱密度函数。

采用MATLAB语言，利用MATLAB语言强大的数据处理和数据可视化能力，通过GUI的对话框模板，使操作更为简便！在一个GUI界面中，同时呈现出不同方法产生出的功率谱。

这里给出了几种不同的方法：BT法，周期图法，平均法以及Welch法。

把几种不同方法所得到的功率谱都呈现在一个界面中，便于对几种不同方法得到的功率谱作对比。

一．题目要求给出一段信号及采样率，利用经典谱估计法估计出信号的功率谱。

二．基本原理及方法经典谱估计的方法，实质上依赖于传统的傅里叶变换法。

它是将数据工作区外的未知数据假设为零，相当于数据加窗，主要方法有BT法，周期图法，平均法以及Welch法。

1. BT法(Blackman-Tukey)● 理论基础：(1)随机序列的维纳-辛钦定理由于随机序列{X(n)}的自相关函数Rx(m)=E[X(n)X(n+m)]定义在离散点m 上，设取样间隔为，则可将随机序列的自相关函数用连续时间函数表示为等式两边取傅里叶变换，则随机序列的功率谱密度(2)谱估计BT法是先估计自相关函数Rx(m)(|m|=0,1,2…,N-1),然后再经过离散傅里叶变换求的功率谱密度的估值。

即其中可有式得到。

2. 周期图法● 理论基础：周期图法是根据各态历经随机过程功率谱的定义来进行谱估计的。

在前面我们已知，各态历经的连续随机过程的功率谱密度满足式中是连续随机过程第i个样本的截取函数的频谱。

对应在随机序列中则有由于随机序列中观测数据仅在的点上存在，则的N点离散傅里叶变换为：因此有随机信号的观测数据的功率谱估计值（称“周期图”）如下：由于上式中的离散傅里叶变换可以用快速傅里叶变换计算，因此就可以估计出功率谱。

3.平均法：● 理论基础：平均法可视为周期图法的改进。

周期图经过平均后会使它的方差减少，达到一致估计的目的，有一个定理：如果是不相关的随机变量，且都有个均值及其方差，则可以证明它们的算术平均的均值为。

功率谱分析在信号处 理中的应用
功率谱分析在信号处理领域具有 广泛的应用，如语音信号分析、 雷达信号处理、通信信号处理等 。通过功率谱分析，可以提取信 号的特征信息，实现信号检测、 识别和分类等任务。
未来发展趋势预测
• 高分辨率功率谱估计：随着信号处理技术的发展，对功率谱估计的分辨率要求 越来越高。未来将继续研究高分辨率的功率谱估计方法，以提高信号处理的精 度和性能。
杂波背景下目标检测
在雷达和声呐应用中，接 收到的信号往往包含杂波 ，即非目标反射的信号。 杂波可能来自地面、海面 、大气等环境因素。
功率谱分析可用于区分目 标回波和杂波。目标和杂 波在功率谱上通常具有不 同的特征，如频率范围、 幅度和形状等。
通过设定合适的阈值和滤 波器，可以在杂波背景下 准确地检测出目标。
定义
随机信号是一种无法用确 定函数描述，但具有一定 统计规律性的信号。
统计规律性
随机信号在大量重复观测 下呈现出一定的统计规律 ，如均值、方差等。
连续性
随机信号通常是时间连续 的，可以用连续时间函数 表示。
随机信号分类
根据信号性质分类
01
非平稳随机信号：统计特性随时间变化的 随机信号。
03
02
平稳随机信号：统计特性不随时间变化的随 机信号。
ARMA模型法
将随机信号建模为自回归滑动平均模型（ARMA），通过求解模型参数得到功率谱估计。 该方法适用于短数据和复杂信号，但模型定阶和参数估计较困难。
不同方法比较与选择
性能比较
现代谱估计方法通常具有更高的分辨率和更低的方差，性能优于经典谱估计方法。其中，MEM和MVM在分辨率 和方差性能方面表现较好，而ARMA模型法在处理短数据和复杂信号时具有优势。

x (n ) e- jΞm
n= 0
2
=
lim
N →∞
1 N
X N (Ξ)
2
由于实际得到的随机信号只能是它的一个样本的
片断, 因此只能用有限长的的样本序列来估计功率谱,
这相当于用一个有限宽度 (N ) 的窗函数 Ξ(n) 去乘样
本序列, 于是有 (用离散频率 K 代替 Ξ) :
参 考 文 献 1 熊沈蜀, 周兆英, 金龙, 陈耘 1 工程图矢量化处理系统 1 清华大学学报 (自然科学版) , 2000, 40 (4) : 35～ 38 2 董海卫, 江早, 王永军 1 基于矢量化的二值工程图符号提取算法 1 计算机辅助设计与图形学学报, 2000, 12 (4)
所示)。
仿真与测试
4 结 语
功率谱估计的实现有许多方法, 也有很多具体的 算法可以参阅。M a tlab 提供的算法函数为我们学习设 计谱估计提供了一条可行的方便途径, 但较为有限。我 们可以在熟悉了估计原理之后, 自己动手编写m 文件 来实现。这对具有一定M a tlab 编程经验的人并不难, 这里就不再赘述了。
2N
1 +
1n= - N x (n) x 3 (n + m )
e- jΞm
N
∑ =
lim
N →∞
2N
1 +
1
x (n ) e- jΞm
n= - N
N
∑ x (n) e- jΞm 3
n= - N
N
∑ =
lim
N →∞
2N
1 +
1
x (n ) e- jΞm
n= - N
2
N- 1
∑ =

现代信号处理作业实验题目：设信号)()8.0cos(25.0)47.0cos()35.0cos()(321n v n n n n x ++++++=θπθπθπ，其中321,,θθθ是[]ππ,-内的独立随机变量，v(n)是单位高斯白噪声。

1.利用周期图法对序列进行功率谱估计。

数据窗采用汉明窗。

2.利用BT 法对序列进行功率谱估计，自相关函数的最大相关长度为M=64,128,256,512采用BARTLETT 窗。

3.利用Welch 法对序列进行功率谱估计，50%重叠，采用汉明窗，L=256,128,64。

4.利用Burg 法对序列进行AR 模型功率谱估计，阶数分别为10，13.要求每个实验都取1024个点，fft 作为谱估计，取50个样本序列的算术平均，画出平均的功率谱图。

实验原理：1）。

周期图法：又称间接法，它把随机信号的N 个观察值x N (n)直接进行傅里叶变换，得到X N (e jw ),然后取其幅值的平方，再除以N ，作为对x （n ）真实功率谱的估计。

2^)(1)(jw e X Nw P N per =， 其中∑-=-=1)()(N n jwn N jwN e n x e X 2）。

BT 法：对于N 个观察值x(0),x(1),。

,x(N-1),令x N (n)=a(n)x(n)。

计算r x （m ）为∑--=-≤+=mN n N Nx N m m n x n xN m r 101),()(1)(,计算其傅里叶变换∑-=--≤=MMm jwm xBT N M e m rm v w P 1 ,)()()(^^，作为观察值的功率谱的估计。

其中v(m)是平滑窗。

3)。

Welch 法：假定观察数据是x(n),n=0,1,2...,N-1，现将其分段，每段长度为M,段与段之间的重叠为M-K,第i 个数据段经加窗后可表示为 1,...,1,0 )()()(-=+=M i iK n x n a n x i M其中K 为一整数，L 为分段数，该数据段的周期图为2)(1)(^w X MU w P i M iper =，其中∑-=-=10)()(M n j w n iM i M e n x w X 。

0 0
1500
50
100
真实谱
10
8
6
4
2
0
150
-2
-1
0
1
2
窗函数，长度为10
2000
20
15
10
5
0
-2
-1
0
1
2
窗函数，长度为20
3000
1000 500
1500 1000
500
2000 1000
0
0
0
100
200
300
0
100
200
例子
bia(Sˆx ())
1
2
S x ()W ( )d S x ()
x(n) sin(2*n) u(n)
高斯白噪声
30 20 10
0 0 10 20 30 40 50
N=50
50
40
30
20
10
0
0
50
100
N=100
250 200 150 100
50 0 0 100 200 300 400 500
(1
)D0
(
2
)d
2
1
2N
S()D0 (1
)D0 ( 2
)d
2
其中：
d0
(n)
1, 0,
| n | N 1 other
F
D0 ()
1 2
Var (Sˆ ( ))
1
2N
2
S ()D0 ( )D0 ( )d
E(Sˆ( ))2
N Var[Sˆx ()] E[Sˆx ()] 0
N=500

10
功率谱估计的参数法（现代谱估计）
AR模型 ARMA模型 有理谱模型 MA模型 线谱模型（谐波）
1 H ( z) A( z ) 1
k a z k k 0 p
参数法谱估计
AR模型 (全极点模型)
MA模型(全零点模型) H ( z ) B( z) bk z k
k 0
13
AR模型参数估计法的功率谱估计：
基本原理： 根据随机采样样本x(0), x(1) … x(N－1) 估计随机时间序列的功
率谱密度： p x(n) ak x(n - k ) u (n)
k 1
式中，u (n) ~ N (0, 2 )
2
H ( z)
1 1 ak z
k 1 p -k
12
第五章 随机信号的功率谱估计
功率谱估计的经典和现代方法
AR模型法的功率谱估计
AR模型法的主要性质 Yule-Walker方程的Levinson – Durbin求解算法 格型滤波器 AR模型参数提取算法 噪声对AR谱估计的影响 ARMA和MA模型法简介
白噪声中正弦波频率的估计
14
③ 对各种阶数下的模型进行比较，应用某种准则估计选择 最好的模型(得阶数p、ak及 2 )。 AR(p)模型的Yule - Walker方程组：
R(0) R(1) R (2) R( p)
R(1) R(0) R(1)
R(2) R(1) R(0)
R x (m)称为取样自相关函数。
5
• 自相关法(B-T法）—直接周期图法
1958年，Blackman 和 Tukey提出。先估计信号的自相关函 数，再求出信号的功率谱密度估计 ：

现代信号处理经典的功率谱估计《现代信号处理》姓名：李建强学号：201512172087专业：电子科学与技术作业内容：在MATLAB平台上对一个特定的平稳随机信号进行经典功率谱估计和现代功率谱估计的比较一、前言功率谱估计是信息学科中的研究热点,在过去的30多年里取得了飞速的发展。

在许多工程应用中，它能给出被分析对象的能量随频率的分布情况。

平滑周期图是一种计算简单的经典方法，它的主要特点是与任何模型参数无关，但估计出来的功率谱很难与信号的真是功率谱相匹配。

与周期图方法不同，现代谱估计主要是针对经典谱估计(周期图和自相关法)的分辨率低和方差性能不好的问题而提出的。

其使用参数化的模型，能够给出比周期图方法高得多的频率分辨率。

其内容极其丰富，涉及的学科和领域也相当广泛，按是否有参数大致可分为参数模型估计和非参数模型估计，前者有AR模型、MA模型、ARMA模型、PRONY指数模型等;后者有最小方差方法、多分量的MUSIC方法等。

二、总体概述本次实验分别使用经典的功率谱估计（如周期图法）与AR模型法对某一特定的平稳随机信号进行其功率谱估计，由图像得到信号的频率。

利用MATLAB平台，直观形象地观察并比较二者估计效果的区别，以便于加深对功率谱估计的理解和掌握。

三、具体的实现步骤1、经典法功率谱估计周期图法又称直接法，它是从随机信号x(n)中截取N长的一段，把它视为能量有限的真实功率谱的估计的一个抽样。

1.1、实现步骤（1）、模拟系统输出参数x(n)=A*sin(2πf1*n)+B*sin(2πf2*n)，包括序列长度N（128或512或1024，加性高斯白噪声（AGWN）功率一定，设置A，B，f1，f2，n的值。

（2）、应用周期图法（不加窗）对信号的功率谱密度进行估计，使用直接法在MATLAB平台上进行编程实现。

（3）、输出相应波形图，进行观察，记录。

1.2 MATLAB源代码实现clear all; %清除工作空间所有之前的变量close all; %关闭之前的所有的figureclc; %清除命令行之前所有的文字n=1:1:128; %设定采样点n=1-128f1=0.2; %设定f1频率的值0.2f2=0.213; %设定f2频率的值0.213A=1; %取定第一个正弦函数的振幅B=1; %取定第一个正弦函数的振幅a=0; %设定相位为0x1=A*sin(2*pi*f1*n+a)+B*sin(2*pi*f2*n+a ); %定义x1函数，不添加高斯白噪声x2=awgn(x1,3); %在x1基础上添加加性高斯白噪声，信噪比为3,定义x2函数temp=0; %定义临时值，并规定初始值为0temp=fft(x2,128); %对x2做快速傅里叶变换pw1=abs(temp).^2/128; %对temp做经典功率估计k=0:length(temp)-1;w=2*pi*k/128;figure(1); %输出x1函数图像plot(w/pi/2,pw1) %输出功率谱函数pw1图像xlabel('信号频率/Hz');ylabel('PSD/傅立叶功率谱估计');title('正弦信号x(n)添加高斯白噪声后的,周期图法功率频谱分析');grid;%------------------------------------------------------------------------- pw2=temp.*conj(temp)/128; %对temp做向量的共轭乘积k=0:length(temp)-1;w=2*pi*k/128;figure(2);plot(w/pi/2,pw2); %输出功率谱函数pw2图像xlabel('信号频率/Hz');ylabel('PSD/傅立叶功率谱估计');title('正弦信号x（n）自相关法功率谱估计');grid;1.3 matlab仿真图形（1）、用直接法，功率谱图像,采样点N=128。

雷达和声呐系统
目标检测
在雷达和声呐系统中，经典功率谱估计常被用于目标检测。通过对接收到的信号进行功率 谱分析，可以判断是否存在目标以及目标的位置和速度等信息。
距离和速度测量
在雷达和声呐系统中，经典功率谱估计还可以用于距离和速度测量。通过对接收到的信号 进行功率谱分析，可以估计出目标与系统之间的距离和相对速度。
信号分类
在雷达和声呐系统中，经典功率谱估计还可以用于信号分类。通过对接收到的信号进行功 率谱分析，可以判断目标的类型，例如区分飞机、船舶或车辆等不同类型目标。
05 经典功率谱估计的改进方 法
基于小波变换的功率谱估计
1
小波变换能够将信号分解成不同频率和时间尺度 的分量，从而更好地揭示信号的内在结构和特征。
然而，这些方法通常需要较长 的数据长度和较为复杂的计算 过程，对于短数据和实时处理 的应用场景具有一定的局限性 。
研究展望
01
随着信号处理技术的发展，经典功率谱估计方法仍有进一步优化的空 间。
02
针对短数据和实时处理的应用场景，研究更为快速、准确的功率谱估 计方法具有重要的实际意义。
03
结合机器学习和人工智能技术，探索基于数据驱动的功率谱估计方法 是一个值得关注的方向。
优点
能够提供较高的频率分辨率和较低的估计误差。
原理
格莱姆-梅尔谱估计利用了信号的模型参数，通过 构造一个模型函数来描述信号的频率响应特性， 并求解该函数的极值问题得到信号的功率谱。
缺点
需要预先设定模型函数的形式和参数，且计算复 杂度较高。
03 经典功率谱估计的优缺点
优点
01
02
03
算法成熟
经典功率谱估计方法经过 多年的研究和发展，已经 相当成熟，具有较高的稳 定性和可靠性。

随机信号功率谱估计1.随机信号功率谱密度定义定义随机信号信号的功率谱()jw x e P 为： ()()m j xe m r P ωω-+∞-∞=∑=m j x e其中()m r x 为随机信号的自相关函数。

功率谱反映了信号的功率在频域随频率ω分布，因此()jw x e P 又称为功率谱密度。

2. 经典谱估计（非参数谱估计）方法简介经典谱估计的方法主要包括两种方法：BT 法和周期图法。

(1) BT 法（间接法）此方法的理论基础是维纳-辛钦定理。

1958年Blackman 和Tukey 给出了这一种方法的具体实现，即由()n N x 估计出自相关函数()m rˆ，然后对()m r ˆ求傅里叶变换得到()n N x 的功率谱，记之为()ωBTP ˆ，并以此作为对()ωP 的估计，即 ()()m j m e m r ˆˆωω--=∑=M MBT P ， 1-≤N M因为这种方法求出的功率谱是通过自相关函数间接得到的，所以称为间接法，又称BT 法或自相关法。

当M 较小时，上式计算量不是很大，因此，该方法是在FFT 问世之前（即周期图法被广泛应用之前）常用的谱估计方法。

(2) 周期图法（直接法）周期图法又称为直接法，它是把随机信号()n x 的N 点观察数据()n N x 视为一个能量有限信号，直接取()n N x 的傅里叶变换，得()jw e x N ，然后再取其幅值的平方，并除于N ，作为对()n x 真实功率谱()jwe p 的估计。

以()ωPERP ˆ表示用周期图法估计的功率谱，则 ()()21ˆωωN PER X NP =3. 信号参数设定模拟信号：()()()()t w t *110cos 5.1t *100cos 2t x ++=ππ其中()t w 为均值为0，方差为1的白噪声。

采样频率：z 1000f s Hnfft=1024时间长度分别取1s ，10s ，100s4. 仿真结果原始信号波形:BT 法：（1）t=1s （2）t=10s （3）t=100s周期法：（1）t=1s （2）t=10s （3）t=100s由图可以看出，在频率50hz和55hz附近处功率谱有两个峰值，说明信号中有50hz和55hz 的周期成分。

随机信号的功率谱估计方法随机信号的功率谱估计方法介绍随机信号是指信号的每个值都是随机的，即在同一时刻下，其取值可以是不同的。

由于随机性导致了随机信号的分布不确定，因此分析随机信号的机理比较复杂。

一个优秀的信号分析方法是估计随机信号的功率谱。

功率谱是一个很有用的统计量，它描述了信号在不同频率上的能量分布。

估计功率谱可以帮助我们了解信号的构成、将信号分解成不同的频率分量、对信号的特征进行定量分析，以及在通信和控制系统中使用。

本文将介绍几种常见的随机信号功率谱估计方法，包括周期图法、自相关函数法、半岭功率谱估计法和最大熵谱估计法。

方法一、周期图法周期图法经常用于信号频谱估计。

当我们有大量采样数据时，可以通过对信号进行傅里叶变换来计算功率谱。

但是，当信号是随机过程时，它的频谱也是一个随机变量，因此我们必须通过使用大量的测量值来确定频谱估计的不确定性。

由此带来的问题是，我们要计算的是随机过程信号的平均功率谱密度函数，而不仅仅是单次测量结果的功率谱。

周期图法通过将数据分成多个重叠的子段，然后计算每个子段的傅立叶变换来估计平均功率谱密度函数。

二、自相关函数法自相关函数法采用的是自相关函数相关的频谱估计方法。

通过对随机信号进行卷积，可以获得信号的自相关函数。

自相关函数是指信号与自身的延迟信号的乘积。

自相关函数可以通过傅立叶变换来计算功率谱密度函数。

这种方法可以用于非平稳和平稳信号，并且在信号较长的情况下效果良好。

三、半岭功率谱估计法半岭功率谱估计法是利用谱曲线的形状确定能量的集中程度。

半岭是谱曲线上右侧的谷底点。

我们可以将信号的谱曲线绘制出来，并计算它到半岭的近似功率谱曲线。

该方法可以适用于处理非平稳信号，需要进行多次计算才能获得准确结果。

四、最大熵谱估计法最大熵谱估计法可以通过最小化误差来估计功率谱密度函数。

该方法通过将信号视为时间序列，然后利用最大熵原理来进行谱估计。

最大熵原理是指在不知道任何关于信号的先验信息的情况下，使用最少的假设来描述数据的过程。

第1篇一、实验目的1. 理解经典功率谱估计的原理和方法；2. 掌握BT法、周期图法、Bartlett法和Welch法等经典功率谱估计方法；3. 通过MATLAB仿真，验证各种方法的性能和特点；4. 分析实验结果，总结经典功率谱估计方法的优缺点。

二、实验原理功率谱估计是信号处理中的一个重要方法，用于分析信号的频率成分。

经典功率谱估计方法主要包括BT法、周期图法、Bartlett法和Welch法等。

1. BT法：先估计自相关函数，然后进行傅里叶变换得到功率谱；2. 周期图法：直接对样本进行傅里叶变换，得到功率谱；3. Bartlett法：将信号分成L段，计算每段的自相关函数，然后进行傅里叶变换得到功率谱；4. Welch法：对信号进行分段，计算每段的自相关函数，然后进行傅里叶变换得到功率谱，并对结果进行加权平均。

三、实验环境1. 操作系统：Windows 10；2. 编程语言：MATLAB；3. 实验数据：随机信号样本。

四、实验步骤1. 生成随机信号样本；2. 使用BT法进行功率谱估计；3. 使用周期图法进行功率谱估计；4. 使用Bartlett法进行功率谱估计；5. 使用Welch法进行功率谱估计；6. 对比分析各种方法的估计结果。

五、实验结果与分析1. BT法：BT法是一种较为精确的功率谱估计方法，其估计结果与真实功率谱较为接近。

但是，BT法需要计算样本的自相关函数，计算量较大。

2. 周期图法：周期图法是一种简单易行的功率谱估计方法，但其估计结果存在较大误差。

当样本长度N较大时，周期图法的估计结果逐渐接近真实功率谱。

3. Bartlett法：Bartlett法在Bartlett窗口的宽度较大时，估计结果较为准确。

但是，当Bartlett窗口的宽度较小时，估计结果误差较大。

4. Welch法：Welch法是一种改进的周期图法，通过分段和加权平均，提高了估计精度。

Welch法在估计精度和计算量之间取得了较好的平衡。

一、随机信号原理分析随机信号的古典法谱估计广义平稳随机过程的功率谱是自相关函数的福利叶变换，它取决于无数个自相关函数值。

但对于许多实际问题，可以利用的数据往往是有限的，所以要准确的计算功率谱是不可能的。

比较合理的目标是设法得出功率谱的一个好的估计，这就是功率谱估计。

功率谱估计有两类大的方法：古典谱估计和现代谱估计。

古典谱估计又有相关法估计，周期图法估计，WOSO法，Bartlett法；现代谱估计有Levinson-Durbin算法和Burg算法。

下面对他们分别作介绍。

1．相关法谱估计这种方法以相关函数为媒介来计算功率谱，所以又叫间接法。

这种方法的具体步骤是：第一步：从无限长随机序列x(n)中截取长度N的有限长序列X N（n）。

第二步：由N长序列求（2M-1）点的自相关函数序列。

即这里，m=-(M-1),。

,-1,0,1,。

,M-1，M N，R^在此处键入公式。

x(m)是双边序列，但是由自相关函数的偶对称性式，只要求出m=0,。

,M-1的，另一半也就知道了。

第三步：由相关函数的傅式变换求功率谱。

即以上过程中经历了两次截断，一次是将x(n)截成N长，称为加数据窗，一次是将截成（2M-1）长，称为加延迟窗。

因此所得的功率谱仅是近似值，也叫谱估计，式中的，代表估值。

一般取M〈〈N，因为只有当M较小时，序列傅式变换的点数才较小，功率谱的计算量才不至于大到难以实现，而且谱估计质量也较好。

因此，在FFT问世之前，相关法是最常用的谱估计方法。

当FFT问世后，情况有所变化。

因为截断后的可视作能量信号，由相关卷积定理可得这就将相关化为线性卷积，而线性卷积又可以用快速卷积来实现。

我们可对上式两边取（2N-1）点DFT，则有于是将时域卷积变为频域乘积。

用快速相关求的完整方案如下：1.对N长的充（N-1）个零，成为（2N-1）长的。

2.求（2N-1）点的FFT，得。

3.求。

由DFT性质，是纯实的，满足共轭偶对称，而一定是实偶的，且以（2N-1）为周期。

功率谱和经典谱估计的应用：
1、功率谱的应用： 功率谱反映了随机信号各频率成分功率能量的分布情况，
可以揭示信号中隐含的周期性及靠得很近的谱峰等有用的信息， 应用及其广泛。例如，在语音信号识别、雷达杂波分析、地震 勘测信号处理、水声信号处理、系统辨识中非线性系统识别、 物理光学中透镜干涉、流体力学的内波分析、太阳黑子活动周 期研究等许多领域，发挥了重要作用。
涡街流量计的信号频率与流体速度成线性比例关系，工 程应用中一般测量该信号的频率，然后根据仪表系数转换算成 实际的流量。因为噪声的原因，数字信号处理必须实现准确的 功率—频率计算。对涡街信号处理的第一步就是直接做功率谱 估计，计算功率谱能量最大的谱线对应的信号频率就是涡街信 号的频率。用这个频率来确定涡街信号的区间范围方便后续进 一步处理。
2、经典谱估计的应用：
经典谱估计法由于假定信号的自相关函数在数据观测区以外等于 零，因此估计出来的功率谱很难与信号的真实功率谱相匹配，是一种低 分辨率的谱估计方法，而现在已有很多质量更好的谱估计方法，所以经 典谱现在主要用于一些要求不高的场合，做一些基础的工作。
（1）涡街流量计
在基于经典谱估计改进方法的涡街流量计中通过经典谱估计的FFT 算法来计算信号频率的区间范围，以待后续进一步的处理。
(2）汽轮机振动信号 当汽轮机产生故障时，其振动信号的频谱能量分布情况会有 所改变，因此对振动信号进行频谱分析是当前常用的汽轮机故障 特征提取方法。周期图法

7.1概述
谱估计方法： 经典方法（非参数法），现代方法（参数法）
• 经典方法：以傅里叶变换为基础， 方法：周期图法 和 Blackman-Tukey(BT)法 (自相关序列估计法)； 适用范围：数据多，对频率分辨率要求不高。
• 现代方法：以随机过程的参数模型为基础， 又称参数方法或模型方法； 最基本的方法： 自回归模型法，线性预测法，最大熵法； 适用范围：数据少，对频率分辨率要求高。 优劣：参数方法较优，利用了“随机过程是如何产生的”信息，
m=−
m=0
7.2 功率谱估计的经典方法------周期图法
• S per (e jw ) 可用 X N (e jw ) 表示为：
N −1
S per (e jw ) =
RN (m)e− jwm =
m=−( N −1)
N −1 [ 1
N m=−( N −1)
N −1− m n=0
xN (n) x*N (n + m)]e− jwm
• RN (m) 的傅里叶变换 S per (e jw ) 为：
+
N −1
S per (e jw ) =
RN (m)e− jwm =
RN (m)e− jwm
m=−
m=−( N −1)
表示 。
• xN (m) 的傅里叶变换 X (e jw ) 为：
+
N −1
X (e jw ) =
xN (m)e− jwm = xN (m)e− jwm
第七章 功率谱估计的经典方法
7.1概述 7.2 功率谱估计的经典方法
7.1概述
功率谱估计， • 是估计平稳随机过程的功率谱， • 根据随机过程的一个取样序列的一段数据, 即有限长数据来估计。 • 假定信号是遍历的，建立在时间平均基础上。

已知信号x(t)=2cos(2πf1t)+3cos(2πf2t)+w(t)。

其中，f1=30HZ，f2=110HZ，w(t)为白噪声，采样频率F S=1000HZ。

求其功率谱1.周期图法，直接法又称周期图法，它是把随机序列x(n)的N个观测数据视为一能量有限的序列，直接计算x(n)的离散傅立叶变换，得X(k)，然后再取其幅值的平方，并除以N，作为序列x(n)真实功率谱的估计。

clc;Fs=1000; %采样频率n=0:1/Fs:1;%产生含有噪声的序列xn=2*cos(2*pi*30*n)+3*cos(2*pi*110*n)+randn(size(n));window=boxcar(length(xn)); %矩形窗nfft=1024;[Pxx,f]=periodogram(xn,window,nfft,Fs); %直接法¨plot(f,10*log10(Pxx));xlabel('频率/kHz');ylabel('相对功率谱密度（dB/Hz）');title('直接法');saveas(gcf, 'zhoujifa','jpg');效果图如图1所示图1 直接法功率谱效果改进的直接法2. Bartlett法clcFs=1000;n=0:1/Fs:1;xn=2*cos(2*pi*30*n)+3*cos(2*pi*110*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);xlabel('频率/kHz');ylabel('相对功率谱密度（dB/Hz）');title('直接法');pause;saveas(gcf, 'zhijie','jpg');figure(2)plot(k,[plot_Pxx plot_Pxx-plot_Pxxc plot_Pxx+plot_Pxxc]); xlabel('频率/kHz');ylabel('相对功率谱密度（dB/Hz）');title('Bartlett法');saveas(gcf, 'bartlett','jpg');相对于直接法，Bartlett法运行结果曲线较为平滑,如图2图3所示：图2直接法图3 Bartlett法3 . Welch法clc;Fs=1000;n=0:1/Fs:1;xn=2*cos(2*pi*30*n)+3*cos(2*pi*110*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)subplot(3,1,1),plot(f,plot_Pxx);subplot(3,1,2),plot(f,plot_Pxx1);subplot(3,1,3),plot(f,plot_Pxx2);xlabel('频率/kHz');ylabel('相对功率谱密度（dB/Hz）');title('Welch法');saveas(gcf, 'Welch','jpg');程序运行结果使信号得以突出，而抑制了噪声,如图4所示：图4 Welch法功率谱效果。

∞
=
Ryy (m) =
p =−∞
∑R
k = −∞ ∞
∑ h( k ) R
xx
∞
xx
( m − k ) = Rxx (m) ∗ h( m)
(m − p) Rhh ( p) = Rxx (m) ∗ Rhh (m)
或
= Rxx (m) ∗ h(m) ∗ h(−m) = Rxy (m) ∗ h(−m)
S yy (e jω ) = S xy (e jω ) H (e− jω )
jω
jω
jω
2
离散随机信号通过线性非移变系统
（4）输入随机过程与输出随机过程的互相关序列Rxy(m)
∞ Rxy ( m) = E [x ( n) y ( n + m) ] = E x ( n) ∑ h( k ) x ( n + m − k ) k = −∞
=
k = −∞
∑ h(k ) E[x(n) x(n + m − k )]
ˆ B =α − E [ α ]
无偏估计， 无偏估计 有偏估计，当观测数据为无穷时B = 0，则称其为渐 渐 B = 0时无偏估计 B ≠ 0 有偏估计 进无偏估计。无偏估计和渐进无偏估计又称为是好估计 进无偏估计 好估计。 好估计
均值 均方值
E[xn ] = mxn = ∫ xpxn ( x, n)dx
∞ −∞
E x = ∫ x 2 pxn ( x, n)dx
2 n −∞
[ ]
2
∞
方差
E xn − mxn
[(
) ]= σ
2 xn
=∫
∞ −∞
(x − m )
xn