基于線陣CCD的風浪環境下水面弱流場速度測量方法

劉 利, 王小青, 陳鵬真, 種勁松

(1. 中國科學院電子學研究所微波成像技術重點實驗室, 北京 100190; 2. 中國科學院大學, 北京 100190)

0 引 言

大量實驗和理論研究表明，水下地形、內波和漩渦等海洋現象都是通過其產生的表面弱流場與表面波的相互作用而被微波遙感手段所觀測到[1]。而對風浪環境下的弱流場進行測量是研究波流相互作用機理的核心步驟。因此，研究風浪環境下表面弱流場測量的方法和裝置對于研究波流相互作用機理具有重要意義。

基于圖像處理技術的傳統流場測量方法，需要在水面散布示蹤粒子，根據示蹤粒子的移動來測量流場[2-3]，這類方法稱為PIV(Particle Image Velocimetry，粒子圖像測速)方法。風浪環境下，風生波會對示蹤粒子的運動產生嚴重影響，用PIV方法測得的是波浪振動速度與表面弱流場的疊加結果，難以將二者區分。而需要測量的表面弱流場是與風生波無關的，因此傳統的流場測量方法無法工作在風浪環境下，難以對流場與風生波相互作用展開研究。

本文提出的基于線陣CCD(Charge-coupled Device，電荷耦合組件)的流場測量方法作為一種非接觸式的光學測量方法，在風浪環境下可以高精度地獲取表面弱流場信息而不受波浪振動的影響。本方法是在基于線陣CCD的波浪測量方法上發展而來。基于線陣CCD的波浪測量方法是近年來在國際上出現的一種高精度、高速的空時二維測量方法[4-7]，其基本原理是在特定光照條件下，水面的光強與水面的傾角有一一對應的映射關系，通過獲取水面光強可以反演出水面斜率，進而反演出水面波高。之前的研究主要集中于使用基于線陣CCD或面陣CCD的光學方法對波浪斜率、波高分布以及波譜、方向譜和能量譜等進行測量[4-11]。目前尚無利用線陣CCD測量水面流場的相關報告。

線陣CCD由于其時間、空間分辨率高，可以對振動頻率很高的微尺度波進行測量，而且可以獲得波浪相位的連續變化。本文所提出的測量方法就是根據表面波浪的相位變化獲得表面波的相速度，其值是波的本征相速度與外界流場的疊加，外界流場包括風生流場及所要測量的表面弱流場。根據經典的波浪理論，波的本征相速度與波浪的振幅高低無關；風速穩定后風生流場也是不變的，因此可以通過表面波相速度的變化計算表面弱流場速度。最后用此方法對水槽實驗獲得的實際CCD數據進行處理，驗證了本方法的正確性。

1 基于線陣CCD的風浪環境下水面流速測量方法

風浪環境下，表面波的相速度可以表示為[12]

c=c0+ud+uc

(1)

式中：c0為本征相速度，ud為風引起的風生流場速度，uc即為所要測量的表面弱流場速度。根據經典波浪理論，波的本征相速度為：

(2)

式中：g為重力加速度，T為表面張力系數，ρ為水的密度，k為表面波波數。由公式(2)可知，波的本征相速度只與波長有關，與波的振幅無關。而且風速穩定后風生流場的速度也是穩定不變的，因此可以通過表面波相速度的變化來測量表面弱流場速度。

表面波的相速度c為波長λ與周期T之比

c=λf

(3)

(4)

式中：波數k=2π/λ，角頻率ω=2πf=2π/T。

因此，為了計算某一波長表面波的相速度，可以有以下兩種方法：

(1) 波頻率法

此方法是計算某一波長的波頻率f，通過公式(3)得到相速度。在計算頻率f時，需要對斜率數據進行時間傅里葉變換，頻率分辨率與所選時間軸長度有關。由于頻率分辨率Δf與時間軸長度ts關系為：Δf=1/ts，因此頻率分辨率若要達到1Hz，時間軸長度須選為1s；頻率分辨率若要達到0.1Hz，時間軸長度須選為10s。因此采用此方法要達到高的頻率分辨就必然會損失時間精度，但此方法的抗噪性能較好。

(2) 波相位法

此方法是計算相位差Δφ，通過公式(4)計算相速度。本方法的優勢是時間分辨率很高，但抗噪性能較差。

由于本研究需要提取較高時間分辨率的流場信息，因此采用第二種計算方法，數據處理中為了減小噪聲對測量結果的影響，對數據進行自適應的二維頻域濾波，對相位差分時間間隔進行了優化，并對不同波長的表面波相速度測量結果進行統計平均得到最終較好的測量結果。

所提出的測量方法主要步驟如圖1所示：

圖1 算法總體流程圖

下面詳細描述圖1中的4個步驟。

(1) 幾何校正及光強斜率轉換

在特定光照條件下，水面的光強與水面的傾角有一一對應的映射關系。在標定映射關系后，就可以通過線陣CCD相機獲取水面的光強變化來反演出水面斜率。由于CCD相機的位置及光學相機拍攝物體的“近大遠小”特征，在標定映射關系前，需要進行幾何校正，將CCD相機拍攝圖像中的像素位置準確投影到空間位置。

(2) 斜率噪聲抑制

由于實驗設備以及測量環境的影響，獲得的數據中不可避免地會混雜有噪聲。噪聲通常較為均勻地分布在時-空二維頻域；而水面波浪由于其特殊的頻散關系，波浪斜率信號在時-空二維頻域分布較為集中。因此利用此原理可以通過二維頻域自適應濾波的方法大幅度提高波浪斜率數據的信噪比。

假設線陣CCD掃描儀測得的斜率表示為s(x,t)，其中x表示空間位置，t表示時間，則斜率的自相關函數可以表示為：

Rs(ξ,τ)=〈s(x,t)s(x+ξ,t+τ)〉

(5)

斜率譜S(kx,ω)為

?Rs(ξ,τ)e-j(kxξ+ωτ)dξdτ

(6)

式中：kx為表面波的波數矢量，ω為表面波的角頻率。根據經典的波浪理論，斜率數據在頻域中主要集中于頻散關系曲線附近，

(7)

式中：ρ為水體密度，T為表面張力系數，g為重力加速度。可用頻域二維濾波器對斜率二維頻譜濾波，濾波范圍可自適應地選擇為其譜寬。

(3) 相位信息提取及修正

假設在時刻t0，沿距離向的斜率分布為s(x,t0)，其空間傅里葉變換為S(k,t0)。傅里葉變換之后可以對不同波數的表面波相位進行計算。對于波長為λ0(波數為k0=2π/λ0)的表面波，假設其傅里葉變換形式為：

S(k0,t0)=SR+jSI

(8)

其相位為：

(9)

由于線陣CCD相機自身存儲空間有限，在其離線采集完一幀數據后需要將數據傳輸回計算機才可以開始下一幀數據的采集。因此相對于幀數據內的時間間隔，幀數據與幀數據間的時間間隔并不穩定。由此導致的結果是相位數據在幀數據間會有小幅度的跳變，這種跳變使得相位差在幀數據切換時刻變化比較劇烈。使用未經修正的相位數據計算出的相速度會出現明顯異常，因此需要對原始的相位數據進行修正。

具體修正方法為：由相位數據得到相位差后，幀數據切換時刻的相位差修正為采用其相鄰數據經過插值后得到的新的相位差。

(4) 表面波相速度計算

計算得到相位數據后，根據公式(4)就可以計算表面波相速度。對于波長為λ0(波數為k0=2π/λ0)的表面波，假設在時刻t0其相位為φ(k0,t0)，在時刻t1=t0+Δt其相位為φ(k0,t1)，則其角頻率為

(10)

在根據公式(4)計算時，由于存在相位噪聲的影響，選擇不同的時間間隔Δt計算得出的相速度精度是不同的。時間間隔過小，相機內部的熱噪聲影響較大，使得最終計算結果中誤差占比過大；加大時間間隔可以減小熱噪聲的影響，但加大時間間隔后由于波浪的形變會導致表面波的相干性降低，相位噪聲增加；此外加大時間間隔還會引起時間分辨率的下降。因此在時間間隔的選擇上需要折中考慮。

以標準差表示精度，標準差越小，精度越高。由公式(4)可知，

(11)

其中，角頻率精度為

(12)

為了選擇合適的時間間隔，2.2節中描述了實驗中所選時間間隔與測量精度的關系。

在本文第2節中，將以實驗數據為例對圖1的算法流程進行詳細的說明。

2 實驗驗證

2.1實驗裝置

實驗于2013年3月在解放軍理工大學分層流實驗水槽中進行。該水槽的長、寬和高分別為12，1.5和1.2m。系統在實驗過程中會產生內波，內波是指發生在密度穩定層化的流體內部的一種波動。本實驗中使用鹽水和淡水產生密度分層，如圖2所示。內波的傳播會在表面產生流場。實驗中風速可變，內波振幅為10cm，傳播方向與水槽長度方向平行，內波傳播過程中遇到水槽壁會發生反射現象。圖2為實驗裝置示意圖。

圖中L表示有效掃描區域。通過特殊的光源以及散光板產生特殊的光照條件，使得水面光強與表面波斜率有近似的一一映射關系。光線經過水面反射后進入圖像采集器件——線陣CCD相機。實驗中CCD相機距離水面的高度為44cm，觀測角度(與垂直方向夾角)為34°。

圖2 水槽內波觀測實驗光學裝置示意圖

光學系統的主要參數如下：有效掃描長度：L=23cm；分辨率：0.38mm(空間分辨率)，1/300s(時間分辨率)；掃描頻率：300Hz；掃描點數：2048。

實驗中采用了6臺CCD相機進行并行數據采集,圖3是實際實驗裝置圖。

圖3 實驗裝置圖

2.2實驗結果及分析

以下給出第1節中所述的4個步驟的具體處理過程及部分結果。

(1) 幾何校正及光強斜率轉換

實驗中得到的時間-空間二維圖像數據如圖4所示。橫坐標x表示距離(單位：cm)，縱坐標t表示時間(單位：s)，灰度值大小表示經過水面反射后進入相機的光強。每幅圖像的尺寸為300pixel×2048pixel，掃描頻率為300Hz，因此一幅圖代表1s時間內掃描區域內波動變化。表面波的運動在圖像中表現為傾斜的條紋，傾斜條紋的斜率即表示表面波的傳播速度。

CCD相機通過記錄掃描區域的光強度從而獲得掃描區域的實時波斜率，如圖5所示。CCD在空間上的掃描范圍約為23cm，每個采樣點的間隔約為0.23m/600=0.38mm,那么CCD可以記錄的表面波的波數范圍理論上為30～8200rad/m。

比較圖4與5，圖4中的傾斜條紋存在一定的幾何畸變(條紋彎曲)，經過幾何校正后這種畸變消失(圖5)；而且圖4中光強分布沿距離方向是不均勻的，但經過光強斜率轉換后的斜率分布比較均勻。

圖4 實驗獲得的CCD圖像數據

圖5 由圖4計算得出的波斜率分布

(2) 斜率噪聲抑制

實驗中測得的斜率s(x,t)和斜率譜S(kx,ω)分別如圖5和6所示。圖6中的黑色曲線為根據頻散關系(公式(7))得出的頻率與波數關系曲線。圖6中波頻集中區域與頻散關系曲線吻合很好，這一方面驗證了實驗所獲取數據的有效性，另一方面也表明了頻率成分主要分布于頻散關系曲線附近，因此可用圖7所示的自適應二維頻域濾波器進行濾波處理，濾波后的斜率分布如圖8所示，與圖5相比斜率分布更加清晰。

某一距離點上波斜率隨時間的變化如圖9所示。

圖6 斜率譜

圖7 頻域二維濾波器

圖8 濾波后的斜率分布

圖9 1s內波斜率隨時間變化

(3) 相位信息提取及修正

以波長λ=4.57cm為例，在觀測時間第28s內相位變化如圖10所示。由于直接計算所得相位值位于[-π,π]區間，在計算相位差時需要對相位進行解卷繞。

圖10 第28s內相位變化

圖11為某波長表面波的相位在一段時間內的變化，可以看出在幀數據與幀數據切換時相位會有小幅度跳變，這個跳變使得相位差變化非常劇烈。

由第1節介紹的方法修正后的相位及相位差如圖12所示，與圖11相比較有明顯的改善。

(4) 表面波相速度計算

在計算得到相位數據后，可由公式(4)計算不同波長的表面波相速度。計算過程中對6臺CCD相機的結果進行了平均處理。

如第1節所述，表面波相速度的測量精度主要與角頻率的精度有關，而角頻率的精度與所選時間間隔有關。圖13為波長分別為2.08和1.52cm的表面波相速度測量精度與所選時間間隔的關系。

由圖13可以明顯看出，在選擇時間間隔Δt≥2.5s時，表面波相速度測量精度優于1cm/s。

(a)

(b)

(a)

(b)

圖14為對波長為4.57cm的表面波采用不同的時間間隔時計算得出的相速度結果。

由圖14可以看出，間隔為2.5s的計算結果要明顯優于間隔為0.5s的計算結果。

計算時采用時間間隔2.5s，經過平滑濾波后得到的波長分別為2.08cm、1.52cm和1.14cm的表面波的相速度如圖15(a)、(b)和(c)所示。

由圖15及圖14(b)可以看出，速度有4次明顯的變化(圖14(b)中A、B、C、D四點)。這是因為在實驗過程中在水槽中有內孤立波經過探測區域，內波的傳播引起了表面波相速度上的變化。圖15及圖14(b)中虛線表示在內波經過之前(前80s)表面波相速度均值。這4次表面波相速度的變化與實驗中觀察到的現象十分吻合，即內波首先逆風向經過探測區域(100s左右)，而后經過3次反射。逆風向經過探測區域時內波傳播使得表面波相速度下降(圖14(b)中A、C)，順風向經過探測區域時內波傳播使得表面波相速度增大(圖14(b)中B、D)。而且由圖15可以看出，雖然不同波長下相速度大小不同，但在內波經過時引起的相速度變化幅值是基本一致的，約為3～4cm/s，整體幅值變化約為6～7cm/s，這種表面波相速度上的變化就是內波傳播引起的表面弱流場導致的。

(a) 波長為2.08cm

(b) 波長為1.52cm

(a) Δt=0.5s

(b) Δt=2.5s

(a) 2.08cm

(b) 1.52cm

(c) 1.14cm

圖16為不同風速下內波經過探測區域前(前80s)不同波長表面波相速度對比。

可以看出，對于波長相同的表面波，風速越大測得的相速度越大。這是由于風速越大，風生流越大，根據公式(1)，表面波相速度越大。同時，可以看出隨著風速增大，相速度曲線在圖中近似于“平移”，這也印證了波的本征相速度與風速無關，亦即與波高無關。

得到如圖15所示的表面波相速度后，將此相速度與前80s(內波經過前)的相速度均值相減，便可得到內波傳播所引起的表面弱流場。由于可以計算不同波長的表面波相速度，將多個波長表面波計算得出的表面弱流場進行平均處理，得到結果如圖17所示。

圖16 不同風速下不同波長表面波相速度(Vwind為風速)

圖17 內波激發的表面弱流場速度(正值表示內波與表面波運動方向相同)

由圖17可以看出，在前80s左右，內波尚未經過觀測區域時，流場速度在零值上下波動且變化幅度很小。對前80s的數據進行統計，其精度達到0.3cm/s，足以應用到對波-流相互作用機理的研究中。

3 結 論

提出了一種基于線陣CCD的水面流場光學測量方法，其最大的優點是可以在風浪環境下高精度地獲取表面弱流場信息而不受波浪振動的影響，克服了傳統PIV流場測量方法不能工作于風浪環境下的缺點。應用此方法對水槽內波觀測實驗數據的處理結果表明，表面弱流場的測量精度達到0.3cm/s。已有學者利用線陣CCD測量表面波浪譜等信息[4-7]，而本文在此基礎上進一步同時測量表面彈流場速度，并且保持實驗中流場方向基本與水槽長度方向平行，因此可以同時獲得風浪環境下的表面弱流場與波浪信息，這在波-流相互作用機理的研究中將會發揮重要作用。

所采用的實驗設備及方法也存在一定的局限性。實驗中對光照條件要求較高，這在實驗室環境中容易實現，但實際測量時測量條件較為惡劣，這對該方法的實用性會造成不利影響。因此今后的研究內容主要是改進該方法使其適用于實際測量以及討論該方法應用于復雜的隨機波浪表面測量可能面臨的問題和對測量精度造成的影響。

參考文獻:

[1]何宜軍. 成像雷達海浪成像機制[J]. 中國科學: D輯, 2000, 30(5): 554-560.

He Yijun. Imaging rader imaging mechanism of ocean waves[J]. Science in Chian, 2000, 30(5): 554-560.

[2]王浩, 曾理江. 二維及三維流場的光學測量方法[J]. 光學技術, 2001, 27(2): 139-142.

Wang Hao, Zeng Lijiang. Optical methods for measuring 2D and 3D flow field[J]. Optical Technology, 2001, 27(2): 139-142.

[3]Willert C E, Gharib M. Digital particle image velocimetry[J]. Experiments in Fluids, 1991, 10(4): 181-193.

[4]Titov V I, Bakhanova V V, Kemarskaja O N, et al. Investigation of sea roughness with complex of optical devices[C]//Proceeding SPIE 7473, Remote Sensing of the Ocean, Sea Ice and Large Water Regions, 2009, 74730T.

[5]Titov V I, Zuikova E M, Luchinin A G, et al. Investigation of surface roughness with optical methods[C]//Proceeding SPIE 7825, Remote swnsing of the Ocear, Sea Ice and Large Water Regions, 2010, 78250G.

[6]Titov V, Bakhanov V, Zuikova E, et al. Remote sensing of water basins using optical range-time images of water surface[C]//Proceeding SPIE 8175, Remote Sensing of the Sea Ice, Coastal Water and Large Water Regions, 2011, 81751D.

[7]Titov V I, Bakhanov V V, Luchinin A G, et al. Remote sensing of sea surface features by optical RTI images[C]//Proceeding SPIE 8888, Remote Sensing ot the Ocean, Sea Ice, Coastal Water and Large Water Regions, 2013, 88880J.

[8]Stockdon H F, Holman R A. Estimation of wave phase speed and nearshore bathymetry from video imagery[J]. Journal of Geophysical Research: Oceans (1978-2012), 2000, 105(C9): 22015-22033.

[9]J?hne B, Klinke J, WAAS S. Imaging of short ocean wind waves: a critical theoretical review[J]. JOSA A, 1994, 11(8): 2197-2209.

[10] Stilwelljr D, Pilon R O. Directional spectra of surface waves from photographs[J]. Journal of Geophysical Research, 1974, 79(9): 1277-1284.

[11] Young I R, Rosenthal W, Ziemer F. A three-dimensional analysis of marine radar images for the determination of ocean wave directionality and surface currents[J]. Journal of Geophysical Research: Oceans (1978-2012), 1985, 90(C1): 1049-1059.

[12] Rozenberg A D, Ritter M J, Melville W K, et al. Free and bound capillary waves as microwave scatterers: Laboratory studies[J]. IEEE Transactions on Geoscience and Remote Sensing, 1999, 37(2): 1052-1065.

作者簡介:

劉利(1989-)，男，山東泰安人，中國科學院電子學研究所碩士研究生。研究方向：信號與信息處理。通訊地址：北京市海淀區北四環西路19號中國科學院電子學研究所一室(100190)。E-mail：kelly07_11@163.com