西南交通大学信号处理期末作业

1、考虑两个谐波信号()x t 和()y t ，其中()cos()c x t A w t φ=+,()cos()c y t B w t =式中A 和c w 为正的常数,φ为均匀分布的随机变量，其概率密度函数为1,02()20,f φπφπ⎧≤≤⎪=⎨⎪⎩其他, 而B 是一个具有零均值和单位方差的标准高斯随机变量，即其分布函数为2()/2),B f b b b =--∞<<∞（1）求()x t 的均值()x u t 、方差2()x t σ、自相关函数()x R τ和自协方差函数()x c τ。

（2）若φ与B 为相互统计独立的随机变量，求()x t 和()x t 的互相关函数()xy R τ与互协方差函数()xy c τ。

解： （1）()x t 的均值()x u t 为：22011(())(cos())cos()sin() 022x c c c u E x t E A w t A w t d A w t ππφφφφππ==+=+=+=⎰方差2()x t σ为：22222222(())(cos ())(1cos(22))(cos(22))2222xc c c A A A A E x t E A w t E w t E w t φσφφφ==+=++=++=自相关函数()x R τ为：22222()(cos()cos(+))(cos()cos(+))(cos(2+22)cos())cos()(cos(2+22))cos()2222x c c c c c c c c c c c c c R E A w t A w t w A E w t w t w A A A A E w t w w w E w t w w τφτφφτφτφτττφτ=++=++=++=++= 自协方差函数()x c τ为：2()()cos()2x x c A c R w τττ==（2）()y t 的均值为：()(())(cos())()cos()0y B B c c u t E y t E B w t E B w t ====,所以()=0E B由互相关函数的定义可知：()(cos()cos())xy c c c R E A w t B w t w τφτ=+-由题意知道φ与B 为相互统计独立的随机变量，所以有()(cos()(cos())(cos()()cos())00cos()0xy c c c c c c c c R E A w t E B w t w AE w t E B w t w A w t w τφτφττ=+-=+-=⨯+⨯-=互协方差函数()xy c τ()()0xy xy c R ττ==2.接收信号由下式给出：cos(),1,2,...,i c i y A i i N ωθω=++=，式中~(0,1)i N ω即i ω是零均值和单位方差的高斯噪声，c ω为载波角频率，而θ是未知的相位。

假设12,,...N ωωω相互独立，求未知相位的最大似然估计^ML θ。

解：由于12,,...N ωωω相互独立，所以1,..N y y 也相互独立并且服从高斯分布，可以得到1,..N y y 与θ的联合概率密度函数分布21[cos()]21(,..|)Ni c i y A i N f y y ωθθ=-+-∑由此，可以得到似然函数211ln(2)[cos()]22Ni c i N L y A i πωθ==---+∑该似然函数对θ求偏导，并令该偏导函数为零，即可得到如下公式： 1[cos()]sin()0Ni c c i L y A i i ωθωθθ=∂=-++=∂∑ 因此，最大似然估计^ML θ 为上述函数的零点值。

则^^^11cos()sin()sin()NNMLML ML cc i c i i A i i y i ωθωθωθ==++=+∑∑该函数为非线性方程，不容易求解，若忽略双倍频率2c ω ，则可简化到如下式子：^1sin()0Nici y i ωθ=+=∑根据三角公式分解得到如下式子：^^11sin cos()cos sin()NNML i c ML i c i i y i y i θωθω===-∑∑由此，可以得到如下公式^11sin tan cos()Nici ML Nici y iy i ωθω===-∑∑所以相位的最大似然估计如下：^11sin()arctan()cos()Nici ML Nici y i y i ωθω===-∑∑3.离散时间的二阶AR 过程由差分方程12()(1)(2)()x n a x n a x n w n =-+-+ 描述，式中()w n 是一零均值、方差为2w σ 的白噪声。

证明()x n 的功率谱为22212122()12(1)cos(2)2cos(4)wx P f a a a a f a f σππ=++---证明：由AR 过程的功率谱公式知222412()1x j fj f P f a ea eωππσ--=--其中2242424121212224422221221212122121221(1)(1)1()()12(1)cos(2)2cos(4)j f j fj f j f j f j f j f f j f j f j f j f a e a e a e a e a e a e a a a e e a a e a a e a e e a a a a f a f ππππππππππππππ--------=----=++-+++-+=++--- 将其带入第一个公式可得：22212122()12(1)cos(2)2cos(4)x P f a a a a f a f ωσππ=++---4、信号的函数表达式为：()()()()sin(2100) 1.5sin(2300)sin(2200)x t t t A t t dn t n t πππ=++++，其中，()A t 为一随时间变化的随机过程，()dn t为经过390-410Hz带通滤波器后的高斯白噪声，()n t为高斯白噪声，采样频率为1kHz，采样时间为2.048s。

分别利用周期图谱、ARMA、Burg最大熵方法估计信号功率谱，其中ARMA方法需要讨论定阶的问题。

解：由题意知采样点数一共为：1000×2.048=2048个数据点。

()A t为一随时间变化的随机过程，由于随机过程有很多类型，如维纳过程、正态随机过程，本文采用了均值为0，方差为1的正态随机过程来作为演示，来代替()A t，高斯白噪声采用强度为2的高斯白噪声代替()dn t。

其中滤波器采用的是契比雪夫数字滤波器。

n t，其带通滤波后为()可得到x(t)如下图所示：1、周期图法matlab中的周期图功率谱法原理是通过计算采样信号的FFT，获得离散点的幅度，再根据幅度与功率之间的关系，转换为离散点的功率，再通过坐标变换将离散点的功率图转换为连续功率谱密度。

Step1:计算采样信号x(n)的DFT，使用FFT方法来计算。

如果此处将复频率处的幅度对称到物理实际频率，得到的就是单边谱，否则就是双边谱Step2:根据正余弦信号功率与幅度的关系以及直流功率与幅度的关系，将幅度转换为离散功率谱。

Step3:对横纵坐标进行转换，横坐标乘以频率分辨率转换为实际连续物理频率，纵坐标除以频率分辨率转换为功率谱密度。

调用MA TLAB中自带的matlab中[Pxx,f]=periodogram(x,window,nfft,fs)函数可得计算结果如下：2、ARMA方法参数模型估计的思想是:①假定研究的过程X(n)是一个输入序列u(n)激励一个线性系统H(z)的输出。

②有已知的X(n)，或其自相关函数来估计H(z)的参数。

③由H(z)的参数来估计X(n)的功率谱。

不论X(n)是确定性信号还是随机信号，u(n)与X(n)之间总有如下输入输出关系：1()()()p qk k k k x n a x n k b u n k ===--+-∑∑()()()k x n h k u n k ∞==-∑对以上两个式子两边分别取Z 变换，并假定b 0=1，可得()()()B z H z A z =其中1()1p kk k A z a z -==+∑,1()1q kk k B z b z -==+∑,0()()k k H z h k z ∞-==∑。

为了保证H(z)是稳定的最小相位系统，A(z)和B(z)的零点都应该在单位圆内。

假定u(n)是一个方差为2σ的白噪声序列，由随机信号通过线性系统的理论可知，输出序列X(n)的功率谱为：222**2()()()()()()()jw jw jwjw x jw jw jwB e B e B e P e A e A e A e σσ==ARMA 阶数确定：本题目采用AIC 准则确定ARMA 的阶数。

分别计算p 、q 从1到20阶数的计算出AIC （p,q ）,如下图所示，当横坐标大概为230左右时，AIC(p,q)取得最小，将此时的p,q 作为带入到模型即可。

ARMA 法谱估计结果：3、Burg 最大熵法Burg 算法的具体实现步骤：步骤1 计算预测误差功率的初始值和前、后向预测误差的初始值，并令m = 1。

210)(1∑==Nn n x NP )()()(00n x n g n f ==步骤2 求反射系数∑∑+=--+=*---+--=Nm n m m Nm n m m m n g n f n g n fK 12121111])1()([21)1()(步骤3 计算前向预测滤波器系数),()()(11i m a K i a i a m m m m -+=*-- 1,...,1-=m imm K m a =)(步骤4 计算预测误差功率12)1(--=m m m P K P步骤5计算滤波器输出)1()()(11-+=--n g K n f n f m m m m )1()()(11-+=--*n g n f K n g m m m m步骤6 令m ←m+1，并重复步骤2至步骤5，直到预测误差功率Pm 不再明显减小。

最后，再利用Levinson 递推关系式估计AR 参数，继而得到功率谱估计。

Burg 最大熵法谱估计结果如下图：5.附件中表sheet1为某地2008年4月28日凌晨12点至2008年5月4日凌晨12点的电力系统负荷数据，采样时间间隔为1小时，利用Kalman 方法预测该地5月5日的电力系统负荷，并给出预测误差（5月5日的实际负荷数据如表sheet2）。

解：卡尔曼滤波是以最小均方误差作为估计的最佳准则，来寻求一套递推估计的算法，其基本思想是：采用信号与噪声的状态空间模型，利用前一时刻地估计值和现在时刻的观测值来更新对状态变量的估计，求得出现时刻的估计值。