本發(fā)明涉及一種基于光場成像的層析PIV測量裝置及方法，屬于流體測量與可視化技術(shù)領(lǐng)域。

背景技術(shù)：
：

流速是描述流場特性的重要參數(shù)之一。在科研領(lǐng)域和工業(yè)領(lǐng)域都涉及到流場流速的測量，如鍋爐燃燒過程中燃料傳輸速度的控制對于燃燒效率、安全穩(wěn)定的燃燒具有重要的意義；微流體的傳熱效率與流場流速有關(guān)等。實現(xiàn)非接觸、瞬態(tài)、三維全場流速的測量是實驗流體力學(xué)中的熱點之一。

隨著計算機技術(shù)、傳感器技術(shù)、光電技術(shù)的快速發(fā)展，PIV技術(shù)已經(jīng)廣泛應(yīng)用于流場流速的測量。該技術(shù)利用CCD在短時間內(nèi)連續(xù)拍攝兩張流場圖像，利用圖像處理技術(shù)來獲取速度矢量分布，該過程無需接觸流場，易獲得流場中多點流速。目前現(xiàn)有的PIV技術(shù)主要歸結(jié)為三類：2D-2C PIV技術(shù)、2D-3C PIV技術(shù)(準(zhǔn)三維PIV技術(shù))、3D-3C PIV技術(shù)。2D-2C PIV技術(shù)是流速測量中最基礎(chǔ)的方法之一，該方法利用約為1mm厚度的激光片光源照明流場，根據(jù)互相關(guān)算法計算兩幅圖像的二維速度。2D-3C PIV技術(shù)利用單個或多個激光片光源照明流場。首先根據(jù)二維PIV技術(shù)計算流場切面的二維速度分量，運用三維流動的連續(xù)性方程或其它技術(shù)計算出切面法線方向的第三個速度分量。該方法是基于2維PIV技術(shù)的方法，又稱為準(zhǔn)三維PIV技術(shù)，但該方法在法線方向上的速度分辨率較低，難以實現(xiàn)體的測量。不同于以上兩種方法，3D-3C PIV是利用具有一定厚度的激光體光源照明某個容積內(nèi)的流場，利用數(shù)學(xué)方法重構(gòu)出示蹤粒子的三維位置坐標(biāo)和強度信息(三維粒子場)，利用三維互相關(guān)技術(shù)計算速度矢量分布。該方法能同時獲得流場的三個速度分量。3D-3C PIV技術(shù)能獲得真正意義上的三維速度矢量分布。層析PIV是3D-3C PIV中的一種，由于該方法具有較高的空間分辨率，適用范圍較廣，已受廣大學(xué)者青睞。但層析PIV技術(shù)往往需要3～6臺相機拍攝流場，大大增加了層析PIV技術(shù)的成本、復(fù)雜程度及重建計算量。每臺CCD需要嚴(yán)格對準(zhǔn)流場，目前的層析PIV技術(shù)的精度非常依賴于CCD的布局及標(biāo)定的準(zhǔn)確性。

技術(shù)實現(xiàn)要素：

本發(fā)明的目的是針對基于傳統(tǒng)相機的層析PIV技術(shù)存在的缺點，提供一種基于光場成像的層析PIV測量裝置及方法，利用光場成像能同時記錄光場的位置和方向信息，使得單臺光場相機代替多臺傳統(tǒng)相機用于層析PIV測量成為可能，且能克服傳統(tǒng)層析PIV的缺點。

上述的目的通過以下技術(shù)方案實現(xiàn)：

一種基于光場成像的層析PIV測量裝置，包括用于產(chǎn)生雙脈沖激光并照亮三維流場的光源設(shè)備和用于拍攝三維流場的成像設(shè)備，所述的光源設(shè)備的發(fā)光端設(shè)置有透鏡組，所述的光源設(shè)備和所述的成像設(shè)備連接時序同步控制器，所述的成像設(shè)備連接計算機；

所述的透鏡組用于將所述的光源設(shè)備產(chǎn)生的雙脈沖激光發(fā)展為一定厚度的體光源和1mm的二維片光源，并照亮流場；

所述的成像設(shè)備，包括主透鏡、CCD探測器以及位于所述的主透鏡與CCD探測器之間的微透鏡陣列。

所述的基于光場成像的層析PIV測量裝置，所述的微透鏡陣列數(shù)為N×N，CCD探測器像面關(guān)于微透鏡的共軛面為主透鏡平面，光線經(jīng)過微透鏡在CCD探測器上所成的像稱為子圖像，每個子圖像包含M×M個像素點，其中N和M都大于2。

用上述基于光場成像的層析PIV測量裝置進行測量的方法，該方法包括如下步驟：

步驟一、坐標(biāo)定義及標(biāo)定板拍攝：流場及標(biāo)定板所在的坐標(biāo)定義為世界坐標(biāo)O-XYZ，光場原圖像中子圖像坐標(biāo)v-o_m-u為光場原圖像中的子圖像所在矩陣的行數(shù)和列數(shù)(v,u)，子圖像中像素坐標(biāo)r-o_p-s為子圖像中像素所在矩陣的行數(shù)和列數(shù)(r,s)，視差圖像中的像素坐標(biāo)r_v-o_v-s_v為視差圖像中像素所在矩陣的行數(shù)和列數(shù)(r_v,s_v)，視差圖像的圖像坐標(biāo)x-o-y為像素在視差圖像中的物理位置(x,y)，原點o位于(r_v,s_v)＝(1,1)的像素中心位置。實驗時，將標(biāo)定板放置于待測流場所在位置，固定光場相機，使標(biāo)定板平面平行于光場相機的探測器面，沿著待測流場深度方向等間隔地移動標(biāo)定板，用光場相機拍攝位于至少三個不同深度位置的標(biāo)定板，獲得光場原圖像；

步驟二、視差圖像的提?。簩⒉襟E一拍攝的標(biāo)定板光場原圖像中的像素進行重新組合，使每個子圖像中位于相同位置的像素組成一幅新圖像，可以獲得M×M幅不同的視差圖像，剔除由子圖像邊緣像素組成的視差圖像，從而得到M′×M′幅視差圖像。記錄標(biāo)定板角點坐標(biāo)(X,Y,Z)與對應(yīng)的每幅視差圖像中的像素坐標(biāo)(r_v,s_v)，利用下式計算角點的圖像坐標(biāo)(x,y)：

式中，p為微透鏡的尺寸。

步驟三、三維流場空間坐標(biāo)與視差圖像中的圖像坐標(biāo)關(guān)系的計算：對每幅視差圖像，利用三階多項式建立由步驟二得到的標(biāo)定板角點坐標(biāo)(X,Y,Z)與圖像坐標(biāo)(x,y)之間的關(guān)系：

式中，已知(X,Y,Z)和(x,y)利用MATLAB擬合工具箱計算出a₁,a₂,…,a₂₀和b₁,b₂,…,b₂₀。對于每一幅視差圖像需要拍攝至少三幅不同深度位置的標(biāo)定板圖像，每一幅視差圖像在不同深度位置計算的a₁,a₂,…,a₂₀和b₁,b₂,…,b₂₀的值不同，然后利用線性插值計算出任意深度位置的a₁,a₂,…,a₂₀和b₁,b₂,…,b₂₀的值，從而確定三維流場空間坐標(biāo)與每幅視差圖像中的圖像坐標(biāo)關(guān)系；

步驟四、權(quán)重矩陣的計算：將流場沿X軸Y軸Z軸劃分成m,n,l個體素，將每個體素的中心坐標(biāo)(X_c,Y_c,Z_c)分別代入步驟三中的式(3)、式(4)得到相應(yīng)的圖像坐標(biāo)，再將圖像坐標(biāo)代入步驟二中的式(2)得到每個體素對應(yīng)的像素坐標(biāo)，利用高斯函數(shù)計算每個體素對每幅視差圖像中的像素的貢獻W(權(quán)重矩陣)，保存權(quán)重矩陣；

步驟五、流場的光場原圖像的獲取與預(yù)處理：將標(biāo)定板移出待測流場，打開雙脈沖激光器，使激光經(jīng)過透鏡組產(chǎn)生具有一定厚度的體光源，當(dāng)激光入射到含有示蹤粒子的待測流動對象時，示蹤粒子會產(chǎn)生散射光，由光場相機連續(xù)兩次拍攝三維流場，獲得光場原圖像，根據(jù)步驟二的方法獲得一系列三維流場的視差圖像的灰度值I_v；

步驟六、層析重建計算三維粒子場：根據(jù)步驟四計算得到的權(quán)重矩陣W及步驟五得到的視差圖像的灰度值I_v，利用MART算法反演三維粒子灰度分布E：

式中，c為體素的序號，c＝1,2,3,…,m×n×l；X_c,Y_c,Z_c為體素中心坐標(biāo)；k為迭代次數(shù)；g為像素的序號，g＝1,2,3,…,N×N×M′×M′；μ為松弛因子。

步驟七、三維互相關(guān)算法計算流速：對步驟六重建所得的三維粒子灰度分布E沿X軸Y軸Z軸分別以x_size,y_size,z_size的尺寸進行網(wǎng)格劃分(查詢窗口)，然后利用下式計算三維流場速度矢量分布：

式中，R為互相關(guān)系數(shù)；x_size,y_size,z_size表示查詢窗口分別在X,Y,Z方向上的尺寸；E_t，E_t+Δt分別為t時刻和t+Δt時刻的三維粒子灰度分布；σ_t,σ_Δ+t分別為t時刻和t+Δt時刻查詢窗口內(nèi)的體素灰度值的標(biāo)準(zhǔn)差；m′,n′,l′為查詢窗口內(nèi)的體素分別在X,Y,Z方向上的坐標(biāo)索引(體素的序號)。式(8)的最大值對應(yīng)的Δm′,Δn′,Δl′為粒子分別在X,Y,Z方向上的位移，已知Δt可以獲得流場的三維速度矢量分布。

所述的基于光場成像的層析PIV測量裝置進行測量的方法，步驟二中所述M×M幅不同視差圖像的提取方法如公式(1)確定：

式中，I(N,N,r,s)為第N行第N列微透鏡子圖像中的第r行第s列的像素的灰度值；r為每個子圖像中的像素所在的行數(shù)，r＝1,2,3,…,M；s為每個子圖像中的像素所在的列數(shù)，s＝1,2,3,…,M。視差圖像里相鄰兩個像素之間的距離為單個微透鏡的尺寸p。

所述的基于光場成像的層析PIV測量裝置進行測量的方法，步驟四中所述利用高斯函數(shù)來計算每個體素對每幅視差圖像中的像素的貢獻(權(quán)重矩陣)，其公式為：

式中，“T”為對矩陣轉(zhuǎn)置；h為每個粒子圖像灰度峰值；d,e,f表征粒子圖像的形狀；W_v(r_v,s_v)為視差圖像中第r_v行第s_v列像素的灰度值，r_v＝1,2,3,…,N，s_v＝1,2,3,…,N；r_v2，s_v2為每個粒子圖像灰度峰值所在的像素坐標(biāo)。

有益效果：

1.本發(fā)明采用單光場相機拍攝流場中示蹤粒子發(fā)出的散射光。相比于傳統(tǒng)相機，光場相機在主透鏡和CCD探測器之間加了一微透鏡陣列。流場本身不會發(fā)光，所以需要往流場中加入示蹤粒子，示蹤粒子被入射光照射時會產(chǎn)生散射光。散射光經(jīng)過鏡頭投影到微透鏡陣列上，再經(jīng)過微透鏡陣列二次成像投影到CCD探測器上，形成一系列子圖像，每個子圖像記錄了光場的位置信息，每個子圖像中位于相同位置的像素記錄了光場的某一視角信息，本發(fā)明成像裝置能同時實現(xiàn)光場方向信息、位置信息及散射光強度的采集。

2.相比于傳統(tǒng)層析PIV裝置，由于光場相機能同時記錄光場的位置和方向信息，采用單光場相機代替多臺傳統(tǒng)相機拍攝流場，從而實現(xiàn)單相機多方向的采集，系統(tǒng)更緊湊，標(biāo)定更方便，能減小多相機同步耦合時帶來的系統(tǒng)誤差。

3.在測量方法上，光場相機中微透鏡陣列數(shù)為N×N，每個子圖像包含M×M個像素點，N和M都大于2。光場相機可以被看成以微透鏡陣列為鏡頭的許多微型相機對物體在不同視角成像的裝置。由于相同方向的光線被微透鏡下相同位置的像素接收，本發(fā)明將每個子圖像中相同位置的像素組成一幅新圖像，從而將光場圖像轉(zhuǎn)換為一系列不同的視差圖像，不同的視差圖像代表對示蹤粒子不同方向的成像，能提取出M×M幅不同的視差圖像。通過標(biāo)定進一步建立三維流場空間坐標(biāo)與每幅視差圖像中像素的圖像坐標(biāo)之間的映射關(guān)系求解權(quán)重矩陣，建立層析重建模型，利用MART算法反演三維粒子灰度分布，利用三維互相關(guān)實現(xiàn)3D-3C速度場測量。

4.對于光場相機的權(quán)重矩陣計算，本發(fā)明對標(biāo)定板圖像中每個子圖像相同位置的像素進行重新組合提取出一系列視差圖像，利用三階多項式建立三維流場空間與這些視差圖像中的像素之間的對應(yīng)關(guān)系，然后利用高斯函數(shù)來計算三維體素對各個視差圖像中的像素的貢獻(權(quán)重矩陣)，該方法計算速度快。

5.本發(fā)明通過建立標(biāo)定板角點坐標(biāo)與視差圖像中像素的圖像坐標(biāo)之間的映射關(guān)系來確定三維流場的每個體素坐標(biāo)與像素之間的關(guān)系，該方法不需要標(biāo)定光場相機的內(nèi)外參數(shù)。

附圖說明

圖1為本發(fā)明所述的基于光場成像的層析PIV測量裝置結(jié)構(gòu)示意圖；

圖中的標(biāo)記含義：1、雙脈沖激光器，2、透鏡組；3、時序同步控制器；4、光場相機；5、計算機；6、待測流場；7、示蹤粒子；8、傳輸線。

圖2為光場相機結(jié)構(gòu)及成像原理圖；圖中的標(biāo)記含義：9、虛擬物面；10、光線；11、主透鏡；12、微透鏡陣列；13、CCD探測器；14、像素；15、子圖像。

圖3(a)為光場原圖像坐標(biāo)定義示意圖；圖中的標(biāo)記含義：16、光場原圖像；橫坐標(biāo)o_mu軸為子圖像所在矩陣的列數(shù)，縱坐標(biāo)o_mv軸為子圖像所在矩陣的行數(shù)；圖3(b)為光場原圖像中子圖像像素坐標(biāo)定義示意圖；圖中的標(biāo)記含義：17、子圖像邊緣像素；橫坐標(biāo)o_ps軸為子圖像中的像素所在矩陣的列數(shù)，縱坐標(biāo)o_pr軸為子圖像中的像素所在矩陣的行數(shù)；圖3(c)為視差圖像坐標(biāo)定義示意圖；圖中的標(biāo)記含義：18、視差圖像，兩條橫坐標(biāo)o_vs_v軸,oy軸分別為視差圖像中像素所在矩陣的列數(shù)及像素在橫向的物理位置，兩條縱坐標(biāo)o_vr_v軸,ox軸分別為視差圖像中像素所在矩陣的行數(shù)及像素在縱向的物理位置。

圖4(a)為三維體標(biāo)定示意圖；圖中的標(biāo)記含義：19、標(biāo)定板；圖4(b)為標(biāo)定板的示意圖；圖中的標(biāo)記含義：20、標(biāo)定板角點。

圖5為體素和查詢窗口劃分示意圖；圖中的標(biāo)記含義：21、體素；22查詢窗口。

具體實施方式

下面結(jié)合具體實施方式，進一步闡明本發(fā)明，應(yīng)理解下述具體實施方式僅用于說明本發(fā)明而不用于限制本發(fā)明的范圍。

本發(fā)明的基于光場成像的層析PIV裝置工作原理為：由于流場本身不會發(fā)光，所以在待測流場中加入跟隨性、散射性較好的示蹤粒子，將具有一定厚度的激光入射到流場中，使示蹤粒子產(chǎn)生散射光。將光場相機設(shè)置在雙幀工作模式下，利用時序同步控制器調(diào)節(jié)激光器的雙脈沖信號的時間間隔，使兩個激光脈沖分別落在光場相機的雙曝光時間上，從而使光場相機連續(xù)兩次拍攝到三維流場圖像。然后利用計算機對采集到的圖像進行視差圖像提取、層析重建、三維互相關(guān)計算，從而得到3D-3C速度場分布。

本發(fā)明裝置能夠?qū)崿F(xiàn)單相機的3D-3C速度場分布測量，且無需標(biāo)定光場相機內(nèi)外參數(shù)。測量裝置如圖1所示，主要包括雙脈沖激光器1、透鏡組2、時序同步控制器3、光場相機4、計算機5、待測流場6、示蹤粒子7、傳輸線8。雙脈沖激光器1是整個測量裝置的光源(波長為532nm，光束直徑D～7mm)。透鏡組2的作用是將雙脈沖激光發(fā)展為具有一定體積(如100mm×100mm×100mm)的激光光束或厚度為1mm的二維片激光光束(如1mm×100mm×100mm)。時序同步控制器3產(chǎn)生多個信號通過傳輸線8分別傳輸給雙脈沖激光器1和光場相機4。光場相機4與雙脈沖激光器1成正交放置。計算機5的作用是圖像存儲、視差圖像提取、層析重建、三維互相關(guān)計算。待測流場6放置于脈沖激光器1和光場相機4正交的交點處，使激光能垂直入射到待測流場6中。示蹤粒子7被均勻地撒入待測流場中。

不同于傳統(tǒng)的層析PIV裝置，本發(fā)明采用單光場相機4代替多臺傳統(tǒng)相機來拍攝流場。光場相機結(jié)構(gòu)如圖2所示，光場相機4在主透鏡11和CCD探測器13之間加了一微透鏡陣列12，使經(jīng)過鏡頭的光線10投影到微透鏡陣列12上，光線經(jīng)過微透鏡陣列12二次投影到CCD探測器13上，產(chǎn)生一幅光場原圖像。光場原圖像里包含一系列子圖像15，每個子圖像15記錄了光場的位置信息，子圖像15內(nèi)每個像素點14記錄了光場的方向信息，能同時實現(xiàn)光場方向和位置信息的采集。這一特性使單光場相機4重建流場中的示蹤粒子7的三維位置和強度信息成為可能。

基于光場成像的層析PIV的測量方法，其步驟包括如下：

步驟一、坐標(biāo)定義及標(biāo)定板拍攝：流場6及標(biāo)定板19所在的坐標(biāo)定義為世界坐標(biāo)O-XYZ。光場原圖像16中子圖像坐標(biāo)v-o_m-u為光場原圖像16中的子圖像15所在矩陣的行數(shù)和列數(shù)(v,u)，子圖像15中像素坐標(biāo)r-o_p-s為子圖像15中像素14所在矩陣的行數(shù)和列數(shù)(r,s)，視差圖像18中的像素坐標(biāo)r_v-o_v-s_v為視差圖像18中像素14所在矩陣的行數(shù)和列數(shù)(r_v,s_v)，視差圖像18的圖像坐標(biāo)x-o-y為像素14在視差圖像18中的物理位置(x,y)，原點o位于(r_v,s_v)＝(1,1)的像素14中心位置，如圖3所示。實驗時，將標(biāo)定板19放置于待測流場6所在位置，固定光場相機4，使標(biāo)定板19平面平行于光場相機4的探測器面13，沿著待測流場深度方向等間隔地移動標(biāo)定板19，用光場相機4拍攝位于至少三個不同深度位置的標(biāo)定板，獲得光場原圖像16，該過程如圖4所示；

步驟二、視差圖像的提取：將步驟一拍攝的標(biāo)定板光場原圖像16中的像素14按式(1)進行重新組合：

式中，I(N,N,r,s)為第N行第N列微透鏡子圖像中的第r行第s列的像素的灰度值；r為每個子圖像中的像素所在的行數(shù)，r＝1,2,3,…,M；s為每個子圖像中的像素所在的列數(shù)，s＝1,2,3,…,M。視差圖像里相鄰兩個像素之間的距離為單個微透鏡的尺寸p。

式(1)使每個子圖像15中位于相同位置的像素組成一幅新圖像，可以獲得M×M幅不同的視差圖像18V，剔除由子圖像邊緣像素17組成的視差圖像，從而獲得M′×M′幅視差圖像。視差圖像里的像素坐標(biāo)(r_v,s_v)與圖像坐標(biāo)(x,y)的關(guān)系按公式(2)計算：

記錄標(biāo)定板角點20坐標(biāo)(X,Y,Z)與對應(yīng)的每幅視差圖像18中的圖像坐標(biāo)(x,y)。

步驟三、三維流場空間坐標(biāo)與視差圖像中的圖像坐標(biāo)關(guān)系的計算：對每幅視差圖像18，利用三階多項式建立由步驟二得到的標(biāo)定板角點20坐標(biāo)(X,Y,Z)與圖像坐標(biāo)(x,y)之間的關(guān)系：

式中，已知(X,Y,Z)和(x,y)利用MATLAB擬合工具箱計算出a₁,a₂,…,a₂₀和b₁,b₂,…,b₂₀。每幅視差圖像18在不同深度位置計算的a₁,a₂,…,a₂₀和b₁,b₂,…,b₂₀的值不同，已知至少三個不同深度位置的a₁,a₂,…,a₂₀和b₁,b₂,…,b₂₀的值，利用線性插值計算出任意深度位置的a₁,a₂,…,a₂₀和b₁,b₂,…,b₂₀的值。從而確定三維流場空間坐標(biāo)與每幅視差圖像18中的圖像坐標(biāo)的函數(shù)關(guān)系。

步驟四、權(quán)重矩陣的計算：

(a)高斯函數(shù)系數(shù)的標(biāo)定：移除標(biāo)定板19，調(diào)整透鏡組2使激光厚度為1mm(稱為片光源)，將片光源入射到含有示蹤粒子7的待測流場區(qū)域，且片光源平面平行于光場相機4的探測器面13。當(dāng)激光入射到含有示蹤粒子7的待測區(qū)域時，示蹤粒子會產(chǎn)生散射光，由光場相機4拍攝至少三幅不同深度位置的示蹤粒子的光場原圖像16，且使每個平面都有約為50個互不重疊的粒子被光場相機記錄，下面以虛擬物面這個深度圖像為例。根據(jù)步驟二中的式(1)獲得位于虛擬物面處粒子的M′×M′幅視差圖像18。對每一幅視差圖像18進行如下處理：識別這50個粒子圖像的灰度峰值p及峰值所在的像素坐標(biāo)14(r_v2,s_v2)，即利用高斯函數(shù)分別對50個粒子進行擬合，高斯函數(shù)公式為：

式中，“T”為對矩陣轉(zhuǎn)置；h為每個粒子圖像灰度峰值；d,e,f表征粒子圖像的形狀；W_v(r_v,s_v)為視差圖像中第r_v行第s_v列像素的灰度值，r_v＝1,2,3,…,N，s_v＝1,2,3,…,N；r_v2，s_v2為每個粒子圖像灰度峰值所在的像素坐標(biāo)。

已知h,r_v,s_v,r_v2,s_v2,W_v(r_v,s_v)的值，可以擬合出50個粒子對應(yīng)的系數(shù)d,e,f，然后求平均值可以確定公式(5)中的系數(shù)d,e,f。其它深度位置的系數(shù)d,e,f按該方法計算，不同深度位置計算所得的d,e,f的值不同，可以利用線性插值計算任意深度位置的d,e,f的值。從而確定了位于任意深度位置每幅視差圖像18對應(yīng)的高斯函數(shù)系數(shù)。

(b)數(shù)據(jù)處理：將流場沿X軸Y軸Z軸劃分成m,n,l個體素21(如m＝n＝1001，l＝101，每個體素的邊長為0.1mm)，將每個體素21的中心坐標(biāo)(X_c,Y_c,Z_c)(X_c,Y_c,Z_c分別表示第c個體素的中心在世界坐標(biāo)系上的坐標(biāo)值，c＝1,2,3,…,m×n×l)的值代入步驟三的式(3)、式(4)中的變量X,Y,Z計算每個體素21的中心坐標(biāo)對應(yīng)的圖像坐標(biāo)，然后利用式(2)將圖像坐標(biāo)轉(zhuǎn)化為像素坐標(biāo)，從而確定了每個體素21對應(yīng)的視差圖像18中的像素坐標(biāo)。將得到的像素坐標(biāo)代入到公式(5)中的變量r_v2,s_v2，從而得到每個體素21對應(yīng)每幅視差圖像18中的所有像素14的灰度值W_v，并轉(zhuǎn)換成一系列一維列向量V′(r,s,c)＝[W_v(1,1),W_v(1,2),W_v(1,3),...,W_v(N,N)]^T，保存形式如下：

式中，W為空間體素21對各幅視差圖像18中的所有像素14的貢獻；V′(M′,M′,m×n×l)為第m×n×l個體素21對由每個子圖像中像素坐標(biāo)為(M′,M′)的像素所組成的視差圖像18的貢獻。

步驟五、流場的光場原圖像的獲取與預(yù)處理：調(diào)整透鏡組2，使激光經(jīng)過透鏡組2產(chǎn)生具有一定厚度的體光源入射到待測流場中。由光場相機4連續(xù)兩次拍攝三維流場，獲得光場原圖像16，根據(jù)步驟二中的式(1)獲得M′×M′幅三維流場的視差圖像18V(r,s)。

步驟六、層析重建計算三維粒子場：

(a)數(shù)據(jù)處理：將步驟五得到的三維流場的每幅視差圖像18的灰度值V(r,s)轉(zhuǎn)化為一維列向量的形式V(r,s)＝[I(1,1,r,s),I(1,2,r,s),I(1,3,r,s),…,I(N,N,r,s)]^T，然后將M′×M′幅不同視差圖像灰度值列向量組成一個新的列向量I_v(x_g,y_g)＝[V(1,1),V(1,2),V(1,3),…,V(M′,M′)]^T。

(b)層析重建：初始化三維粒子灰度分布E¹＝[1,2,3,…,m×n×l]′；已知權(quán)重矩陣W和視差圖像18灰度值I_v的值，利用MART算法反演三維粒子灰度分布E：

式中，c為體素的序號，c＝1,2,3,…,m×n×l；X_c,Y_c,Z_c為體素中心坐標(biāo)；k為迭代次數(shù)；g為像素的序號，g＝1,2,3,…,N×N×M′×M′；μ為松弛因子。

步驟七、三維互相關(guān)算法計算流速：對步驟六重建所得的三維粒子灰度分布E沿X軸Y軸Z軸分別以x_size,y_size,z_size(如x_size＝y(tǒng)_size＝z_size＝16個體素)的尺寸進行網(wǎng)格劃分(查詢窗口)22。然后利用下式計算三維流場速度矢量分布：

式中，R為互相關(guān)系數(shù)；x_size,y_size,z_size表示查詢窗口分別在X,Y,Z方向上的尺寸；E_t，E_t+Δt分別為t時刻和t+Δt時刻的三維粒子灰度分布；σ_t,σ_Δ+t分別為t時刻和t+Δt時刻查詢窗口內(nèi)的體素灰度值的標(biāo)準(zhǔn)差；m′,n′,l′為查詢窗口內(nèi)的體素分別在X,Y,Z方向上的坐標(biāo)索引(體素的序號)。式(8)的最大值對應(yīng)的Δm′,Δn′,Δl′為粒子分別在X,Y,Z方向上的位移，已知Δt可以獲得流場的三維速度矢量分布。