專利名稱：用于粒子圖像測速中的多網(wǎng)格處理方法
技術(shù)領(lǐng)域：
本發(fā)明涉及圖像處理，特別是涉及一種用于粒子圖像測速中的多網(wǎng)格 處理方法。
背景技術(shù)：
自然界圖像的運(yùn)動形式千差萬別，流體運(yùn)動是一種典型的非剛體運(yùn)
動，流體運(yùn)動圖像計(jì)算與分析即粒子圖像測速(Particle Image Velocimitry，簡稱PIV)是一種新型的非接觸式的測量技術(shù)。通常PIV是 在流體中投入示蹤粒子，在激光片光的照射下，在與片光垂直方向用攝像 機(jī)拍攝隨流體運(yùn)動的粒子圖像，再對圖像進(jìn)行處理分析計(jì)算，最終得到流 場橫切面上二維速度場的一種測量方法。作為一種全流場、無接觸、無擾 動、高精度的流動可視化方法，PIV適用于湍流、非定常流等時(shí)變復(fù)雜流 場的測量。如今，PIV已是一門跨學(xué)科的交叉綜合技術(shù)，其結(jié)合了激光技 術(shù)、視頻圖像處理技術(shù)、計(jì)算機(jī)技術(shù)、流體力學(xué)和近代光學(xué)技術(shù)的最新成 果，廣泛應(yīng)用于實(shí)驗(yàn)流體力學(xué)、生物醫(yī)學(xué)、工業(yè)制造等許多領(lǐng)域。PIV將整 個(gè)運(yùn)動矢量場作為一個(gè)整體連續(xù)函數(shù)進(jìn)行估計(jì)，其分辨率可以達(dá)到圖像中 的每一個(gè)像素，為微觀、精確地進(jìn)行局部圖像運(yùn)動計(jì)算和分析提供了可能， 同時(shí)克服了傳統(tǒng)相關(guān)方法中的一些固有缺陷。PIV系統(tǒng)在實(shí)驗(yàn)流體力學(xué)、 空氣動力學(xué)、制造業(yè)、醫(yī)學(xué)、工業(yè)、飛機(jī)制造、水利水電、氣象學(xué)等領(lǐng)域 都有著重要的科學(xué)和經(jīng)濟(jì)價(jià)值。
許多圖像分析問題需要解決帶有邊值問題的橢圓偏微分方程，并最終 通過求解線性或非線性方程組實(shí)現(xiàn)如基于變分方法的圖像非線性擴(kuò)散、 主動輪廓計(jì)算、圖像恢復(fù)等。雖然光流計(jì)算無論在計(jì)算精度、可靠性等方 面都取得了較大的飛躍，但其數(shù)值方法和手段并沒有重大突破。
多網(wǎng)格思想提出于上世紀(jì)60年代，70年代后被應(yīng)用于實(shí)際科學(xué)計(jì)算 中。如今多網(wǎng)格計(jì)算已經(jīng)廣泛應(yīng)用于物理學(xué)、力學(xué)、電磁場等科學(xué)計(jì)算，工 程數(shù)學(xué)和有限元數(shù)值計(jì)算中廣泛使用的多網(wǎng)格方法。網(wǎng)格是偏微分方程數(shù) 值解法的基礎(chǔ)，網(wǎng)格體系的好壞直接影響計(jì)算結(jié)果的精度。多網(wǎng)格方法的 研究經(jīng)歷了從結(jié)構(gòu)化到非結(jié)構(gòu)化，從單一網(wǎng)格到混合網(wǎng)格的過程，并不斷
3出現(xiàn)新的針對不同情況的網(wǎng)格生成技術(shù)。

發(fā)明內(nèi)容
本發(fā)明的主要目的就是針對現(xiàn)有技術(shù)的不足，提供一種多網(wǎng)格處理方 法，能有效加快粒子圖像測速的圖像處理過程。 為實(shí)現(xiàn)上述目的，本發(fā)明采用以下技術(shù)方案
一種用于粒子圖像測速中的多網(wǎng)格處理方法,所述粒子圖像測速包括
求解線性方程組AU-F以實(shí)現(xiàn)對運(yùn)動場的最佳逼近的步驟，U為待求解， A為系數(shù)矩陣，其至少包含運(yùn)動矢量場參數(shù)，F(xiàn)為向量，所述方法包括以
下步驟
A、 設(shè)置第一網(wǎng)格，獲得方程組AU-F對應(yīng)在第一網(wǎng)格h下的方程組 AhUh = Fh;
B、 迭代求解方程組AhUh-Fh，根據(jù)所得實(shí)際解tJh得到方程誤差 Rh =Fh-AhtJh，
C、 設(shè)置尺寸較第一網(wǎng)格更粗的第二網(wǎng)格H，將在第一網(wǎng)格h下對方程 組AV =Rh的求解尺度映射為在第二網(wǎng)格H下對方程組AV1 =RH的求解；
D、 迭代求解在第二網(wǎng)格H下的高頻計(jì)算誤差eH;
E、 對第二網(wǎng)格H下的高頻計(jì)算誤差eH進(jìn)行尺度映射，得到在第一網(wǎng) 格h下的低頻計(jì)算誤差e、
F、 用第一網(wǎng)格h下的低頻計(jì)算誤差eh更新實(shí)際解tI|;e,v=eh+Ch。 優(yōu)選地，還包括
在更新實(shí)際解后，再在第一網(wǎng)格h下迭代求解AhUh = Fh ，來消除t^ew的 高頻計(jì)算誤差。
優(yōu)選地，循環(huán)多次進(jìn)行所述處理，其中第二網(wǎng)格的尺度逐次遞增。 本發(fā)明有益的技術(shù)效果是
圖像分析問題需要解決帶有邊值問題的橢圓偏微分方程，并最終通過 求解線性或非線性方程組實(shí)現(xiàn)，如基于變分方法的圖像非線性擴(kuò)散、主 動輪廓計(jì)算、圖像恢復(fù)等，在大尺寸網(wǎng)格上線性方程的高頻誤差很容易消 除，相反即使經(jīng)過大量迭代低頻誤差卻很難消除，但在小尺度網(wǎng)格上這種 低頻誤差只需經(jīng)少量迭代就能明顯消除。因此，在粒子圖像測速的圖像處 理過程中，采用多網(wǎng)格方法，通過不同網(wǎng)格的迭代進(jìn)行相關(guān)線性方程組的 求解處理，可以有效加快線性方程組迭代的收斂過程，大大加快圖像處理 速度，同時(shí)不會降低處理的可靠性和精度。采用網(wǎng)格多尺度(有可能尺度因子小于2)逐級逼近的方法，可以降低 方程模型誤差帶來的影響，改善計(jì)算精度。每次尺度變化都會進(jìn)行一次基 于當(dāng)前估計(jì)結(jié)果的補(bǔ)償更新。將圖像尺度的變化看作是不同尺寸網(wǎng)格的變 化，本質(zhì)上多尺度逼近計(jì)算可以看作是多網(wǎng)格算法的一種特例。


圖1為粒子圖像測速處理系統(tǒng)的結(jié)構(gòu)框圖2為用水平集函數(shù)來劃分內(nèi)部光滑連續(xù)而邊界不連續(xù)的場的示意
圖3為一種實(shí)施例中利用多網(wǎng)格進(jìn)行加速處理的流程圖； 圖4a和圖4b分別為V型和W型多網(wǎng)格結(jié)構(gòu)的示意圖； 圖5為一種優(yōu)選實(shí)施例的多網(wǎng)格結(jié)構(gòu)的示意圖。 本發(fā)明的特征及優(yōu)點(diǎn)將通過實(shí)施例結(jié)合附圖進(jìn)行詳細(xì)說明。
具體實(shí)施例方式
首先介紹可應(yīng)用本發(fā)明實(shí)施例的多網(wǎng)格處理方法的一種粒子圖像測 速方法。
一個(gè)完整的粒子圖像測速處理系統(tǒng)關(guān)鍵部分如圖1所示，主要包括數(shù) 據(jù)采集模塊、圖像預(yù)處理模塊、運(yùn)動矢量場(光流場)計(jì)算模塊和后續(xù)數(shù)據(jù) 分析模塊，后續(xù)數(shù)據(jù)分析模塊包括旋度計(jì)算模塊和散度計(jì)算模塊。其中的 運(yùn)動矢量場計(jì)算模塊是整個(gè)PIV系統(tǒng)的核心模塊。
在計(jì)算2維運(yùn)動矢量場、散度和旋度場時(shí)，采用了基于光流運(yùn)動方程 的計(jì)算思想，以整體函數(shù)的形式計(jì)算運(yùn)動矢量。為了能夠抑止噪聲同時(shí)又 能保持矢量場、散度和旋度場的局部結(jié)構(gòu)信息，還結(jié)合了水平集的思想。
2維運(yùn)動矢量場是指在大小為NXM像素范圍內(nèi)(如通常熟知的PAL制 768 X 576或640 X 480個(gè)像素)能夠描述每一像素位置(NXM內(nèi))下該像素當(dāng)
前時(shí)刻的位移矢量，該矢量為一浮點(diǎn)數(shù)二維矢量，表示為v=|;Vx,Vy;|，
其中、，vy為x軸和y軸方向的兩個(gè)速度分量。為了統(tǒng)一描述，用大寫粗體
表示矢量，小寫細(xì)體表示標(biāo)量。速度場V實(shí)際是坐標(biāo)位置變量X氣X,y]的函數(shù)。
流體運(yùn)動圖像計(jì)算出了獲得瞬時(shí)運(yùn)動矢量場v外，還需計(jì)算v的散度場《和 旋度場n，如(1)(2)式所示。
<formula>formula see original document page 5</formula><formula>formula see original document page 6</formula>
由于實(shí)際的流體運(yùn)動可能會有多矢量流交界的現(xiàn)象，即運(yùn)動邊界和多 運(yùn)動問題，為了解決運(yùn)動邊界和多運(yùn)動問題，同時(shí)為了提高流體運(yùn)動特征 描述的精度，將水平集的策略和思想集成于原有的光流計(jì)算框架。
如圖2所示，在某一區(qū)域Q內(nèi)運(yùn)動矢量、散度和旋度有可能不連續(xù)， 為此可以將Q劃分為兩個(gè)子區(qū)域^,Q,， Q,區(qū)域內(nèi)對應(yīng)的運(yùn)動場、散度場和 旋度場為[Hn,]; Q,區(qū)域內(nèi)對應(yīng)的運(yùn)動場、散度場和旋度場為[V2,m2]。 假設(shè)在區(qū)域q內(nèi)[Hn,光滑連續(xù)，在區(qū)域Q,內(nèi)[V^,T^光滑連續(xù)，但在Q^2 的交界處運(yùn)動場、散度場和旋度場并不連續(xù)。為了保證計(jì)算時(shí)[V^,t!,和 [V",，ri,]的連續(xù)性，同時(shí)又不破壞^與Q,之間的結(jié)構(gòu)信息以及便于計(jì)算，構(gòu) 建水平集函數(shù)(()和指示函數(shù)H ， c))為定義在圖像大小范圍內(nèi)(NXM)的連續(xù)函 數(shù)。H(x)-l,x20and H(x)=0,x<0。
用水平集函數(shù)小來劃分^,Q:， po的區(qū)域?qū)儆赒,，而"0的區(qū)域?qū)儆赹。 由于采用了水平集的思想，因此在計(jì)算[Hr^和[V",、]時(shí)同時(shí)要不斷計(jì)
算更新(t)o
為了估計(jì)出[Hn,]和[v"^]，構(gòu)造以下能量函數(shù)。
E二E,+XE2 (3) 其中E為全局能量泛函，EpE2分別為以下式(4)、 (5)形式的數(shù)據(jù)約束
能量和光滑約束能量。
<formula>formula see original document page 6</formula>
這里I(X,t)代表粒子圖像序列。最終估計(jì)的光流場可以描述為<formula>formula see original document page 7</formula>通過全局能量泛函的最小化來實(shí)現(xiàn)光流場的最佳逼近。
能量泛函(3)的最小化可以通過求解線性方程組(6)實(shí)現(xiàn)，。
AU = F (6) 其中U-[VpV2,^，^，TipTi2,小]，A為9X9的系數(shù)矩陣，F(xiàn)為9X1的向
量表示線性方程的奇次項(xiàng)。
線性方程組(6)的獲得和計(jì)算可以通過經(jīng)典的變分計(jì)算方法實(shí)現(xiàn)，這 里不詳細(xì)闡述。
對于上述的PIV方法，可以采用多網(wǎng)格加速策略與之相結(jié)合，具體而 言，即通過多網(wǎng)格方法來處理線性方程組(6)，加速PIV計(jì)算。
多網(wǎng)格計(jì)算的一個(gè)中心任務(wù)是求解各種各樣來自實(shí)際問題的偏微分 方程及線性方程組。由于許多圖像分析問題需要解決帶有邊值問題的橢圓 偏微分方程，并最終通過求解線性或非線性方程組實(shí)現(xiàn)，如基于變分方 法的圖像非線性擴(kuò)散、主動輪廓計(jì)算、圖像恢復(fù)等，在大尺寸網(wǎng)格上線性 方程的高頻誤差很容易消除，相反即使經(jīng)過大量迭代低頻誤差卻很難消除， 但在小尺度網(wǎng)格上這種低頻誤差只需經(jīng)少量迭代就能明顯消除。因此多網(wǎng) 格策略就是通過不同網(wǎng)格的迭代計(jì)算來加速收斂過程。對于圖像分析而言， 網(wǎng)格的定義比較簡單，原始圖像的每一個(gè)像素表示一個(gè)網(wǎng)格，即最細(xì)的網(wǎng) 格，計(jì)算都是基于這樣的離散網(wǎng)格點(diǎn)上的值運(yùn)算完成的。
請參考圖3， 一種實(shí)施例的多網(wǎng)格處理方法具體包括如下步驟 在網(wǎng)格尺寸為h,xhy的第一網(wǎng)格h下，可以將式(6) AU-F寫為式(7):
AhUh=Fh (7) Ah，Uh分別表示在網(wǎng)格尺寸為h^hy (不同網(wǎng)格尺寸對應(yīng)不同的圖像
分辨率)下的方程系數(shù)矩陣和解。
假設(shè)方程(7)經(jīng)過n,次迭代后得到的解為tJh ，并令u-為理想的真實(shí)解， 則運(yùn)動矢量的計(jì)算誤差為
eh=Uh-!Jh (8)
由于計(jì)算結(jié)果的高頻誤差經(jīng)過少量迭代就可以基本消除，因此e-在這 里主要包含計(jì)算誤差的低頻能量。如果能知道eh則可以對tjh進(jìn)行補(bǔ)償達(dá)到 減少計(jì)算誤差的目的。雖然直接得到eh并不容易，但可以通過式(9)和式(10) 進(jìn)行間接求解。通過已知的Ch可以計(jì)算出方程(7)的誤差R、有式(9)。
Rh=Fh-AhtIh
(9)
由于Ah是線性的因此有
AV=Rh (10)這里的eh主要是在網(wǎng)格h下的低頻計(jì)算誤差。根據(jù)前面的分析，知道僅從方程(10)很難得到低頻誤差eh 。多網(wǎng)格計(jì)算的巧妙之處就在于eh不是在原先網(wǎng)格h下計(jì)算得到的，而是在更粗的網(wǎng)格上計(jì)算得到的。
對于較粗網(wǎng)格H的尺寸為H,xHy，并有H^^、,H,hy，在粗網(wǎng)格上可
以將(10)寫為
AHeH=RH (11)其中AH，RH為V，Rh經(jīng)過尺度映射后得到的新方程系數(shù)，映射關(guān)系如
下
AH 二 p"h〖Ab) ， RH = p"h〖R"
Ph—h為映射函數(shù)，eH反映的是網(wǎng)格H下的高頻計(jì)算誤差，同時(shí)又是細(xì)網(wǎng)格h下的低頻計(jì)算誤差。由于eH和eh不在一個(gè)尺度上，因此還需要將eh映射為eH，如下
在網(wǎng)格H下的高頻計(jì)算誤差eH只需通過很少的迭代就能得到，通過映
射(□)最終可以得到eh ，有了 eh就可以通過(12)來更新Ob 。
tJL=eh+tJh (12)
由于尺度之間的映射也會帶來一定的高頻誤差，因此最后還可以在網(wǎng)格h下再進(jìn)行&次迭代來消除tJL的高頻誤差。
可將多網(wǎng)格結(jié)構(gòu)與圖像計(jì)算中的多尺度思想相結(jié)合應(yīng)用于PIV計(jì)算。采用多尺度(有可能尺度因子小于2)逐級逼近的策略，可以降低方程模型誤差帶來的影響，改善計(jì)算精度。每次尺度變化都會進(jìn)行一次基于當(dāng)前估計(jì)結(jié)果的補(bǔ)償更新。將圖像尺度的變化看作是不同尺寸網(wǎng)格的變化，本質(zhì)上多尺度逼近計(jì)算可以看作是多網(wǎng)格算法的一種特例。
因此，根據(jù)實(shí)際計(jì)算需要也可以采用更多尺寸的網(wǎng)格。多網(wǎng)格算法在實(shí)際中通常分為三類V型多網(wǎng)格計(jì)算、W型多網(wǎng)格計(jì)算和完全型多網(wǎng)格計(jì)算，圖4a和圖4b給出了 V型和W型多網(wǎng)格計(jì)算的示意圖。
優(yōu)選地，采用如圖5所示的多網(wǎng)格結(jié)構(gòu)。圖5中的圓點(diǎn)代表不同尺度的網(wǎng)格(從下至上網(wǎng)格尺寸依次增大)，上升箭頭表示從小尺度網(wǎng)格向大尺
度網(wǎng)格過渡，向下箭頭則相反。每一完成一次v型多網(wǎng)格計(jì)算后，便將式
(7)中的方程系數(shù)A和F進(jìn)行更新，進(jìn)行多次循環(huán)處理。通過多網(wǎng)格計(jì)算可以有效加快線性方程組迭代的收斂速度，同時(shí)不會降低計(jì)算的可靠性和精度。
以上內(nèi)容是結(jié)合具體的優(yōu)選實(shí)施方式對本發(fā)明所作的進(jìn)一步詳細(xì)說明，不能認(rèn)定本發(fā)明的具體實(shí)施只局限于這些說明。對于本發(fā)明所屬技術(shù)領(lǐng)域的普通技術(shù)人員來說，在不脫離本發(fā)明構(gòu)思的前提下，還可以做出若干簡單推演或替換，都應(yīng)當(dāng)視為屬于本發(fā)明的保護(hù)范圍。
權(quán)利要求
1.一種用于粒子圖像測速中的多網(wǎng)格處理方法，所述粒子圖像測速包括求解線性方程組AU＝F以實(shí)現(xiàn)對運(yùn)動場的最佳逼近的步驟，U為待求解，A為系數(shù)矩陣，其至少包含運(yùn)動矢量場參數(shù)，F(xiàn)為向量，其特征在于，所述方法包括以下步驟A、設(shè)置第一網(wǎng)格，獲得方程組AU＝F對應(yīng)在第一網(wǎng)格h下的方程組AhUh＝FhB、迭代求解方程組AhUh＝Fh，根據(jù)所得實(shí)際解<overscore>U</overscore>h得到方程誤差Rh＝Fh-Ah<overscore>U</overscore>h，C、設(shè)置尺寸較第一網(wǎng)格更粗的第二網(wǎng)格H，將在第一網(wǎng)格h下對方程組Aheh＝Rh的求解尺度映射為在第二網(wǎng)格H下對方程組AHeH＝RH的求解；D、迭代求解在第二網(wǎng)格H下的高頻計(jì)算誤差eH；E、對第二網(wǎng)格H下的高頻計(jì)算誤差eH進(jìn)行尺度映射，得到在第一網(wǎng)格h下的低頻計(jì)算誤差eh；F、用第一網(wǎng)格h下的低頻計(jì)算誤差eh更新實(shí)際解<maths id="math0001" num="0001" ><math><![CDATA[ <mrow><msubsup> <mover><mi>U</mi><mo>&OverBar;</mo> </mover> <mi>new</mi> <mi>h</mi></msubsup><mo>=</mo><msup> <mi>e</mi> <mi>h</mi></msup><mo>+</mo><msup> <mover><mi>U</mi><mo>&OverBar;</mo> </mover> <mi>h</mi></msup><mo>.</mo> </mrow>]]></math> id="icf0001" file="A2009101094290002C1.tif" wi="27" he="4" top= "134" left = "145" img-content="drawing" img-format="tif" orientation="portrait" inline="yes"/></maths>
2. 如權(quán)利要求l所述的多網(wǎng)格處理方法，其特征在于，還包括 在更新實(shí)際解后，再在第一網(wǎng)格h下迭代求解VlJh-Fh，來消除C^的高頻計(jì)算誤差。
3. 如權(quán)利要求1或2所述的粒子圖像測速處理方法，其特征在于,循 環(huán)多次進(jìn)行所述處理，其中第二網(wǎng)格的尺度逐次遞增。
全文摘要
本發(fā)明公開了一種用于粒子圖像測速中的多網(wǎng)格處理方法，包括以下步驟設(shè)置第一網(wǎng)格，獲得待求方程組對應(yīng)在第一網(wǎng)格下(對應(yīng)高分辨圖像)的方程組；迭代求解方程組，根據(jù)所得實(shí)際解得到方程誤差；設(shè)置尺寸較第一網(wǎng)格更粗的第二網(wǎng)格，將在第一網(wǎng)格下對方程組的求解尺度映射為在第二網(wǎng)格下對方程組的求解；迭代求解在第二網(wǎng)格下的高頻計(jì)算誤差；對第二網(wǎng)格下的高頻計(jì)算誤差進(jìn)行尺度映射，得到在第一網(wǎng)格下的低頻計(jì)算誤差；用第一網(wǎng)格下的低頻計(jì)算誤差更新實(shí)際解。本發(fā)明的方法能大大加快粒子圖像測速的處理過程，同時(shí)不會降低處理的可靠性和精度。
文檔編號G01P5/20GK101629965SQ200910109429
公開日2010年1月20日 申請日期2009年8月18日 優(yōu)先權(quán)日2009年8月18日
發(fā)明者盧宗慶, 廖慶敏 申請人:清華大學(xué)深圳研究生院