基于Python和GIS信息的南海氣象雷達專題圖自動制作

韋凱華， 梁偉杰,2， 韋海寧， 袁宇東,2

(1.廣東省氣象臺，廣東 廣州 510080；2.中國南海石油聯合服務總公司氣象服務公司，廣東 廣州 510080；3.廣西水利電力勘測設計研究院有限責任公司，廣西 南寧 530023)

0 引言

南海海域是我國重要的戰略縱深，地理位置十分重要。隨著近年來海洋主權意識的不斷增強，在南海進行生產作業和駐守的軍民也在增加。受地理環境和氣象環境條件影響，南海災害性天氣頻發，海洋氣象條件復雜，季風和熱帶云團導致的強對流天氣以及路徑曲折、強度多變的臺風[1-3]，使得在南海航行及進行漁業作業的船舶常受惡劣天氣的威脅[4]。天氣雷達作為監測近海臺風及強對流天氣的有效工具，其相關產品對預報員分析判斷未來天氣動向，以及為公眾和行業提供防御決策指導至關重要。但現行的雷達產品圖在精細化地理信息(如南海島礁)的繪制方面的不足限制了其使用價值的進一步挖掘。本文在此前工作的基礎上[5-7]利用Python編程語言及其相關類庫重新制作南海海域氣象雷達可視化專題圖，并實現自動化生成，以期進一步推進相關資料應用，充分發揮氣象防災減災第一道防線的作用。

Python作為一種跨平臺的計算機程序設計語言，結合了解釋性、編譯性、互動性等特點，并且是一種面向對象的高級語言。其設計之初定位自動化腳本(Shell)語言的編寫，隨著相關類庫的不斷更新和發展，逐漸被用于大型項目的開發。其提供的便捷的數據分析、擬合等功能，并被應用于金融分析、數據挖掘等多個領域[8]。本研究選NumPy、Pandas進行數據解碼，利用強大的繪圖類庫Matplotlib完成專題圖繪制。

1 數據預處理

1.1 氣象雷達基數據資料的讀取

南海海域以及周邊的海南、廣東、廣西三省沿海共布設有S波段(SA、SC)氣象雷達12部，不同型號的雷達其基數據格式不盡相同。中國氣象局于2020年度完成了氣象雷達數據格式標準化工作，實現標準格式單軌運行。本文選取位于南海中部三沙市永興島的西沙SC型雷達作為代表，一方面考慮西沙附近島礁密集，對于精確地理信息與氣象信息結合需求相對較大；另一方面，SC型雷達在上述業務變更前后，基數據格式的改動較大。標準格式雷達基數據的數據結構如圖1所示。

圖1 標準格式基數據整體結構

其中，N表示第N個仰角；M表示第M個徑向；K表示第K個數據類型；I表示徑向數據長度，數據讀取過程中涉及到上述多個結構字段的循環。本文首先選用NumPy類庫中的numpy.dtype函數構建類似結構體的字典類型。并采用NumPy類庫的frombuffer函數循環遍歷讀取雷達基數據的各個數據塊并經解碼后存入對應numpy.ndarray類型的三維數組(仰角、方位角和徑向)中。其關鍵代碼如圖2所示。

圖2 標準格式基數據讀取及解碼關鍵代碼

由于Numpy底層采用Fortran和C++ BOOST進行數據處理和運算的工作，本文采用numpy.frombuffer函數可以在確保運行效率的同時充分發揮Python易讀性強的特點。數據讀取和解碼完成后，根據各個徑向的方位角信息，以方位角升序的方式對變量三維數組進行排序，形成排序后的反射率因子(ref)和多普勒速度(vel)三維數組數據，用于后續進一步繪圖或計算處理。

1.2 GIS信息的預制作

從國家地理信息共享平臺(全國地理信息資源目錄服務系統)獲取1：100萬全國地理信息基礎數據，由于獲取的是分幅數據，故需要通過ArcGIS軟件按照所需要的南海區域地理信息范圍合并相應的幾幅，再進行裁剪，最后得到一組以SouthChinaSea命名的地理信息數據。對于主要島礁名稱的文字標注，需要讀取其中的庫表關系獲取名稱和位置信息。為簡化步驟，提高讀取效率，用ArcGIS軟件中的要素轉點功能，將.shp文件中原有的面元素和線元素轉換為點元素，同時運用軟件中的Excel轉換工具，將已轉換為點元素的.shp文件中包含各元素的屬性表轉換為.xls文件并另存，每行中就包含了元素所屬類別、中英文名稱、經緯度坐標等。后續即可通過Python程序讀取該.xls文件實現地理信息的標注。

2 氣象雷達可視化專題圖自動制作程序設計

本研究利用Python語言設計程序實現氣象雷達可視化專題圖自動制作，程序框圖及數據處理流程如圖3所示。

圖3 程序框圖及數據處理流程

氣象雷達站觀測到的雷達回波信息數據上傳至國省兩級氣象信息中心數據庫后，再通過推送下發或數據接口訪問的方式將基數據提供給氣象臺等業務單位使用[9]。設計的程序主要分為3大功能模塊，即時間調度與實時資料檢索匹配模塊、雷達數據讀取計算處理模塊、雷達數據與GIS信息繪圖模塊。經程序處理后自動制作生成的氣象雷達可視化專題圖可推送到本地或網絡指定位置，供氣象部門內及部門外各渠道調用。

2.1 實時資料的檢索與匹配

無論通過數據接口訪問或是本地推送(或FTP共享)的方式獲取氣象雷達基數據，要實現實時的專題服務圖自動制作，首先需完成最新資料檢索下載。以本地(或FTP共享)方式訪問為例，讀取配置文件settings.ini中的資料路徑fileFolder，列出該目錄下后綴名為.bz2的文件。標準雷達基數據文件名的形式為Z_RADR_I_Ziiii_YYYYMMDDhhmmss_O_DOR_TT_CAP_FMT.bin.bz2，其中iiii為4位數字的雷達區站號；TT為雷達型號；YYYYMMDDhhmmss為基數據產生時間。本文利用正則表達式提取雷達站號和數據時間信息。通過提取后的雷達站號確定西沙雷達基數據存放的次級目錄，進入到目錄中后對列表中的所有文件根據文件名時間進行逆序排列，排序后的第一個文件名即為最新的資料。考慮到可能出現的雷達數據傳輸故障，設定為機器當前時間與最新雷達數據文件名時間的延遲小于15分鐘以內時，才進行相應的數據處理工作。具體代碼如圖4所示。

圖4 實時資料檢索與匹配關鍵代碼

2.2 氣象雷達基數據資料的可視化專題圖繪制

在讀取雷達基數據資料的基礎上進行可視化專題圖的繪制。自定義一個帶參數的draw_radar_ppi函數，傳入的參數為要繪制的專題圖數據類型。先利用循環按投影比例繪制雷達距離圈和方位角射線，此后根據傳入的參數確定繪制數據類型及其對應的色標(色標預先以數組形式定義在mycolors_變量名.npy文件中，通過自定義的get_cmap函數傳遞)。以雷達所在坐標經緯度為中心，按方位和徑向計算雷達極坐標投影至地理坐標的經緯度位置。本文充分利用NumPy的向量運算的特性，運用Pyproj包中的Geod函數，將極坐標通過WGS84坐標系投影轉換為等經緯網上的笛卡爾坐標，得到雷達數據在平面上的lon(經度)與lat(緯度)值，再利用pcolormesh函數按排序處理后的三維數組數據值對應的色標繪制圖形。pcolormesh函數可用于均勻/非均勻色標的繪制，通過對BoundaryNorm類型的指定實現，其構造函數為BoundaryNorm(boundaries, ncolors, clip)。各參數含義如下：levels為傳入的色標閾值；ncolors為色標的色階數；clip為控制邊界值是否截斷進行繪制，選用True表示超出閾值上下邊界的均截斷為閾值。

再確定子坐標軸，分別以文本形式繪制相關的雷達基本信息和觀測信息，以矩形填色方式繪制色標與閾值，并從基數據中提取時間信息，用datetime類完成儒略歷、地方時轉換，用作保存專題圖的文件名使用，完成整幅圖的基本繪制。關鍵代碼如圖5所示。

圖5 專題圖繪制關鍵代碼

2.3 GIS信息的繪制

運用Matplotlib類庫中的Basemap庫指定數據繪圖的投影方式以及經緯度范圍，并通過readshapefile函數依次讀入國省市及縣界的.shp地理信息，本設計中，使用的.shp地理信息基礎數據包括了南海所有島礁的信息以及南海九段線，符合地理信息規范化制圖的規定。drawmapboundary函數使用指定顏色繪制底圖，此處指定為海洋顏色，fillcontinents函數將在底圖上用指定顏色填充繪制陸地地形。此外，南海島礁名稱的標注，通過Pandas類庫中的read_excel函數，按指定格式讀取第2.2節中預先制作好的SouthChinaSea.xls中按行列存儲的地理信息，提取其中的名稱、經緯度值，最后通過循環的方式逐條通過text函數將島礁名稱按前述繪圖投影方式的經緯度以標注文字的形式繪制在底圖上。

3 專題圖效果展示與應用

原程序制作的專題圖如圖6(a)所示。設計的Python程序自動運行后制作的南海氣象雷達專題圖如圖6(b)、圖6(c)所示。其中圖6(b)為繪制的雷達反射率因子(ref)資料。相較于原有VB(Visual Basic)程序繪制的專題圖(圖6(a))，可以發現本文設計程序制圖的若干優點。其一是增加了精確的地理信息，三沙市范圍內西沙群島、中沙群島的島礁位置及名稱均清晰標注在圖內，同時在海南島范圍內將地理信息精確繪制至縣級(原來僅為市級)，并增加標注沿海重要港口及雷達探測距離圈范圍，幫助終端用戶能清楚獲悉惡劣天氣的影響區域、范圍和過去及未來一段時間的移動趨勢，特別是在臺風臨近登陸時，上述功能可以輔助確定精確到縣級的可能登陸地段，以便更有針對性地采取防御措施。其二是改進了雷達數據處理的邏輯，原有的VB專題圖制作程序設計直接采取徑向數據方位角取整的方法處理和繪制圖像，造成了個別徑向的缺失，如圖6(a)中150°-180°方位角的部分徑向缺失影響了雷達回波圖像的完整性。本文的程序由于在繪制圖像前先進行了排序，并考慮到可能缺失徑向的插補，避免了繪制缺失徑向問題的出現。其三是高精度圖像繪制，由于雷達進行了技術升級與數據格式標準化統一工作，分辨率從原來的1.2 km提升了4倍至300 m，相應程序也進行了修改，使得繪制的圖像更精細。同時本文的程序設計可以通過改變圖像dpi的方式調整輸出圖像的分辨率(默認dpi為100)。其四是在原有基礎上新增了徑向速度(vel)資料的專題圖，如圖6(c)所示。

圖6 各程序制作的專題圖

所選個例為2020年第22號臺風“環高”過程，通過制作的專題圖，預報員可以分析計算得到專題圖代表的雷達探測0.5°仰角上，該臺風底層中心附近的最大風速為51 m/s(15級上限、16級下限)，考慮到臺風中心距離雷達約150 km，根據測高公式算得回波離地高度約2.5 km，故結合預報經驗公式折算至地面估計的風應在14到15級(46-48 m/s)左右。本臺預報員最終制作發布的臺風路徑和強度在未來18 h內均維持在45-50 m/s，同時繼續發布了臺風緊急警報，指出“受冷空氣和臺風‘環高’的共同影響，南海中西部海面：旋轉風10到11級陣風13級，其中臺風中心經過的附近海面14級陣風15級”，預報結論與觀測實況較為一致，應用專題圖起到了較好的服務效果。

此外程序自動生成的專題圖可以按配置文件settings.ini中的目標路徑輸出到指定目錄下或FTP服務器，除可以生成網站及App端常用的.jpg或.png格式外，還可以生成矢量圖.eps及.svg格式文件，為研究應用印刷出版提供支持。

4 總結

本文基于Python語言，結合精細的GIS信息實現了南海氣象雷達專題服務圖的自動制作，支持雷達基數據自動檢索下載、解碼讀取、圖像投影繪制等功能。程序流程簡明清晰，從數據檢索到成圖分發耗時在5 s以內，滿足業務時效需求。制作的專題服務圖供氣象部門內預報分析網站及對內/外App等多渠道應用，特別是為近海臺風、強對流等天氣的氣象防災減災工作提供了有力的資料支撐。添加的精確地理信息對于天氣影響區域的研判準確性較此前大為提升，且可以通過修改相關預制GIS文件的方式增加或突出顯示，在后續也可以按需增加近海航線、漁場等信息[10]，為更具體的應用場景服務。此外，程序設計是基于標準格式的雷達基數據處理，可擴展應用到南海海域及周邊任意一部業務雷達進行專題服務圖的制作，具有一定的推廣應用價值。