Signals, Systems and Inference, Chapter 11: Wiener Filtering (mit.edu)

基本

在图像处理的时候, 遇到了这个维纳滤波, 其推导的公式不是很理解, 于是上网查了查, 并做个简单的总结.

符号 说明
\(x[k]]\) 观测信号\(x\)的第k个元素
\(\hat{y}\) 为\(y\)的一个估计
\(v\) 噪声信号
\(e[k]\) 误差, 为\(e[k]=\hat{y}[k] - y[k]\)
\(R_{xy}(i, j)\) 相关系数: \(E\{x[n-i]y[n-j]\}\)
\(S_{xy}\) 傅里叶变换: \(\mathcal{F} \{R_{xy}\}\)

基本的假设:

\(y[k]\)服从wide-sense stationary (WSS), 即

  1. \(E[y[k]] = E[y[0]]\);
  2. \(R_{yy}(i, j) = R_{yy}[j-i]\).

维纳滤波可以应用于很多场景, 但是这里只讨论下面的去噪的情形:

\[x[n] = y[n] + v[n],
\]

且假设\(y, v\)之间相互独立, \(E[v]=0\).

我们的目标是找到一个滤波\(h\), 得到一个估计

\[\hat{y}[n] = h\star x [n],
\]

使得下式最小

\[E[e^2[n]].
\]

维纳滤波需要分情况讨论, 这里只关注离散的情形, 包括

  1. non-causal: \(\hat{y}[n] = \sum_{k=-\infty}^{+\infty} h[k]y[n-k]\);
  2. FIR: \(\hat{y}[n] = \sum_{k=0}^{N-1} h[k]y[n-k]\).

causal的情况这里就不写了.

滤波的推导

自然地, 寻找驻点:

\[\begin{array}{ll}
\frac{\partial E[e^2[n]]}{\partial h[m]}
&=2E[e[n]\frac{\partial e[n]}{\partial h[m]}] \\
&=2E[e[n]x[n-m]] \\
&=2(h\star R_{xx}[m] -R_{yx}[m]) \\
&=0.
\end{array}
\]

于是必须满足

\[\tag{1}
h \star R_{xx}[m] = R_{yx}[m].
\]

相应的在频率域内存在(假设DFT存在, 参考文献用的z变换, 这个不是特别了解):

\[\tag{2}
H[u] S_{xx}[u] = S_{yx}[u].
\]

在FIR情况下, 可以通过(1)推导出一个线性方程组从而求解, non-causal下可用(2)推导出结果.

注: \(H[u]\)用了[]是为了保持一致, 在non-causal情况下\(H(z)\)可能更加妥当.

故最优解为:

\[\tag{*}
H = S_{yx} / S_{xx}.
\]

特别的情况

进一步, 由于

\[S_{xx}
= \mathcal{F}\{R_{xx}\}
= \mathcal{F}\{E[(y[n]+v[n])(y[n-m]+v[n-m])]\}
=S_{yy} + S_{vv},
\]

最后一步成立的原因是\(E{v}=0\), 且\(v, y\)相互独立.

同时

\[S_{yy} = S_{yx}.
\]

所以:

\[\tag{3}
H = \frac{S_{yy}}{S_{yy} + S_{vv}}.
\]

特别的例子

在图像数字处理中, 给出的这样的情形(FIR):

\[x[n] = g \star y [n] + v[n],
\]

记\(r[n] = g \star y [n]\),

\[\begin{aligned}
S_{rr}[u]
&= \mathcal{F}\{R_{rr}[m]\}[u] \\
&= \mathcal{F}\{E[r[n] r[n-m]]\}[u] \\
&= E[r[n] \mathcal{F}\{r[n-m]\}][u] \\
&= E[r[n] G[-u] Y[-u] e^{-j2\pi nu/N}] \\
&= G[-u] E[r[n]Y[-u] e^{-j2\pi nu/N}] \\
&= G[-u] \mathcal{F}\{R_{ry}[m]\}[u] \\
&= G[-u] \mathcal{F}\{E[r[n+m]y[n]]\}[u] \\
&= G[-u] E[\mathcal{F}\{r[n+m]\}[u]y[n]] \\
&= G[-u] E[G[u]Y[u]e^{-j2\pi nu/N}y[n]] \\
&= G[-u]G[u] \mathcal{F}\{R_{yy}[m]\}[u] \\
&= G[-u]G[u] S_{yy}[u]
\end{aligned}
\]

由于

\[S_{yx} = S_{yr} = G[-u]S_{yy}[u],
\]

证明和上面是类似的,

\[S_{xx} = S_{rr} + S_{vv},
\]

于是

\[H[u] = \frac{1}{G[u]}\frac{1}{1 + S_{vv}[u] / (G[-u]G[u]S_{yy}[u])}.
\]

当\(g\)为实的时候, 有\(G[-u]G[u]=|G[u]|^2\), 在数字图像处理书中, 给出的公式中:

\[S_{yy}[u] = |F[u]|^2, S_{vv}[u] = |V[u]|^2,
\]

个人觉得是这里的期望\(E\)用

\[R_{yy}[m] = \sum_{n=0}^{N-1} y[n]y[n-m],
\]

代替了,

所以

\[\mathcal{F}\{R_{yy}[m]\} = F^*(u)F(u) = |F(u)|^2,
\]

这里假设\(y[n] \in \mathbb{R}\).

\(S_{vv}[u]\)也是类似的.

最新文章

  1. mongoose数据库连接和操作
  2. Android—实现自定义相机倒计时拍照
  3. sys.sysprocesses视图的使用小结
  4. C语言实现单链表-03版
  5. 【原】macbook不睡眠的排查与解决
  6. IE6 BUG 汇总
  7. BZOJ 1071组队
  8. Linux 修改hostname 文件
  9. SPOJ CNTPRIME 13015 Counting Primes (水题,区间更新,求区间的素数个数)
  10. Error 1937.An error occurred during the installation of assembly...
  11. WildFly 9.0.2+mod_cluster-1.3.1 集群配置
  12. BZOJ 1668: [Usaco2006 Oct]Cow Pie Treasures 馅饼里的财富
  13. 基于FPGA的cordic算法的verilog初步实现
  14. Rails 执行 rails server 报错 Could not find a JavaScript runtime
  15. (办公)json报错的解决问题的小经验.
  16. pycharm安装配置
  17. Confluence 6 配置手动备份
  18. maven项目update报错
  19. Iowait的成因、对系统影响及对策--systemtap
  20. 转:如何捕获winform程序全局异常?

热门文章

  1. A Child's History of England.39
  2. day14搭建博客系统项目
  3. java9 模块化 jigsaw
  4. Spark(八)【利用广播小表实现join避免Shuffle】
  5. [php安全]原生类的利用
  6. 【编程思想】【设计模式】【创建模式creational】Borg/Monostate
  7. 【Java多线程】Java 原子操作类API(以AtomicInteger为例)
  8. 【Java 调优】Java性能优化
  9. 30个类手写Spring核心原理之MVC映射功能(4)
  10. 4、Linux下安装tomcat