现代谱估计计算机仿真实验报告.doc

资源描述

现代谱估计计算机仿真实验报告.doc

《现代谱估计计算机仿真实验报告.doc》由会员分享，可在线阅读，更多相关《现代谱估计计算机仿真实验报告.doc（21页珍藏版）》请在冰豆网上搜索。

现代谱估计计算机仿真实验报告.doc

现代谱估计计算机仿真实验报告

胡敏

在许多工程应用中，利用观测到的一组样本数据估计并分析一个平稳随机信号的功率谱密度是十分重要的。

例如，在雷达信号处理中，由回波信号的功率谱密度、谱峰的宽度、高度和位置，可以确定目标的位距离和运动速度；在阵列信号处理中，空间功率谱描述了信号功率随空间角度的分布情况。

在许多信号处理应用中，谐波过程经常会遇到，它对应的功率谱为线谱，谐波过程的功率谱估计就是要确定谐波的个数，频率和功率（合称谐波恢复）。

为了更好的学习现代信号处理中该部分的内容，我们做了相应的计算机仿真实验。

1实验目的

1、深入理解现代谱估计和谐波恢复的基本理论，包括ARMA模型，ARMA谱估计，ARMA模型识别，Pisarenko谐波分解，信号子空间和噪声子空间，旋转不变技术（ESPRIT）；

2、熟悉与上述谱估计和谐波恢复理论相关的数学方法以及各自的特点，包括最小二乘估计（LS），奇异值分解（SVD），总体最小二乘估计（TLS），特征值分解和广义特征值分解；

3、体会ARMA功率谱估计中的Cadzow谱估计子和Kaveh谱估计子，ARMA模型的识别方法，Pisarenko谐波恢复方法，ARMA建模谐波恢复方法，MUSIC方法进行谐波恢复，两种ESPRIT方法（LS-ESPRIT和TLS-ESPRIT进行谐波恢复；

2实验原理

2.1ARMA谱估计

相当多的平稳随机过程都可以通过用白噪声激励线性时不变系统来产生，而线性系统又可以用线性差分方程进行描述，这种差分模型就是自回归—滑动平均（ARMA）模型。

而且，任何一个有理式的功率谱密度都可以用一个ARMA随机过程的功率谱密度精确逼近。

ARMA随机过程定义为服足下列线性差分方程的离散随机过程：

（1）

式中是一离散白噪声；式

（1）所示的差分方程称为ARMA模型，系统和分别称为自回归（AR）参数和滑动平均（MA）参数，而和分别叫做AR阶数和MA阶数。

ARMA过程可以简记为ARMA，使用移位算子可以把它写作如下形式：

（2）

其中：

若，则平稳ARMA过程的功率谱密度为：

（3）

用（3）式进行谱估计必须事先辨识ARMA模型和激励噪声的方差，而MA参数的估计需要求解非线性方程。

为了避开非线性运算，Cadzow和Kaveh分别提出了一种线性谱估计子：

1、Cadzow谱估计子

（4）

（5）

（6）

（7）

其中为的协方差函数。

因此用Cadzow谱估计子只需要确定AR阶数和AR参数就能进行ARMA谱估计。

2、Kaveh谱估计子

（8）

（9）

Kaveh谱估计子需要确定AR阶数，AR参数和MA阶数来进行ARMA谱估计。

3、AR阶数和参数的确定

对于一个ARMA过程，可以推出其相关函数满足如下方程：

（10）

其中为系统的冲激响应，根据其定义可以得到：

（11）

由（10）式和（11）式就能推导得到著名的修正Yule-Walker方程，简称MYW方程：

，（12）

若已知AR阶数，就能通过求解个MYW方程来求解AR参数：

（13）

其中：

该方程可以用最小二乘法估计的值：

（14）

而实际问题中，AR阶数往往是未知的，此时可用奇异值分解法确定AR阶数，总体最小二乘估计AR参数，合称SVD-TLS算法。

考虑超定方程：

（15）

其中：

，，。

若，就可以通过对奇异值分解：

（16）

中包含个奇异值，将其归一化：

，（17）

选择一个接近于零的数作为阈值，把大于此值得最大整数作为有效秩，它就是AR阶数。

根据总体最小二乘方法可以得到矩阵：

（18）

其中：

是酉矩阵第列的一个加窗段，定义为：

（19）

由的可以估计AR参数为：

，（20）

4、MA阶数和参数的确定

在AR定阶和参数估计的SVD-TLS算法中，取，令，构造超定矩阵：

。

计算其SVD，计算比值：

（21）

式中是对应的第个奇异值，若大于某个给定的阈值，则接受。

在推导Kaveh估计子的过程中可以得到一组非线性方程组，使用Newton-Raphson算法求解该方程组可以得到MA参数，其过程如下：

（1）、令MA参数的初始值为：

，，，

（2）、计算迭代的拟合误差函数：

，（22）

式中可由（9）式求解。

（3）计算矩阵：

（4）更新MA参数估计向量：

（23）

（5）判断MA参数估计向量是否收敛，若收敛，则迭代停止，若发散，令返回

（2）计算，直到MA参数估计收敛为止。

2.2Pisarenko谐波分解理论

考虑由个无重复频率的正弦波组成的过程：

（24）

是在内均匀分布的随机数，则其个频率由特征多项式：

，，（25）

的对共轭复根决定。

一般在加性白噪声中观测该过程，所以观测过程是一个特殊的ARMA过程，且AR参数和MA参数完全相等。

在使用和统计独立的假设下，可以得到一个重要的法方程：

（26）

其中：

这表明是观测过程的自相关矩阵的特征值，其对应的特征向量正好是特征多项式（25）的系数，Pisarenko谐波分解法的思想就是找出最小的特征值并将其对应的特征向量代入（25）式以求得对共轭复根，再由下式确定频率：

（27）

2.3ARMA建模法谐波恢复

2.2中分析的ARMA过程不满足MYW方程的条件，但可以推导出其服从的法方程和MYW方程的形式是一致的，所以谐波恢复的ARMA建模算法如下：

（1）利用观测数据计算样本自相关矩阵：

其中：

，；

（2）用SVD-TLS算法确定AR阶数和系数向量的总体最小二乘估计；

（3）计算特征多项式（25）的共轭根对（，有（27）式确定频率。

2.4MUSIC方法

考虑白噪声中的个谐波信号

，（28）

，引入以下向量：

则有：

（29）

对进行特征值分解得到：

（30）

式中是无加性噪声时信号自相关矩阵的特征值。

所以的特征值由个信号特征值和个噪声特征值组成，与此对应把特征向量矩阵的列向量分为两部分，即

分别称为和分别由信号特征向量和噪声特征向量组成；由和分别张成的空间分别叫做信号子空间和噪声子空间。

可以推导出：

（31）

所以用MUSIC方法进行谐波恢复的思想是：

以很小的步长对进行搜索，寻找

（32）

或（33）

的个极大值点，其对应的频率就是所求的谐波频率。

2.5旋转不变技术ESPRIT

2.4中描述的问题，引入向量：

于是有：

（34）

（35）

酉矩阵称为旋转算符，在上面描述的过程是平移的结果，可以看作是最简单的旋转。

向量和的互相关矩阵为：

（36）

其中：

因为平移保持了和信号子空间的不变性，所以可以构造矩阵束：

其中是的最小特征值。

对矩阵束进行广义特征值分解，其特征值矩阵与矩阵有如下关系：

（37）

基本ESPRIT算法（LS-ESPRIT算法）的思想就是用矩阵束的广义特征值矩阵的前个特征值来估计谐波频率，计算公式使用（27）式。

TLS-ESPRIT算法的思想是：

对进行奇异值分解，确定有效秩，并存储与个主奇异值对应的，和；求的广义特征值分解，得到个广义特征值，它们给出了个谐波频率。

3实验内容

仿真的观测数据由下式给出：

其中是一列高斯白噪声，。

进行如下各项实验：

1、取AR阶数分别为4和6，用最小二乘法估计AR参数，然后使用Cadzow谱估计子进行功率谱估计，并试根据该谱确定谐波频率；

2、假定AR阶数未知，用SVD-TLS方法确定AR阶数和参数，然后使用Cadzow谱估计子进行功率谱估计，并试根据该谱确定谐波频率；

3、用ZHANG方法确定MA阶数，用Kaveh谱估计子进行功率谱估计，用Newton-Raphson方法估计MA参数，结合2确定ARMA模型，计算ARMA功率谱，并试根据该谱确定谐波频率；

4、用Pisarenko谐波恢复方法确定谐波频率；

5、用ARMA建模谐波恢复方法确定谐波频率；

6、用MUSIC方法确定谐波频率；

7、用LS-ESPRIT方法确定谐波频率；

8、用TLS-ESPRIT方法确定谐波频率。

4实验过程

1、编写基于最小二乘法和Cadzow谱估计子的计算机仿真程序lsestimate.m（代码见附录1），独立运行程序5次，表1给出了频率的估计数据。

从表中数据可以看出第三次程序运行效果比较理想，图1所示的是这次仿真得到的功率谱估计结果。

表1LS法+Cadzow估计子得到的频率估计数据

AR阶数

MA阶数

第1次

第2次

第3次

第4次

第5次

0.2031

0.2891

0.2031

0.3125

0.2031

0.1953

0.2031

0.2656

0.2031

0.1953

0.2031

0.1953

0.2031

0.2578

0.2031

0.1953

0.2031

0.1484

0.2031

0.1953

0.2031

0.1953

0.2

展开阅读全文