一種基于地轉(zhuǎn)風(fēng)模型的平流層風(fēng)場(chǎng)插值方法

胡康，周江華，*，張曉軍

（1.中國科學(xué)院空天信息創(chuàng)新研究院，北京 100094； 2.中國科學(xué)院大學(xué)光電學(xué)院，北京 100049）

高空氣球作為平流層科學(xué)實(shí)驗(yàn)的主要載體，具有有效載荷大、發(fā)放成本低等顯著優(yōu)點(diǎn)。但是作為一種無動(dòng)力飛行器，高空氣球在駐空階段只能通過調(diào)節(jié)浮力和重力來控制垂直方向的運(yùn)動(dòng)，而無法主動(dòng)調(diào)節(jié)飛行的水平方向，因此，平流層風(fēng)場(chǎng)水平方向風(fēng)速分量的數(shù)據(jù)精度和密度將直接影響高空氣球的飛行軌跡［1］。目前，高空氣球飛行仿真和實(shí)際實(shí)驗(yàn)使用的全球風(fēng)場(chǎng)數(shù)據(jù)主要來自于歐洲中期天氣預(yù)報(bào)中心（ECWMF）的ERA系列全球氣象數(shù)據(jù)模型，其中最新的ERA5模型提供了全球氣象再分析資料和未來5 d內(nèi)的氣象預(yù)報(bào)數(shù)據(jù)，每一時(shí)刻的數(shù)據(jù)水平分辨率為0.25°×0.25°（經(jīng)緯度），包涵從地面到距地面80 km處的37個(gè)等壓面風(fēng)場(chǎng)模型［2］。然而在使用ERA5數(shù)據(jù)對(duì)高空氣球進(jìn)行飛行仿真和實(shí)際操控的過程中，通過控制氣球高度處于合適的風(fēng)速帶來達(dá)成區(qū)域駐留目標(biāo)，其有效駐留控制半徑在100 km以上，控制精度差，對(duì)高空氣球?qū)嶒?yàn)的軌跡預(yù)測(cè)、設(shè)備回收和空域申請(qǐng)等工作都帶來很大程度的影響，因此，需要研究一種根據(jù)現(xiàn)有數(shù)據(jù)生成大范圍、高精度、高密度風(fēng)場(chǎng)數(shù)據(jù)的插值方法來提高高空氣球的控制精度。

現(xiàn)有的風(fēng)場(chǎng)插值方法可以分成2類：直接插值類方法和融合插值類方法。直接插值類方法是對(duì)風(fēng)場(chǎng)的預(yù)報(bào)數(shù)據(jù)或者觀測(cè)數(shù)據(jù)采用各種插值方法進(jìn)行插值，從而獲得最終結(jié)果。Jardin和Erzberger［3］利用全局多項(xiàng)式插值方法和風(fēng)場(chǎng)預(yù)報(bào)對(duì)風(fēng)場(chǎng)格點(diǎn)數(shù)據(jù)進(jìn)行了插值。Tao等［4］研究了Hermite插值在風(fēng)場(chǎng)光譜表示中的誤差傳遞和數(shù)值模型。湯新民等［5］研究了3種不同變異模型下的克里格插值法在風(fēng)場(chǎng)插值中的應(yīng)用效果。董志南等［6］對(duì)地面區(qū)域內(nèi)的實(shí)測(cè)風(fēng)場(chǎng)數(shù)據(jù)進(jìn)行交叉驗(yàn)證，比較了反距離權(quán)重插值、全局多項(xiàng)式插值、局部多項(xiàng)式插值、徑向基函數(shù)插值、最近鄰插值及普通克里金插值6種方法。上述直接插值類方法［3-6］是對(duì)風(fēng)場(chǎng)的預(yù)報(bào)數(shù)據(jù)或者觀測(cè)數(shù)據(jù)直接進(jìn)行插值，具有計(jì)算速度快、數(shù)據(jù)質(zhì)量要求低等優(yōu)點(diǎn)。這些方法主要面對(duì)地面及近地風(fēng)場(chǎng)，觀測(cè)的時(shí)間和空間尺度范圍小且采樣難度低，通過密集的風(fēng)場(chǎng)實(shí)測(cè)數(shù)據(jù)，補(bǔ)插出與實(shí)際風(fēng)場(chǎng)高度擬合的插值函數(shù)。但是平流層風(fēng)場(chǎng)時(shí)間和空間觀測(cè)尺度大、采樣成本高、數(shù)據(jù)稀疏且缺乏連續(xù)性特征［7］，直接插值類方法難以給出有效的插值結(jié)果。

融合插值類方法則先對(duì)風(fēng)場(chǎng)數(shù)據(jù)進(jìn)行融合，在融合過程中分析風(fēng)場(chǎng)特征并設(shè)計(jì)插值函數(shù)，從而對(duì)融合結(jié)果進(jìn)行插值，獲得最終結(jié)果。Gandin［8］采用最優(yōu)插值法進(jìn)行了風(fēng)場(chǎng)融合插值，以最小方差作為約束條件，對(duì)觀測(cè)風(fēng)場(chǎng)和預(yù)報(bào)風(fēng)場(chǎng)進(jìn)行數(shù)據(jù)融合，從而獲得對(duì)風(fēng)場(chǎng)插值結(jié)果的最優(yōu)估計(jì)。朱成陣等［9］通過風(fēng)場(chǎng)的實(shí)測(cè)數(shù)據(jù)和氣象預(yù)報(bào)數(shù)據(jù)，用內(nèi)插外推法預(yù)測(cè)特定位置的風(fēng)場(chǎng)信息。上述融合插值類方法［8-9］將小范圍的實(shí)測(cè)數(shù)據(jù)作為觀測(cè)場(chǎng)，將風(fēng)場(chǎng)預(yù)報(bào)數(shù)據(jù)作為背景場(chǎng)，進(jìn)行數(shù)據(jù)融合，分析風(fēng)場(chǎng)特征，再根據(jù)特征進(jìn)一步設(shè)計(jì)風(fēng)場(chǎng)插值模型。但是統(tǒng)計(jì)數(shù)據(jù)表明，風(fēng)場(chǎng)在時(shí)間和空間維度上有顯著的差異性［10］，只采用小范圍的實(shí)測(cè)數(shù)據(jù)進(jìn)行數(shù)據(jù)融合獲得的插值方法面對(duì)更大的區(qū)域難以具有良好的泛化性能。同時(shí)，計(jì)算高精度的風(fēng)場(chǎng)預(yù)報(bào)作為模型的背景場(chǎng)需要巨大的算力和時(shí)間開銷。以國內(nèi)全球／區(qū)域同化預(yù)測(cè)系統(tǒng)（global／regional assimilation prediction system，GRAPES）中的全球版本為例，作為序貫?zāi)Ｐ?，?duì)某一時(shí)刻0.5°×0.5°（經(jīng)緯度）網(wǎng)格31層模型的積分計(jì)算需要200個(gè)以上的CPU節(jié)點(diǎn)計(jì)算10 d［11］，面對(duì)需要更高數(shù)據(jù)密度的高空氣球?qū)嶒?yàn)，傳統(tǒng)的融合插值類方法很難符合其時(shí)效性要求。

為了解決平流層風(fēng)場(chǎng)插值面臨的插值精度低、計(jì)算時(shí)間長和算力開銷大等問題，進(jìn)行速度快、精度高的插值計(jì)算，本文提出一種基于地轉(zhuǎn)風(fēng)（based on geostrophic wind，BGW）模型的平流層水平方向風(fēng)場(chǎng)插值方法，用2019年7月ERA5全球氣象再分析數(shù)據(jù)對(duì)平流層的全球風(fēng)場(chǎng)進(jìn)行插值計(jì)算實(shí)驗(yàn)驗(yàn)證。相比于直接插值類方法，BGW 插值對(duì)平流層水平分量的插值效果有明顯提升；相比于融合插值類方法，BGW 插值大幅度降低了計(jì)算復(fù)雜度，有效提高了平流層風(fēng)場(chǎng)插值的計(jì)算速度。

1 基于地轉(zhuǎn)風(fēng)模型的風(fēng)場(chǎng)插值

本文提出的BGW 風(fēng)場(chǎng)插值方法是一種改進(jìn)的融合插值類方法，主要包含兩部分：①為了提升計(jì)算效率，用重力位勢(shì)數(shù)據(jù)計(jì)算地轉(zhuǎn)風(fēng)模型作為插值背景場(chǎng)；②為了減小背景場(chǎng)誤差和觀測(cè)誤差對(duì)插值結(jié)果的影響，將已有的風(fēng)速數(shù)據(jù)作為觀測(cè)場(chǎng)，用迭代的二元線性回歸對(duì)背景場(chǎng)數(shù)據(jù)進(jìn)行修正，獲得最終的插值結(jié)果。

1.1 地轉(zhuǎn)風(fēng)模型改進(jìn)

一般的融合插值類方法需要使用氣象預(yù)報(bào)模型計(jì)算的目標(biāo)精度風(fēng)場(chǎng)預(yù)報(bào)數(shù)據(jù)作為背景場(chǎng)。氣象預(yù)報(bào)模型一般是基于大氣運(yùn)動(dòng)方程組的數(shù)值模型，能夠高精度地預(yù)測(cè)風(fēng)速的東向、北向及垂直方向分量，對(duì)中短期風(fēng)場(chǎng)都有著良好的預(yù)測(cè)效果。然而數(shù)值預(yù)報(bào)模型的實(shí)際計(jì)算過程需要大量的歷史數(shù)據(jù)、算力和時(shí)間作為支撐。因此，當(dāng)前的融合類插值方法需要大量的時(shí)間來進(jìn)行計(jì)算，或依賴于外部輸入的目標(biāo)密度風(fēng)場(chǎng)預(yù)報(bào)數(shù)據(jù)，難以在實(shí)驗(yàn)現(xiàn)場(chǎng)根據(jù)需求進(jìn)行實(shí)時(shí)的離線計(jì)算。為了解決這一問題，本文提出一種用地轉(zhuǎn)風(fēng)模型替代風(fēng)場(chǎng)預(yù)報(bào)作為背景場(chǎng)的方法。

地轉(zhuǎn)風(fēng)模型由白貝羅于1857年提出，是大氣運(yùn)動(dòng)方程組在大尺度風(fēng)場(chǎng)條件下推導(dǎo)而出，當(dāng)大氣運(yùn)動(dòng)滿足科里奧利力和氣壓梯度力平衡，即地轉(zhuǎn)平衡條件時(shí)的理論風(fēng)［12］，在局地直角坐標(biāo)系下，其公式如下：

式中：Ω 為地球的旋轉(zhuǎn)角速度，值為7.29×10-5rad／s；φ為緯度。

地轉(zhuǎn)風(fēng)雖然是根據(jù)氣壓分布計(jì)算出的估計(jì)值，并非實(shí)際存在的風(fēng)，但是根據(jù)實(shí)驗(yàn)驗(yàn)證，以地轉(zhuǎn)平衡作為條件計(jì)算的中高緯度自由大氣風(fēng)場(chǎng)與真實(shí)風(fēng)場(chǎng)十分接近［13］，因此用地轉(zhuǎn)風(fēng)模型作為風(fēng)場(chǎng)插值的背景場(chǎng)能夠在提高計(jì)算速度的同時(shí)，準(zhǔn)確地顯示中高緯度地區(qū)的風(fēng)速分布特征。

在緯度逼近赤道的過程中，由于科里奧利頻率f的值逐漸趨近于0，導(dǎo)致在低緯度地區(qū)地轉(zhuǎn)風(fēng)計(jì)算值趨于無窮，因此對(duì)f的計(jì)算公式進(jìn)行改進(jìn)。對(duì)于一般意義上的赤道附近區(qū)域，即南北緯10°之間的區(qū)域，選取邊界對(duì)科里奧利頻率公式進(jìn)行截?cái)?，截?cái)噙吔绲倪x擇范圍為（0°，10°］，越接近0則赤道附近區(qū)域風(fēng)場(chǎng)誤差越大，越接近10則初值場(chǎng)和準(zhǔn)地轉(zhuǎn)運(yùn)動(dòng)偏離越多。綜合考慮選取中值5°為截?cái)噙吔?，將北?°和南緯5°之間的f值固定為2Ωsin 5°，改進(jìn)后的公式如下：

經(jīng)過實(shí)驗(yàn)驗(yàn)證，相比于使用原始的地轉(zhuǎn)風(fēng)模型，導(dǎo)致最終插值精度較差，用式（3）的BGW 插值在低緯度地區(qū)也能獲得良好的插值效果。

1.2 基于二維卷積的地轉(zhuǎn)風(fēng)快速計(jì)算方法

式（1）的地轉(zhuǎn)風(fēng)模型需要大氣密度ρ作為輸入，而在實(shí)際實(shí)驗(yàn)中無法實(shí)時(shí)獲得大氣密度數(shù)據(jù)，因此對(duì)式（1）進(jìn)行改進(jìn)。另外，為了進(jìn)一步提升背景場(chǎng)的計(jì)算速度，結(jié)合GPU的矩陣運(yùn)算優(yōu)勢(shì)改進(jìn)了背景場(chǎng)計(jì)算方法。

1.2.1 基于重力位勢(shì)的地轉(zhuǎn)風(fēng)模型

由于氣壓p是局地直角坐標(biāo)系（x，y，z，t）的單值函數(shù)，可以將坐標(biāo)系中的高度z用氣壓p替代，從而對(duì)于任意坐標(biāo)系（x，y，z，t）下的函數(shù)F（x，y，z，t）有

式（9）把式（1）中的大氣密度和氣壓這2個(gè)變量轉(zhuǎn)換為只有重力位勢(shì)這一個(gè)變量，在降低計(jì)算復(fù)雜度的同時(shí)，解決了地轉(zhuǎn)風(fēng)模型需要大氣密度參與計(jì)算的問題。

1.2.2 基于二維卷積的優(yōu)化計(jì)算

式中：“*”表示矩陣的二維卷積運(yùn)算；A為提取的卷積核。

由于輸入是全球重力位勢(shì)數(shù)據(jù)的平面展開投影，而卷積計(jì)算對(duì)矩陣邊緣有計(jì)算損失，考慮到投影過程是將地球沿180°經(jīng)線東西向分開，南北兩極將緯度90°處展開，在進(jìn)行卷積計(jì)算時(shí)將重力位勢(shì)矩陣左右進(jìn)行循環(huán)擴(kuò)展，上下進(jìn)行復(fù)制擴(kuò)展，這種改進(jìn)有效抑制了卷積計(jì)算的邊緣損失。

相比于傳統(tǒng)模型沿等壓線根據(jù)大氣密度進(jìn)行逐格計(jì)算來求解地轉(zhuǎn)風(fēng)［14］，BGW 插值通過提取模型卷積核，用二維卷積進(jìn)行快速的矩陣計(jì)算，有效提升了地轉(zhuǎn)風(fēng)模型的計(jì)算速度。

1.3 迭代二元線性回歸確定插值函數(shù)

BGW 插值方法將地轉(zhuǎn)風(fēng)模型計(jì)算的目標(biāo)密度風(fēng)場(chǎng)格點(diǎn)數(shù)據(jù)作為背景場(chǎng)，用部分格點(diǎn)的實(shí)際風(fēng)速觀測(cè)數(shù)據(jù)作為觀測(cè)場(chǎng)，用觀測(cè)場(chǎng)對(duì)背景場(chǎng)進(jìn)行修正以獲得最終插值結(jié)果。每一個(gè)格點(diǎn)上的插值結(jié)果是由背景場(chǎng)值和修正值用迭代二元線性回歸確定，修正值由一定范圍內(nèi)N個(gè)格點(diǎn)上的觀測(cè)值與背景場(chǎng)值的偏差加權(quán)獲得，最終獲得如下線性回歸方程式（12）和其約束條件式（13）：

在地理信息領(lǐng)域，傳統(tǒng)的距離相關(guān)權(quán)重矩陣一般采用反距離權(quán)重矩陣。反距離權(quán)重受到觀測(cè)點(diǎn)和待插值點(diǎn)的距離及冪次參數(shù)控制，需要根據(jù)數(shù)據(jù)密度、搜索范圍和插值效果不斷調(diào)整冪次參數(shù)［15］。在權(quán)重wi中引入搜索半徑參數(shù)r并進(jìn)行重新設(shè)計(jì)，使得權(quán)重可以根據(jù)搜索范圍進(jìn)行自適應(yīng)調(diào)整，解決了對(duì)不同搜索范圍都要重新設(shè)計(jì)權(quán)重的問題，其公式如下：

式中：d為觀測(cè)點(diǎn)到待插值點(diǎn)的距離。同時(shí)為了加快線性回歸的收斂速度，對(duì)wi進(jìn)行歸一化處理。

為了進(jìn)一步提高插值精度，在完成一輪線性回歸計(jì)算后，將插值結(jié)果矩陣作為背景場(chǎng)數(shù)據(jù)再一次輸入回歸模型進(jìn)行迭代，重復(fù)迭代過程，直到插值精度收斂。經(jīng)過驗(yàn)證，在兩輪迭代時(shí)模型的插值精度已經(jīng)收斂，取兩輪迭代的模型作為最終的插值函數(shù)。

2 實(shí)驗(yàn)驗(yàn)證和對(duì)比分析

2.1 數(shù)據(jù)集

平流層風(fēng)場(chǎng)插值方法的目的是為了獲得盡可能貼近真實(shí)風(fēng)場(chǎng)情況的插值結(jié)果，ECWMF的ERA系列氣象再分析模型使用了全球氣象資料融合系統(tǒng)和完善的氣象數(shù)據(jù)數(shù)據(jù)庫，從1970年至今，對(duì)全球各地來自地面雷達(dá)測(cè)站、船載氣象測(cè)站、海上浮標(biāo)、探空氣球、飛機(jī)和衛(wèi)星的多源觀測(cè)資料進(jìn)行質(zhì)量評(píng)價(jià)和數(shù)據(jù)融合，從而生成高度接近真實(shí)風(fēng)場(chǎng)的完整數(shù)據(jù)集［2］。以往的高空氣球飛行試驗(yàn)也以ERA系列模型作為參考，歷史數(shù)據(jù)也驗(yàn)證了ERA系列模型的可靠性。因此，選擇ERA系列最新的ERA5模型2019年7月的氣象再分析數(shù)據(jù)進(jìn)行插值計(jì)算和實(shí)驗(yàn)驗(yàn)證。由于高空氣球的駐空階段一般駐留在平流層的底部區(qū)域，距地面20 km左右，位于50 hPa大氣壓等壓面附近，選用50 hPa等壓面模型的0.5°×0.5°（經(jīng)緯度）網(wǎng)格的重力位勢(shì)數(shù)據(jù)作為輸入，用經(jīng)緯度網(wǎng)格0.25°×0.25°的風(fēng)速數(shù)據(jù)進(jìn)行驗(yàn)證。

2.2 評(píng)價(jià)指標(biāo)

通過實(shí)驗(yàn)對(duì)BGW 插值方法中，背景場(chǎng)對(duì)實(shí)際風(fēng)場(chǎng)特征的估計(jì)能力、插值模型的優(yōu)劣程度和插值精度進(jìn)行評(píng)價(jià)，分別使用相關(guān)系數(shù)矩陣、R2檢驗(yàn)和平均絕對(duì)誤差（mean absolute error，MAE）作為評(píng)價(jià)指標(biāo)。

2.2.1 相關(guān)系數(shù)矩陣

背景場(chǎng)對(duì)實(shí)際風(fēng)速分布特征的表達(dá)能力用背景場(chǎng)數(shù)據(jù)與再分析資料的相關(guān)系數(shù)矩陣表示，矩陣X和Y的相關(guān)系數(shù)矩陣ρXY的公式為

2.3 BGW 插值背景場(chǎng)風(fēng)場(chǎng)特征估計(jì)能力實(shí)驗(yàn)

為驗(yàn)證BGW 插值的背景場(chǎng)模型對(duì)實(shí)際風(fēng)場(chǎng)特征的估計(jì)能力，實(shí)驗(yàn)選擇ERA 模型2019年7月第一個(gè)星期的50 hPa等壓面氣象數(shù)據(jù)進(jìn)行實(shí)驗(yàn)，對(duì)一周內(nèi)每日UTC時(shí)間8：00的平流層數(shù)據(jù)，計(jì)算其BGW 插值的背景場(chǎng)和再分析風(fēng)場(chǎng)的相關(guān)系數(shù)矩陣，計(jì)算時(shí)將風(fēng)場(chǎng)數(shù)據(jù)沿緯度線剪開，每一緯度的數(shù)據(jù)作為一維特征。由于每一個(gè)完整的相關(guān)系數(shù)矩陣包涵721個(gè)相關(guān)系數(shù)，難以直接分析，取相關(guān)系數(shù)矩陣的中位數(shù)進(jìn)行分析，結(jié)果保留6位有效數(shù)字，數(shù)據(jù)如表1所示。

表1 相關(guān)系數(shù)矩陣中位數(shù)Table 1 Median of correlation coefficient matrix

從表1中可以看出，西風(fēng)、南風(fēng)和合風(fēng)的背景場(chǎng)與再分析風(fēng)場(chǎng)的相關(guān)系數(shù)中位數(shù)都在0.4以上，屬于中等相關(guān)和強(qiáng)相關(guān)范疇，而且合風(fēng)的相關(guān)系數(shù)都在0.6以上，屬于強(qiáng)相關(guān)，說明背景場(chǎng)數(shù)據(jù)能夠一定程度上表達(dá)實(shí)際風(fēng)場(chǎng)的風(fēng)速分布特征。

牛沖含有優(yōu)質(zhì)蛋白質(zhì)及豐富的維生素 B1、維生素 B2、維生素 C、維生素E，以及雄性激素，具有補(bǔ)腎益精，壯腰膝的功效。常食對(duì)于改善性功能有一定作用。

由于表1只對(duì)風(fēng)速的數(shù)值特征進(jìn)行分析，為了更直觀地分析風(fēng)場(chǎng)的幾何特征，以ERA5模型2019年7月24日UTC時(shí)間8：00的50 hPa等壓面氣象數(shù)據(jù)為例，用數(shù)據(jù)中經(jīng)緯度網(wǎng)格0.5°×0.5°的重力位勢(shì)再分析數(shù)據(jù)作為輸入計(jì)算經(jīng)緯度網(wǎng)格0.25°×0.25°的風(fēng)場(chǎng)水平方向合風(fēng)速的西風(fēng)分量u和南風(fēng)分量v背景場(chǎng)數(shù)據(jù)，風(fēng)速數(shù)據(jù)用熱圖表示，西風(fēng)分量和南風(fēng)分量背景場(chǎng)風(fēng)速熱圖分別如圖1（a）和圖1（b）所示。以ERA5氣象再分析模型的0.5°精度西風(fēng)u和南風(fēng)v數(shù)據(jù)作為實(shí)際風(fēng)場(chǎng)進(jìn)行對(duì)比，其風(fēng)速熱圖分別如圖1（c）和圖1（d）所示。

圖1中每幅子圖的橫坐標(biāo)為經(jīng)度，縱坐標(biāo)為緯度。在同一幅西風(fēng)分量圖中，顏色越偏向黃色表示西風(fēng)越大，顏色越偏向藍(lán)色表示東風(fēng)越大；在同一幅南風(fēng)分量圖中，顏色越偏向黃色表示南風(fēng)越大，顏色越偏向藍(lán)色表示北風(fēng)越大。

從圖1（a）中可以發(fā)現(xiàn)，模型計(jì)算的背景場(chǎng)西風(fēng)數(shù)據(jù)在南緯60°附近有一個(gè)沿緯度圈環(huán)繞地球一周、風(fēng)速40 m／s以上的大風(fēng)速西風(fēng)帶，而其他區(qū)域的西風(fēng)風(fēng)速在沿緯度方向分布較為均勻。從圖1（b）可以發(fā)現(xiàn)，背景場(chǎng)在南緯60°緯線圈上西經(jīng)60°和西經(jīng)140°附近有2處南風(fēng)超過20 m／s的集中分布區(qū)域。而從實(shí)際風(fēng)場(chǎng)圖1（c）和圖1（d）中可以看出，背景場(chǎng)估計(jì)的西風(fēng)分量和南風(fēng)分量的大風(fēng)速帶和風(fēng)速集中分布區(qū)域的主要特征與觀測(cè)場(chǎng)基本一致，結(jié)合相關(guān)系數(shù)分析結(jié)果說明，用改進(jìn)的地轉(zhuǎn)風(fēng)模型作為插值背景場(chǎng)能夠替代風(fēng)場(chǎng)預(yù)報(bào)數(shù)據(jù)，對(duì)風(fēng)速分布特征進(jìn)行估計(jì)。

圖1 風(fēng)速背景場(chǎng)和實(shí)際風(fēng)場(chǎng)數(shù)據(jù)Fig.1 Background field and actual field of wind speed data

2.4 平流層風(fēng)場(chǎng)插值實(shí)驗(yàn)

用搜索半徑分別為30 km、60 km和90 km的BGW 插值方法與直接插值類方法中使用最多的最近鄰域插值法（nearest neighbor，NN）和反距離權(quán)重插值法（inverse distance weighted，IDW）進(jìn)行對(duì)比，IDW 的冪次參數(shù)取2，比較幾種方法在2019年7月每天24個(gè)50 hPa全球風(fēng)場(chǎng)數(shù)據(jù)上的模型插值性能，取R2檢驗(yàn)結(jié)果的均值作為評(píng)價(jià)指標(biāo)，結(jié)果保留6位有效數(shù)字，數(shù)據(jù)如表2所示。

表2 不同插值方法的R2 檢驗(yàn)結(jié)果Table 2 R2 square test results of different interpolation methods

對(duì)于數(shù)據(jù)網(wǎng)格間距為0.5°（經(jīng)緯度）的輸入數(shù)據(jù)，待插值點(diǎn)與觀測(cè)點(diǎn)的最遠(yuǎn)平均距離約為20 km，從表2的R2檢驗(yàn)結(jié)果可以發(fā)現(xiàn)，BGW 插值的搜索半徑為30 km時(shí)模型性能最好，該結(jié)果驗(yàn)證了1.3節(jié)中搜索半徑的最佳參數(shù)范圍應(yīng)該在待插值點(diǎn)與觀測(cè)點(diǎn)最遠(yuǎn)距離的一倍到兩倍之間的結(jié)論。之后的實(shí)驗(yàn)BGW 插值默認(rèn)搜索半徑選擇30 km。

表3 不同插值方法的平均絕對(duì)誤差Table 3 Mean absolute error of different interpolation methods

從表3中可以發(fā)現(xiàn)，BGW 插值在西風(fēng)、南風(fēng)和合風(fēng)上的MAE都要小于NN和IDW 插值。在西風(fēng)上，BGW 插值誤差是 NN 插值誤差的42.32%，是IDW 插值誤差的20.77%；在南風(fēng)上，BGW 插值誤差是NN插值誤差的48.41%，是IDW 插值誤差的85.68%；在合風(fēng)上，BGW 插值誤差是NN插值誤差的59.26%，是IDW 插值誤差的29.16%。

同時(shí)對(duì)BGW 插值方法的計(jì)算速度進(jìn)行實(shí)驗(yàn)，對(duì)于單等壓面風(fēng)場(chǎng)的平面插值問題，用Windows系統(tǒng)下的Python代碼進(jìn)行算法實(shí)現(xiàn)，用英特爾7代i7CPU進(jìn)行計(jì)算的平均耗時(shí)為1.37 s，用英偉達(dá)GTX1050GPU進(jìn)行同樣工作的平均耗時(shí)為1.26 s。而且在實(shí)際的高空氣球飛行仿真和外場(chǎng)實(shí)驗(yàn)中，需要對(duì)風(fēng)場(chǎng)進(jìn)行持續(xù)、實(shí)時(shí)的插值計(jì)算，GPU在矩陣計(jì)算上的速度優(yōu)勢(shì)在面對(duì)大量矩陣計(jì)算時(shí)還會(huì)進(jìn)一步擴(kuò)大。另外，一般的融合類插值背景場(chǎng)采用氣象預(yù)報(bào)模型，無法對(duì)每一層等壓面風(fēng)場(chǎng)數(shù)據(jù)解耦進(jìn)行單獨(dú)計(jì)算，而完整的預(yù)測(cè)模型需要大量的算力設(shè)備支持和天級(jí)的運(yùn)算時(shí)間，BGW 插值有效解決了融合插值類方法背景場(chǎng)計(jì)算時(shí)間長、算力開銷大的問題。

實(shí)驗(yàn)結(jié)果表明，相比于直接插值類方法，BGW 插值對(duì)一維的重力位勢(shì)矢量進(jìn)行插值，再用觀測(cè)風(fēng)速對(duì)背景場(chǎng)進(jìn)行修正，誤差絕對(duì)值要小于對(duì)三維的觀測(cè)風(fēng)速矢量進(jìn)行直接插值，有效減小了插值方法的誤差。而相比于融合插值類方法，BGW 插值是用當(dāng)前時(shí)刻的重力位勢(shì)數(shù)據(jù)結(jié)合地轉(zhuǎn)風(fēng)模型來計(jì)算風(fēng)速的估計(jì)值作為背景場(chǎng)，與直接引入預(yù)報(bào)數(shù)據(jù)的方法相比，BGW 插值無需引入外部目標(biāo)精度的氣象預(yù)報(bào)數(shù)據(jù)，可以只通過已有的數(shù)據(jù)進(jìn)行離線運(yùn)算；與其他自主計(jì)算風(fēng)場(chǎng)預(yù)報(bào)的方法相比，由于預(yù)報(bào)模型是一個(gè)高度耦合的數(shù)值系統(tǒng)，無法對(duì)風(fēng)速的水平方向分量單獨(dú)進(jìn)行預(yù)測(cè)，導(dǎo)致其模型的運(yùn)算時(shí)間和算力消耗極大，而BGW 插值能只對(duì)單一等壓面的背景場(chǎng)進(jìn)行計(jì)算，且計(jì)算速度在秒級(jí)，有效提高了模型的計(jì)算效率。

3 結(jié) 論

面對(duì)稀疏的平流層風(fēng)場(chǎng)數(shù)據(jù)，基于地轉(zhuǎn)風(fēng)模型和迭代二元線性回歸提出一種針對(duì)全球平流層風(fēng)場(chǎng)水平分量數(shù)據(jù)的BGW 插值方法，有效提升了計(jì)算速度和精度。

1）提出了一種適用于全緯度地區(qū)地轉(zhuǎn)風(fēng)模型計(jì)算的科里奧利頻率改進(jìn)公式，同時(shí)推導(dǎo)出一種適用于高空氣球?qū)嶒?yàn)且能有效表達(dá)真實(shí)風(fēng)場(chǎng)特征的地轉(zhuǎn)風(fēng)快速計(jì)算方法。

2）通過迭代的二元線性回歸和改進(jìn)的自適應(yīng)觀測(cè)場(chǎng)偏差權(quán)重矩陣減小了觀測(cè)場(chǎng)誤差對(duì)插值結(jié)果的影響，并解決了面對(duì)不同密度的風(fēng)場(chǎng)數(shù)據(jù)都要重新設(shè)計(jì)權(quán)重矩陣的問題。

3）實(shí)驗(yàn)結(jié)果表明，與2種直接插值類方法相比，NN插值的平均絕對(duì)誤差約為0.21 m／s，IDW插值的平均絕對(duì)誤差約為0.43 m／s，BGW 插值將風(fēng)速平均絕對(duì)誤差減小到約0.13 m／s，有效提升了插值精度；與融合插值類方法相比，一般使用的背景場(chǎng)模型需要10 d以上的計(jì)算時(shí)間，BGW 插值對(duì)單一等壓面的風(fēng)場(chǎng)模型的插值計(jì)算速度則提升到1.26 s，有效提升了插值速度。

為使提出的平流層風(fēng)場(chǎng)插值方法能對(duì)高空氣球仿真和實(shí)驗(yàn)提供進(jìn)一步支持，下一階段工作應(yīng)將目前的風(fēng)場(chǎng)插值方法結(jié)合大氣位溫方程，研究時(shí)間維度上的風(fēng)場(chǎng)插值方法。