一種基于網格重心的等值線追蹤算法

王俊驕 陳 晴 張姍姍

(1.杭州市氣象局，浙江 杭州 310008;2.浙江省氣象信息網絡中心，浙江 杭州 310017)

0 引 言

在氣象要素圖形顯示中經常用到等值線圖，它能較直觀地給出氣象要素的分布情況。不管是網格化的數值預報產品分析，還是在離散站點數據網格化應用中，等值線圖的應用都非常廣泛。等值線圖是以相等數值點的連線表示連續分布且逐漸變化數量特征的一種圖型，等值線生成常用算法有標準網格法和三角網格法兩種[1]。而標準網格法利用區域內的高程離散點,建立逼近的曲面模型,然后通過內插得到一組規則的矩形網格數據,再在矩形網格的基礎上線性插值追蹤得到各條等值線。

目前，國內對標準網格等值線追蹤的算法很多，如按一定方向直接對每個網格各邊的等值點進行等值線追蹤[2]、通過對角線轉軸法進行追蹤[3]、利用數據關聯表生成矩形網格等值線[4]等等，對追蹤出來的等值線再用一些插值擬合方法進行平滑處理,如Bezier、B樣條擬合等[5]，或為了使等值線圖更直觀,還需要對等值線進行顏色填充[7]。

以上等值線追蹤算法的設計均旨在提高整體的搜索效率，提升速度，在搜索過程中只考慮小網格的進邊和出邊，未分析等值線在小網格內的走向，在小網格尺寸很大時尤其明顯。因此，本文提出一種基于網格重心的等值線追蹤算法，在等值線追蹤過程中，引入重心點分析每個網格，通過增加等值點序列中點的數量，不僅可更真實地反應等值線整體逐漸變化的特征，而且為后續的等值線光滑提供更多的控制點，同時由于本文算法在判斷出邊時僅需通過重心判斷，避免了網格出口邊判斷復雜的問題。

1 等值線的原理及算法

對于標準網格的等值線生成，首先需獲得網格數據，在氣象應用中，網格數據一般通過內插方法將離散數據轉換成網格數據，也可從氣象應用中常用的Micaps四類格式、Grib格式和Netcdf文件獲取。

在整個大網格中，分析每個小網格獲取等值點，再判斷走向，將它們依次用線連接起來，最后光滑等值線，從而構成等值線。本文重點討論分析網格數據等值線追蹤[8]。

1.1 等值點的判斷

在幾何中，對于由點A和點B構成的一條邊AB，其中XA，YA分別表示點A的橫坐標和縱坐標，XB，YB分別表示點B的橫坐標和縱坐標，VA，VB分別代表點A和點B的高程值。判斷在該邊上是否存在滿足高程值為W的等值點R的條件如下：

(VA-W)(VB-W)<0

(1)

若邊AB滿足(1)式，則通過內插可知等值點R的橫坐標XR和縱坐標YR：

(2)

假設規則網格是由m×n個網格數據點構成,用Di,j表示縱向第i行，橫向第j列的網格點，則用Vi,j表示該點的高程值，Di,jDi,j+1表示由點Di,j和Di,j+1所構成的橫向邊，Di,jDi+1,j表示由點Di,j和Di+1,j所構成的縱向邊，其中1≤i≤m，1≤j≤n。

定義1對于任意一條邊Di,jDi+1,j，若滿足(1)式，則邊Di,jDi+1,j定義為等值邊，而等值點R的橫向位置x和縱向位置y：

(3)

同理對于任意一條橫邊Di,jDi,j+1為等值邊，可得到等值點相應的位置。

定理1由等值點的滿足條件(1)式可知，在任意一個小網格中，若某邊為等值邊，即存在等值點，則該網格的另外3條邊中至少有1條邊同時為等值邊。

另外，在網格數據等值線追蹤過程中，須對每次搜索到的等值點進行判別，分析該點是否為終點，以停止搜索。而根據等值線特性可知，在搜索某高程值等值線時，終點需滿足以下條件：該點為邊界上的點或者該點為搜索起始點。

1.2 傳統網格追蹤

一條等值線的傳統追蹤算法描述如下：首先，依照自上而下，自左向右的順序搜索每個小網格，將第一個找到等值邊上的等值點作為搜索起始點，以該等值邊為進邊分析該網格，計算該網格上另外一個等值點和相應的等值邊，連接兩個等值點，即完成第一個網格的追蹤。然后，以新的等值邊作為進邊，得到下一個網格，依此類推，直至等值點為終點。

然而，在傳統追蹤過程中，當網格內不止存在一條出邊時(見圖1)，就需判斷等值線的走向。不僅降低了追蹤的效率，而且可能由于不同的判斷方法導致產生不同的等值線。當網格尺寸較大時，雖然能加快等值線搜索速度，但無法更真實地反映網格內部和整體的變化特征。

圖1 傳統追蹤示意圖(括號中代表該點高程值)

2 本文算法

2.1 基于小網格重心的追蹤

本文在傳統的追蹤算法基礎上，提出利用網格的重心來控制等值線的走向，通過重心引導，避免了等值線走向問題，同時使等值線在網格內的走向更符合等值線的實際特征(見圖2中b，c，d)。

圖2 等值點追蹤時小網格內部的走向

在一個小網格中，已知四個點A,B,C,D，括號內為該點的高程值。則重心點E的橫坐標Ex、縱坐標Ey和高程值VE為：

基于小網格重心的追蹤算法可描述為：當某一邊進入小網格時，分析以該網格重心與進入邊構成的三角形。由于在該三角形中另外兩條邊中有且僅有一條為等值邊，則找到該等值邊及相應的等值點。繼續分析下一個三角形，直到找到的等值點在該網格的四條邊上，即結束該小網格的搜索。一般該追蹤算法可出現圖2中b，c，d 3種情況。

以圖2中的c為例，追蹤高程值為5的等值線。首先得到網格ABCD的重心E，以邊AB作為進入邊，則分析三角形ABE得到新的等值邊AE，由于AE不是網格的邊，則繼續分析三角形AED，CDE，直到找到網格的出邊CD，才完成整個小網格的追蹤。在該網格的追蹤中，共增加了2個新點，將這2個點存入到整條等值線點的序列中，可為后續的光滑，增加新控制點。

2.2 算法

等值線追蹤到的所有等值點坐標按順序連接就形成了等值線。然而該線條為一條折線，故需對這些等值點進行曲線平滑處理，將處理后的坐標連接起來即可繪出連續平滑曲線。曲線平滑常用的方法有Bezier方法、B樣條方法、三次樣條方法及最小二乘法，同時對等值線進行填充顏色，邊界裁剪等后續加工處理。

為避免在追蹤過程中，對某些等值點重復追蹤，造成算法冗余，在追蹤順序時，本文先依序搜索網格最外面4條邊，根據等值線要么自封閉,要么邊界封閉的特性，先將邊界封閉的等值線追蹤完畢，再追蹤自封閉的等值線。因此，本文的算法流程圖如下。

圖3 本文算法的流程框圖

3 本文算法的應用

為驗證本文算法的可用性，選取2016年4月5日08時歐洲中心細網格的小時降水預報，在數值預報中繪制高程值為8的等值線，為能清楚地展示網格內部等值線追蹤的效果，截取了中心10×10網格數據,用于展示等值線追逐效果，并通過傳統和本文的算法進行比較(結果見圖4和圖5)。

圖4 傳統追蹤的示范圖

圖5 本文基于重心追蹤的示范圖

由圖4和圖5可知：采用傳統的追蹤算法比基于重心追蹤的等值點少，所以等值線會相對比較粗糙，基于重心的追蹤法更能反映等值線在小網格內的走向。傳統的追蹤方法因為需對網格出口邊判斷，從而導致整體生成的等值線會不同，而等于重心追蹤在小網格中不涉及出口邊判斷問題，確保了等值線的唯一性，可以應用到實際業務中。

4 結 語

本文在分析傳統矩形網格等值線追蹤算法的基礎上,引入小網格重心的算法加以改進,實現尋點與追蹤同步。與傳統方法比較,本文的方法有以下改進：1)可更真實地反應等值線整體逐漸變化的特征；2)在等值線點序列中增加了節點，為后續的等值線光滑提供更多的控制點；3)在追蹤過程中，不涉及出口邊判斷問題，確保了等值線的唯一性；4)可利用同樣原理擴展到三角網格的等值線追蹤中。

本文算法通過編碼，已成功應用到氣象數值預報產品、Micaps四類格式的資料處理中，并應用到實際業務中。

[1] 林毅,金燁，馬登哲等.正規網格等值線的虛路徑掃描算法[J].計算機工程及應用，2001,37(13):92-94.

[2] 余衛東,李湘閣,王靖.氣象場等值線自動繪制[J].氣象科技，2002,30(4):222-225.

[3] Dayhoff M O. A contour map program for X-ray crystallogra-phy[J]. Co mmunications of the Association for Computing Machinery,1963,6(10):620-62.

[4] 孫科峰,孫根正,李潔.一種新的矩形網格生成等值線算法[J].東華大學學報(自然科學版)，2005,31(4):66-69.

[5] 韓麗娜,張紅祥,張群會.Bezier曲線修改的一種分割算法[J].計算機工程與科學，2006(7):77-79.

[6] 韓家新,王家華.一種等值線填充的連通區域搜索算法[J].計算機應用與軟件，2003(10)：5-6.

[7] 陳正旭,封秀燕,王亞云.多島嶼地圖上繪制氣象要素等值線色塊的自適應方法[J].氣象科技，2009,37(3):356-359.

[8] 王亞強.等值線相關算法類庫的開發與應用[J].氣象科技，2010(4):478-482.

[9] 張紅祥,車鵬飛.一種多段Bezier曲線光順擬合方法[J].科技信息,2008(13):48-49.

[10] 孫桂茹,馬亮,路登平等.等值線生成與圖形填充算法[J].天津大學學報,2000,33(6):816-818.