Wiener Filter

因为最近看文章接触了维纳滤波,所以这里写一下Weiner Filter的一些简单理解和推导。

基本定义

维纳滤波是一种在含噪声的时序信号把信号提取出来的滤波器,其基本框图如下:

简单的维纳滤波其实就是通过一个FIR滤波器,去除噪声的过程。在这里,hhh的作用也可以理解为: 通过训练集的数据对信号和噪声的建模,然后通过前几个点的信息,预测当前时刻的噪声信号所占的比例,然后去除掉,剩下的就是预测的时序信号了。维纳滤波作为一种使用很广泛的滤波器,其变化的形式也有很多种,可以是单输入输出的,也可以是多输入输出的。hhh所表示的变换也可以写成非线性;hhh可以是有限长的FIR滤波器,也可以是无限长的IIR滤波器。要取决于当前你所解决的问题。但是维纳滤波的基本思想还是一致的。通过滤波(矩阵或者其他模型的形式)来从信号和噪声的混合中提取信号。所以维纳滤波的核心,就是计算这个滤波器(矩阵hhh或者模型的参数)。也就是解Wiener-Hopf方程。

本文用比较简单的单输入输出,且只考虑有限长滤波(即认为当前时刻的信号只和前有限个时间点的信号相关)。

公式推导

首先,对于图1中的滤波器:

y(n)=x(n)∗h(n)=(s(n)+v(n))∗h(n)(1)y(n) = x(n) * h(n) = (s(n)+v(n))*h(n) \tag{1}y(n)=x(n)∗h(n)=(s(n)+v(n))∗h(n)(1)

其中∗*∗表示卷积,x(n)x(n)x(n)表示输入信号, y(n)y(n)y(n)表示输出信号, s(n)s(n)s(n)表示输入信号中,有用的信号部分;v(n)v(n)v(n)表示输入信号中的噪声部分。

维纳滤波的目标是,保证输出y(n)y(n)y(n)和真实信号s(n)s(n)s(n)的差别最小,由于y(n)y(n)y(n)和s(n)s(n)s(n)是时序信号,所以要保证两者的均方误差最小,所以有:

E{e2(n)}=E{(y(n)−s(n))2}=E{(x(n)∗h(n)−s(n))2}(2)E\{e^2(n)\} = E\{(y(n)-s(n))^2\} = E\{(x(n)*h(n)-s(n))^2\} \tag{2} E{e2(n)}=E{(y(n)−s(n))2}=E{(x(n)∗h(n)−s(n))2}(2)

即求使得Eq(2)最小的hhh。所以E{e2}E\{e^2\}E{e2}对hhh求偏导。有:

∂E{e2(n)}∂h=2E{e(n)∗∂e(n)∂h}=0(3)\frac{\partial{E\{e^2(n)\}}}{\partial{h}} = 2E\{e(n) * \frac{\partial{e(n)}}{\partial{h}}\} = 0 \tag{3} ∂h∂E{e2(n)}​=2E{e(n)∗∂h∂e(n)​}=0(3)

∂E{e2(n)}∂h=2E{[∑m=0N−1h(m)x(n−m)−s(n)]x(n−j)},j=0,1,...,N−1(4)\frac{\partial{E\{e^2(n)\}}}{\partial{h}} = 2E\{[\sum_{m=0}^{N-1}{h(m)x(n-m) - s(n)}]x(n-j)\}, j = 0, 1, ... , N-1 \tag{4} ∂h∂E{e2(n)}​=2E{[m=0∑N−1​h(m)x(n−m)−s(n)]x(n−j)},j=0,1,...,N−1(4)

∂E{e2(n)}∂h=2∑m=1N−1h(m)E{x(n−j)x(n−m)}−2E{s(n)x(n−j)}=0,j=0,1,...,N−1(5)\frac{\partial{E\{e^2(n)\}}}{\partial{h}} = 2\sum_{m=1}^{N-1}{h(m)}E\{x(n-j)x(n-m)\} - 2E\{s(n)x(n-j)\} = 0, j = 0, 1, ..., N-1 \tag{5} ∂h∂E{e2(n)}​=2m=1∑N−1​h(m)E{x(n−j)x(n−m)}−2E{s(n)x(n−j)}=0,j=0,1,...,N−1(5)

我们设xxx和sss的相关系数为RxsR_{xs}Rxs​,则有:

Rxs(j)=∑m=0N−1h(m)Rxx(j−m),j=0,1,...,N−1(6)R_{xs}(j)=\sum_{m=0}^{N-1}{h(m)R_{xx}(j-m)}, j=0,1,...,N-1 \tag{6}Rxs​(j)=m=0∑N−1​h(m)Rxx​(j−m),j=0,1,...,N−1(6)

其中,Rxx(j−m)R_{xx}(j-m)Rxx​(j−m)表示x(n−j)x(n-j)x(n−j)和x(n−m)x(n-m)x(n−m)的相关系数,这里mmm是固定的,jjj是变化的。且m>=0m>=0m>=0,Rxs(j)R_{xs}(j)Rxs​(j)表示x(n−j)x(n-j)x(n−j)和s(n)s(n)s(n)的相关系数。上述公式中,nnn表示的是时序信号中的时间点。

然后,根据Eq(6),可以得到NNN个线性方程:
{Rxs(0)=h(0)Rxx(0)+h(1)Rxx(1)+...+h(N−1)Rxx(N−1)Rxs(1)=h(1)Rxx(1)+h(0)Rxx(0)+...+h(N−1)Rxx(N−2)...Rxs(N−1)=h(N−1)Rxx(N−1)+h(N−2)Rxx(N−2)+...+h(0)Rxx(0)(7)\begin{cases} R_{xs}(0)=h(0)R_{xx}(0)+h(1)R_{xx}(1)+...+h(N-1)R_{xx}(N-1)\\ R_{xs}(1)=h(1)R_{xx}(1)+h(0)R_{xx}(0)+...+h(N-1)R_{xx}(N-2)\\ ...\\ R_{xs}(N-1)=h(N-1)R_{xx}(N-1)+h(N-2)R_{xx}(N-2)+...+h(0)R_{xx}(0)\\ \end{cases} \tag{7} ⎩⎪⎪⎪⎨⎪⎪⎪⎧​Rxs​(0)=h(0)Rxx​(0)+h(1)Rxx​(1)+...+h(N−1)Rxx​(N−1)Rxs​(1)=h(1)Rxx​(1)+h(0)Rxx​(0)+...+h(N−1)Rxx​(N−2)...Rxs​(N−1)=h(N−1)Rxx​(N−1)+h(N−2)Rxx​(N−2)+...+h(0)Rxx​(0)​(7)
写成矩阵形式,有:

RxxH=Rxs(8)\displaystyle \boldsymbol{R_{xx}H}=\boldsymbol{R_{xs}} \tag{8}Rxx​H=Rxs​(8)

其中, H=[h(0),h(1),...,h(N−1)]T\displaystyle \boldsymbol{H} = [h(0), h(1),...,h(N-1)]^TH=[h(0),h(1),...,h(N−1)]T是需要求的滤波器参数

Rxs=[Rxs(0),Rxs(1),...,Rxs(N−1)]T\displaystyle \boldsymbol{R_{xs}} = [R_{xs}(0),R_{xs}(1), ..., R_{xs}(N-1)]^TRxs​=[Rxs​(0),Rxs​(1),...,Rxs​(N−1)]T是xxx和sss的相关系数
Rxx=[Rxx(0)Rxx(1)...Rxx(N−1)Rxx(1)Rxx(0)...Rxx(N−2)⋮⋮...⋮Rxx(N−1)Rxx(N−2)...Rxx(0)](9)\displaystyle \boldsymbol{R_{xx}} = \begin{bmatrix} R_{xx}(0)&R_{xx}(1)&...&R_{xx}(N-1)\\ R_{xx}(1)&R_{xx}(0)&...&R_{xx}(N-2)\\ {\vdots}&{\vdots}&...&{\vdots}&\\ R_{xx}(N-1)&R_{xx}(N-2)&...&R_{xx}(0)\\ \end{bmatrix} \tag{9} Rxx​=⎣⎢⎢⎢⎡​Rxx​(0)Rxx​(1)⋮Rxx​(N−1)​Rxx​(1)Rxx​(0)⋮Rxx​(N−2)​............​Rxx​(N−1)Rxx​(N−2)⋮Rxx​(0)​​⎦⎥⎥⎥⎤​(9)

所以根据Eq(8)可以求得:

H=Rxx−1Rxs(10)\displaystyle \boldsymbol{H} = \boldsymbol{R_{xx}^{-1}R_{xs}} \tag{10}H=Rxx−1​Rxs​(10)

此时,信号的均方误差最小,根据Eq(2),可得:

E{e2(n)}=E{(s(n)−∑m=0N−1h(m)x(n−m))2}(11)E\{e^2(n)\} = E\{(s(n)-\sum_{m=0}^{N-1}h(m)x(n-m))^2\} \tag{11}E{e2(n)}=E{(s(n)−m=0∑N−1​h(m)x(n−m))2}(11)

E{e2(n)}=E{s2(n)−2s(n)∑m=0N−1h(m)x(n−m)+∑m=0N−1∑r=0N−1h(m)x(n−m)h(r)x(n−r)}E\{e^2(n)\} = E\{s^2(n) - 2s(n)\sum_{m=0}^{N-1}h(m)x(n-m)+\sum_{m=0}^{N-1}\sum_{r=0}^{N-1}{h(m)x(n-m)h(r)x(n-r)}\}E{e2(n)}=E{s2(n)−2s(n)m=0∑N−1​h(m)x(n−m)+m=0∑N−1​r=0∑N−1​h(m)x(n−m)h(r)x(n−r)}

E{e2(n)}=Rss(0)−2∑m=0N−1h(m)Rxs(m)+∑m=0N−1h(m)∑r=0N−1h(r)Rxx(m−r)E\{e^2(n)\}=R_{ss}(0)-2\sum_{m=0}^{N-1}{h(m)R_{xs}(m)+\sum_{m=0}^{N-1}{h(m)}\sum_{r=0}^{N-1}{h(r)R_{xx}(m-r)}}E{e2(n)}=Rss​(0)−2m=0∑N−1​h(m)Rxs​(m)+m=0∑N−1​h(m)r=0∑N−1​h(r)Rxx​(m−r)

根据Eq(5),可得:

E{e2(n)}=Rss(0)−∑m=0N−1h(m)Rxs(m)(12)E\{e^2(n)\} = R_{ss}(0) - \sum_{m=0}^{N-1}{h(m)R_{xs}(m)} \tag{12}E{e2(n)}=Rss​(0)−m=0∑N−1​h(m)Rxs​(m)(12)

假设信号s(n)s(n)s(n)和噪声v(n)v(n)v(n)互相独立,那么有:

Rsv=Rvs=0R_{sv}= R_{vs} = 0Rsv​=Rvs​=0

Rxs=Rss+Rvs=RssR_{xs} = R_{ss} + R_{vs} = R_{ss}Rxs​=Rss​+Rvs​=Rss​

Rxx=Rss+Rsv+Rvs+Rvv=Rss+RvvR_{xx} = R_{ss}+R_{sv}+R_{vs}+R_{vv} = R_{ss}+R_{vv}Rxx​=Rss​+Rsv​+Rvs​+Rvv​=Rss​+Rvv​

则,根据Eq(12),有:

E{e2(n)}=Rss(0)−∑m=0N−1h(m)Rss(m)(14)E\{e^2(n)\} = R_{ss}(0) - \sum_{m=0}^{N-1}{h(m)R_{ss}(m)} \tag{14}E{e2(n)}=Rss​(0)−m=0∑N−1​h(m)Rss​(m)(14)

至此,最简单的维纳滤波的基本公式推导完成,如果涉及到多输入多输出的维纳滤波,会更加复杂,这里不做推导。

维纳滤波(Wiener Filter)相关推荐

  1. 数字图像处理Python语言实现-图像滤波-维纳滤波(Wiener Filter)

    维纳滤波(Wiener Filter) 1.前言 维纳滤波器(Wiener Filter)是最早用于图像复原经典滤波之一,目前被广泛用于信号滤波降噪和图像预处理中.维纳滤波器的目的是使用相关信号作为输 ...

  2. 维纳滤波原理(Wiener Filter)

    维纳滤波原理(Wiener Filter) - 知乎 自适应滤波:维纳滤波器--FIR及IIR设计 - 桂. - 博客园

  3. 数字图像处理实验(13):PROJECT 05-04,Parametric Wiener Filter

    实验要求: Objective: To understand the high performance of the parametric Wiener Filter in image restora ...

  4. Wiener Filter

    假设分别有两个WSS process:$x[n]$,$y[n]$,这两个process之间存在某种关系,并且我们也了解这种关系.现在我们手头上有process $x[n]$,目的是要设计一个LTI系统 ...

  5. 【图像处理】参数维纳滤波(Parametric Wiener Filter)

    实验要求   (a) 编写一个给图像中添加高斯噪声的程序,程序的输入参数为噪声的均值与方差.   (b) 编写程序实现公式(5.6-11)所示的污损滤波:   (c) 如图5.26(b)所示,对图像5 ...

  6. 参数维纳滤波(Parametric Wiener Filter)

     摘   要:本实验主要使用维纳滤波法(又名为最小均方误差滤波)实现图像复原与重建.首先我们通过对一幅图像加入运动污损滤波和高斯噪声,然后从噪声中提取出原始图像信号.在各种估计方法中,维纳滤波是一 ...

  7. Wiener Filter维纳滤波器halcon算子,持续更新

    目录 gen_psf_defocus gen_psf_motion simulate_defocus simulate_motion wiener_filter wiener_filter_ni ge ...

  8. 维纳滤波器(Wiener Filter)在图像处理中的应用(一)

    维纳滤波器是一种自适应的滤波器,在数字信号处理中有着广泛的应用.ispforfun会在从今天开始定期给大家带来维纳滤波器在图像处理中应用.本节讲诉维纳滤波器在图像去噪中的简单应用. 让我们从Matla ...

  9. 维纳滤波器推导以及MATLAB代码(Wiener Filter)

    维纳滤波器 1. 背景及术语介绍 随机信号或者随机过程:它是是普遍存在的.一方面,任何确定性信号经过测量后往往会引入随机性误差:另一方面,任何信号本身都存在随机干扰. 噪声:按照功率谱密度(power ...

  10. 【youcans 的 OpenCV 例程 200 篇】107. 退化图像的维纳滤波

    欢迎关注 『youcans 的 OpenCV 例程 200 篇』 系列,持续更新中 欢迎关注 『youcans 的 OpenCV学习课』 系列,持续更新中 [youcans 的 OpenCV 例程 2 ...

最新文章

  1. 实现Singleton 模式——七种实现方式
  2. 设计模式之---Factory
  3. Manjaro下带供电的USB Hub提示error -71
  4. 十八、彻底掌握金融量化交易库Tushare
  5. linux系统上传代码到gitlab服务器
  6. part-time job
  7. 更改路由器的外网IP
  8. CodeForces 877E DFS序+线段树
  9. Visio帮你轻松画出3D效果示意图
  10. 【代码笔记】多线程游戏开发——伏魔记:第一步——开始游戏界面实现(一)...
  11. 免安装版VSCode配置(便携模式)
  12. 远程控制工具——Centos7上向日葵安装使用(xy)
  13. Volatility
  14. HTML|颜色的设置方法
  15. 开源的C++静态分析工具
  16. java.lang.IllegalStateException: Failed to convert message:‘‘ to outbound message.
  17. Javaweb支付宝支付
  18. AR涂涂乐⭐四、 获取截图、赋值给物体,将数据传递给shader
  19. 新办的卡为什么显示无服务器,为什么插入卡后显示无服务,有时有有时又没有?...
  20. 7.Docker容器使用辅助工具汇总

热门文章

  1. 苹果计算机重装系统步骤,苹果台式电脑重装系统教程,适合imac恢复出厂设置...
  2. MAC - 必备软件安装与使用
  3. mysql数据库的基本操作
  4. FX3SA三菱PLC使用软件GX Works2编写程序(梯形图等)
  5. mapminmax 用法
  6. python爬虫豆瓣网的模拟登录实现
  7. 对话改写论文笔记(2021年初 )
  8. 浙大计算机专业硕士专业代码,浙江大学海洋学院电子信息(专硕)专业代码
  9. Matlab基础(5)——符号运算
  10. 在pycharm中查看opencv版本