基于图像结构张量的粗到细光流估计
Abstract
提出了一种基于图像局部结构张量的粗到细光流方法。将所提出的结构张量恒定假设和灰色值恒定假设与鲁棒惩罚函数相结合,投影了光流模型的数据项,并规划了一个各向同性非线性函数来平滑项。采用从粗到细的多尺度扭曲策略,建立了一种线性迭代方案。实验结果表明,该方法具有较高的精度和较好的鲁棒性。
1. Introduction
众所周知,变分光流策略通常需要构成一个包含两部分的能量函数:数据项和平滑项,其中也包括一个加权因子。数据项通常包括图像数据的各种恒定假设,这决定了不同图像区域的光流估计。帕彭伯格[4]对所提出的恒定假设进行了详细的比较。为了使Brox[5]在光照变化条件下具有更好的性能,我们结合了经典的亮度恒定假设和梯度恒定假设,由Bruhn[6]、Zimmer[7]和Zhang[8]提出了该思想的一些改进策略。
平滑数据由不同的多种扩散策略组成,控制了不同图像区域的光流扩散程度。平滑项的第一个概念是由霍恩和申克提出的,他们在假设流动在任何区域都是连续的情况下使用均匀扩散。由于均匀扩散往往会影响运动边缘,因此提出了一些各向同性[9-10]和各向异性[11-12]扩散策略,可以降低运动边缘的光流扩散,增加运动区域的光流扩散,以保持运动边界。
为了应对大位移光流估计,Anandan[13]和Black[14]采用了从粗到细的方法来处理通常对图像噪声敏感的线性数据项的缺陷。Brox[15]结合了粗到细和区域匹配的方法来估计光流,其他相关研究由Amiaz[16]和Nir[17]提出。
利用图像局部结构张量的几何流动法,可以获得比经典方法更好的性能。本文的其余部分组织如下。在第Ⅱ节中,我们给出了光流估计的计算模型。详细的实验项目在Ⅲ部分提供,简要总结。
2. Method
A. The data term
对于灰度序列,设I(X)为第一帧点X的灰度值,而下一帧对应点的灰度值为I(X+W)。假设相应点具有相同的灰度值,灰色值恒定假设如下:
I
x
u
+
I
y
v
+
I
t
=
0
(2)
I_xu+I_yv+I_t=0 \tag{2}
Ixu+Iyv+It=0(2)由于灰色值恒定假设对光照变化非常敏感,Tretiak[18]可能首先提出了梯度常量假设,并证明了它不受光照变化的影响:
∇
I
(
X
+
w
)
−
∇
I
(
X
)
=
0
(3)
\nabla I(X+w)-\nabla I(X)=0 \tag{3}
∇I(X+w)−∇I(X)=0(3)众所周知,灰值恒定假设在光照变化中是可靠的,将灰值恒定假设和梯度恒定假设相结合的数据项可以得到准确的光学估计,但梯度恒定假设与平移运动不同。虽然存在三阶恒定假设,如哈森恒定假设和拉普拉斯恒定假设,它们比梯度恒定假设有更多的方向信息,但高阶微分对图像噪声更敏感。
首先由福斯特纳[19]和Bigun[20]首先提出的图像局部结构张量可以用后面的张量矩阵来表示。
T
=
G
σ
∗
(
∇
I
)
T
(
∇
I
)
=
[
I
x
2
I
x
I
y
I
y
I
x
I
y
2
]
(4)
T=G_{\sigma} *(\nabla I)^{T}(\nabla I)=\left[\begin{array}{cc} I_{x}^{2} & I_{x} I_{y} \\ I_{y} I_{x} & I_{y}^{2} \end{array}\right] \tag{4}
T=Gσ∗(∇I)T(∇I)=[Ix2IyIxIxIyIy2](4)在传统的光流策略[21]中,通常用于投影平滑项的局部结构张量。众所周知,局部结构张量包含灰度值变化的方向和大小,我们可以假设图像局部结构张量在短时间内是恒定,因此我们可以得到结构张量常量假设,可以写成:
T
(
X
+
W
)
−
T
(
X
)
=
0
(5)
T(X+W)-T(X)=0 \tag{5}
T(X+W)−T(X)=0(5)利用泰勒展开式,得到了等式的线性公式(5)如下:
{
I
x
x
u
+
I
x
y
v
+
I
x
t
=
0
(
I
x
x
I
y
+
I
x
I
y
x
)
u
+
(
I
x
y
I
y
+
I
x
I
y
y
)
v
+
I
x
t
I
y
+
I
x
I
y
t
=
0
I
y
x
u
+
I
y
y
v
+
I
y
t
=
0
(6)
\left\{\begin{array}{c} I_{x x} u+I_{x y} v+I_{x t}=0 \\ \left(I_{x x} I_{y}+I_{x} I_{y x}\right) u+\left(I_{x y} I_{y}+I_{x} I_{y y}\right) v+I_{x t} I_{y}+I_{x} I_{y t}=0 \\ I_{y x} u+I_{y y} v+I_{y t}=0 \end{array}\right. \tag{6}
⎩⎨⎧Ixxu+Ixyv+Ixt=0(IxxIy+IxIyx)u+(IxyIy+IxIyy)v+IxtIy+IxIyt=0Iyxu+Iyyv+Iyt=0(6)不难发现等式(6)适用于各种运动,因为它包含了比梯度常量假设更多的信息。我们结合灰色值恒定假设和结构张量恒定假设,创建了一个准确和稳健的数据项如下:
E
D
a
t
a
=
∬
Ω
{
ψ
[
(
I
(
X
+
w
)
−
I
(
X
)
)
2
]
+
ψ
[
(
T
(
X
+
w
,
t
+
1
)
−
T
(
X
,
t
)
)
2
]
}
d
x
d
y
(7)
\begin{aligned} E_{Data}=& \iint_{\Omega}\left\{\psi\left[(I(X+w)-I(X))^{2}\right]+\right.\\ &\left.\psi\left[(T(X+w, t+1)-T(X, t))^{2}\right]\right\} dxdy \end{aligned} \tag{7}
EData=∬Ω{ψ[(I(X+w)−I(X))2]+ψ[(T(X+w,t+1)−T(X,t))2]}dxdy(7)
B. The smoothing term
为了设计平滑项,我们首先引入了Horn和Schunck[1]提出的经典平滑策略如下:
E
S
m
o
o
t
h
=
∬
Ω
(
∣
∇
u
∣
2
+
∣
∇
v
∣
2
)
d
x
d
y
(8)
E_{Smooth}=\iint_{\Omega}\left(|\nabla u|^{2}+|\nabla v|^{2}\right)dxdy\tag{8}
ESmooth=∬Ω(∣∇u∣2+∣∇v∣2)dxdy(8)这是一种线性齐次平滑的平滑策略,我们采用非线性平滑策略,在平滑项中引入惩罚函数
ψ
(
s
2
)
ψ(s^2)
ψ(s2),以保持光流的不连续性。
E
Smooth
=
∬
Ω
ψ
(
∣
∇
u
∣
2
+
∣
∇
v
∣
2
)
d
x
d
y
(9)
E_{\text{Smooth }}=\iint_{\Omega} \psi\left(|\nabla u|^{2}+|\nabla v|^{2}\right)dxdy\tag{9}
ESmooth =∬Ωψ(∣∇u∣2+∣∇v∣2)dxdy(9)在传统的光流方法中,平滑项的加权因子通常是一个固定的值,它使图像的任何区域具有相同的平滑程度。为了很好地保持图像的边缘,我们使用一个与图像梯度相关的函数作为自适应加权因子,如下[8]:
λ
=
A
+
B
(
1
/
σ
2
π
)
e
−
(
∇
I
)
2
/
2
σ
2
(10)
\lambda=A+B(1 / \sigma \sqrt{2 \pi}) e^{-(\nabla I)^{2} / 2 \sigma^{2}}\tag{10}
λ=A+B(1/σ2π)e−(∇I)2/2σ2(10)其中A和B为固定值,σ为高斯标准差,∇I为图像梯度。利用自适应加权因子,平滑项可以实现各向同性的非线性扩散,从而更好地保持图像的边缘。
C. The coarse-to-fine method for optical flow estimation
通过数据项和平滑项的设计,所提出的光流模型可以写成如下的能量函数:
E
=
∬
Ω
{
[
ψ
[
(
I
(
X
+
w
)
−
I
(
X
)
)
2
]
+
ψ
[
(
T
(
X
+
w
,
t
+
1
)
−
T
(
X
,
t
)
)
2
]
+
λ
⋅
ψ
(
∣
∇
u
∣
2
+
∣
∇
v
∣
2
)
}
d
x
d
y
E=\iint_{\Omega}\{[\psi[(I(X+w)-I(X))^{2}]+\psi[(T(X+w,t+1)-T(X, t))^{2}]+\lambda \cdot \psi(|\nabla u|^{2}+|\nabla v|^{2})\}dxdy
E=∬Ω{[ψ[(I(X+w)−I(X))2]+ψ[(T(X+w,t+1)−T(X,t))2]+λ⋅ψ(∣∇u∣2+∣∇v∣2)}dxdy
等式对应的欧拉-拉格朗日方程(11)通过最小化光流向量(u,v)可以写成如下:
{
ψ
D
1
′
⋅
I
z
I
x
+
ψ
D
2
′
⋅
(
I
x
x
I
x
z
+
I
x
y
I
y
z
+
(
I
x
x
I
y
+
I
x
I
y
x
)
(
I
x
z
I
y
+
I
x
I
y
z
)
)
−
λ
div
(
ψ
S
′
⋅
∇
u
)
=
0
ψ
D
1
′
⋅
I
z
I
y
+
ψ
D
2
′
⋅
(
I
x
y
I
x
z
+
I
y
y
I
y
z
+
(
I
x
y
I
y
+
I
x
I
y
y
)
(
I
x
z
I
y
+
I
x
I
y
z
)
)
−
λ
div
(
ψ
S
′
⋅
∇
v
)
=
0
(12)
\left\{\begin{array}{c} \psi_{D 1}^{\prime} \cdot I_{z} I_{x}+\psi_{D 2}^{\prime} \cdot\left(I_{x x} I_{x z}+I_{x y} I_{y z}+\left(I_{x x} I_{y}+\right.\right. \\ \left.\left.I_{x} I_{y x}\right)\left(I_{x z} I_{y}+I_{x} I_{y z}\right)\right)-\lambda \operatorname{div}\left(\psi_{S}^{\prime} \cdot \nabla u\right)=0 \\ \psi_{D 1}^{\prime} \cdot I_{z} I_{y}+\psi_{D 2}^{\prime} \cdot\left(I_{x y} I_{x z}+I_{y y} I_{y z}+\left(I_{x y} I_{y}+\right.\right. \\ \left.\left.I_{x} I_{y y}\right)\left(I_{x z} I_{y}+I_{x} I_{y z}\right)\right)-\lambda \operatorname{div}\left(\psi_{S}^{\prime} \cdot \nabla v\right)=0 \end{array}\right. \tag{12}
⎩⎪⎪⎨⎪⎪⎧ψD1′⋅IzIx+ψD2′⋅(IxxIxz+IxyIyz+(IxxIy+IxIyx)(IxzIy+IxIyz))−λdiv(ψS′⋅∇u)=0ψD1′⋅IzIy+ψD2′⋅(IxyIxz+IyyIyz+(IxyIy+IxIyy)(IxzIy+IxIyz))−λdiv(ψS′⋅∇v)=0(12)其中
ψ
′
(
s
2
)
ψ'(s^2)
ψ′(s2)是
ψ
(
s
2
)
ψ(s^2)
ψ(s2)相对于
s
2
s^2
s2的一阶导数,我们定义:
ψ
D
1
′
=
ψ
′
[
(
I
(
X
+
w
)
−
I
(
X
)
)
2
]
ψ
D
2
′
=
ψ
′
[
(
T
(
X
+
w
,
t
+
1
)
−
T
(
X
,
t
)
)
2
]
ψ
S
′
=
ψ
′
(
∣
∇
u
∣
2
+
∣
∇
v
∣
2
)
,
I
x
=
∂
x
I
(
X
+
w
)
,
I
y
=
∂
y
I
(
X
+
w
)
,
I
z
=
I
(
X
+
w
)
−
I
(
X
)
,
I
x
x
=
∂
x
x
I
(
X
+
w
)
,
I
x
y
=
∂
x
y
I
(
X
+
w
)
I
y
y
=
∂
y
y
I
(
X
+
w
)
,
I
x
z
=
∂
x
I
z
,
I
y
z
=
∂
y
I
z
.
\begin{array}{l} \psi_{D 1}^{\prime}=\psi^{\prime}\left[(I(X+w)-I(X))^{2}\right] \\ \psi_{D 2}^{\prime}=\psi^{\prime}\left[(T(X+w, t+1)-T(X, t))^{2}\right] \\ \psi_{S}^{\prime}=\psi^{\prime}\left(|\nabla u|^{2}+|\nabla v|^{2}\right), \quad I_{x}=\partial_{x} I(X+w), \\ I_{y}=\partial_{y} I(X+w), \quad I_{z}=I(X+w)-I(X), \\ I_{x x}=\partial_{x x} I(X+w), \quad I_{x y}=\partial_{x y} I(X+w) \\ I_{y y}=\partial_{y y} I(X+w), \quad I_{x z}=\partial_{x} I_{z}, \quad I_{y z}=\partial_{y} I_{z} . \end{array}
ψD1′=ψ′[(I(X+w)−I(X))2]ψD2′=ψ′[(T(X+w,t+1)−T(X,t))2]ψS′=ψ′(∣∇u∣2+∣∇v∣2),Ix=∂xI(X+w),Iy=∂yI(X+w),Iz=I(X+w)−I(X),Ixx=∂xxI(X+w),Ixy=∂xyI(X+w)Iyy=∂yyI(X+w),Ixz=∂xIz,Iyz=∂yIz.采用粗到细的多尺度翘曲策略,假设金字塔采样层为n。对于k层(1≤k≤n)中的迭代方案,初始化为
w
k
=
(
u
k
,
v
k
)
w^k=(u^k,v^k)
wk=(uk,vk),通过计算以下方程可以得到下一层的光流结果
w
k
+
1
=
(
u
k
+
1
,
v
k
+
1
)
w^{k+1}=\left(u^{k+1}, v^{k+1}\right)
wk+1=(uk+1,vk+1)
为了消除图像导数中的非线性,我们使用一阶泰勒展开式使
I
∗
k
+
1
I_*^{k+1}
I∗k+1呈线性:
{
I
z
k
+
1
≈
I
z
k
+
I
x
k
d
u
+
I
y
k
d
v
I
x
z
k
+
1
≈
I
x
z
k
+
I
x
x
k
d
u
+
I
x
y
k
d
v
I
y
z
k
+
1
≈
I
y
z
k
+
I
y
x
k
d
u
+
I
y
y
k
d
v
(14)
\left\{\begin{array}{l} I_{z}^{k+1} \approx I_{z}^{k}+I_{x}^{k} du+I_{y}^{k}dv \\ I_{x z}{ }^{k+1} \approx I_{xz}^{k}+I_{x x}^{k}du+I_{x y}^{k}dv \\ I_{y z}^{k+1} \approx I_{y z}^{k}+I_{y x}^{k}du+I_{yy}^{k}dv \end{array}\right.\tag{14}
⎩⎨⎧Izk+1≈Izk+Ixkdu+IykdvIxzk+1≈Ixzk+Ixxkdu+IxykdvIyzk+1≈Iyzk+Iyxkdu+Iyykdv(14)其中
w
k
+
1
=
(
u
k
+
1
,
v
k
+
1
)
≈
(
u
k
+
d
u
k
,
v
k
+
d
v
k
)
w^{k+1}=\left(u^{k+1}, v^{k+1}\right) \approx\left(u^{k}+d u^{k}, v^{k}+d v^{k}\right)
wk+1=(uk+1,vk+1)≈(uk+duk,vk+dvk)。未知的
u
k
+
1
,
v
k
+
1
u^{k+1},v^{k+1}
uk+1,vk+1在最后一层
u
k
,
v
k
u^k,v^k
uk,vk的结果和当前层的未知增量
d
u
k
,
d
v
k
du^k,dv^k
duk,dvk中被分裂。
为了消除非线性函数中的非线性,应实现内、外不动点迭代方案,让 d w k , l = ( d u k , l , d v k , l ) dw^{k, l}=\left(du^{k, l}, dv^{k, l}\right) dwk,l=(duk,l,dvk,l)表示步骤 l l l中的迭代变量, ψ ∗ ∗ k , l ψ_{**}^{k,l} ψ∗∗k,l表示等式(13)中的鲁棒惩罚函数和扩散系数在迭代k,l时。在内层迭代中,函数 ψ ∗ ∗ k ψ_{**}^{k} ψ∗∗k关于固定点 d w k , l = ( d u k , l , d v k , l ) dw^{k, l}=\left(du^{k, l}, dv^{k, l}\right) dwk,l=(duk,l,dvk,l)是不变的,它将通过内部迭代进行更新。然后,第k层的补充光流 d w k , l + 1 = ( d u k , l + 1 , d v k , l + 1 ) dw^{k, l+1}=\left(du^{k, l+1}, dv^{k, l+1}\right) dwk,l+1=(duk,l+1,dvk,l+1),可以用下面的线性方程求解。
在从粗到细的多尺度翘曲策略中,光流估计的精度和效率都与降采样因子η和η∈(0.5,1)有关。虽然光流结果的精度随着下采样因子η值的增加而提高,但时间消耗也增加了。为了估计具有可接受时间消耗的光流,这里我们设定η=0.9,内部迭代为500,外部迭代为3。
3. EXPERIMENTS
A. Measurement index
本文给出了平均角度误差(AAE)和平均端点误差(AEE)来评价该方法的质量和效率
式中
(
u
E
,
v
E
)
(u_E,v_E)
(uE,vE)表示估计的流结果,
(
u
G
,
v
G
)
(u_G,v_G)
(uG,vG)表示真实值,N为图像像素的总数。