现代数字信号处理笔记

现代数字信号处理

Chapter2：维纳滤波卡尔曼滤波

2.3 信号的自相关序列$R_{xx}=0.8^{|m|}, m=0,\pm1,\pm2,\cdots$。观察信号：$x(n)=s(n)+v(n)$。其中$v(n)$为零均值方差为0.45的白噪声,$s(n),v(n)$统计独立。设计长度N=3的FIR滤波器$x(n)$。使其输出$\hat s(n)$与$s(n)$ 的差的均方差值最小。

解法：

(1).由 $x(n)=s(n)+v(n)$推出$x(n+m)=s(n+m)+v(n+m)$

$x(n) \cdot x(n+m) = [s(n)+v(n)] \cdot [s(n+m)+v(n+m)]=\underline {s(n)\cdot s(n+m)}+\underline{v(n) \cdot v(n+m)} + s(n) \cdot v(n+m)+v(n) \cdot s(n+m)$

其中下划线的表示自相关，所以也就是推出

$R_{xx}(m)=R_{ss}(m)+R_{vv}(m)+0 \\ \qquad \quad= R_{ss}(m)+0.45\sigma(m) \\ \qquad \quad=0.8^{|m|}+0.45\sigma(m)$

计算出:

$R_{xx}(0)=1+0.45=1.45 \\ R_{xx}(1)=0.8 \\ R_{xx}(2)=0.64$

计算标准方程互相关$P=R \cdot h ​$ 推出$h_{opt} =R^{-1} \cdot P​$

计算公式：

$$\begin{bmatrix} R_{sx}(0) \\ R_{sx}(1) \\ \vdots \\ \vdots \\ R_{sx}(N-1) \end{bmatrix}=\begin{bmatrix} R_{xx}(0) & R_{xx}(1) & \cdots & R_{xx}(N-1) \\ R_{xx}(1) & R_{xx}(0) & \cdots & R_{xx}(N-2) \\ \vdots & R_{xx}(1) & \cdots & \vdots \\ \vdots & \vdots & \cdots & \vdots \\ R_{xx}(N-1) & R_{xx}(N-2) & \cdots & R_{xx}(0) \\ \end{bmatrix} \begin{bmatrix} h(0) \\ h(1) \\ \vdots \\ \vdots \\ h(N-1) \end{bmatrix}$$
其中， $P=R_{sx}(n)=E[s(n)x(n)]=R_{ss}(n)+R_{sv}(n)=\underline {R_{ss}(n)+0}$
$$\begin{vmatrix} R_{xx}(0) & R_{xx}(1) & R_{xx}(2) \\ R_{xx}(1) & R_{xx}(0) & R_{xx}(1) \\ R_{xx}(2) & R_{xx}(1) & R_{xx}(0) \\ \end{vmatrix} \begin{vmatrix}h(0)\\h(1)\\h(2)\end{vmatrix}= \begin{vmatrix}R_{ss}(0) \\ R_{ss}(1) \\ R_{ss}(2) \\ \end{vmatrix}$$

WF & KF 计算

已知信号功率谱为$S_{ss}(z)=\frac {0.36}{(1-0.8z^{-1})(1-0.8z)}$，测量时，混入了加性噪声$v(n)$，及测量信号为$x(n)=s(n)+v(n)$。其中$v(n)$方差为1均值为0的白噪声，与$s(n)$不相关。设计一因果IIR WF对$x(n)$处理以得到最佳估计。

解法：

一、给出已知条件$S_{sx}(z),S_{xx}(z)$

可以由以下几种方式得出：

1. 隐含的关系：$S_{xx}(z)=S_{ss}(z)+S_{vv}(z)$ 和 $S_{sx}(z)=S_{ss}(z)+0$ 推导

2. 给出$R_{ss}(n)$和$R_{vv}(n)$求Z变换

3. 信号模型： $s(n)=as(n-1)+\underline {w(n)} \quad |a|<1,|c|<1$

   测量模型：$x(n)=cs(n)+\underline {v(n)}$

   $w(n)$：白噪声，方差Q

   $v(n)$：测量白噪声，方差R

二、谱分解

求出$S_{ss}(z),S_{sx}(z),S_{xx}(z)$ 的功率谱

$$S_{ss}(z)= \frac{Q}{(1-az^(-1))(1-az)} \\ S_{sx}(z)=S_{s(cs+v)}(z)=CS_{ss}(z)+S_{sv}(z) \\ S_{xx}(z)=C^2S_{ss}+S_{vv}(z)= \frac{c^2Q}{(1-az^(-1))(1-az)} + R$$
根据谱分解定理：
$$S_{xx}=\sigma_{\varepsilon} ^2 B(z)B(z^{-1})=\sigma_\varepsilon ^2 \frac {1-f \cdot z^{-1} }{1-az^{-1}} \frac {1-fz}{1-az} \qquad |f|<1$$
解出 $f$和$a$之后，求出 $B(z)$和 $B(z^{-1})$

三、求因果部分

$$\frac {S_{sx}(z)}{B(z^{-1})} = \frac{cQ}{(1-az^(-1))(1-az)}/{\frac {1-fz}{1-az}} \\ = \frac {A}{1-az^{-1}} + \frac {Bz}{1-fz}$$

四、求$H_c(z)$

求$H_c(z)$之前我们需要先求出$\frac {S_{sx}(z)}{B(z^{-1})}$。也就是求出其中的因果部分。

得到$[\frac {S_{sx}(z)}{B(z^{-1})}]_+$之后，

$$H_c(z) = \frac {1}{\sigma_\varepsilon ^2 B(z)} [\frac {S_{sx}(z)}{B(z^{-1})}]_+$$
通过 $H_c(z) \quad \underrightarrow{反变换} \quad h(n)$。

五、计算最小均方差

$$\zeta_{min}(n) = \frac {1}{2 \pi j} \int\limits {[S_{ss}(z)-H_c(z)S_{sx}(z)]z^{-1}} dz$$

2-5 因果IIR WF计算

(1) 根据 a、c、R、Q，得： $Q=P-\frac {PRa^2}{R+c^2P}$

(2)求维纳增益: $G=\frac {CP}{R+c^2P}$

(3) 滤波器系数：$f=a(1-cG)$

(4)求$H_c(z)$的表达式：$H_c(z)=\frac {G}{1-fz^{-1}}$

其中，a就是信号模型中的a，c就是测量模型，Q是$w(n)$的方差，R是$v(n)$的方法

互补维纳滤波器用于飞机盲目着陆系统：

飞机着陆时以较慢的恒定速度沿着一个无线电波速下降。为了自动对准跑道，通常需要为系统提供二个信号。

① 一个是由无线电波束提供的信号（它与飞机航向的滑离跑道方向的大小成比例），有飞机信号与高频噪声叠加

②另一个信号由飞机本身通过自身方位的测量得到，在下降过程中，由于风向的改变在信号中引入低频噪声

Chapter3 自适应滤波器

性能曲面：

$$\xi(n)=E[d^2(n)]+W^TRW-2P^TW$$
性能曲面主轴：就是 $R$的特征向量

1. 先求特征值 $(R-\lambda I)Q_n=0$。
   • 利用$QQ^T=I$正交求出其值。

   • 性能曲面沿主轴的二阶导数：就是特征值的2倍

     例3 P43

     已知: $R=\begin{vmatrix} 3 & 1 \\ 1 & 3\end{vmatrix} \qquad P=\begin{vmatrix}2 \\ 5 \end{vmatrix} \qquad E[d^2(n)]=10$

     （1）求出性能曲面表达式

     （2）求最佳权矢量

     （3）求性能曲面的主轴

     （4）求性能曲面沿主轴的二阶导数

     全矢量的计算公式：

     迭代器形式： $W(n+1)=(I-2uR)W(n)+2uP$

     任意时刻的表达式： $W(n)=W^*+[I-2uR]^n[W(0)-W^*]$ 或者 $\xi(n)=\xi_{min}+\sum\limits^{L}_{k=0}[v^{’}_k]^2 \lambda_k(1-2u\lambda_k)^{2n}$。

     chapter4 功率谱估计

     取样自相关 $\hat{R}_{xx}(m)=\frac {1}{N} \sum\limits_{n=0}^{N-1-|m|}x(n)x(n+m) \quad |m| \leq N-1$。

     Yule-Walker方程
     

     $$\begin{cases} \delta^2_k=R(0)+\sum\limits_{i=1}^ka_{k,i}R(i) \\ D_k = \sum\limits_{i=0}^ka_{k,i}R_{(k+1-i)}, \quad a_{k,0}=1 \\ v_{k+1}=\frac {D_k}{\delta_k^2} \\ \delta^2_{k+1} = (1-v_{k+1}^2)\delta^2_k \\ a_{k+1,i}=a_{k,i}-v_{k+1} \cdot a_{k,k+1-i}, \quad i = 1,2,\cdots, k,k+1 \end{cases}$$

     AR模型的激励

     $H_c(z)=\frac {1}{1+\sum\limits_{k=1}^2a_k \cdot z^{-k}}$

     自相关值得估计与模型参数计算估计

     例：已知$x(n)=(1,1,1,1,1)$，求$x(n)$的AR(2)模型参数

     （1） 求出自相关函数

     求出$\hat R(0)=1，\hat R(1)=0.8，\hat R(2)=0.6$

     方法一：

     (1)、0接AR模型

     $\delta_0^2 = R(0)$ 和 $D_0=R(1)$

     (2)、1阶AR模型
     

     $$v_1 = \frac {D_0}{\delta^2_0} \\ \delta_1^2=(1+v_1^2)\cdot \delta_0^2 \\ a_{1,i} = a_{0,i} - v_1 \cdot a_{0,1-i}, \quad i = 1,2,\cdots, k \\ D_1 = \sum\limits_{i=0}^{1} a_{1,i} \cdot R(1+1-i), \quad a_{1,0}=1$$
     (3)、二阶AR模型