西南交大现代信号处理作业

现代信号处理作业1.(5″)证明下面定理：任何一个无偏估计子方差的下界叫作Cramer-Rao下界定理：令1(,,)Nxxx为一样本向量，(|)fx是x的条件密度，若ˆ是的一个无偏估计子，且(|)/fx存在，则221ˆˆvar()()[ln(|)]EEfx式中ˆln(|)()()fxK。其中()K是的某个不包含x的正函数。2.(10″)Wiener滤波是信号处理中最常用和基础的波形估计工具之一，对其在自己研究领域的应用情况进行一个简单综述。3.(5″)二阶滑动平均过程由2()()1(1)2(2),{()~(0,)}xnwnbwnbwnwnN定义，式中2(0,)N表示正态分布，其均值为零、方差为2。求x(n)的功率谱。4.(20″)信号的函数表达式为：()sin(2100)sin(2300)()sin(2200)()()xtttAttdntnt，其中，A(t)为一随时间变化的随机过程，dn(t)为经过390-410Hz带通滤波器后的高斯白噪声，n(t)为高斯白噪声，采样频率为1kHz，采样时间为2.048s。(1)利用现代信号处理知识进行信号的谱估计；(2)利用现代信号处理知识进行信号的频率提取；(3)分别利用Wiener滤波和Kalman滤波进行去噪；(4)利用Wigner-Ville分布分析信号的时频特征。5.(10″)附件中表sheet1为某地2008年4月28日凌晨12点至2008年5月4日凌晨12点的电力系统负荷数据，采样时间间隔为1小时，利用ARMA方法预测该地5月5日的电力系统负荷，并给出预测误差（5月5日的实际负荷数据如表sheet2）。1、定理：令1(,,)Nxxx为一样本向量，(|)fx是x的条件密度，若参数估计ˆ是真实参数的一个无偏估计子，且(|)/fx、22(|)/fx存在，则ˆ的均方误差所能达到的下界（称为Cramer-Rao下界）等于Fisher信息的导数，即：221ˆˆvar()()[ln(|)]EEfx（1-1）不等式中等号成立的充分必要条件是：ˆln(|)()()fxK（1-2）其中()K是的某个正函数，与样本1(,,)Nxxx无关。证明：由假设条件知，ˆ{}E或ˆ{}0E，因此有1ˆˆ{}()()0NEfxdxdx（1-3）对上式两边求关于的偏导，得ˆˆˆ{}()()[()()]0Efxdxfxdx即有ˆ()()()0fxdxfxdx（1-4）另一方面，由复合函数的求导法，又有()[ln()]()fxfxfx（1-5）由于()fx是x的条件概率密度，故()1fxdx（1-6）将式（1-5）和式（1-6）带入式（1-4），得ˆ[ln()]()()1fxfxdx或改写作ˆ[ln()()][()()]1fxfxfxdx（1-7）由Cauchy-Schwartz不等式知，对于任意两个复函数()fx和()gx，恒有不定式：222()()()()fxgxdxfxdxgx成立，并且当且仅当()()fxcgx，等号成立，将Cauchy-Schwartz不等式应用于式（1-7），则有22ˆ[ln()]()()()1fxfxdxfxdx或等价为221ˆ()()[ln()]()fxdxfxfxdx（1-8）由Cauchy-Schwartz不等式等号成立的条件知：当且仅当ˆln()()()()()fxfxKfx，即式（1-2）成立时，不等式（1-8）才取等号。注意到ˆ{}E，故有22ˆˆˆvar(){()}()()Efxdx（1-9）另由公式{()}()()Egxgxfxdx知22{ln()}ln()()Efxfxfxdx（1-10）将式（1-9）和式（1-10）代入式（1-8），直接得到不等式（1-1）。根据前面的分析，不等式等号成立的充分必要条件是式（1-2）成立3、解：对()xn取延迟形式：12()()(1)(2)xnwnbwnbwn于是有：1212[()()]{[()(1)(2)][()(1)(2)]}ExnxnEwnbwnbwnwnbwnbwn展开上式得到：122111222122[()()][()()()(1)()(2)(1)()(1)(1)(1)(2)(2)()(2)(1)(2)(2)]ExnxnEwnwnbwnwnbwnwnbwnwnbwnwnbbwnwnbwnwnbbwnwnbwnwn即：221211221122()(1)()()(1)(2)()(1)(2)x对上式进行傅里叶变换，则有：22221211221122()(1)()()()()()()()jwjwjwjwx从而得到：2222222221211221122222222212112222222121122()(1)()()(1)()()()(1)2()cos2cos2jwjwjwjwxjwjwjwjwPwbbbbbebebbbebebbbbbeebeebbbbbwbw即22222121122()(1)2()cos2cos2xPbbbbbb4、music算法原理：music(MultipleSignalClassification)算法是针对多元天线阵测向问题提出的。假定M元的均匀线阵，阵元间距为d，信号的工作波长为λ。空间信号源共有D个，各信号不相关，各阵元的噪声(),m1,2,,Mmnt互不相关，噪声和信号(),1,2,,kStkD也不相关。因此，第m个阵元的输出为(1)1()()()kDjmmkmkxtStent（4-1）式中2sinkkd，k为第ｋ个信号源的方向。将式（4-1）写成矩阵形式：()()()XtAStNt（4-2）式中：12[(),(),,()]DAaaa、(1)()[1,,,]kkjTjMTTkaee。求各阵元输出的相关矩阵，有：2()()HHREXtXtAPAI（4-3）()()HPEStSt（4-4）式中：2———噪声的方差。对式（4-3）的相关矩阵R作特征分解，其各特征值及其相对应的特征向量分别为：λ1≥λ2≥…≥λＤ≥λＤ+1≥…≥λＭ（4-5）ｖ1ｖ2…ｖＤｖＤ+1…ｖＭ（4-6）据式式（4-3），可得以下结论：（1）R的最小特征值等于2，重数为(MD)，即λＤ+1=…=λＭ=2。据此，空间信号源的个数D可由下式得出：D=M-（R最小特征值的重数）。最小重数为1，因此，M阵元可测向的信号源数目的最大值为max1DM。（2）各特征向量相互正交。这些向量为矩阵Ｒ列空间的基，由于最小特征值为噪声的贡献，因此与最小特征值对应的那些特征向量所张成的子空间也是噪声的贡献，称之噪声子空间，记为N。这样Ｒ的列空间被划分成两个子空间,即信号子空间S和噪声子空间N：1,,NDMspanvv（4-7）1,,SDspanvv（4-8）由于各特征向量相互正交，故有：SN。在信号源所在方向上，诸方向向量(),1,,kakD，均处于信号子空间S中，故()kNa。构造矩阵：1[,,]NDMEvv（4-9）显然有()0,1,,NkEakDmusic算法就是根据式（4-9）来求空间谱()MUP，有221()()MUHNPEa（4-10）谱峰所对应θ值就是信号源方向的估值。维纳滤波算法原理：维纳(Wiener)是用来解决从噪声中提取信号的一种过滤（或滤波）方法。这种线性滤波问题，可以看做是一种估计问题或一种线性估计问题。一个线性系统，如果它的单位样本响应为()hn，当输入一个随机信号()xn，且：()()()xnsnvn（4-1）其中：()xn表示信号，()vn表示噪声，则输出()yn为：y()()()mnhnxnm（4-2）我们希望()xn通过线性系统()hn后得到的()yn尽量接近于()sn，因此称()yn为()sn的估计值，用()sn表示，即：y()()nsn（4-3）则维纳滤波器的输入—输出关系可用下面图1表示。图4-1维纳滤波器的输入—输出关系实际上，式（4-2）所示的卷积形式可以理解为从当前和过去的观察值(),(1),,()xnxnxnm，来估计信号的当前值()sn。因此，用()hn进行过滤问题实际上是一种统计估计问题。一般地，从当前的和过去的观察值(),(1),,()xnxnxnm，估计当前的信号值y()()nsn成为过滤或滤波；从过去的观察值，估计当前的或者将来的信号值y()()(N0)nsnN称为外推或预测；从过去的观察值，估计过去的信号值y()()(N0)nsnN称为平滑或内插。因此维纳滤波器又常常被称为最佳线性过滤与预测或线性最优估计。这里所谓的最佳与最优是以最小均方误差为准。如果我们分别以()sn与()sn表示信号的真实值与估计值，而用()en表示他们之间的误差，即：()()()ensnsn（4-4）显然()en可能是正值，也可能是负值，并且它是一个随机变量。因此，用它的均方误差来表达误差是合理的，所谓均方误差最小即它的平方的统计期望最小。卡尔曼滤波算法原理：卡尔曼滤波是基于状态空间方法的一套递推滤波算法，在状态空间方法中，引入了状态变量的概念。卡尔曼滤波的模型包括状态空间模型和观测模型。状态模型是反映状态变化规律的模型，通过状态方程来描写相邻时刻的状态转移变化规律；观测模型反映了实际观测量与状态变量之间的关系。Kalman滤波问题就是联合观测信息及状态转移规律来得到系统状态的最优估计。假设动态系统的状态空间模型为(1)()()XtXtWt（4-1）()()()YtHXtVt（4-2）其中，X(t)为系统在时刻t的状态；Y(t)为对状态的观测值；W(t)为系统噪声，方差阵为Q；V(t)为观测噪声，方差阵为R；为状态转移矩阵；H为观测矩阵；为系统噪声驱动矩阵。卡尔曼滤波的计算流程为：计算状态估计值：ˆˆ(1|1)(1|)(1)(1)XttXttKtt（4-3）计算状态一步预测：ˆˆ(1|)(|)XttXtt（4-4）计算新息：ˆ(1)(1)(1|)tYtHXtt（4-5）计算卡尔曼滤波增益：1(1)(1|)[(1|)]TTKtPttHHPttHR（4-6）计算一步预测均方误差：(1|)(|)TTPttPttQ（4-7）计算一步预测估计均方误差：(1|1)[(1)](1|)nPttIKtHPtt（4-8）下面给出卡尔曼滤波的系统模型框图：Matlab图形：图1原始型号时域波形图图2music方法对信号进行谱估计图3信号频率提取图4两种滤波方式时域结果图图5Wiener滤波去噪信号频域波形图6Kalman滤波去噪信号频域波形图7信号Wigner-Ville分布图8信号伪Wigner-Ville分布Matlab程序：cl