基于回歸分析的表面流場粒子匹配算法

陳 紅，周國梁，閆 靜，嵇 陽，晏成明

(1. 河海大學水利水電學院，江蘇 南京 210098； 2. 廣東水利電力職業技術學院，廣東 廣州 510635)

我國河流上興建了大量水利、水運、交通等涉河工程，受河流復雜水沙運動作用，大中型涉河工程關鍵設計參數多通過物理模型試驗研究確定。其中，流速是河流物理模型試驗中最基本的參量[1-2]，試驗工況確定、水流結構分析、模型驗證或工程方案比選均依賴流速數據分析，流速測量準確性直接關系到試驗成果水平。目前常用流速測量儀器有畢托管[1]、微型旋槳流速儀[3]、超聲波多普勒流速儀(acoustic Doppler velocimetry，ADV)[3]、大尺度流場圖像測速儀( large-scale particle image velocimetry，LSPIV)[3]、熱線熱膜流速儀(hot wire/hot film anemometry, HWFA)[4]、激光多普勒流速儀 (laser Doppler velocimetry/anemometry, LDV/LDA)[5]、粒子圖像測速儀(particle image velocimetry，PIV)[6-10]等。物理模型試驗水體含雜質較多，易阻塞畢托管或損壞HWFA敏感元件，因此畢托管和HWFA不適用于物理模型試驗。PIV和LDV測量范圍有限，布置煩瑣，易受雜質干擾，也無法適用于物理模型試驗[11-13]。微型旋槳流速儀和ADV常用于物理模型流速測量，微型旋槳流速儀操作簡便、價格低廉、測量精度較高，然而，對流態干擾較大，存在測量盲區，流速較小時測量誤差大。ADV測量精度高，可同時測量三維點流速和流向，但超聲波在空氣和水體中傳播速度不同，測量時應將發射接收探頭置于水面下才能進行準確測量，對流態有一定干擾，而且其價格昂貴，難以大規模應用。LSPIV通過測量水流表面示蹤粒子計算對應水流質點流速，具有測量范圍大、測量效率高、不干擾流場形態等顯著優點，目前已廣泛應用于物理模型試驗。

LSPIV系統主要包括示蹤粒子、圖像采集和圖像處理等[14]，王興奎等[15]基于圖像采集卡研發了LSPIV系統，建立了最鄰近算法實現了粒子匹配，并測量了三峽工程物理模型表面流場，分析了工程前后流場變化。唐洪武等[16]改進了粒子匹配算法，采用Hopfield神經網絡對流場數據進行校核和修正，并應用到快速自航船模船尾流場及航道斷面流場的實時測量。陳誠[17]對示蹤粒子跟隨性和圖像處理技術進行了深入研究，并應用于黃河模型試驗中。LSPIV最核心技術是示蹤粒子匹配，常規采用最鄰近算法，但當粒子密度較大時，粒子匹配誤差會大幅增加，降低流場測量精度，影響工程前后流速變化分析。

本文針對粒子密度大時粒子匹配誤差大等問題，改進了最鄰近匹配算法，提出多幀回歸粒子匹配算法，在此基礎上開發了改進型大尺度模型表面流場圖像測試系統(ILSPIV)，并將該系統應用于八卦洲河工模型，同時與旋槳流速儀測量結果進行了對比。

1 流場圖像特征

在模型水流表面播散跟隨性能好的示蹤粒子跟隨水流運動，利用CCD攝像頭記錄一定時間間隔的流場圖像(如圖1所示)，經圖像處理提取示蹤粒子像，通過粒子匹配，獲取t1時刻示蹤粒子經過Δt時間后在t2時刻的圖像，得到Δt時間內粒子位移，從而求得示蹤粒子的速度，即對應水流點的速度。

圖1 流場粒子圖

圖1流場圖像中示蹤粒子為白色，RGB灰度值較大。在圖像中選取5個粒子像，依次讀取中心點、邊緣點及相鄰背景點灰度值，RGB灰度值如表1所示?；叶戎当砻鞑煌恢檬聚櫫Ｗ覴GB灰度值存在一定差異，粒子中心點及邊緣點處灰度值大于背景點灰度值，且隨著背景點灰度值增大，粒子中心點及邊緣點處灰度值也相應增大，粒子中心點、邊緣點、背景點灰度值變化趨勢表現一致，示蹤粒子中心點灰度值與對應背景點灰度值存在一定相關關系。

表1 示蹤粒子RGB灰度值

2 粒子像提取

粒子像提取是將示蹤粒子從圖像中分離出來，獲取示蹤粒子中心坐標。傳統粒子提取算法主要采用閾值分割，分割算法如下式:

(1)

式中：gi(x,y)為第i幀圖像在坐標(x,y)處灰度值；ki(x,y)為閾值分割后圖像灰度值；T為閾值。

圖2 閾值分割結果對比

閾值分割處理前景和背景差異較小的圖像時分割效果較差。對圖1直接進行閾值分割，處理結果如圖2(a)所示。幀差法的基本原理是首先對相鄰兩幀圖像做差分運算，再通過閾值變換消除背景、分離運動目標信息，該方法能提高背景和前景差異，有利于閾值分割，適用于運動目標和背景差異較小的圖像檢測，具有算法簡單、運算速度快等優點。幀差法計算方法如下:

pi(x,y)=gi+1(x,y)-gi(x,y)

(2)

式中：pi(x,y)為幀差法處理后圖像灰度值。

幀差法能消除背景影響，即使背景存在因外部條件導致的灰度差異，只要該差異小于背景與示蹤粒子間的差異，仍可采用幀差法和閾值分割組合(幀差法閾值分割)提取示蹤粒子。幀差法閾值分割先采用式(2)計算出pi(x,y)，然后將gi(x,y)=pi(x,y)代入式(1)得出ki(x,y)。對圖1進行幀差法閾值分割，結果如圖2(b)所示，幀差法閾值分割能有效消除模糊背景，示蹤粒子提取準確度更高。

選取幀差法閾值分割處理圖像，對第i幀和i+1幀圖像進行幀差法閾值分割得到ki(x,y)，對第i+1幀和i+2幀圖像進行幀差法閾值分割得到ki+1(x,y),對ki(x,y)和ki+1(x,y)進行和運算，得到hi(x,y)：

(3)

hi(x,y)是第i幀、i+1幀、i+2幀粒子像全部疊加在一起得到的粒子圖像，背景灰度值為0，粒子灰度值為255。沿x、y方向依次讀取坐標灰度值，當灰度值為255時，即可判別該坐標為粒子像坐標，記為(xm,yn),m=0,1,…,l1;n=0,1,…,l2;l1、l2分別為提取出的x、y方向上粒子像的最大數。

3 粒子匹配

示蹤粒子像包含了第i幀、i+1幀、i+2幀示蹤粒子像，粒子濃度較大時，無法通過最鄰近匹配算法進行粒子同源匹配。根據流體運動連續性，多幀同源示蹤粒子像遵循以下原則：①同源示蹤粒子運動位移量較??；②同源示蹤粒子運動軌跡連續光滑，不存在突變，表現為多幀同源示蹤粒子像運動軌跡近似直線，回歸分析時相關系數最大。本文改進了最鄰近匹配算法，提出了多幀回歸粒子匹配算法，算法的步驟如下：

步驟1搜索第i+1時刻同源疑似點。原則①表明示蹤粒子運動范圍是以初始時刻坐標為圓心，最大運動位移量d0(與運動速度相關，可由斷面平均流速估計)為半徑的圓，根據該原則初步確定同源點，稱為第i+1時刻同源疑似點。從示蹤粒子像坐標提取起始點(xm,yn)，為提高圖像處理速度，設立正方形查詢窗口，邊長2r，在該正方形范圍內提取點(xm+k,yn+h)，其中k=-r, -r+1,…,r，h=-r, -r+1,…,r，兩個點的距離d1為

(4)

如果d1

步驟2搜索第i+2時刻同源疑似點。以點(xm+k,yn+h)為中心，在邊長為2r的正方形查詢窗口內提取點(xm+k+g,yn+h+f)，其中g=-r, -r+1,…,r，f=-r, -r+1,…,r。計算這兩個點的距離d2，如果d2

(5)

總平方和為

(6)

假設上述三點擬合直線方程為

y=ax+b

(7)

(8)

(9)

相關系數Rm為

(10)

圖3 同源疑似示蹤粒子像

依次對其他同源疑似點與起始點(xm,yn)進行回歸分析，相關系數最大的作為同源點，完成全流場粒子匹配。

圖3為根據原則①處理獲取的同源示蹤粒子像，黑色點坐標為(x0,y0)，其同源疑似點有點(x0+1,y0-1)、(x0+1,y0)和(x0,y0+1)；點(x0+1,y0-1)的同源疑似點有(x0+2,y0-2)、(x0+2,y0-1)和(x0+1,y0);點(x0+1,y0)的同源疑似點有(x0+1,y0-1)、(x0,y0+1)和(x0+2,y0-1);點(x0,y0+1)的同源疑似點有(x0-1,y0+2)和(x0+1,y0)，則對應相關系數最高的點為點(x0,y0)、(x0+1,y0-1)和(x0+2,y0-2)，其擬合直線為

y=x

(11)

擬合直線相關系數為1，依次類推可得，黑色點與紅色點屬于同源粒子。

基于回歸分析的粒子匹配算法開發了改進型大尺度模型表面流場圖像測試系統(ILSPIV)，對圖1表面流場粒子圖像進行處理，得到流速矢量如圖4所示，總體上矢量圖錯誤較少，矢量與粒子軌跡重合度高。

圖4 表面流場矢量圖

4 實際應用

長江南京河段八卦洲左汊河道沿程阻力遠大于右汊，導致左汊逐漸淤積、萎縮，既影響左汊沿岸生產生活，又影響南京長江河段整體河勢穩定。因此，需要采取有效工程措施穩定八卦洲汊道河勢，提高左汊分流比。為對比研究工程效果，采用物理模型研究了工程前后洲頭流速分布等變化情況。根據試驗要求，確定八卦洲汊道整治物理模型模擬范圍，上邊界為大勝關，下邊界至九鄉河，模擬河段總長度約43 km。選定模型平面比尺λL= 480，垂直比尺λH= 120，流速比尺λv= 10.95。

試驗工況流量為15 290 m3/s，選定測流斷面位于分流口，分別采用旋槳流速儀和ILSPIV對流場進行對比測試。旋槳流速儀測量數據為斷面垂線平均流速，ILSPIV測量流速為表面流速，根據有關理論兩者流速關系為

vc=0.85vb

(12)

式中：vc為斷面垂線平均流速；vb為表面流速。

依據式(12)將ILSPIV測量的表面流速數據換算成斷面垂線平均流速，圖5給出了兩種測試技術獲得的工程前后斷面垂線平均流速。總體上，兩種測試技術測量的斷面垂線平均流速分布形態相似，局部流速點存在一定偏差，相對誤差小于2.7%。ILSPIV測量的工程前后流速變化范圍為-0.25～0.22 m/s，旋槳流速儀測量的工程前后流速變化范圍為-0.18～0.22 m/s，表明工程前后兩種技術測量結果相對變化量基本一致，對工程影響評估無差異。

圖5 兩種測試技術獲得的工程前后斷面垂線平均流速對比

5 結 語

粒子匹配是影響LSPIV測量精度最重要的因素之一，為提高粒子匹配精度，改進了最鄰近匹配算法，提出了多幀回歸粒子匹配算法，從動態圖像中提取連續幀，分組疊加再做減運算，采用一元線性回歸，以相關系數為評估參數，確保了高粒子濃度條件下粒子的精確匹配。運用多幀回歸粒子匹配算法開發了ILSPIV，將ILSPIV應用到八卦洲河工模型試驗，ILSPIV和旋槳流速儀測量的斷面垂線平均流速分布形態一致，相對偏差小于2.7%，工程前后兩種技術測量結果相對變化量偏差小，表明多幀回歸粒子匹配算法能提升高濃度粒子條件下流速測量精度，更有利于工程前后流速變化分析和評估。