李佳田 李顯凱 李應(yīng)蕓 李果家 錢堂慧

(昆明理工大學(xué)國(guó)土資源工程學(xué)院 昆明 650093)

?

總變差一范數(shù)光流場(chǎng)計(jì)算綜述①

李佳田②李顯凱③李應(yīng)蕓 李果家 錢堂慧

(昆明理工大學(xué)國(guó)土資源工程學(xué)院 昆明 650093)

重點(diǎn)回顧與總結(jié)了計(jì)算光流場(chǎng)的總變差線性一范數(shù)方法，簡(jiǎn)稱總變差一范數(shù)(TV-L1)方法。通過一致的形式化描述，將亮度守恒模型與梯度守恒模型過渡至全局平滑模型，并著重描述了模型的數(shù)值計(jì)算過程，透析計(jì)算原理與算法底層。最后給出了光流場(chǎng)可視化及質(zhì)量評(píng)估方法，并對(duì)總變差一范數(shù)方法中的技術(shù)難點(diǎn)與下一步工作進(jìn)行總結(jié)。

光流場(chǎng)， 總變差， 數(shù)值計(jì)算， 算法， 綜述

0 引 言

光流場(chǎng)目標(biāo)的運(yùn)動(dòng)速度在圖像上的投影[1，2]，其在目標(biāo)提取、圖像重建與圖像匹配等方面已有初步應(yīng)用[3-7]。光流場(chǎng)計(jì)算(或光流場(chǎng)估計(jì))實(shí)際上是依據(jù)圖像亮度模式的變化對(duì)物體的運(yùn)動(dòng)進(jìn)行估計(jì)，因此，解決該問題的一個(gè)切入點(diǎn)為：假設(shè)運(yùn)動(dòng)前后像素點(diǎn)保持亮度一致。然而，亮度一致性假設(shè)僅提供了一個(gè)約束方程來限定光流的兩個(gè)分量，不能獲取唯一解[8]。為此，Horn和Schunck(HS)、Lucas和Kanade(LK)等分別通過附加全局平滑約束[8]和多點(diǎn)鄰域約束[9]的方法以適定光流的解，形成了以HS算法為代表的全局法和LK算法為代表的局部法，其中全局法因能獲得稠密光流場(chǎng)、易于集成不同模型進(jìn)行改進(jìn)與獲取的光流精度高等顯著優(yōu)點(diǎn)而更受青睞。全局法通過數(shù)據(jù)約束和先驗(yàn)約束下的總誤差構(gòu)成的能量泛函來獲取全局最優(yōu)解，通常將此類方法稱為總變差方法?？傋儾罘椒ǖ膬?yōu)化技術(shù)主要分為三個(gè)方面：一是對(duì)數(shù)據(jù)模型進(jìn)行改進(jìn)，如利用光學(xué)不變測(cè)度或多模型結(jié)合的方法構(gòu)建數(shù)據(jù)約束；二是對(duì)平滑模型進(jìn)行改進(jìn)，如采用高階先驗(yàn)約束、加權(quán)的先驗(yàn)約束或各向異性的先驗(yàn)約束等；三是總體約束誤差的表達(dá)，早期的方法采用L2范數(shù)來整合單點(diǎn)約束誤差，L2范數(shù)要求誤差滿足高斯分布，這在大部分真實(shí)場(chǎng)景中不成立。因此，提出了不同的穩(wěn)健函數(shù)來優(yōu)化誤差表達(dá)，其中L1范數(shù)因效果較好而被廣泛使用，此時(shí)，認(rèn)為數(shù)據(jù)項(xiàng)與平滑項(xiàng)以線性組合方式作用于光流計(jì)算，所以稱這種方法為總變差線性一范數(shù)(total variation linearly1 norm，TV-L1)，簡(jiǎn)稱為總變差一范數(shù)，它是目前最主要的光流場(chǎng)計(jì)算方法。

通過一致的形式化描述，本文總結(jié)了總變差一范數(shù)計(jì)算方法的基本原理與模型演變，以及透視光流計(jì)算中非線性約束下的凸優(yōu)化問題。在數(shù)值計(jì)算方法上，通過變分法最小化能量泛函把問題轉(zhuǎn)變?yōu)榉蔷€性方程組的求解，考慮到計(jì)算效率及避免陷入局部最小值，歸納了固定點(diǎn)迭代和金字塔結(jié)合的求解策略，使得問題最終簡(jiǎn)化為大型稀疏線性方程組的求解。

1 基本原理與模型演變

光流場(chǎng)計(jì)算產(chǎn)生了多種基本假設(shè)模型，這些模型各有優(yōu)劣，如何選取合適的模型進(jìn)行集成以獲取高精度光流是近年來的研究熱點(diǎn)?？傋儾钍遣煌Ｐ瓦M(jìn)行集成、優(yōu)化的強(qiáng)有力工具，形式如下：

EGlobal=EData+αEPrior

(1)

式中EData表示數(shù)據(jù)項(xiàng)，用于描述原始圖像數(shù)據(jù)隨時(shí)間的變化不滿足特定恒等關(guān)系時(shí)的誤差，由數(shù)據(jù)約束模型(如亮度守恒模型、梯度守恒模型)給出；EPrior表示先驗(yàn)項(xiàng)，用于描述光流場(chǎng)變化不滿足先驗(yàn)恒等關(guān)系時(shí)的誤差，由平滑約束模型給出；α>0，表示正則化參數(shù)。

1.1 數(shù)據(jù)項(xiàng)

1.1.1 亮度守恒假設(shè)

亮度守恒假設(shè)是光流估計(jì)的基本假設(shè)，假定目標(biāo)點(diǎn)發(fā)生位移之后像素強(qiáng)度或顏色值保持不變，其隱含了三個(gè)假設(shè)條件，即目標(biāo)反射特性(如：假設(shè)目標(biāo)表面為理想散射表面)、場(chǎng)景狀態(tài)(如：假設(shè)場(chǎng)景光照不發(fā)生變化)與圖像采集傳輸(如：假設(shè)相機(jī)成像過程不產(chǎn)生暈影)。如果以f(x, y, t)表示像點(diǎn)(x, y)在t時(shí)刻的像素值，以(u(x, y, t), v(x, y, t))表示光流，則亮度守恒假設(shè)可表示為

f(x, y, t)=f(x+u, y+v, t+1)

(2)

線性化式(2)，即對(duì)其進(jìn)行一階泰勒級(jí)數(shù)展開，則有等式

(3)

fxu+fyv+ft=0

(4)

1.1.2 懲罰函數(shù)選擇

式(2)和(4)給出了單點(diǎn)約束誤差，關(guān)于圖像總體誤差的表達(dá)，Horn和Schunck采用L2范數(shù)[8]，則以亮度守恒假設(shè)為約束的數(shù)據(jù)項(xiàng)表示如下：

EData=∑x,y(fxu+fyv+ft)2

(5)

從概率分布的角度來看，采用L2范數(shù)意味著假設(shè)約束誤差滿足高斯獨(dú)立同分布，然而在多數(shù)真實(shí)場(chǎng)景中，尤其是在邊界發(fā)生遮擋的區(qū)域，此假設(shè)并不成立。Black和Anandan[10]使用Lorentzian懲罰函數(shù)有效地解決了這一問題。Lorentzian懲罰函數(shù)雖然具有很強(qiáng)的穩(wěn)健性，但因?qū)儆诜峭购瘮?shù)，使得后續(xù)優(yōu)化過程變得復(fù)雜，且容易陷入局部最小值。因此，Brox[11]和Wedel[12]提出采用L1范數(shù)克服上述缺陷，并用一個(gè)基于Charbonnier函數(shù)改寫的可微函數(shù)來近似L1范數(shù)

(6)

式中E是關(guān)于單點(diǎn)約束誤差Ex, y的向量；‖·‖1表示L1范數(shù)；是一個(gè)很小的正常數(shù)。

1.1.3 光學(xué)不變特征

亮度守恒假設(shè)是一種理想假設(shè)，在自然場(chǎng)景中由于成像設(shè)備噪聲、光照變化、表面反射以及遮擋等原因往往會(huì)使得亮度值發(fā)生變化，尤其在某些應(yīng)用情況下這種變化十分明顯，極大違背了亮度守恒假設(shè)。早期的算法[13]主要通過對(duì)圖像進(jìn)行濾波預(yù)處理來消除亮度變化的影響，但改善效果不是很明顯，后來有學(xué)者提出基于光學(xué)不變特征構(gòu)建數(shù)據(jù)項(xiàng)，其中梯度特征、紋理特征在算法中被廣泛應(yīng)用。此外，也有學(xué)者采用更復(fù)雜的特征，如：Liu采用尺度不變特征變換(SIFT)特征[14]、Mileva等采用基于彩色信息的光度不變量[15]。

(1) 梯度特征

圖像梯度是指像素值在時(shí)間、空間上的導(dǎo)數(shù)，由于數(shù)字圖像是離散化的，每一個(gè)像點(diǎn)梯度的求解與周圍像點(diǎn)有著密切關(guān)系，當(dāng)發(fā)生均勻光照變化時(shí)，像點(diǎn)強(qiáng)度值發(fā)生變化而梯度值不會(huì)改變，為此Brox等[11]提出了梯度守恒假設(shè)，即

(7)

(2) 紋理特征

圖像可以看作是結(jié)構(gòu)部分和紋理部分的線性組合，光照變化發(fā)生后，產(chǎn)生的負(fù)面影響主要集中在結(jié)構(gòu)部分，而紋理部分幾乎不受影響[12]。因此，以圖像紋理特征作為光流計(jì)算的輸入可避免光照變化的影響。紋理特征的提取通過紋理結(jié)構(gòu)分解(structure texture decomposition，STD)完成，STD基于總變差ROF(rudin osher fatemi)去噪模型實(shí)現(xiàn)[16]，對(duì)于圖像f(x)，其結(jié)構(gòu)部分fs(x)的求解模型為

(8)

其中Ω表示所求圖像范圍，θ為常數(shù)，模型(8)的求解采用原始對(duì)偶算法[17]。圖像紋理部分fτ(x)為原始圖像與結(jié)構(gòu)部分之差，即fτ(α, x)=f(x)-αfs(x)，其中α∈(0,1)，經(jīng)驗(yàn)證，α取值為0.95時(shí)，平均端點(diǎn)誤差(Average End-point Error，AEE)最小[12]。

1.1.4 多通道圖像

在構(gòu)建數(shù)據(jù)項(xiàng)時(shí)，如果僅基于單通道圖像進(jìn)行考慮，則數(shù)據(jù)項(xiàng)無法滿足多通道圖像的光流場(chǎng)計(jì)算，一個(gè)較簡(jiǎn)單的解決辦法就是以各通道的總體誤差之和來表達(dá)數(shù)據(jù)項(xiàng)：

(9)

1.2 先驗(yàn)項(xiàng)

由于每個(gè)像點(diǎn)的光流均包含了兩個(gè)分量，而式(2)或(4)僅給出一個(gè)約束方程來限定這兩個(gè)分量，即該問題是不適定的，這種現(xiàn)象被稱為“孔徑問題”[8]。因此，必須增加一個(gè)先驗(yàn)項(xiàng)來正則化該病態(tài)問題，先驗(yàn)項(xiàng)一般由平滑函數(shù)構(gòu)成。

1.2.1 一階先驗(yàn)約束

Horn和Schunck根據(jù)有限區(qū)域內(nèi)的光流具有相同變化規(guī)律的先驗(yàn)知識(shí)，提出全局平滑假設(shè)[8]，并用光流場(chǎng)的一階導(dǎo)數(shù)來表示，

(10)

約束式(10)主要應(yīng)用于2幀圖像的光流估計(jì)，當(dāng)應(yīng)用情況變?yōu)槎鄮蛘麄€(gè)圖像序列的光流場(chǎng)估計(jì)時(shí)，需要在式(10)中增加時(shí)間平滑約束，并使用單向[19]或雙向[20]時(shí)間校正技術(shù)。單向時(shí)間校正技術(shù)假設(shè)光流場(chǎng)隨時(shí)間是漸變的，即在先驗(yàn)約束中增加時(shí)間平滑項(xiàng)?u/?t和?v/?t。雙向時(shí)間校正技術(shù)需要選擇中間幀作為參考進(jìn)行雙向迭代校正，因此，只能在圖像幀數(shù)已知的情況下使用。

1.2.2 懲罰函數(shù)選擇

先驗(yàn)項(xiàng)(式(10))采用了L2范數(shù)來表示先驗(yàn)約束誤差，如1.1.2中所述，真實(shí)場(chǎng)景難以滿足L2范數(shù)的使用條件。Anandan[10]提出在先驗(yàn)項(xiàng)中，通過引進(jìn)穩(wěn)健懲罰函數(shù)來對(duì)該問題進(jìn)行修正，并在算法中對(duì)數(shù)據(jù)項(xiàng)和先驗(yàn)項(xiàng)運(yùn)用了一致的懲罰函數(shù)。實(shí)際使用中，數(shù)據(jù)項(xiàng)和先驗(yàn)項(xiàng)可不必采取一致的懲罰函數(shù)。在多樣的懲罰函數(shù)中，L1范數(shù)的懲罰函數(shù)因其顯著優(yōu)點(diǎn)而被廣泛使用[11,12]。

在先驗(yàn)項(xiàng)中，懲罰函數(shù)有兩種作用方式：一是對(duì)每一項(xiàng)先驗(yàn)約束誤差分別利用懲罰函數(shù)，再將其總和作為先驗(yàn)項(xiàng)；二是先求取所有約束誤差的平方和或絕對(duì)值之和，再對(duì)其利用懲罰函數(shù)作為先驗(yàn)項(xiàng)。第一種方式有助于獲取更加平滑的光流場(chǎng)，Black和Anandan等[10]就采用了這種作用方式；第二種方式因能有效濾除異質(zhì)點(diǎn)而應(yīng)用更為廣泛，Brox[11]、Wedel[12]與Bruhn[19]均采用了這種方式。

1.2.3 帶權(quán)先驗(yàn)約束

另一種比較常用的先驗(yàn)約束是帶權(quán)的一階約束，即基于圖像梯度對(duì)L2范數(shù)的懲罰函數(shù)進(jìn)行加權(quán)：Eprior=∑x, yw(

(11)

其中w是關(guān)于像素梯度的權(quán)重函數(shù)，在邊緣處像素梯度較大，對(duì)應(yīng)權(quán)重w(f)較小，減弱了對(duì)光流的懲罰力度，有利于保留邊緣信息。當(dāng)然，權(quán)重函數(shù)也可不必基于圖像梯度，如Seitz[21]采用分段的權(quán)重函數(shù)。

1.2.4 各向異性的先驗(yàn)約束

式(11)給出了一個(gè)各向同性的先驗(yàn)項(xiàng)，其缺點(diǎn)是無法避免先驗(yàn)約束與數(shù)據(jù)約束產(chǎn)生的相互干擾問題，應(yīng)削弱先驗(yàn)約束在數(shù)據(jù)約束方向的作用力度，增強(qiáng)在數(shù)據(jù)約束正交方向的作用力度[1]。為此，Nagel[22]和Welberger[20]提出各向異性的平滑模型，使得先驗(yàn)約束在像素梯度方向的作用力度小于像素梯度正交方向；Sun等[23]通過構(gòu)建可控隨機(jī)場(chǎng)，根據(jù)光流分布及圖像結(jié)構(gòu)自適應(yīng)調(diào)整約束方向；Zimmer等[18]提出不再基于像素梯度衡量方向信息，根據(jù)數(shù)據(jù)約束來考慮作用方向。

1.2.5 高階先驗(yàn)約束

在弱紋理區(qū)域，光流估計(jì)存在分區(qū)域恒定現(xiàn)象[1]，為此Trobin等[24]提出采用二階先驗(yàn)項(xiàng)來懲罰分區(qū)恒定，以歐式向量范數(shù)表示如下：

‖·‖=

(12)

高階的先驗(yàn)約束能產(chǎn)生相對(duì)較平滑的光流場(chǎng)，但在邊緣處容易產(chǎn)生過度平滑現(xiàn)象，為此，Nir等[21，25]提出采用仿射平滑先驗(yàn)約束獲取更精確的光流場(chǎng)。首先參數(shù)化光流，然后對(duì)于每一像點(diǎn)不再直接求解光流矢量(u(x, y, t), v(x, y, t))，而是求解6個(gè)仿射參數(shù)，通過所求參數(shù)來表示光流：

(13)

(14)

其中(x0, y0)是圖像的中點(diǎn)，通過式(13)和(14)可取代數(shù)據(jù)項(xiàng)直接進(jìn)行光流求解。

1.3 總變差一范數(shù)模型

1.3.1 復(fù)合數(shù)據(jù)項(xiàng)構(gòu)建

考慮到光照影響及特征提取的簡(jiǎn)易性，梯度特征是構(gòu)建數(shù)據(jù)項(xiàng)的較好選擇。然而，由于梯度守恒假設(shè)包含了方向信息使得約束式(7)僅在物體發(fā)生平移運(yùn)動(dòng)時(shí)成立[26]，且在后續(xù)優(yōu)化過程中線性化能量函數(shù)時(shí)圖像需滿足二次可微條件，在多運(yùn)動(dòng)模式并存或一些更復(fù)雜的運(yùn)動(dòng)模式中，上述限制將會(huì)嚴(yán)重降低算法性能[1]。這些缺陷在亮度守恒假設(shè)中是完全可以避免的，因此，將基于亮度守恒模型和梯度守恒模型構(gòu)建復(fù)合數(shù)據(jù)項(xiàng)。

線性化的亮度守恒約束式(4)建立在兩個(gè)基本假設(shè)之上：圖像像素強(qiáng)度是連續(xù)的；物體的運(yùn)動(dòng)是微小的。然而，自然場(chǎng)景中目標(biāo)的運(yùn)動(dòng)難以嚴(yán)格滿足這兩個(gè)條件，因此，在建模階段依然沿用原始的非線性守恒約束式(2)和式(7)[11]。根據(jù)1.1.2中的描述，廣義Charbonnier懲罰函數(shù)表示如下：

(15)

定義x=(x, y, t)T， w=(u, v, 1)τ，則復(fù)合數(shù)據(jù)項(xiàng)表示為

ED(u, v)=∫Ωψ(|f(x+w)-f(x)|2

(16)

其中γ是兩個(gè)假設(shè)之間的權(quán)重，用于調(diào)節(jié)模型間的相對(duì)作用程度。

1.3.2 先驗(yàn)項(xiàng)選擇

在1.2中總結(jié)了5種不同類型的先驗(yàn)項(xiàng)，從作用方式看，這些模型主要分為兩大類：全局平滑模型和有向平滑模型[27]。其中有向平滑模型對(duì)噪聲比較敏感，且計(jì)算較為復(fù)雜，應(yīng)用不如全局平滑模型廣泛?；诖?，選用一階先驗(yàn)約束，并采用L1范數(shù)的懲罰函數(shù)式(15)來保持邊緣信息，用光流場(chǎng)的空間梯度算子表示[11]，約束式(10)改寫為

EP(u,v)=∫Ωψ(|u|2+|v|2)dx

(17)

式(16)和(17)分別給出了數(shù)據(jù)項(xiàng)和先驗(yàn)項(xiàng)，根據(jù)式(1)可得出總變差一范數(shù)光流場(chǎng)計(jì)算模型[14]

EG(u, v)=ED+αEP

(18)

現(xiàn)在的目標(biāo)是需要找到使得式(18)所示泛函EG(u, v)最小的函數(shù)u和v，即問題轉(zhuǎn)變?yōu)镋G(u, v)的最優(yōu)化問題，通常通過連續(xù)優(yōu)化方法實(shí)現(xiàn)。

2 優(yōu)化方法

在光流場(chǎng)計(jì)算中常用的連續(xù)優(yōu)化技術(shù)有兩種[26]：梯度下降法； 全變分法或極值法。近年也提出了一些比較新的方法，如Trobin[24]和Wedel[12]等采用的解耦算法、Seitz和Baker[21]提出的基于線性規(guī)劃的凸優(yōu)算法，但這些方法仍不成熟。梯度下降算法中以最速下降法[26]最為簡(jiǎn)明，但其需要大量迭代計(jì)算才能達(dá)到收斂，為此提出了Hessian矩陣的方法[28]，當(dāng)圖像含有N個(gè)像素點(diǎn)時(shí)，該方法需要求解2N×2N的矩陣及其逆矩陣，因此，梯度下降優(yōu)化方法僅適用于以LK算法為代表的一類求取未知量較少的算法[9，29]。對(duì)于能量函數(shù)式(18)的優(yōu)化，采用全變分法比較合適，根據(jù)變分原理，如果式(18)存在最小值，則必滿足歐拉—拉格朗日方程[11]：

(19)

另一種與全變分法類似的方法是極值法，以κ表示每個(gè)像點(diǎn)的光流(u,v)，極值法通過?EG/?κ=0來獲取約束方程組，該方法與全變分法的區(qū)別僅在于能量函數(shù)的逐像素化是在式(19)之前還是之后進(jìn)行[26]。選用全變分法進(jìn)行優(yōu)化，由于EG(u,v)是高度非線性的，最小化過程比較繁瑣，因此，為了提高可讀性，定義下列縮寫表示：

(20)

其中用含z的變量表示特定的差值便于同時(shí)間導(dǎo)數(shù)進(jìn)行區(qū)分。

用上述變量表示式(19)的結(jié)果如下：

(21)

3 模型的數(shù)值計(jì)算方法

通過變分最小化后，每個(gè)像素點(diǎn)將產(chǎn)生如式(21)所示兩個(gè)非線性方程，假設(shè)圖像像素點(diǎn)個(gè)數(shù)為N，則結(jié)果為一個(gè)含有2N個(gè)歐拉-拉格朗日方程式的大型非線性方程組。此時(shí)，光流場(chǎng)計(jì)算簡(jiǎn)化為該方程組的求解，通過雙重不動(dòng)點(diǎn)迭代策略[11]消除方程組非線性后可把問題進(jìn)一步簡(jiǎn)化為稀疏線性方程組的求解，運(yùn)用一般的線性方程組解法[26](如雅克比迭代、高斯-賽德爾迭代、超松弛迭代等)即可解出光流。同時(shí)，為了避免陷入局部最小值及加快迭代收斂速率，采用不動(dòng)點(diǎn)迭代和多尺度方法相結(jié)合的求解策略[11，30]。

3.1 多尺度運(yùn)動(dòng)估計(jì)

實(shí)現(xiàn)多尺度運(yùn)動(dòng)估計(jì)主要依靠由粗到細(xì)的金字塔結(jié)構(gòu)實(shí)現(xiàn)，首先通過下采樣建金字塔來構(gòu)造不同的尺度空間，然后光流估計(jì)從粗尺度開始進(jìn)行，在獲取粗略解后，將粗略解映射到精細(xì)尺度下做進(jìn)一步計(jì)算。尺度越粗糙，所含像點(diǎn)越少，收斂速度越快，以粗略的解作為精細(xì)尺度解的初始值有利于減少迭代次數(shù)，提高算法的效率[1]；且粗尺度下的圖像含有較少的高頻分量，能有效減少局部最小值的數(shù)量；在粗尺度下，光流較小，因此對(duì)亮度守恒、梯度守恒進(jìn)行線性化是合理的[26]。

3.1.1 高斯金字塔構(gòu)建

高斯金字塔構(gòu)建主要包括高斯濾波和降采樣[27]。首先采用高斯濾波作用于原始圖像以消除噪聲，降低異質(zhì)點(diǎn)影響，然后通過降采樣為金字塔每一層獲取數(shù)據(jù)。

(1) 高斯濾波

高斯濾波是一種線性平滑濾波，被廣泛用于圖像去噪，其實(shí)質(zhì)是對(duì)圖像進(jìn)行加權(quán)平均的過程[31]。定義濾波器gfilter，根據(jù)正態(tài)分布3原則可知位于當(dāng)前像點(diǎn)3σ個(gè)像素外的像點(diǎn)基本不會(huì)對(duì)當(dāng)前像點(diǎn)產(chǎn)生影響，則可設(shè)置濾波器寬度k=3σ×2+1，得到下式：

(22)

其中i表示濾波器中元素的位置，di表示其對(duì)應(yīng)點(diǎn)距離當(dāng)前像點(diǎn)的距離，σ通常設(shè)置為1。實(shí)現(xiàn)高斯濾波的方式一般有兩種，一是采用濾波模板與圖像做卷積，另一種是通過傅里葉變換實(shí)現(xiàn)。σ為1時(shí)，模板大小為7×7，卷積運(yùn)算量較小。由于高斯函數(shù)有可分離性，可采用分離濾波器來加快速度，即把二維窗口卷積分離為兩次一維卷積運(yùn)算。

(2) 降采樣

降采樣過程是建立高斯金字塔的關(guān)鍵步驟，整個(gè)過程有兩方面需要注意：

① 金字塔采樣因子和層數(shù)。金字塔采樣因子ratio是指當(dāng)前層與下一層圖像行或列的比值，一般采用隔行隔列采樣，即ratio=0.5。金字塔層數(shù)直接關(guān)系著光流計(jì)算精度及計(jì)算效率，通常由一個(gè)給定的最小寬度來控制。

② 下采樣插值。下采樣插值主要用于構(gòu)建金字塔過程中不是間隔固定行列采樣時(shí)非整數(shù)坐標(biāo)取值問題。主要的插值方法有雙線性插值[7]和雙三次插值[32]，后者獲得的圖像具有更好的平滑性。

3.1.2 光流逐層估計(jì)

逐層估計(jì)包含了金字塔不同層級(jí)間的增量校正過程及同一層級(jí)間的光流迭代更新過程，其計(jì)算過程將在3.2節(jié)給出，本節(jié)著重關(guān)注逐層估計(jì)中的Warping處理。Warping處理是指通過已解出的光流在第2幀上封裝出一幅新圖像[11]，其主要作用于兩個(gè)過程：一是在低分辨率層光流場(chǎng)向高分辨率層映射；二是內(nèi)層迭代(詳見3.3節(jié))結(jié)束后通過Warping處理縮小與第1幀圖像的殘差。其中，對(duì)于非整數(shù)坐標(biāo)取值問題采用雙線性或雙三次插值解決，光流(u，v)越界問題通過下式處理：

(23)

3.2 外層不動(dòng)點(diǎn)迭代

(24)

然后把未知量ul+1，vl+1分離為前一步的迭代結(jié)果ul，vl和未知增量dul，dvl。設(shè)

(25)

3.3 內(nèi)層不動(dòng)點(diǎn)迭代

(26)

3.4 梯度、散度離散近似

采用有限差分近似像素梯度與光流梯度。設(shè)定當(dāng)前點(diǎn)位置(x, y)，對(duì)于寬M、高N的圖像，光流增量du水平方向梯度的前向差分如下：

(27)

豎直方向梯度，

(28)

梯度的散度采用后向差分結(jié)果表達(dá)，目的是通過在前向差分結(jié)果上進(jìn)行后向差分等價(jià)實(shí)現(xiàn)二階中心差分：

(29)

式(26)經(jīng)過差分近似后可得到兩個(gè)線性方程式，則圖像上所有像點(diǎn)對(duì)應(yīng)的方程式就構(gòu)成了一個(gè)大型的線性方程組。其中，待求解的列向量即為光流場(chǎng)增量，為加快收斂速度，采用超松弛迭代法[33]求解。

總變差一范數(shù)光流算法流程圖如圖1所示。

4 光流場(chǎng)可視化與質(zhì)量評(píng)價(jià)

4.1 可視化

通過總變差一范數(shù)算法，可以獲得一個(gè)致密的光流場(chǎng)，如何在圖像上以直觀的方式對(duì)光流場(chǎng)進(jìn)行顯示也是一個(gè)重要的問題。目前光流場(chǎng)顯示主要有三種方法：灰度圖、二維矢量圖與彩色蘊(yùn)含圖[26]。

(1) 灰度圖。對(duì)于流場(chǎng)的兩個(gè)分量場(chǎng)，由于其值有正有負(fù)，因此將光流絕對(duì)值化后再進(jìn)行歸一化處理，就可以用兩個(gè)灰度圖來進(jìn)行顯示，但是這種顯示方案在絕對(duì)值化時(shí)舍去了流場(chǎng)的方向信息，因此只能用于運(yùn)動(dòng)物體檢測(cè)，無法提供更多運(yùn)動(dòng)相關(guān)信息。

(2) 二維矢量圖。二維矢量圖是比較直觀顯示方案，主要通過有向線段來表示光流矢量。由于每個(gè)像點(diǎn)具有一個(gè)光流矢量，為了觀察方便，通常采用稀疏流場(chǎng)或放大的流場(chǎng)顯示。

圖1 總變差一范數(shù)光流算法流程圖

(3) 彩色蘊(yùn)含圖。由于人眼對(duì)彩色的分辨能力比對(duì)灰度的強(qiáng)，因此基于彩色信息來表達(dá)光流有利于肉眼觀察。彩色蘊(yùn)含圖是一種比較常用的顯示方法，它一般基于面向視覺感知的顏色模型(如HIS、HSV等)來進(jìn)行顯示。其中，最常用的是HIS模型，它采用三個(gè)參數(shù)—色調(diào)、飽和度、亮度來描述顏色。當(dāng)飽和度、亮度固定時(shí)，把最敏感的紅、綠、藍(lán)三種色調(diào)間的連續(xù)變化用離散圓域來表示，構(gòu)成孟塞爾色輪。由于光流場(chǎng)是一個(gè)二維向量場(chǎng)，因此，可以將光流值與色輪上的顏色值一一對(duì)應(yīng)，按相應(yīng)顏色值進(jìn)行顯示即可得到一個(gè)彩色表示的光流場(chǎng)。

4.2 質(zhì)量評(píng)價(jià)

(30)

AAE定義為

(31)

5 總結(jié)與展望

光流場(chǎng)算法一直致力于對(duì)孔徑問題、成像噪聲、運(yùn)動(dòng)不連續(xù)性、無紋理區(qū)域、光照變化、大位移光流、非剛體運(yùn)動(dòng)、小目標(biāo)等方面進(jìn)行改善，并取得了顯著成果，但要在實(shí)際中深入應(yīng)用，還需著重解決以下幾個(gè)方面的問題：

(1) 設(shè)計(jì)合適的優(yōu)化算法。目前廣泛使用的優(yōu)化方法都屬于連續(xù)優(yōu)化算法，此類算法通過數(shù)值迭代獲取光流解，但每次迭代只能實(shí)現(xiàn)微小的更新，使得算法效率不高。許多學(xué)者提出了離散優(yōu)化算法(如圖割法[35]、置信度傳播[36]等)。離散優(yōu)化算法是在一個(gè)事先定義的離散的光流空間里搜索最優(yōu)解，有利于可提高收斂速度，但此類方法尚屬于起步階段，存在著許多問題需要解決[1]。

(2) 復(fù)雜光照變化的影響。由于受陰影遮擋、拍攝角度與成像傳輸設(shè)備的影響，光流算法對(duì)亮度變化異常敏感，無法進(jìn)行精確光流估計(jì)。為此，對(duì)圖像進(jìn)行濾波預(yù)處理、尋求光學(xué)不變特征、多模型結(jié)合等方法被用來改善算法在光照變化下的估計(jì)效果，此舉在一定程度上降低了光照變化帶來的不利影響，但由于真實(shí)場(chǎng)景中光照變化的復(fù)雜性，目前還沒有一種算法能夠自適應(yīng)地解決復(fù)雜光照變化問題[1]，該問題嚴(yán)重限制著光流技術(shù)在實(shí)際中的應(yīng)用。

(3) 大位移運(yùn)動(dòng)估計(jì)。雖然采用了由粗到細(xì)的金字塔法來解決大位移光流的獲取，但是所能估計(jì)的位移尺度依然有限，Alvarez等[37]通過對(duì)基本光流場(chǎng)進(jìn)行多方面改進(jìn)后能解決10個(gè)像素大小的位移估計(jì)問題，但如何獲取更大位移的光流依然是值得考慮的問題。由于傳統(tǒng)的金字塔法提供的改進(jìn)空間有限，可考慮融合其他圖像特征作為引導(dǎo)信息來降低最優(yōu)化的搜索范圍及難度[38]。其次是小目標(biāo)的大位移問題，為實(shí)現(xiàn)多尺度的方法，需要構(gòu)造許多不同的尺度空間，這樣將會(huì)造成低分辨率下小目標(biāo)的丟失而無法估計(jì)其光流，融合通道表示下描述符匹配的方法[39]可以實(shí)現(xiàn)對(duì)小目標(biāo)的捕捉，但是該方法不能有效處理結(jié)構(gòu)復(fù)雜的運(yùn)動(dòng)圖像，因此如何更好地對(duì)小目標(biāo)的運(yùn)動(dòng)進(jìn)行估計(jì)依然值得深入研究。

(4) 算法實(shí)時(shí)性。當(dāng)光流技術(shù)應(yīng)用到機(jī)器人導(dǎo)航、視頻運(yùn)動(dòng)目標(biāo)跟蹤這些場(chǎng)景下時(shí)，必須考慮算法實(shí)時(shí)計(jì)算性能。為了求解非線性系統(tǒng)及獲取大位移光流，需要采用不動(dòng)點(diǎn)迭代、多尺度運(yùn)動(dòng)估計(jì)等求解策略，這就會(huì)極大地增大算法的時(shí)間復(fù)雜度，且整個(gè)估計(jì)過程需要逐層、逐點(diǎn)進(jìn)行迭代計(jì)算，隨著圖像增大，大量浮點(diǎn)運(yùn)算將會(huì)給CPU帶來沉重負(fù)荷，因而設(shè)計(jì)結(jié)合GPU的并行處理架構(gòu)或采用多核處理器將會(huì)為實(shí)時(shí)系統(tǒng)的構(gòu)建提供有力支持[1，38]。

[1] 宋爽, 楊健, 王涌天. 全局光流場(chǎng)估計(jì)技術(shù)及展望. 計(jì)算機(jī)輔助設(shè)計(jì)與圖形學(xué)學(xué)報(bào), 2014, 26(5): 841-850

[2] 涂志剛, 謝偉, 熊淑芬等. 一種高精度的TV-L1光流算法. 武漢大學(xué)學(xué)報(bào)·信息科學(xué)版, 2012, 37(4): 496-499

[3] 范沖, 龔健雅, 朱建軍等. ALOS-PRISM遙感影像超分辨率重建. 遙感學(xué)報(bào), 2009, 13(1): 75-82

[4] 陳王麗, 孫濤, 陳喆等. 利用光流配準(zhǔn)進(jìn)行嫦娥一號(hào)CCD多視影像超分辨率重建. 武漢大學(xué)學(xué)報(bào)·信息科學(xué)版, 2014, 39(9): 1103-1108

[5] 閆利, 鞏翼龍, 張毅等. 光流動(dòng)態(tài)紋理在土地利用/覆蓋變化檢測(cè)研究中的應(yīng)用. 光譜學(xué)與光譜分析, 2014, 34(11): 3056-3061

[6] 張正鵬, 江萬壽, 張靖. 光流特征聚類的車載全景序列影像匹配方法. 測(cè)繪學(xué)報(bào), 2014,43(12): 1266-1273

[7] 王昊京, 王建立, 王鳴浩等. 采用雙線性插值收縮的圖像修復(fù)方法. 光學(xué)精密工程, 2010, 18(5): 1234-1241

[8] Horn B K P, Schunck B G. Determining optical flow.ArtificialIntelligence, 1981, 17(1):185-203

[9] Lucas B D, Kanade T. An iterative image registration technique with an application to stereo vision. In: Proceedings of the 7th International Joint Conference on Artificial Intelligence, San Francisco, USA, 1981, 81: 674-679

[10] Black M J, Anandan P. The Robust Estimation of Multiple Motions: Parametric and Piecewise-Smooth Flow Fields.ComputerVisionandImageUnderstanding, 1996, 63(1): 75-104

[11] Brox T, Bruhn A, Papenberg N, et al. High accuracy optical flow estimation based on a theory for warping. In: Proceedings of the 8th European Conference on Computer Vision, Heidelberg: Springer, 2004, 3024: 25-36

[12] Wedel A, Pock T, Zach C, et al. An improved algorithm for TV-L1 optical flow. In: Lecture Notes in Computer Science. Heidelberg: Springer, 2008, 5064: 23-45

[13] Anandan P. A computational framework and an algorithm for the measurement of visual motion.InternationalJournalofComputerVision, 1989, 2(3): 283-310

[14] Liu C, Yune J, Torralba A, et al. SIFT flow: dense correspondence across difference scenes. In： Proceedings of the 10th European Conference on Computer Vision, 2008, 5304:28-42

[15] Mileva Y, Bruhn A, Weickert J. Illumination robust variational optical flow with photometric invariants. In: Proceedings of the DAGM Symposium on Pattern Recognition. Heidelberg: Springer-Verlag, 2007. 152-159

[16] Rudin Li, Osher S, Fatemi E. Nonlinear total variation based noise removal algorithms.PhysicaD, 1992, 60(1-4): 259-268

[17] 李秀智, 譚君, 賈松敏等. 基于結(jié)構(gòu)紋理分解和多重網(wǎng)格的光流估計(jì)算法. 上海交通大學(xué)學(xué)報(bào), 2014, 48(7): 959-964, 970

[18] Zimmer H, Bruhn A, Weickert J, et al. Complementary optic flow. In: Proceedings of the 7th international workshop on energy minimization methods in computer vision and pattern recognition, Bonn, Germany, 2009, 5681: 207-220

[19] Bruhn A, Weickert J, Schnorr C. Combine local and global optical flow methods.InternationalJournalofComputerVision, 2005, 61(3): 211-231

[20] Werlberger M, Trobin W, Pock T, et al. Anisotropic Huber-L1 optical flow. In: Proceedings of the British machine vision conference, Manchester, UK, 2009. 1-11

[21] Seitz S M, Baker S. Filter flow. In: Proceedings of the 12th IEEE International Conference on Computer Vision. Los Alamitos, USA, 2009. 143-150

[22] Nagel H, Enkelmann W. An investigation of smoothness constraints for the estimation of displacement vector fields from image sequences.IEEETransactionsonPatternAnalysisandMachineIntelligence, 1986, 8(5): 565-593

[23] Sun D, Roth S, Lewis J, et al. Learning optical flow. In: Proceedings of the 10th European conference on computer vision, Heidelberg: Springer, 2008, 5304: 83-97

[24] Trobin W, Pock T, Cremers D, et al. An unbiased second-order prior for high-accuracy motion estimation. In: Proceedings of the 30th DAGM Symposium Munich, Heidelberg: Springer, 2008, 5096: 396-405

[25] Nir T, Bruckstein A M, Kimmel R. Over-parameterized variational optical flow.InternationalJournalofComputerVision, 2008, 76(2): 205-216

[26] Baker S, Scharstein D, Lewis J P, et al. A Database and Evaluation Methodology for Optical Flow.InternationalJournalofComputervision, 2011, 92(1): 1-31

[27] 項(xiàng)學(xué)智, 開湘龍, 張磊等. 一種變分偏微分多模型光流求解方法. 儀器儀表學(xué)報(bào), 2014, 35(1): 109-116

[28] Baker S, Matthews I. Lucas-Kanade 20 years on: a unifying framework.InternationalJournalofComputerVision, 2004, 56(4): 221-255

[29] Le Besnerais G, Champagnat F. Dense optical flow by iterative local window registration. In: Proceedings of the International Conference on Image Processing, 2005, 1:137-140

[30] Lei C, Yang Y H. Optical flow estimation on coarse-to-fine region-trees using discrete optimization. In: Proceedings of the 12th IEEE International Conference on Computer Vision. Los Alamitos: IEEE Computer Society Press, 2009, 1562-1569

[31] 李秀智, 尹曉琳, 賈松敏等. 改進(jìn)的TV-L1平滑光流估計(jì). 光學(xué)學(xué)報(bào), 2013, 33(10): 188-194

[32] 龐志勇, 譚洪舟, 陳弟虎. 一種改進(jìn)的低成本自適應(yīng)雙三次插值算法及VLSI實(shí)現(xiàn). 自動(dòng)化學(xué)報(bào), 2013, 39(4): 407-417

[33] 李慶揚(yáng), 王能超, 易大義. 數(shù)值分析. 北京：清華大學(xué)出版社, 2006

[34] Baker S, Scharstein D, Lewis J P, et al. The middlebury computer vision pages [BE/OL]. http://vision.middlebury.edu/flow/

[35] 路子赟, 唐土生, 高雋等. 帶有局部控制因子的圖割光流估計(jì). 中國(guó)圖像圖形學(xué)報(bào), 2012, 17(7): 846-852

[36] 郁理, 郭立, 袁紅星. 基于分級(jí)置信度傳播的立體匹配新方法. 中國(guó)圖像圖形學(xué)報(bào), 2011, 16(1): 103-109

[37] Alvarez L, Weickert J, Sanchez J. Reliable estimation of dense optical flow field with large displacements.InternationalJournalofComputerVision, 2000, 39(1): 41-56

[38] 李秀智, 賈松敏, 尹曉琳等. 視覺光流矢量場(chǎng)估計(jì)算法綜述. 北京工業(yè)大學(xué)學(xué)報(bào), 2013, 39(11): 1638-1643

[39] Laura S, Deqing S, Erik G, et al. Optical flow estimation with channel constancy. In: Proceedings of the 13th European Conference, Switzerland: Springer, 2014, 8689: 423-438

A review on TV-L1optical flow field calculation

Li Jiatian, Li Xiankai, Li Yingyun, Li Guojia, Qian Tanghui

(Faculty of Land Resource Engineering, Kunming University of Science and Technology, Kunming 650093)

The method for optical flow field calculation using the Total Variation Linearly 1 norm, called the TV-L1method for short, is reviewed and summarized emphatically. By means of the coherent formal description, the bright conservation model and the gradient conservation model are transferred to a global smooth model, and the models’ numerical computing progress, analytical calculation principle and physical level of algorithm are emphatically described. At the end, the method for visualization of optical flow fields and the method for quality assessment optical flow fields are given, and the technical difficulties of the TV-L1and its future work are concluded.

optical flow field, total variation, numerical calculation, algorithm, review

10.3772/j.issn.1002-0470.2016.04.004

①國(guó)家自然科學(xué)基金(41561082， 41161061)資助項(xiàng)目。

，E-mail: 1132468376@qq.com(

2016-01-26)

②男，1975年生，博士，碩士生導(dǎo)師，副教授；研究方向：視覺測(cè)量與機(jī)器場(chǎng)景理解；E-mail： ljtwcx@163.com