数据项是变分光流算法中能量泛函的重要组成部分之一,主要包含各种常值守恒假设,如灰度守恒假设、梯度守恒假设、Hessian矩阵守恒假设、拉普拉斯守恒假设等,由这些守恒假设构成的约束条件是光流计算过程中决定运动模型的主要因素。而这些守恒假设的适用性和准确性等依赖于许多先验知识,它们和所处理的图像序列的获取条件(如光照变化、噪声大小等)及图像序列的运动形式是密切相关的。下面就这几种常用的守恒假设做一简要介绍。
灰度守恒假设即是前面所提到的Horn变分方法中使用的假设约束,其形式如式(2-1),设X= (x,y) T ,w= (u,ν) T 将式(2-1)写为
用泰勒公式将式(3-33)展开并忽略高阶项,即线性化此假设便得著名的光流基本约束方程,即
此公式不仅应用在Horn方法中,在经典的Lucas局部光流一致方法等有效算法中也被广泛使用。在许多情况下,用灰度守恒假设计算能得到较准确的光流,但是在全局或者局部光照发生变化的情况下,灰度守恒假设将不再有效,这时需要寻找对亮度变化不敏感的守恒假设。
众所周知,光照变化会对图像序列的灰度产生影响,它会使灰度值发生总体偏移,也可能使灰度尺寸进行重新定义,对于前者,它不会对图像序列中的梯度因素产生任何影响,对于后者,它仅仅改变梯度矢量的长度,而不会改变梯度矢量的方向。基于此可定义如下梯度守恒假设
式中,∇= (∂ x ,∂ y ) T ,表示图像空间梯度,可见式(3-35)为矢量形式,用泰勒展开线性化此公式,便得到如下联立的方程组
上式可用矢量的形式表示为
Hessian矩阵守恒假设是关于图像中像素灰度高阶导的一种守恒假设,它能包含图像的更多信息,并且不受光照变化的影响,对此假设进行公式化便得到如下形式
其中H ss 表示求Hessian矩阵,即
那么式(3-39)等价于如下联立的方程组
对式(3-40)用泰勒公式展开可得
同样上式可用矢量的形式表示为
Hessian矩阵守恒假设虽然能够避免光照变化的不利影响,但与灰度守恒假设不同,梯度守恒假设和Hessian矩阵守恒假设这些基于一阶或高阶导的假设约束,并不适合各种运动形式。因为梯度守恒假设和Hessian矩阵守恒假设都包含方向信息,所以它们比较适应于计算平移运动、扩散运动及速度较缓慢的旋转运动等运动形式,并且它们比灰度守恒假设包含更多的约束条件,这就更增强了运动估计的鲁棒性。但另一方面,当图像序列中含有快速旋转运动的成分时,这两种守恒假设将不再适合。
下面介绍一种能够避免旋转运动不利影响的守恒假设,即Laplacian守恒假设,也就是取Hessian矩阵中对角线元素的和,认为在运动中它保持不变,可以看出Laplacian守恒假设对于光照变化也不受影响。其公式化的形式为
其中
示空间Laplacian算子,即
,这里直接给出其泰勒公式展开后的形式
前面介绍了几种比较常见的用于光流计算的守恒假设约束,下面给出几种其他的旋转不变的守恒假设形式,为简便起见,这里不再讨论其线性化形式。我们知道,图像序列任意像素点的空间梯度是一个矢量值,它可以分解为模的信息和方向信息,不考虑方向信息,认为梯度的模是旋转不变的,这个量可以得到另一个守恒假设,即
同理,对空间Hessian矩阵求模,可以得到一个守恒假设
另外,空间Hessian矩阵的值,认为它在运动过程中不变,又可以得到一个守恒假设
当然,还可以制定出 3 阶或者更高阶数的守恒假设,但是需要注意到,更高的阶数会对噪声更加敏感。另外,在图像的某些区域,其像素灰度的导数为零,增加阶数也不会获得更多的有用信息,所以本书中对更高阶数的守恒假设不予讨论。
平滑项是变分光流算法中能量泛函的又一重要组成部分,主要包含各种平滑策略,它使变分光流算法能够求得唯一解。各种平滑策略与应用于图像的各种扩散思想密切相关。许多学者借鉴各种图像扩散的算法策略,形成了各种应用于变分光流算法中平滑项的平滑策略,使得处理光流计算问题,特别是处理边界光流计算问题更加灵活和多变。本节主要介绍平滑项设计的物理学背景,为下一节不同平滑策略的介绍做铺垫。
晶体中存在的杂质,因为其浓度分布不均匀,杂质将从浓度高的区域流向浓度较低区域,这种迁移过程在物理学上称为扩散。与扩散相类似,当介质中的温度分布不均匀时,那么热量将从温度较高的区域流向温度较低区域,这种迁移过程称之为热传导。
用函数u(x,y,z,t)表示杂质浓度随空间和时间变化的程度,用梯度∇u来描述空间的不均匀性,那么由梯度产生的作用力−∇u将推动杂质做定向迁移,这里的负号表示作用力的方向指向u值减小的方向。那么在这一作用力的推动下,介质中将产生流密度,如果流密度的方向与梯度的方向一致,即流密度为单位时间通过与梯度矢量垂直的单位面积的杂质。那么称这种介质为各项同性(Isotropic)介质。用公式表示为
式(3-48)中的系数a称为传导系数。传导系数a为一标量,它并不一定只是一个固定的常值,这只是最简单的情况,它还可以是一个标量函数,函数值依赖于空间的位置,即a(x,y,z) ,或者函数值依赖于空间位置和u本身,即a(x,y,z,u) 。
如图 3-1 所示,讨论包围某一定点p的闭合曲面Ω的总流量F。
由高斯定理得
图 3-1 通过闭合曲面的流量示意图
式(3-50)中,V表示由闭合曲面Ω所包围的体积,由上节可知f表示流密度,可见F的物理意义是在单位时间内,通过曲面Ω流出的杂质的总量。于是有下式成立
式(3-51)的等价表达式为
由前面的讨论知道,各项同性介质中,f= −a∇u,代入式(3-52)可得
式(3-53)在任意正交坐标系(x,y)中,可展开表示为
不难看出,式(3-54)右边可以看成在两个相互垂直的方向上的一维扩散之和,并且两个相互垂直方向上的传导系数相同,故称为“各项同性扩散”。而“各项同性扩散”根据传导系数a的不同又可分为如下几种情况。
(1)传导系数a为一个常数,式(3-54)称为线性扩散。
(2)传导系数a为一个依赖于空间位置的标量函数,即a(x,y,z) ,式(3-54)称为拟线性扩散。
(3)传导系数a为一个既依赖于空间位置又依赖u本身,即a(x,y,z,u) ,式(3-54)称为各项同性的非线性扩散。
本节前一部分介绍了各项同性扩散的各种情况,但是自然界中的晶体中的杂质扩散,其流密度并不一定与梯度方向一致,这样一类介质称为各项异性(Anisotropic)介质。各项异性介质中,流密度可这样表达
其中,D为 2× 2(二维)或3× 3(三维)矩阵,它描述了流密度与梯度方向的不一致性,称为扩散张量。将式(3-55)代入式(3-52)得
本书中,扩散张量D取图像处理中常用的 2× 2对称矩阵进行讨论。
将式(3-56)在任意正交坐标系(x,y)中进行展开得到如下形式
式(3-58)并不能看做两个相互垂直的方向上的一维扩散之和,不过可以利用矩阵本征分解定理进行变换,式(3-58)表示为在两个本征矢量构成的局部坐标系上的形式。由本征分解定理得
在两个本征矢特征向量构成的局部坐标系(ν 1 ,ν 2 )中,展开∇u得
再由本征矢量的正交归一化性质得
可见,张量扩散也可以看做是沿着两个相互垂直的本征矢量方向上的一维扩散之和。但是由式(3-61)可以看出,这两个方向上的传导系数不同,这与式(3-54)表示的各项同性扩散有本质的区别,故称为“各项异性扩散”。
上面介绍的物理学中各种扩散思想,在图像处理中得到了广泛的应用和发展。不难想象,把整幅 2 维图像当做一个 2 维晶体,用图像中各个像素点的灰度值I(x,y,t)代表这个位置的杂质浓度u(x,y,t),图像的梯度∇I即表示杂质的扩散作用力∇u,那么就可以把各种扩散思想应用到图像处理上来。
以线性扩散为例,假设传导系数为常值 1,图像线性扩散方程为
可以想象,原始图像经时间t不断演化,灰度值从高的区域逐渐流向灰度值低的区域,实际上,式(3-62)也可通过Fourier变换的方法得到它的解
其中,
表示中心位于坐标原点,宽度为
的二维Gaussian函数。由此可见,图像的线性扩散,等价于传统图像处理中对图像采用高斯滤波器进行滤波。特别的,当时间趋于无穷大时,线性扩散方程的稳态解
趋于原始图像的平均灰度值,可见,随着时间的增长,图像将变得越来越模糊,最后以图像灰度变为平均值而告终。
在Horn变分光流算法中的能量泛函,经最速下降法得到一个扩散反应方程,如式(3-26)所示。而扩散反应方程中的扩散项与图像线性扩散形式一致,只不过是图像灰度I的扩散变成了光流w的扩散,如式(3-32)所示。其中传导系数g为常数 1,导致光流在图像中任何区域的扩散程度都是一样的,这使得光流在运动边缘处很模糊。
值得注意的是:扩散反应方程中的演化参数τ与图像扩散中的扩散时间t不能混淆。计算光流稳态值,是扩散反应方程[式(3-26)]到达稳态,能量泛函取极值时的光流值,即∂ τ w=Δw−(1/ α) ∇I(∇I T w+I t )=0。并不是像图像扩散中时间趋于无穷大时的稳态,实际应用中,图像扩散也不是一直扩散下去,而是规定一个扩散时间范围。可见,光流扩散实际上是继承了图像扩散中的扩散趋势,一直进行到能量泛函取得极值时的光流状态。
图像线性扩散,模糊了包含图像重要信息的边缘,这在图像处理中往往是不能接受的。如果在图像扩散的过程中,传导系数能依赖图像的局部特性,即在图像平坦区域,传导系数自动增大,而在图像边缘区域,传导系数自动减小,这样就能达到既能消除噪声又能保护边缘的效果。Perona和Malik首先提出了如下的扩散PDE方程
式(3-64)为P-M方程。其中传导系数为g,一般情况下它是图像边缘检测(Edge Detection)中用到的边缘函数,比如式(3-65)为一常用的边缘函数
式中,r是自变量,k是选定的常数,它可以控制g的下降速率。原则上,g可以是任何单调递减的非负函数,它依赖于图像梯度∇I的模的大小,可见P-M方程为各项同性的非线性扩散。它将传统的图像滤波与边缘检测统一起来,达到既能滤除噪声又能保护边缘的目的。有关P-M方程的具体行为分析,不再赘述。
受到各向同性非线性扩散图像滤波的启发,Weickert等人给出了基于各向同性非线性扩散、图像驱动平滑项的具体形式
式中,g为单调递减非负函数。这里的图像驱动表明传导系数g只与图像数据有关,与光流数据无关,因此认为运动边缘是图像边缘的子集,这样就可以保护运动边缘的光流不被过于平滑和模糊。可以证明平滑项[式(3-66)]对应的反应扩散方程中的扩散项为
由式(3-67)和式(3-64)对比可以看出,基于各向同性、图像驱动的平滑项对应的光流扩散过程与图像中各向同性非线性扩散相一致。使用平滑项[式(3-67)]对数据项进行约束,在计算光流时能够在图像边缘,分别减弱光流u方向和ν方向的扩散程度,从而保护了光流边缘。
顾名思义,各向异性、图像驱动的平滑项与各向同性、图像驱动的平滑项的区别为:在图像边缘,不是各个方向都减小扩散,而是在不同的方向施加不同的作用力。一个基本的方法是在垂直图像边缘的地方减小扩散,在平行图像边缘的地方增加扩散。最早把Horn变分光流算法中的全局平滑约束替换成各向异性、图像驱动平滑项的算法是由Nagel提出的,基于此方法可得到这样一个光滑项
式中
,其中E是单位矩阵。可以证明平滑项[式(3-68)]对应的反应扩散方程中的扩散项为
下面说明扩散张量矩阵D(∇ I)的设计过程,首先假设扩散张量矩阵D(∇ I)的特征向量已知,分别为
其方向分别为平行于图像梯度方向和垂直于梯度方向。那么根据矩阵本征分解定理
,只需要设计扩散张量矩阵的特征值λ
1
、λ
2
就可以了。扩散张量矩阵D(∇ I)的两个特征值设计为
可以看出,在图像中灰度梯度|∇ I|→ 0 的平坦区域,λ 1 → 0 .5 ,λ 2 → 0.5,上述扩散变为各项同性的扩散。在灰度边缘区域|∇ I|→∞ ,λ 1 → 0,λ 2 →1 ,这样就增大了ν 2 方向,即垂直梯度方向上的扩散,并且减小了ν 1 即平行梯度方向上的扩散。
在实际中,图像中灰度边缘往往大于运动边缘,如果图像中存在纹理比较强的区域,那么使用图像驱动往往造成过分割现象,为了仅仅在运动边缘减小过度平滑,引入光流驱动的光滑项。Schnörr和Weickert最早提出了如下各项同性、光流驱动的平滑项
其中,ψ(s 2 )是一个可微的递增函数,为了保证其收敛性,要求函数ψ(s 2 )是凸性的,例如
可以证明平滑项[式(3-71)]对应的反应扩散方程中的扩散项为
从式(3-73)可以看出,传导系数
是个标量函数,并且依赖于光流模值
的大小,这表明这个模型是各项同性的并且是光流驱动的。由式(3-72)定义的ψ函数可得
因此,这个非线性扩散关于自变量是递减的,所以在光流的模值
+
大的区域(即运动边缘),扩散效果被削弱,从而使得运动边缘得到保护。
要设计各项异性、光流驱动的平滑项,必然要使在光流扩散过程中,扩散张量与光流值相关,以消除图像驱动带来的过分割现象,并且还要在运动边缘处对各个方向施加不同的影响,一个基本的目标是在垂直于光流边缘(运动边缘)的地方减小扩散,在平行光流运动边缘的地方增加扩散。
首先设标量函数ψ(s
2
)为式(3-72)所示,A为n×n阶矩阵,其正交归一化的特征向量为
,对应的特征值为
,把标量函数枞的扩展为关于矩阵A的函数为
同理,对于ψ(s 2 )的导函数扩展为
定义矩阵A对角线元素之和为
通过以上定义,Weickert给出了如下基于各向异性、光流驱动的平滑项
式(3-78)中,令
式中,J是一个对称半正定 2×2 阶矩阵,称为结构张量。含有两个正交特征矢量ν 1 、ν 2 ,其对应的特征值为μ 1 、μ 2 ,这两个特征值能反映光流在这两个特征矢量ν 1 、ν 2 方向上的变化。Zenzo曾经在有关多通道图像的边缘分析中引入此概念并进行了详细介绍,在图像中,这个张量能反映图像局部结构信息的变化,这里则反映了光流的局部变化信息。使用如式(3-78)所示的平滑项,经最速下降法能得到所期待的扩散项,即
注意到,式(3-80)中的扩散矩阵,其特征值ψ'(u 1 ) 、ψ'(u 2 )是不同的,从而实现光流的各项异性扩散。