粒子拖尾長度對高頻DPIV測量誤差的影響

陳啟剛，李丹勛，鐘 強(qiáng)，王興奎

（清華大學(xué)水沙科學(xué)與水利水電工程國家重點實驗室，北京 100084）

粒子拖尾長度對高頻DPIV測量誤差的影響

陳啟剛，李丹勛，鐘 強(qiáng)，王興奎

（清華大學(xué)水沙科學(xué)與水利水電工程國家重點實驗室，北京 100084）

基于連續(xù)激光源的高采樣頻率DPIV系統(tǒng)在進(jìn)行曝光時可能會產(chǎn)生粒子拖尾現(xiàn)象。明渠湍流試驗和數(shù)值模擬試驗結(jié)果顯示，粒子拖尾長度的增加會導(dǎo)致平均流速減小和湍流強(qiáng)度增大。分析表明，導(dǎo)致平均流速誤差的原因之一是診斷窗口邊緣局部粒子圖像的缺失，且圖像缺失程度會隨著粒子拖尾長度增加而增大，從而使平均速度隨粒子拖尾長度增加而減??；粒子拖尾的另一個直接后果是增大了粒子圖像的有效粒徑，從而導(dǎo)致湍流強(qiáng)度增大。通過增大第二個診斷窗口的尺寸，進(jìn)行相關(guān)函數(shù)標(biāo)準(zhǔn)化和相關(guān)系數(shù)補(bǔ)償?shù)确椒?，可以消除粒子拖尾長度引起的誤差。

PIV；拖尾長度；連續(xù)激光源；誤差；平均速度；湍流強(qiáng)度

0 引 言

PIV是一種無干擾全場測速技術(shù)，可以快速獲取瞬時二維或三維流速場。自20世紀(jì)80年代以來，隨著對示蹤粒子、光學(xué)設(shè)備、成像系統(tǒng)、速度提取算法等研究的不斷深入，傳統(tǒng)的PIV技術(shù)日趨成熟。

目前，標(biāo)準(zhǔn)的商用DPIV系統(tǒng)是以兩臺交錯工作的脈沖激光為光源，數(shù)字?jǐn)z像機(jī)按雙幀單曝光模式拍攝示蹤粒子圖片，采用移動、變形或多重迭代的診斷窗口，基于快速傅里葉變換（FFT）對兩個診斷窗口進(jìn)行互相關(guān)計算得到速度場，通過相關(guān)系數(shù)峰值亞像素擬合可以提高測量精度［1－3］。文獻(xiàn)中對基于脈沖光源PIV 系統(tǒng)的誤差分析全面而深入［4－11］。

在標(biāo)準(zhǔn)的商用DPIV系統(tǒng)中，粒子曝光時間等于脈沖持續(xù)時間，可以設(shè)置到ns量級，但受限于激光的脈沖頻率，其采樣頻率一般低于50Hz［12］。為了滿足對高頻湍流研究的需要，基于高頻率脈沖激光的TRPIV系統(tǒng)近年來得到了較快發(fā)展，但其價格昂貴且對同步設(shè)備要求極高，尚無法在明渠湍流研究領(lǐng)域推廣普及。提高DPIV采樣頻率的另一個重要途徑是采用連續(xù)激光代替脈沖激光［12－14］。在基于連續(xù)激光的DPIV系統(tǒng)中，采樣頻率由高速攝像機(jī)控制，但由于連續(xù)激光的強(qiáng)度低于脈沖激光，相機(jī)曝光時間（相機(jī)快門從打開到關(guān)閉的時間間隔）一般要長于脈沖激光源的脈沖持續(xù)時間。Sante等［15］在將基于連續(xù)激光的TR－PIV應(yīng)用于旋轉(zhuǎn)槽道流的測量時，其相機(jī)曝光時間為200～300μs。

隨著曝光時間的增加，示蹤粒子在成像過程中會出現(xiàn)拖尾現(xiàn)象。這一現(xiàn)象早就引起了研究者的注意，例如Prasad［2］和 Adrian［3］在有關(guān)PIV 的綜述中對此就有所提及，最近Sante等［13］在將連續(xù)激光器應(yīng)用于實際流動測量時，也注意到粒子會隨著曝光時間的增加而產(chǎn)生拖尾現(xiàn)象。但是，目前還未見有系統(tǒng)分析粒子拖尾長度對DPIV測量誤差影響的研究成果。

為深入研究粒子拖尾現(xiàn)象對DPIV測量精度造成的影響，首先應(yīng)用基于連續(xù)激光源的高頻DPIV對同一明渠湍流時均速度和湍流強(qiáng)度進(jìn)行測量，對比不同曝光時間（拖尾長度）下測量結(jié)果的差異；然后，采用蒙特－卡羅隨機(jī)模擬方法生成兩幀粒子圖像，分析拖尾長度對平均速度的影響；最后對粒子拖尾導(dǎo)致DPIV測量誤差的原因進(jìn)行分析，并提出消除這種誤差的方法。

1 明渠湍流試驗

1.1 試驗條件

構(gòu)建了基于連續(xù)激光源的高頻DPIV系統(tǒng)，在不同曝光時間條件下，對同一明渠均勻湍流進(jìn)行測量，分析拖尾長度對時均流速和湍流強(qiáng)度測量精度的影響。

試驗水槽長11m、寬0.25m、高0.20m，可調(diào)坡降范圍為0～0.01，觀測段距上游進(jìn)水口8m。PIV硬件系統(tǒng)由2W連續(xù)激光源、工業(yè)用高速攝像機(jī)和光學(xué)透鏡構(gòu)成，相機(jī)CCD大小為640pixel×480pixel，采樣頻率604Hz。流場計算采用基于FFT的互相關(guān)算法；采用二重網(wǎng)格迭代模式進(jìn)行計算，最大診斷窗口尺寸為64pixel×16pixel，最小診斷窗口尺寸為16pixel×4pixel；采用式（1）進(jìn)行三點高斯亞像素擬合以確定相關(guān)系數(shù)峰值的位置。

在相同的水流條件下進(jìn)行了4組試驗，相機(jī)分辨率Rs＝54.2μm／pixel，水深2cm，坡降0.4%，垂線平均流速49.7cm／s，摩阻流速u＊＝2.6cm／s，雷諾數(shù)Re＝8590。4組試驗的曝光時間Te分別為300、600、1000、1500μs。每組試驗均連續(xù)采集40000幀圖片進(jìn)行統(tǒng)計計算。其中，在曝光時間為300μs時，基本上無拖尾現(xiàn)象，可以作為對比的標(biāo)準(zhǔn)。

1.2 試驗結(jié)果

定義拖尾長度L、平均流速Um、平均流速誤差Rb、湍流強(qiáng)度等參數(shù)如下：

其中，Ui為單次采樣得到的瞬時速度值，Ur為相同水流條件下300μs曝光時間（基本沒有拖尾影響）的平均流速。

圖1點繪了由（4）式計算的平均流速的誤差。從圖中可以看出，隨著曝光時間的增加，平均流速減小，即平均流速隨粒子拖尾長度的增加而減??；同時，近壁區(qū)速度梯度的存在進(jìn)一步加劇了平均流速減小的幅度。

圖1 平均流速誤差沿水深的變化Fig.1 Mean－velocity error along water depth

圖4點繪了縱向湍流強(qiáng)度隨曝光時間的變化，從圖中可以看出，隨著曝光時間的增加，縱向湍流強(qiáng)度有增大的趨勢，這反映了粒子拖尾長度的影響。

圖2 平均流速誤差隨速度梯度的變化Fig.2 Mean－velocity error with respect to velocity gradient

圖3 平均流速誤差隨拖尾長度的變化Fig.3 Mean－velocity error with respect to tailing length

圖4 縱向湍流強(qiáng)度隨曝光時間的變化Fig.4 Streamwise turbulent intensity with respect to exposure time

2 數(shù)值模擬試驗

2.1 圖像生成

模擬均勻流場的情況。采用蒙特－卡羅模擬方法在400pixel×400pixel的圖片上隨機(jī)生成了一系列示蹤粒子作為第一幀粒子圖像（圖5（a）），其中每個示蹤粒子的大小為1個像素點；將第一幀圖像中的每個像素點從左往右平行移動相同的距離，即為第二幀粒子圖像，即每個診斷窗口均具有相同的縱向流速。為了模擬粒子的拖尾現(xiàn)象，只需在圖5（a）中每個粒子的右邊補(bǔ)充相同數(shù)量的粒子即可，圖5（b）給出了圖5（a）的粒子拖尾2個像素點后對應(yīng)的圖像。

2.2 數(shù)值模擬結(jié)果

共分析了拖尾長度為0～9的7種情況。模擬得到各種拖尾長度的粒子圖像后，采用64pixel×64pixel的診斷窗口對粒子圖像進(jìn)行相關(guān)計算，將每個流場中得到的流速矢量求平均后即得到各拖尾長度條件下的平均流速。選擇粒子無拖尾時的平均流速為參考流速Ur，用（4）式分別計算各拖尾長度條件下的平均流速與Ur的誤差，計算結(jié)果如圖6。

圖5 模擬粒子圖像Fig.5 Simulated particle images

圖6 平均流速誤差的數(shù)值模擬結(jié)果Fig.6 Numerical simulation result of mean－velocity error

圖6的結(jié)果表明，平均流速測量誤差隨著拖尾長度的增加而逐漸增大，這與圖3中由試驗得出的結(jié)果是吻合的。雖然圖3和圖6中二者的變化趨勢是相同的，但在相同拖尾長度時，平均流速誤差的數(shù)值有差別，這可能是由于湍流脈動及試驗顆粒的三維運動特性等原因造成的，而基于簡化的二維模擬圖片進(jìn)行PIV誤差分析時一般會低估實際誤差［9］。

3 誤差分析

3.1 平均流速誤差的原因

當(dāng)粒子不發(fā)生拖尾時，經(jīng)過8像素位移以后，第一個診斷窗口中離右邊框8像素范圍內(nèi)的粒子將跑出第二個診斷窗口（如圖7）。當(dāng)粒子發(fā)生拖尾時，以拖尾2個像素為例，相當(dāng)于粒子水平方向的粒徑為原來的3倍，此時除了離右邊框8像素范圍內(nèi)的粒子將跑出第二個診斷窗口以外，離右邊壁8～10像素距離范圍內(nèi)的粒子圖像將部分跑出第二個診斷窗口，導(dǎo)致粒子圖像局部缺失（如圖8），這是導(dǎo)致平均流速產(chǎn)生誤差的原因之一。

為了驗證以上分析，采用2.1節(jié)中已經(jīng)生成的粒子拖尾2像素的圖像進(jìn)行計算，在計算之前，首先將第一幀診斷窗口中距離右邊框10像素范圍的粒子全部剔除，然后再與第二幀圖片中相同位置的診斷窗口進(jìn)行基于FFT的互相關(guān)運算。經(jīng)過以上處理以后，可以保證粒子在拖尾長度不超過2像素時，第一幀診斷窗口中距離右邊框大于8像素范圍內(nèi)的所有粒子圖像均能完整的出現(xiàn)在位移8像素后的第二幀診斷窗口中，不致產(chǎn)生粒子圖像的局部缺失。圖9點繪了經(jīng)過以上處理后得到的平均流速誤差隨粒子拖尾長度的變化趨勢。

對比圖6與圖9可以發(fā)現(xiàn)，剔除距離右邊框10像素范圍內(nèi)的示蹤粒子后，當(dāng)拖尾長度不超過2像素時，粒子拖尾不會導(dǎo)致平均流速誤差。當(dāng)拖尾長度超過2像素后，平均流速之所以再次產(chǎn)生誤差，是因為距離右邊框大于10像素的部分粒子又會產(chǎn)生圖像的局部缺失。實際上，即使粒子不發(fā)生拖尾現(xiàn)象，當(dāng)粒子圖像的粒徑大于1個像素時，也會因為診斷窗口右邊緣的粒子圖像局部缺失而導(dǎo)致平均流速誤差。

如圖10所示，當(dāng)粒子圖像均保持完整時，在將第一幀診斷窗口相對第二幀診斷窗口進(jìn)行移動相關(guān)的過程中，由于實際位移點左右兩側(cè)的有效圖像信息完全對稱，F(xiàn)FT計算得到的最大相關(guān)點及其左右兩個次最大相關(guān)點分布在以實際位移點為峰值的高斯曲線上，因此，按（1）式進(jìn)行亞像素擬合不會使相關(guān)系數(shù)峰值偏離實際位移。當(dāng)粒子圖像存在局部缺失現(xiàn)象時，實際位移點左側(cè)的圖像信息量不變，而實際位移點右側(cè)的有效圖像信息量將會減少，這使得FFT計算得到的實際位移點右側(cè)的相關(guān)系數(shù)偏小，亞像素擬合得到的相關(guān)系數(shù)峰值位置將偏向?qū)嶋H位移點的左側(cè)，這是導(dǎo)致平均流速減小的原因。當(dāng)粒子拖尾長度增加時，第二幀診斷窗口邊緣缺失的圖像信息量會增大，這是使得平均流速誤差隨拖尾長度增大而增大的原因。

圖7 無拖尾時粒子位移前后示意圖Fig.7 Schematic diagrams of particles moving without tailing

圖8 拖尾2像素時粒子位移前后示意圖Fig.8 Schematic diagrams of particles moving with two pixels'tailing

圖9 剔除粒子后拖尾長度與平均流速誤差的關(guān)系Fig.9 Relationship between mean－velocity error and tailing length after eliminating particles near the boundary

圖10 基于FFT計算結(jié)果的高斯亞像素擬合示意圖Fig.10 Schematic diagrams of Gaussian sub－pixel fit based on the result of FFT calculation

3.2 湍流強(qiáng)度誤差產(chǎn)生的原因

根據(jù) Adrain［16］，Willert［17］，Westerweel［6］等的研究，當(dāng)粒子粒徑大于2個像素后，計算得到的湍流強(qiáng)度會隨著粒子半徑的增大而增大。這是由于當(dāng)粒子粒徑增大時，F(xiàn)FT計算得到的相關(guān)系數(shù)峰值會變寬，在進(jìn)行亞像素擬合時，相關(guān)系數(shù)峰值位置的波動幅度會增大，從而導(dǎo)致湍流強(qiáng)度增大。當(dāng)示蹤粒子出現(xiàn)拖尾時，等價于增大了粒子的有效粒徑，因此，試驗得到的湍流強(qiáng)度會呈現(xiàn)隨拖尾長度的增大而增加的趨勢。

3.3 消除誤差的方法

上述分析表明，當(dāng)粒子圖像出現(xiàn)拖尾現(xiàn)象時，平均速度產(chǎn)生誤差是由于第二個診斷窗口右邊框附近的粒子圖像局部缺失引起的；而拖尾導(dǎo)致湍流強(qiáng)度增大的原因在于增大了粒子的有效粒徑。Anandarajah等［10］的研究表明，在進(jìn)行互相關(guān)計算時，可以增大第二個診斷窗口的尺寸，使得第一個診斷窗口中的所有粒子均能完整的出現(xiàn)在第二個診斷窗口中，從而消除粒子圖像局部缺失引起的平均流速誤差。Huang等［7］的研究則發(fā)現(xiàn)，對相關(guān)函數(shù)進(jìn)行標(biāo)準(zhǔn)化和對相關(guān)系數(shù)峰值附件的像素點處的相關(guān)系數(shù)進(jìn)行補(bǔ)償后，湍流強(qiáng)度將不再隨粒徑變化。

基于以上研究成果，可以通過以下方法消除粒子拖尾引起的DPIV測量誤差。首先，采用大小不同的診斷窗口I1和I2（I2＞I1），按（6）式對互相關(guān)函數(shù)進(jìn)行標(biāo)準(zhǔn)化后進(jìn)行相關(guān)運算：

其中A是第一和第二幀診斷窗口在移動相關(guān)過程中的重疊區(qū)域。計算得到各整像素點處的相關(guān)系數(shù)后，設(shè)R（xi，yi）為像素點相關(guān)系數(shù)峰值，R（xk，yk）∈｛R（xi＋1，yi），R（xi－1，yi），R（xi，yi＋1，R（xi，yi－1）｝是R（xi，yi）周圍用于亞像素擬合的四個相鄰像素點的相關(guān)系數(shù)。采用（7）式對R（xk，yk）進(jìn)行補(bǔ)償后，再按照（1）式進(jìn)行亞像素擬合，即可消除粒子拖尾引起的DPIV測量誤差。

其中，k（0＜k＜1）為補(bǔ)償系數(shù)。

4 結(jié) 論

對湍流等復(fù)雜流動的研究推動了基于連續(xù)激光源的高頻DPIV系統(tǒng)的研發(fā)和推廣，這種系統(tǒng)在對粒子圖像進(jìn)行曝光時可能會產(chǎn)生粒子拖尾現(xiàn)象。明渠湍流試驗和數(shù)值模擬試驗表明，隨著粒子拖尾長度的增加，計算得到的平均速度會減小，而湍流強(qiáng)度則會增大。

由于診斷窗口邊緣附近的粒子圖像局部缺失，使得FFT計算的實際位移點及其右側(cè)的相關(guān)系數(shù)偏小，最終導(dǎo)致高斯擬合后的相關(guān)系數(shù)峰值位置偏向左側(cè)，從而使得計算得到的平均流速偏小，而粒子拖尾本身不會直接導(dǎo)致平均速度誤差。當(dāng)粒子拖尾長度增加時，圖像局部缺失程度增加，F(xiàn)FT計算結(jié)果的不對稱性加大，使得平均流速隨拖尾長度的增大而減小。

粒子發(fā)生拖尾時，實際上相當(dāng)于增大了粒子的有效粒徑，等效粒徑的增大是導(dǎo)致計算湍流強(qiáng)度增大的原因。

通過增大第二個診斷窗口，對相關(guān)函數(shù)進(jìn)行標(biāo)準(zhǔn)化以及對相關(guān)系數(shù)峰值附近像素點處的相關(guān)系數(shù)進(jìn)行補(bǔ)償，可以避免平均速度誤差的產(chǎn)生，也可使得湍流強(qiáng)度不隨粒子粒徑的增大而增加，從而消除拖尾長度對PIV測量誤差的影響。因此，將基于連續(xù)激光源的高頻DPIV系統(tǒng)應(yīng)用于實際流動的測量是可行的。

［1］WILLERT C E，GHARIB M.Digital particle image velocimetry［J］.Experiments in Fluids.1991，10（4）：181－193.

［2］PRASAD A K.Particle image velocimetry［J］.Current Science，2000，79（1）：51－60.

［3］ADRIAN R J.Twenty years of particle image velocimetry［J］.Experiments in Fluids，2005，39（2）：159.

［4］GUEZENNEC Y G，KIRITSIS N.Statistical investigation of errors in particle image velocimetry［J］.Experiments in Fluids，1990，10（2－3）：138－146.

［5］KEANE R D，ADRIAN R J.Theory of cross－correlation analysis of PIV images［J］.Applied Scientific Research，1992，49（3）：191－215.

［6］WESTERWEEL J.Fundamentals of digital particle image velocimetry［J］.Measurement Science ＆ Technology，1997，8（12）：1379－1392.

［7］HUANG H，DABIRI D，GHARIB M.On errors of digital particle image velocimetry［J］.Measurement Science ＆ Technology，1997，8（12）：1427－1440.

［8］FORLITI D J，STRYKOWSKI P J，DEBATIN K.Bias and precision errors of digital particle image velocimetry［J］.Experiments in Fluids，2000，28（5）：436－447.

［9］WESTERWEEL J.Theoretical analysis of the measurement precision in particle image velocimetry［J］.Experiments in Fluids，2000，29（Supplement 1）：S3－S12.

［10］ANANDARAJAH K，HARGRAVE G K，HALLIWELL N A.Digital particle image velocimetry：particle image error（PIE）［J］.Journal of Physics：Conference Series，2006，45（1）：174－185.

［11］CARR Z R，AHMED K A，F(xiàn)ORLITI D J.Spatially correlated precision error in digital particle image velocimetry measurement of turbulent flow［J］.Experiments in Fluids，2009，47（1）：95－106.

［12］王龍，李丹勛，王興奎.高幀頻明渠湍流粒子圖像測速系統(tǒng)的研制與應(yīng)用［J］.水利學(xué)報，2008，39（7）：781－787.

［13］DI Sante A，THEUNISSEN R，VAN den Braembussche R A.A new facility for time－resolved PIV measurements in rotating channels［J］.Experiments in Fluids，2008，44（2）：179－188.

［14］余俊，萬津津，施鎏鎏，等.基于連續(xù)式激光光源的TR－PIV測試技術(shù)［J］.上海交通大學(xué)學(xué)報，2009，（8）：1254－1257.

［15］SANTE A D，CASTANO J G，BRAEMBUSSCHE R A V.Time resolved PIV in a rotating diverging channel：13th Int Symp on Application of Laser Techniques to Fluid Mechanics［Z］.Lisbon：2006.

［16］ADRIAN R J.Particle－imaging techniques for experimental fluids mechanics［J］.Annual Review of Fluid Mechanics，1991，23（1）：261－304.

［17］WILLERT C.The fully digital evaluation of photographic PIV recordings［J］.Applied Scientific Research，1996，56（2－3）：79－102.

陳啟剛（1987－），漢族、四川宜賓人，博士研究生，主要從事明渠湍流和泥沙運動試驗研究。通訊地址：北京市海淀區(qū)清華大學(xué)水利系泥沙研究所（100084）；聯(lián)系電話：15120001870；E－mail：cqg09＠ mails.tsinghua.edu.cn

Bias error induced by particle tailing in DPIV measurement

CHEN Qi－gang，LI Dan－xun，ZHONG Qiang，WANG Xing－kui
（State Key Laboratory of Hydroscience and Engineering，Tsinghua University，Beijing 100084，China）

The presence of particle tailing affects the performance of DPIV systems that use continuous wave laser as light source.Open－channel flow experiments and numerical stimulations show that the measured mean velocity of the flow decreases and the turbulent intensity increases when the exposure time prolongs.Underestimation of the mean velocity is due to the fact that the existence of particle tailing aggravates the partial lost of particle images near the boundary of the interrogation windows.Overestimation of the turbulent intensity is caused by the enlargement of the seeding particles due to particle tailing.An algorithm was proposed to eliminate or reduce these errors by enlargement of the second interrogation window，normalization of correlation function and complementation for integer peak neighbors.

particle image velocimetry；particle tailing；continuous wave laser；error；mean velocity；turbulent intensity

O352

A

1672－9897（2011）06－0077－06

2010－11－24；

2011－08－12

國家自然科學(xué)基金（50779023）