下載及裁剪美國(guó)GFS全球數(shù)值天氣預(yù)報(bào)資料的方法

史小康，程文聰，忤建勛，李 森

(1.北京航空氣象研究所，北京100085；2.西北核技術(shù)研究院，烏魯木齊 841700)

0 引言

美國(guó)NCEP機(jī)構(gòu)制作的GFS全球數(shù)值天氣預(yù)報(bào)產(chǎn)品在中國(guó)氣象業(yè)務(wù)中有著非常廣泛的應(yīng)用，除了可提供預(yù)報(bào)信息，還可作為背景場(chǎng)驅(qū)動(dòng)有著同樣廣泛應(yīng)用的區(qū)域數(shù)值天氣模式WRF開(kāi)展短期和短時(shí)預(yù)報(bào)等。隨著GFS預(yù)報(bào)模型空間分辨率的不斷提高，GFS全球數(shù)值天氣預(yù)報(bào)產(chǎn)品的文件大小急速增加。相比過(guò)去1.0°分辨力的數(shù)據(jù)，目前0.25°的grib2格式的預(yù)報(bào)產(chǎn)品的數(shù)據(jù)大小增加了10多倍，單個(gè)文件達(dá)到330 Mb左右，這對(duì)數(shù)據(jù)傳輸和存儲(chǔ)帶來(lái)較大壓力。

然而，及時(shí)下載和存儲(chǔ)GFS的0.25°grib2格式的預(yù)報(bào)產(chǎn)品有著非常重要的業(yè)務(wù)價(jià)值。文章重點(diǎn)分析了當(dāng)前下載GFS數(shù)據(jù)存在的困難和問(wèn)題，給出了解決方法。另外，通過(guò)裁剪關(guān)心區(qū)域和減少無(wú)用變量的方法，縮小了原始grib2格式文件大小，并可給WRF模式提供背景場(chǎng)。

1 方法分析

1.1 現(xiàn)有方法和問(wèn)題

1)在Windows環(huán)境下，訪問(wèn)數(shù)據(jù)網(wǎng)站https：//ftp.ncep.noaa.gov/data/nccf/com/gfs/prod/，點(diǎn)擊下載目標(biāo)文件。0.25°分辨力預(yù)報(bào)結(jié)果的文件命名格式為：gfs.t**z.pgrb2.0p25.f###，其中，“**”取值為00、06、12或18，代表1 d中的4次預(yù)報(bào)；“###”為預(yù)報(bào)時(shí)效。

但由于不明原因，當(dāng)前數(shù)據(jù)下載速度非常慢(測(cè)試環(huán)境帶寬為200 Mb/s)，并且由于數(shù)據(jù)文件比較大，時(shí)常下載停止或掉線，從而導(dǎo)致下載失敗。此現(xiàn)象在其他測(cè)試環(huán)境中也會(huì)出現(xiàn)。

2)在Linux環(huán)境下(測(cè)試環(huán)境帶寬為200 Mb/s)，數(shù)據(jù)下載相對(duì)穩(wěn)定一些。經(jīng)測(cè)試，每日清晨是下載速度最快的時(shí)段。但即使在此時(shí)段，下載12 h預(yù)報(bào)結(jié)果(比如，取3 h間隔，則共5個(gè)數(shù)據(jù)文件)也常需要數(shù)個(gè)小時(shí)，使得GFS資料的可用性降低。

3)訪問(wèn)數(shù)據(jù)網(wǎng)站https：//nomads.ncep.noaa.gov/cgi-bin/filter_gfs_0p25.pl，針對(duì)GFS預(yù)報(bào)產(chǎn)品，可以選取變量、層次和關(guān)心區(qū)域，這樣下載的文件比較小。在相同網(wǎng)速下，這種方法耗時(shí)明顯減少，且成功率非常高。但從2020-06-24起，該網(wǎng)站停止提供此服務(wù)，且尚未恢復(fù)。

1.2 初步解決方法

1)申請(qǐng)1臺(tái)騰訊云計(jì)算服務(wù)器(Linux環(huán)境)，最低配置即可。以騰訊云服務(wù)器為平臺(tái)，利用其IP，可以較為穩(wěn)定的訪問(wèn)數(shù)據(jù)網(wǎng)站https：//ftp.ncep.noaa.gov/data/nccf/com/gfs/prod/。

2)編寫(xiě)腳本，實(shí)現(xiàn)數(shù)據(jù)下載。

3)如需將騰訊云上的下載數(shù)據(jù)轉(zhuǎn)移至本地，數(shù)據(jù)轉(zhuǎn)移速度與購(gòu)買(mǎi)的滕訊云帶寬有關(guān)。為加快數(shù)據(jù)傳輸速度，減少存儲(chǔ)空間，利用python軟件及相關(guān)程序包編寫(xiě)程序，生成僅包含關(guān)心區(qū)域和變量的grib2文件，提供給WRF模式作為背景場(chǎng)。

2 grib2文件格式和WRF模式背景場(chǎng)變量

2.1 grib2文件格式

grib是世界氣象組織(WMO)開(kāi)發(fā)的一種用于交換和存儲(chǔ)規(guī)則分布數(shù)據(jù)的二進(jìn)制文件格式?，F(xiàn)行的grib碼版本有g(shù)rib1和grib2兩種格式，grib2相比grib1有較大優(yōu)點(diǎn)，比如支持多種壓縮方式、表示多維數(shù)據(jù)等。其中，jasper和jpeg是兩種常用且高效的壓縮方式，大大減小了文件大小。

1份grib2文件一般情況下可包含9部分[1]：0段為指示段；1段為標(biāo)志段；2段為本地使用段；3段為網(wǎng)格定義段；4段為產(chǎn)品定義段；5段為數(shù)據(jù)表示段；6段為位圖段；7段為數(shù)據(jù)段；8段為結(jié)束標(biāo)志碼7777。

由于文章僅對(duì)grib2文件進(jìn)行區(qū)域裁剪，不作其他修改，故涉及修改的地方主要在3段和5段。其中，在3段中重新定義區(qū)域，在5段中重新寫(xiě)入數(shù)據(jù)。

2.2 WRF模式背景場(chǎng)變量

grib2格式的GFS資料可作為WRF模式背景場(chǎng)，驅(qū)動(dòng)WRF模式開(kāi)展更高時(shí)空分辨率的數(shù)值模擬。GFS資料包含的變量特別多，但在WRF模式中，通過(guò)編譯安裝grib2讀寫(xiě)庫(kù)，再借助Vtable變量表的形式，可實(shí)現(xiàn)指定變量的讀取[2]。

使用包含以下物理量的grib2格式的GFS資料，可驅(qū)動(dòng)WRF模式正常運(yùn)行：各等壓面層的位勢(shì)高度、溫度、相對(duì)濕度、U風(fēng)速、V風(fēng)速；各層的土壤溫度和土壤濕度；累計(jì)雪深的水相當(dāng)量、雪深、2 m氣溫、2 m相對(duì)濕度、10 mU風(fēng)速、10 mV風(fēng)速、地表氣壓、海陸分布標(biāo)志、海平面氣壓和地形高度。

3 具體實(shí)現(xiàn)

3.1 數(shù)據(jù)下載腳本

核心內(nèi)容是使用Linux系統(tǒng)自帶的wget命令實(shí)現(xiàn)數(shù)據(jù)下載，具體為：

/usr/bin/wget ${filepath}/${filename}

其中，${filepath}為包含數(shù)據(jù)的網(wǎng)址；${filename}為被下載的GFS數(shù)據(jù)文件名。

經(jīng)測(cè)試，利用騰訊云計(jì)算服務(wù)器的數(shù)據(jù)下載速度比使用之前的Linux環(huán)境快5倍左右。

3.2 基于python、ecCodes和pygrib軟件的grib2文件裁剪

1)軟件安裝

參考https：//www.cnblogs.com/sunshine8/p/10348204.html[3]中方法，依次安裝anaconda3(V4.4.0)、ecCodes支持庫(kù)、ecCodes(V 2.18.0)和pygrib(V 2.0.5)等軟件。針對(duì)pygrib軟件安裝過(guò)程中存在的問(wèn)題，文章給出以下解決步驟：

①修改setup.cfg文檔內(nèi)容，主要是去掉下列行的注釋，并修改

tar -xzvf v2.0.5rel.tar.gz

cp setup.cfg.template setup.cfg

packages_to_install=["pygrib"，"ncepgrib2"]

grib_api_dir=/usr/local

jasper_dir=/usr

openjpeg_dir=/usr

grib_api_libname=grib_api_lib

②創(chuàng)建下面的2個(gè)鏈接，因?yàn)閜ygrib需要訪問(wèn)grib_api*庫(kù)

cd /usr/local/lib

ln -sf libeccodes.so libgrib_api_lib.so

ln -sf libeccodes_f90.so libgrib_api_lib_f90.so

③安裝

cd pygrib-2.0.5rel/

python3 setup.py build

python3 setup.py install

2)grib2文件的讀寫(xiě)和區(qū)域裁剪

編寫(xiě)python程序，利用pygrib庫(kù)實(shí)現(xiàn)grib2文件的讀寫(xiě)和區(qū)域裁剪。

①導(dǎo)入庫(kù)

import pygrib

from ncepgrib2 import Grib2Decode，Grib2Encode

②打開(kāi)讀寫(xiě)文件

infile=pygrib.open(infile_str) #被讀的且完整的GFS文件

outfile=open(outfile_str，′wb′) #存儲(chǔ)結(jié)果的GFS文件

③設(shè)置WRF模式運(yùn)行需要的變量

#將要讀取的變量包含在varnames列表中，這些變量可供WRF模式正常運(yùn)行使用。

varnames=[′Geopotential Height′，′Temperature′，′Relative humidity′，′U component of wind′，′V component of wind′，′Soil Temperature′，′Volumetric soil moisture content′，′Water equivalent of accumulated snow depth (deprecated)′，′Snow depth′，′2 metre temperature′，′2 metre relative humidity′，′10 metre U wind component′，′10 metre V wind component′，′Surface pressure′，′Land-sea mask′，′Pressure reduced to MSL′，′MSLP (Eta model reduction)′，′Sea ice area fraction′，′Land -sea coverage (nearest neighbor) [land=1，sea=0]′，′Orography′]

④讀取每一個(gè)gribmessage單元，并進(jìn)行區(qū)域裁剪

for varname in varnames：

grbs=infile.select(name=varname)

for i in range(len(grbs))：

grbmsg=grbs[i]

#獲得關(guān)心區(qū)域的數(shù)據(jù)datsubset和原文件的編碼信息

datsubset，latsubset，lonsubset=grbmsg.data(lat1=lat1，lat2=lat2，lon1=lon1，lon2=lon2)

grb=Grib2Decode(grbmsg.tostring()，gribmsg=True)

#修改編碼信息

grb.grid_definition_info[1]=number_data_unpack #關(guān)心區(qū)域的格點(diǎn)數(shù)

grb.grid_definition_template[7]=points_x #東西方向格點(diǎn)數(shù)

grb.grid_definition_template[8]=points_y #南北方向格點(diǎn)數(shù)

grb.grid_definition_template[11]=lat2*1000000 #起始緯度

grb.grid_definition_template[12]=lon1*1000000 #起始經(jīng)度

grb.grid_definition_template[14]=lat1*1000000 #終點(diǎn)緯度

grb.grid_definition_template[15]=lon2*1000000 #終點(diǎn)經(jīng)度

grb.latitude_first_gridpoint=lat2*1.0

grb.latitude_last_gridpoint=lat1*1.0

grb.longitude_first_gridpoint=lon1*1.0

grb.longitude_last_gridpoint=lon2*1.0

grb.number_of_data_points_to_unpack=number_data_unpack

grb.points_in_x_direction=points_x

grb.points_in_y_direction=points_y

#存儲(chǔ)datsubset和新的編碼信息

grbo=Grib2Encode(grb.discipline_code，grb.identification_section)

grbo.addgrid(grb.grid_definition_info，grb.grid_definition_template)

grbo.addfield(grb.product_definition_template_number，grb.product_definition_template，grb.data_representation_template_number，grb.data_representation_template，datsubset) #寫(xiě)數(shù)據(jù)

grbo.end()

outfile.write(grbo.msg) #寫(xiě)輸出文件

通過(guò)以上步驟，可實(shí)現(xiàn)GFS預(yù)報(bào)資料的區(qū)域裁剪。其中，lat1

3.3 結(jié)果分析

通過(guò)以上方法，根據(jù)關(guān)心區(qū)域的格點(diǎn)數(shù)與全球的格點(diǎn)數(shù)的比值，以及存儲(chǔ)變量個(gè)數(shù)與全部變量個(gè)數(shù)的比值，可以大概算出文件最終的大小。

以生成區(qū)域——lat：17°E～57°E；lon：70°E～120°E的WRF模式背景場(chǎng)為例，實(shí)測(cè)表明采用文章所述方法可將1個(gè)GFS文件的大小從330 Mb縮小至7 Mb，在很大程度上節(jié)省了存儲(chǔ)空間并加快了數(shù)據(jù)傳遞速度。

通過(guò)編寫(xiě)作業(yè)腳本，利用Linux系統(tǒng)自帶的crontab命令設(shè)置定時(shí)作業(yè)，可實(shí)現(xiàn)GFS文件的批量自動(dòng)下載和裁剪。經(jīng)測(cè)試，測(cè)試平臺(tái)可實(shí)現(xiàn)多個(gè)數(shù)據(jù)文件的并行下載，但基于python文件裁剪功能的進(jìn)程數(shù)量不宜同時(shí)超過(guò)4個(gè)，因此需在作業(yè)提交時(shí)給予合理時(shí)間安排。

3.4 pygrib在處理grib文件中的其他應(yīng)用

grib格式是存放氣象數(shù)值預(yù)報(bào)產(chǎn)品的一種重要數(shù)據(jù)格式，開(kāi)展grib格式文件的讀寫(xiě)是一件非常重要的能力。在以往研究中，?；贔ortran或C程序進(jìn)行開(kāi)發(fā)，難度較大，文章基于python環(huán)境，利用pygrib庫(kù)可較容易地實(shí)現(xiàn)數(shù)據(jù)讀寫(xiě)與裁剪。

但此過(guò)程的實(shí)現(xiàn)，也用到了pygrib庫(kù)的其他函數(shù)。比如，grbs.seek()函數(shù)可獲得grib文件中所有變量的文字介紹；grbs.select()函數(shù)可實(shí)現(xiàn)對(duì)指定變量或指定位置的內(nèi)容的讀取；grbs.latlons()函數(shù)可獲得詳細(xì)的經(jīng)緯度信息等。

此外，在pygrib庫(kù)的test子目錄下，有大量的python示例程序可供參考；在utils子目錄下，有多個(gè)可執(zhí)行程序可供使用，比如grib1和grib2的互轉(zhuǎn)、grib文件的拆分和合并等，可加以關(guān)注和使用，這些都將方便用戶對(duì)grib格式文件的使用。

4 結(jié)束語(yǔ)

針對(duì)當(dāng)前美國(guó)GFS全球數(shù)值天氣預(yù)報(bào)資料下載速度較慢和文件所需存儲(chǔ)空間較大的問(wèn)題，文章提出借助騰訊云計(jì)算服務(wù)器實(shí)現(xiàn)數(shù)據(jù)的穩(wěn)定、快速下載；基于python、ecCodes和pygrib軟件實(shí)現(xiàn)grib2格式文件的區(qū)域裁剪，減少無(wú)用變量，縮小文件大小，并可給WRF模式提供背景場(chǎng)，具有較強(qiáng)的業(yè)務(wù)應(yīng)用價(jià)值。