基于經驗模態分解的粒子圖像測速流場重構研究

王曉宇,賀靜,柯漢兵,林梅

(1.西安交通大學能源與動力工程學院,710049,西安;2.武漢第二船舶設計研究所,430205,武漢)

粒子圖像測速(PIV)技術是一種誕生于20世紀80年代的無擾動流場測試技術,它克服了傳統接觸式流場測量中單點測量的局限性,可以對存在變化流速及各類渦結構的復雜流動實現全流場瞬時非接觸測量[1-2]。在PIV流場圖像后處理中常采用光學楊氏條紋法、自相關法或互相關法,逐點處理PIV圖像,從而獲得流場速度數據[3-4]。由于在實驗中可能存在示蹤粒子分布不均勻、PIV圖像變形或扭曲、PIV圖像背景噪聲、成像系統變形等因素,獲得的流速數據中會存在干擾正常流場分析的錯誤數據。當錯誤數據造成的誤差超過2%時,便會影響正常的流場分析研判[5]。

針對PIV流場數據中錯誤數據的消除方法,國內諸多學者都進行了相關研究。王平讓總結了錯誤數據產生的原因,并從PIV圖像濾波處理的角度提出了消除錯誤數據的方法[6]。王燦星等建立了基于互相關技術的連續PIV圖像之間的粒子對應方法,提出了識別錯誤數據的方法和消除準則[7]。高琪等基于速度場的本征正交分解后處理技術,通過迭代方法有效地實現了錯誤數據的識別與修復[8]。吳龍華等在對人腦判別PIV錯誤數據方式進行模擬的基礎上,建立了可以進行準確錯誤數據識別的多證據推理Hopfield神經網絡模型[9]。王晉軍等引入不可壓縮連續性方程對PIV測得的流場數據進行后處理,使得修正后的流場嚴格滿足差分形式的連續性方程,實現了流場數據的準確修正[10]。阮曉東等采用Delaunay網格技術建立了基于流體連續性假設的自動檢測PIV數據中錯誤數據的算法,實現了錯誤數據的識別與替換[11]。張工利用Visual C++開發了通用的PIV圖像后處理軟件,實現了PIV圖像的增強等功能,可以從PIV圖像中消除部分錯誤數據[12]。王元等廣泛研究了各類提高流場測量及數據后處理精度的粒子追蹤算法,減少了錯誤數據的產生[13-16]。上述方法對PIV流場數據中錯誤數據的剔除和修正都有較好的效果,但也存在處理流程復雜、效率不高、引入非物理偏差等不足。

本研究提出了一種基于經驗模態分解(EMD)的PIV流場數據后處理方法,將數據分解得到的波形異常的本征模態分量(IMF)進行濾波處理并與其余波形平滑的正常本征模態分量進行反向疊加,重構流場數據,使錯誤數據得到有效剔除。當分解得到的本征模態分量較多且錯誤數據集中在某一個本征模態分量時,直接將與原始流場數據相關的本征模態分量反向疊加重構流場數據,摒棄不相關的本征模態分量,實現錯誤數據的直接剔除。本方法已經在多個模擬或實驗獲取的變化較為緩和的平穩流場數據處理中得到了良好的使用,本文利用該方法對2個典型的人為添加錯誤矢量的模擬標準流場進行數據處理,同時對某實驗的原始流場數據進行修正,以上述3個流場為例,說明并驗證該方法的可行性。

1 基于經驗模態分解的流場數據重構原理

1.1 經驗模態分解概述

經驗模態分解是Huang等于1998年提出的一種無需基函數的自適應性信號數據處理新方法,在各工程領域得到了廣泛、高效的應用[17-18]。該方法利用數據自身的尺度特征進行數據分解,特別適用于非平穩、非線性數據的平穩化處理。經驗模態分解處理信號的基本步驟如下。

(1)確定隨機信號x(t)在整個時間尺度上的所有極大值點和極小值點,并用三次樣條差值曲線分別連結所有極大值點和極小值點,形成上包絡線a(t)和下包絡線b(t),求取上、下包絡線的平均值

(1)

(2)求取原始隨機信號x(t)與m(t)之間的差值

h(t)=x(t)-m(t)

(2)

(3)如果h(t)的極值點與過零點的個數最多相差一個且h(t)的上、下包絡線的均值為0,則h(t)為第一個本征模態分量c1。求取x(t)與c1的差值

r(t)=x(t)-c1

(3)

將r(t)視為最新的隨機信號,重新開始上述分解過程(1)～(3)求取剩余的本征模態分量,直至r(t)變為單調信號。如果h(t)不能滿足上述條件,則將h(t)視為最新的隨機信號,重新開始上述分解過程(1)～(3)。

關于經驗模態分解更加具體的數學化描述及證明過程請參見文獻[19-20]。經過經驗模態分解后,一組數據可以被分解為若干個表征不同數據尺度和數據特征的本征模態分量,它們是組成原始信號的不同尺度的一系列分量,是對原始信號特征的一系列表達。理論上連續的兩個本征模態分量的標準差在0.2～0.3之間時,則分解得到的本征模態分量既能滿足線性條件,穩定性良好,又能夠較好地反映原始隨機信號的性質。經驗模態分解具有正交性和完備性的特征,因此將分解得到的本征模態分量反向疊加即可得到原始數據。在處理現代流動測試技術所獲得的流場數據時,數據量的大小并不會直接影響所分解出的本征模態分量數量,而原始流場的復雜程度則對本征模態分量的數量影響較大。原始數據經過經驗模態分解后一般可以得到5～15個本征模態分量,其中本征模態分量在5～10個時可以認為較少,超過10個時可以認為分解得到的本征模態分量較多。

1.2 相關系數概述

在數據處理中常用相關系數來表征兩組數據之間的相關程度,在本研究中通過求解不同的本征模態分量與原始流場數據的相關系數,可以判斷該分量與原始流場數據的相關程度。相關系數的劃分并沒有固定的標準,而是根據所處理數據的具體特征來決定的。表1所示的是一種相關性的劃分方式,本研究根據這種相關系數大小的劃分方式將本征模態分量與原始流場數據的關系歸類為不相關、微相關、實相關、高度相關4類,其中微相關、實相關、高度相關統稱為相關。

表1 相關性劃分

1.3 流場數據處理流程

圖1是本研究所提出的基于經驗模態分解的PIV流場數據后處理方法流程圖,其中實線表示的流程是通用方法,虛線表示的流程是在原始流場較為復雜、原始流場數據經過經驗模態分解后得到的本征模態分量較多且錯誤數據集中在某一個本征模態分量時對通用方法的一種簡化。

圖1 基于EMD的PIV流場數據后處理流程圖

錯誤數據在數值大小和方向上與正確數據有明顯的差異,在原始流場數據沒有經過經驗模態分解時,錯誤數據摻雜在正確數據中難以識別。當用經驗模態分解將原始流場數據分解為若干個本征模態分量時,錯誤數據就會造成本征模態分量曲線出現明顯的局部不平滑。由于錯誤數據造成本征模態分量曲線局部的不平滑十分明顯,因此常用的濾波方法就可以將錯誤數據點在本征模態分量中濾除而不誤刪正確數據。根據經驗模態分解的正交性和完備性,波形異常的本征模態分量經濾波處理后與波形平滑正常的本征模態分量反向疊加就可以重構剔除了錯誤數據的流場數據。上述的原始流場數據處理思路可表示為圖1中實線所示的流程。

原始流場數據經過經驗模態分解后得到的本征模態分量越多,每個本征模態分量攜帶的原始流場數據就越少,所占原始流場數據的比例就越小。在有些實驗中,粒子圖像測速技術所獲得的流場數據相對比較復雜且錯誤數據比較明顯,這時分解得到的本征模態分量較多且錯誤數據會集中在某一個本征模態分量。摻雜了錯誤數據的本征模態分量波形異常且與原始流場數據的相關性較低(呈不相關),因此在該流程中可以通過求解各本征模態分量與原始流場數據的相關性大小來判斷錯誤數據所在的本征模態分量。為了提高處理速度和效率,可以按照圖1中虛線所示的處理流程,省略相關性低(呈不相關)的本征模態分量的濾波過程,直接將該分量刪除并將其余相關性較高的本征模態分量反向疊加直接重構流場數據。隨著錯誤數據的刪除,也會損失該本征模態分量中的正確數據,但在分解得到的本征模態分量較多的情況下,損失的小部分正確數據不會對整體流場數據造成明顯的影響,流場的特性仍然可以清楚地表達。

2 算例驗證與應用分析

2.1 模擬流場中的驗證與分析

圖2是用計算機模擬的流速為1 m/s的單向均勻流場,在4 000個原始流場數據中加入了明顯的錯誤數據,人為添加的誤差約為1.7%,因此速度矢量場中出現了數值大小和方向明顯錯誤的速度矢量。圖3是原始流場數據按圖1中實線流程處理后得到的速度矢量場。圖4和圖5是對該原始流場x、y方向全部的4 000個流場數據進行經驗模態分解后得到的本征模態分量,分別為6個和8個,不符合使用簡化流程的要求,因此采用的是通用方法。由圖可知x、y方向的錯誤數據都集中在了第1個本征模態分量中,因此該本征模態分量的曲線中出現了不平滑的跳變點,也正是錯誤數據存在的位置。對x、y方向上錯誤數據集中的第1個本征模態分量進行濾波處理,即可把錯誤數據存在的跳變點濾除而不誤刪正確數據。將x、y方向上濾波處理后的第1個本征模態分量與其他本征模態分量反向疊加,即可重構流場數據。圖2與圖3對比表明,實線表示的處理流程能有效實現該單向均勻流場錯誤數據的剔除和速度矢量場的修正,修正并重構后的流場數據與正確流場數據的誤差約為0.002%。

圖2 原始單向流場數據得到的速度矢量場

圖3 處理后的單向流場數據得到的速度矢量場

圖4 原始單向流場x方向數據的本征模態分量

圖5 原始單向流場y方向數據的本征模態分量

為了驗證所提出的方法在處理相對復雜的流場數據時的有效性,本研究模擬了速度大小不定的單向旋轉流場,并在原始流場數據中加入了一些明顯的錯誤數據,人為添加的誤差約為3.3%。圖6和圖7是對該原始流場x、y方向的流場數據進行經驗模態分解后得到的本征模態分量與原始流場數據的相關系數,本征模態分量的數量和特點符合上文提出的虛線表示的簡化處理流程的要求,因此處理該流場數據時采用的是簡化流程。圖8和圖9分別是由原始流場數據和按圖1中虛線表示的簡化流程處理后的流場數據得到的速度矢量場。對比表明,圖1中虛線表示的簡化處理流程能夠有效地實現該速度大小不定的單向旋轉流場中錯誤數據的剔除和速度矢量場的修正,修正并重構后的流場數據與正確流場數據的誤差約為0.18%。

圖6 旋轉流場x方向流場數據的各本征模態分量與原始流場數據的相關系數

圖7 旋轉流場y方向流場數據的各本征模態分量與原始流場數據的相關系數

圖8 原始旋轉流場數據得到的速度矢量場

圖9 處理后的旋轉流場數據得到的速度矢量場

2.2 實驗流場中的應用與分析

為了驗證本文所提出的基于經驗模態分解的PIV流場數據后處理方法在實驗流場數據處理中的準確性與實用性,本研究對低速風洞中某一T型通道的PIV流場原始數據進行了處理與分析。關于該T型通道PIV流動測試實驗的具體細節請參見文獻[21-22]。

正如上文所述,實驗流場的復雜性一般相對較高,流場數據經過經驗模態分解之后得到的本征模態分量較多。該實驗流場x、y方向的原始流場數據經過經驗模態分解之后分別得到了10個本征模態分量和12個本征模態分量,可以用圖1中的簡化流程進行數據處理。圖10和圖11分別表示了x方向流場數據、y方向流場數據的本征模態分量與原始流場數據之間的相關程度。由圖可知,x方向流場數據的第1個本征模態分量和y方向流場數據的第6個本征模態分量與原始流場數據的整體相關性很小(呈不相關),是錯誤數據集中的具體表現。該流場數據可以按照圖1中虛線表示的流程進行處理,直接刪除這兩個本征模態分量,將x、y方向上剩余的9個、11個本征模態分量反向疊加重構流場數據。圖12、圖13分別表示了該流場數據處理前后的速度矢量場,可知原始流場數據中位于速度矢量場上部的錯誤數據(圖12中區域1)得到了準確的剔除,一些明顯的錯誤矢量得到了有效的消除,速度矢量場更加平滑,品質得到了顯著改善,速度矢量場上部的流動特征(圖13中區域1＇)也更加趨近于文獻[22]報道的各工況下的真實流動情況,但在個別不存在錯誤矢量的區域出現了流場特征的輕微變化。這是隨著簡化流程中呈不相關的本征模態分量的直接刪除,造成少部分正確數據損失所導致的。如果對局部流場特性的準確性要求較高,可以不采用簡化流程,而采用相對較為煩瑣的通用方法來進行處理。

圖10 實驗流場x方向流場數據的各本征模態分量與原始流場數據的相關系數

圖11 實驗流場y方向流場數據的各本征模態分量與原始流場數據的相關系數

圖12 原始實驗流場數據得到的速度矢量場

圖13 處理后的實驗流場數據得到的速度矢量場

3 結 論

本研究提出了一種基于經驗模態分解的PIV流場數據后處理方法,對存在錯誤數據的本征模態分量進行濾波處理或直接刪除,將其余本征模態分量反向疊加,實現了錯誤數據的剔除與流場數據的重構。利用該方法處理了計算機模擬流場和某實驗流場的原始數據,結果表明錯誤數據得到了有效剔除,重構的流場數據可更加真實地反映流場特征。

目前任何一種PIV流場數據后處理方法都不能完全將錯誤數據剔除而不影響正確數據,本研究所提出的方法也是如此,但本方法具備以下兩個優點:首先,該方法最大限度地從原始流場數據中剔除了錯誤數據,并非通過各類算法在計算時“避開”錯誤數據,具備處理的徹底性;其次,該方法簡單易行,沒有復雜的錯誤數據識別流程與標準,可以快捷簡便地提高流場分析的準確性。在PIV技術中,錯誤數據超過2%就已經會影響流場分析,證明實驗中存在的誤差太大,此時解決誤差大的最好辦法是重新進行實驗,而不是消除已有的錯誤數據。就這一點而言,本方法對誤差在2%以內的流場數據是有效的,但本方法也有自身的缺點,當所處理的流場變化非常劇烈、雜亂(例如射流撞擊壁面)時,本方法容易把部分真實反映流場的數據視為錯誤數據,因此本方法只適用于處理變化較為緩和的平穩流場。

受該方法的啟發,也可以將原始的PIV圖像視為一個二維數組,利用經驗模態分解的方法通過上述流程來對其進行處理。這時處理的對象是PIV原始粒子圖像,各個本征模態分量也是以圖像或數組的形式來呈現的,最終通過本征模態分量的反向疊加得到處理后的PIV圖像,并從中提取流場數據。此外,本方法在三維流場數據中的應用及可行性還有待驗證,在接下來的研究工作中,上述兩部分內容將是研究的重點。