基于海洋調查實測資料的中尺度渦旋識別結果的驗證及邊界擬合技術

彭漢幫，王金虎，吳克亮，楊春輝

( 1.解放軍92192部隊，浙江 寧波 315122)

1 引言

中尺度渦旋是大洋中廣泛存在的一種海洋現象，海洋中的中尺度渦旋對海洋中的水文物理性質有著重要影響[1–2]，在海洋中扮演著重要的角色。西北太平洋海域每天分布著幾十個中尺度渦旋，這些中尺度渦旋水平尺度為幾十千米到百千米量級[3–4]。無論從海洋動力學等理論研究方面還是從國防軍事等實際應用方面，西北太平洋中尺度渦旋都具有很高的研究價值。渦旋的識別是中尺度渦旋研究的重要內容之一，也因此成為國內外學者研究的重點[3,5]。

目前中尺度渦旋的識別方法可主要歸為4類：基于物理參數的 O?W（Okubo–Weiss）法、基于流動幾何學的WA（Winding-Angle）法、海表面高度異常（Sea Level Anomaly，SLA）等值線法和O?W與SLA綜合識別法（下文簡稱綜合法）。O?W法通過定義W參數（Okubo–Weiss參數）對渦旋進行判斷，在進行渦旋識別時，通常采用W的一個閾值圈定渦旋的形態和位置[6–7]。WA法是一種基于流場中瞬時流線的幾何狀態識別中尺度渦的渦旋自動識別算法，其主要思想是通過計算并挑選閉合的流線來達到識別渦旋的目的，該算法由Sadarjoen和Post[8]首次提出。SLA法根據SLA的絕對值大于某個設定的閾值來界定冷暖渦旋的形態[9–10]。綜合法首先使用O-W方法識別出渦旋存在的區域，再利用SLA的極值找出渦旋的中心，最后利用SLA等值線找出渦旋的邊界[11]。

為了定量評判中尺度渦旋識別方法的優劣，Chaigneau等[10]定義了識別準確率（the Success of Detection Rate，SDR）和識別冗余率（the Excess of Detection Rate，EDR）兩個參數。SDR表征某種渦旋識別算法對海區渦旋的成功檢測率，EDR表征該算法檢測渦旋的冗余比例。因此SDR越高，EDR越小，說明該算法對渦旋識別的準確性越高。Yi等[11]利用SDR和EDR分別對O?W法、WA法和綜合法進行檢測，結果顯示綜合法SDR最高且EDR最小。綜合法基于海區的SLA資料對渦旋進行識別，渦旋邊界用SLA的閉合等值線表示、渦旋中心點用SLA>區域極值點表示[11]。然而SLA只是中尺度渦旋引起的海表面高度異常，并不能反映渦旋的全部性質。因此，綜合法識別出的渦旋邊界形態是否代表真實渦旋的水平形態、識別的中心點與實際的渦旋中心點是否一致等仍是研究者需要解決的問題。目前包括綜合法在內的各類渦旋識別方法，其識別結果均缺少實測資料的檢驗。因此本文利用渦旋的實測資料，通過診斷渦旋的邊界和中心點來診斷綜合法的識別效果。

得益于中尺度渦旋識別方法的發展和衛星觀測資料的日益增多，一大批針對海洋中尺度渦旋的統計特征分析也相繼出現[12–15]。其中在對中尺度渦旋的水平尺度進行統計分析時，Chen等[13]采用圓半徑表征渦旋的水平尺度；Wang等[12]和楊光[14]將渦旋邊界東西向最長距離與南北向最長距離平均的一半作為渦旋的半徑；江偉等[15]將渦旋邊界各個方向離中尺度渦旋核心位置的距離平均值作為中尺度渦的半徑。由此可見，為了描述中尺度渦旋，多數學者利用圓來表征渦旋的水平形態，并利用圓半徑表征渦旋的水平尺度。然而江偉等[15]指出，實際海洋中水平形態接近圓形的中尺度渦旋占極少數，更多的是接近橢圓形，但對于該說法并未給出相關的證明。為了評判圓和橢圓對渦旋水平形態表征的優劣，本文基于綜合識別法，分別利用圓和橢圓對其識別的邊界進行擬合，通過設置一定的比對參數評判二者的擬合效果。

文章的第二部分介紹本文使用的數據和渦旋綜合識別法；第三部分為利用實測數據診斷識別渦旋的邊界、中心點；第四部分為渦旋邊界擬合效果的比對；最后是總結和討論。

2 數據和方法介紹

2.1 數據介紹

本文所用的海表面高度計資料來自于哥白尼海洋環境監測局（Copernicus Marine and Environment Monitoring Service，CMEMS）。該衛星觀測資料由 T/P、ERS、Sentinel-3A、Jason-3、HY-2A、Saral/Altika 等多顆衛星高分辨率海面高度異常數據整合而成，可以準實時地為世界上主要的海氣預報中心提供服務。其時間分辨率為日平均，即每天提供1個全場衛星觀測數據，空間分辨率為（1/4）°×（1/4）°。

渦旋實測數據來源于2019年一次遠海綜合調查測量航次中對一個冷渦進行的觀測，其中心位置位

于 18.2°N，128.5°E 附近（圖1）。走航觀 測航線 為A→B→C→D→E→F→G→H→I→J→A，并在其中兩條斷面上開展13個大面站觀測。走航時采用300 kHz ADCP對上80 m層海流進行觀測，最上層觀測深度為10 m，垂向采樣間隔為2 m；采用38 kHz ADCP對上900 m層海流進行觀測，最上層觀測深度為38 m，垂向采樣間隔為16 m；選取300 kHz ADCP的80 m以淺和38 kHz ADCP的80 m以深的海流數據用于本文的研究。對航線中橫跨渦旋東西方向的FG斷面、 橫跨南北方向的AJ斷面，利用下放CTD對溫鹽剖面進行測量，測量深度為2 000 m；溫鹽觀測站點中，W1至W7相鄰站點距離約為33 km，W8至W13相鄰站點距離為28～75 km。

2.2 渦旋識別綜合法

本文借鑒Yi等[11]發展的渦旋識別綜合法，綜合法按照確定渦旋位置→確定渦旋中心點→確定渦旋邊界的先后順序，對海洋中的中尺度渦旋進行識別。具體過程如下：（1）利用O?W參數W確定渦旋的位置，W

2.3 圓與橢圓擬合法

本文擬合的對象是閉合的渦旋邊界。采用的圓擬合法首先分別以閉合邊界上的橫坐標和縱坐標平均值作為圓心的橫坐標和縱坐標，再取閉合邊界上的點到圓心的平均距離作為圓半徑，最后繪制圓形如圖2a所示。

圖2 渦旋邊界的圓擬合（a）和橢圓擬合（b）結果Fig.2 Circle (a) and ellipse (b) fitting results of the eddy boundary

本文使用的橢圓擬合法來自Gander等[16]開發的最小二乘法橢圓擬合法。最小二乘法橢圓擬合是較常用的橢圓擬合方法[17–18]，其主要技術是尋找參數集合，從而將數據點與橢圓邊界之間的距離度量最小化，常見的距離度量有幾何距離和代數距離[19]。無論是幾何距離還是代數距離，其目的就是求如橢圓一般方程式（1）的解，從而根據方程系數A、B、C、D、E和F的值確定渦旋長軸、短軸、焦點和傾斜角等參數。根據獲取的渦旋參數，繪制橢圓如圖2b所示。

3 綜合法渦旋識別結果的診斷

3.1 渦旋識別邊界的診斷

本文利用實測流場診斷渦旋識別的邊界。由2.2節內容可知，綜合法識別出的渦旋邊界是SLA閉合等值線，因此欲診斷識別的渦旋邊界，應選取上層的實測海流數據。本文選取200 m以淺中15 m、50 m、100 m和200 m為代表的標準層實測海流數據，用于診斷渦旋邊界形態。如圖3所示，實測流場為15 m、50 m、100 m和200 m層平均流場，紅實線為識別渦旋的邊界。從圖3可以看出，該冷渦主體部分（17°N以北）大致呈傾斜橢圓狀，實測流場圍繞該橢圓逆時針旋轉，與識別的渦旋邊界形態較為吻合。因此，識別的渦旋邊界形態基本可以反映實測渦旋水平形態。

圖3 渦旋邊界（紅線）與實測流場（流場為15 m、50 m、100 m和200 m層平均流場）Fig.3 The eddy boundary (red line) and observed currents average of 15 m, 50 m, 100 m and 200 m

為進一步診斷識別的渦旋邊界形態與實測渦旋水平形態的吻合情況，如圖4所示，分別將15 m、50 m、100 m和200 m層流場與渦旋邊界投影到一個平面上，從各層流場中選取位于渦旋邊界附近的實測海流點（為保證診斷數據的質量，實測海流點與渦旋邊界點間的水平距離控制在2 km范圍內）。通過計算每個實測海流點的流矢量與渦旋邊界切線的夾角α診斷識別渦旋的邊界形態與實測渦旋的形態吻合度，α越小說明識別的渦旋邊界形態越能反映實測渦旋的水平形態。結果顯示15 m層α在0.6°～30°范圍內，平均值約為 12.8°；50 m 層α在 0.4°～32°范圍內，平均值約為 13.8°；100 m 層α在 1.3°～36°范圍內，平均值約為 14.6°；200 m 層α在 2.2°～38°范圍內，平均值約為14.9°。因此，識別的渦旋邊界形態基本可以反映實測渦旋水平形態。

圖4 渦旋邊界及其附近15 m (a)、50 m (b)、100 m (c)和200 m (d)層實測流場Fig.4 The eddy boundary and the nearby observed currents of 15 m (a), 50 m (b), 100 m (c) and 200 m (d)

3.2 渦旋中心位置的診斷

本文利用實測海流和溫度數據反演渦旋的中心位置，并用來診斷識別的渦旋中心位置。由圖1可知，航次任務對渦旋進行觀測時，根據綜合法識別結果，選擇分別位于中心點的北面和南面兩個緯向斷面（FG和HI斷面）進行走航式海流觀測。由圖3可知，該冷渦主體部分大致呈橢圓形，且長軸呈西北–東南走向。由于橢圓長軸兩端曲率大，因而渦旋在此處流向變化大和流速小。根據此特點，利用FG和HI斷面的海流數據（200 m以淺平均海流）可以反演出該橢圓渦旋的長軸線（本文記為渦旋中心線，如圖5a所示），并初步判斷實際渦旋中心位置大致位于該中心線上。同時，根據冷渦中心溫度較低的特點，通過計算渦旋內每個大面觀測站剖面200 m以淺平均溫度，可以進一步縮小渦旋中心位置的分布區域（本文記為中心區域），如圖5b所示。由圖5b可知在該冷渦區域內所有觀測站點200 m以淺平均溫度中，W4、W5和W10這3個站點是溫度最低的3個點。同時從圖5c溫度剖面圖可以看出，W4、W5和W10站點的溫度剖面處于低溫一側。由于溫度由冷渦中心向外呈增加趨勢，可以判斷渦旋中心不會出現在站點 W4至W3和W11、W5至W6以及W10至 W9區域內，而大概率出現在由站點W4、W5和W10圍成的區域內。根據以上推斷，利用經過同一平面內不在同一條直線上的3點可確定一個圓的定理，可基于W4、W5和W10 3個站點反演出冷渦的中心區域（如圖5b中的紅圓所示）。

結合實測資料反演的渦旋中心線和中心區域，通過求二者的“交集”可以進一步縮小實際渦旋中心的位置范圍。如圖5d所示，渦旋中心線與中心區域圓相交于弦AB，故可將渦旋中心位置鎖定在弦AB上。通過距離計算，識別中心點C至弦AB的最小距離為C到D的距離（約8 km），最大距離為C到B的距離（約25 km）。因此，綜合法識別的中心點與實測資料反演的渦旋中心點距離范圍在8～25 km之間，由于本文用到的衛星高度計資料空間分辨率為（1/4）°×（1/4）°（約28 km），因此可以說兩中心點位置較為吻合。

4 渦旋邊界擬合的準確率

本文從兩方面比對圓與橢圓對識別渦旋邊界的擬合效果。第一個是識別邊界區域和擬合邊界區域重疊和不重疊部分的面積占識別邊界區域面積的比例，即重疊率（Overlap Rate, OR）和冗余率（Redundancy Rate, RR），表達式分別為式（2）和式（3）；第二個是擬合結果與SLA反演的地轉流吻合度，即地轉流矢量與擬合邊界切線的夾角α大小。顯然，OR越大，RR、α越小，擬合效果越好。

圖5 基于實測數據的渦旋中心位置診斷過程Fig.5 Diagnosis processes of the eddy center based on observed data

為敘述方便，分別將渦旋識別邊界、圓擬合結果和橢圓擬合結果圍成的區域簡稱邊界區域、圓區域和橢圓區域。式中，A為邊界區域面積；AO為圓區域或橢圓區域與邊界區域重疊部分的面積；AR為圓區域或橢圓區域中與邊界區域不重疊部分的面積，如圖6所示。

圖7a為西北太平洋海域某個渦旋分別進行圓和橢圓擬合的結果。通過計算可知，該渦旋邊界區域的面積A為6.8（本文近似的將地理坐標中1°×1°圍成的面積記為1），圓區域中與其重疊的面積AO為5.3、與其不重疊的面積AR為2.5，因此該渦旋邊界圓擬合的結果OR為0.78、RR為0.36；橢圓區域中與邊界區域重疊的面積AO為6.5，與其不重疊的面積AR為0.6，因此橢圓擬合的結果OR為0.95、RR為0.09。

圖6 邊界區域與圓區域疊加Fig.6 The eddy boundary overlaps with the circle

與3.1節中利用實測海流點診斷渦旋邊界方法相似，從地轉流中選取與擬合邊界距離最近的地轉流點（如圖7b所示），通過計算地轉流矢量與擬合邊界切線的夾角α診斷邊界擬合效果。結果顯示，圓與橢圓平均α分別約為33°和10°。因此，就圖7a所示的渦旋識別的邊界而言，綜合考慮OR、RR和α，橢圓擬合效果比圓更佳。

為進一步比對圓與橢圓對渦旋邊界的擬合效果，本文隨機選取西北太平洋某區域，對該區域進行渦旋識別并分別進行圓與橢圓擬合（如圖8a所示），計算每個渦旋的OR和RR。同時從地轉流中選取與擬合邊界距離最近的地轉流點（如圖8b所示），計算每個渦旋的平均α。OR、RR和α結果如表1所示，從表中可知： （1）該區域9個渦旋中，有8個渦旋橢圓擬合OR大于圓擬合OR、橢圓擬合RR有7個小于圓擬合RR；（2）圓擬合 OR平均為 0.83、RR平均為 0.26，橢圓擬合OR平均為0.88、RR平均為0.17；（3）有9個渦旋橢圓擬合α小于圓擬合α，平均小8.3°；（4）橢圓擬合OR方差與圓擬合OR方差相等，橢圓擬合RR、α方差分別小于圓擬合RR、α方差。因此，綜合考慮OR、RR和α的結果，可以得出的結論是利用橢圓比利用圓擬合渦旋邊界效果更佳。

5 結論

本文在渦旋OW與SLA綜合識別法的基礎上，利用實測資料對識別的渦旋形態和中心點進行比對驗證。結果表明綜合法識別的渦旋與實測結果吻合度較高。

圖7 單個渦旋邊界擬合結果（a）（矢量箭頭為地轉流）及擬合邊界附近地轉流分布（b）Fig.7 Fitting results of the eddy boundary (a) (vectors indicate the geostrophic current) and the distribution of geostrophic current nearby the fitting boundary (b)

圖8 多個渦旋擬合結果（a）及其擬合邊界附近地轉流分布（b）（Nxx為渦旋編號）Fig.8 Fitting results of several eddy boundaries (a) and the distribution of geostrophic current nearby the fitting boundaries (b) (Nxx denote the eddy number)

表1 N01–N09號渦旋圓與橢圓擬合結果（OR、RR和α）Table 1 The results (OR, RR and α) of circle and ellipse fitting for number N01 to N09 eddies

本文通過診斷渦旋邊界附近的實測海流矢量與渦旋邊界夾角α來評判渦旋邊界與實測渦旋水平形態的吻合度。通過計算15 m、50 m、100 m和200 m層實測海流矢量與渦旋邊界的夾角α，結果顯示α平均在12.8°～14.9°范圍內，因此可以判斷綜合法識別的邊界形態能夠較好地反映實際渦旋的水平形態。

本文利用實測海流和溫度剖面資料反演實際渦旋中心位置，并與識別的中心位置進行比對，結果顯示渦旋識別的中心位置與實測渦旋中心位置基本吻合。比對過程如下：（1）根據渦旋產生的旋轉流，利用渦旋中心線處流向變化大、流速小的特點確定渦旋中心線；（2）利用冷渦中心溫度較低的特性，分別對渦旋內每個溫鹽觀測站點200 m以淺溫度剖面求平均值，選取平均溫度較低的3個站位點并以此圈定渦旋中心區域；（3）渦旋中心線和中心區域共同確定的渦旋中心與渦旋識別的中心點距離范圍在8～25 km，小于衛星高度計資料的空間分辨率（28 km）。

本文隨機選取西北太平洋某個海域，比對圓與橢圓對識別渦旋的邊界擬合效果。利用擬合參數面積重疊率、面積冗余率和地轉流矢量與擬合邊界夾角診斷擬合效果，結果顯示橢圓OR比圓OR平均大0.05、橢圓RR比圓RR平均小0.09、橢圓α比圓α平均小8.3°。可以看出橢圓擬合效果要明顯高于圓擬合效果。

本文在研究中存在以下幾點不足：（1）僅利用1個渦旋的實測資料進行識別結果的診斷，得出的結論不具備足夠的說服力；（2）由于渦旋的觀測時間跨度為3 d，而渦旋識別利用的SLA數據的日期為渦旋觀測期間的某一日。考慮到西北太平洋中的渦旋在以每天約10 km的量級向西傳播[3]，本文在利用實測數據對識別渦旋邊界和中心點進行比對驗證時，會存在一定的誤差。下一步將結合更多的渦旋同步觀測資料進一步驗證渦旋識別結果，并利用橢圓識別開展相關的研究。