本發(fā)明涉及氣象數(shù)據(jù)質量控制技術領域，特別是一種區(qū)域自動氣象站小時雨量數(shù)據(jù)質量控制方法。

背景技術：

近幾年，隨著我國地面氣象觀測現(xiàn)代化建設的推進，全國自動氣象站數(shù)量已經(jīng)超過5萬個，但這也加大了維護的難度。自動氣象站的雨量數(shù)據(jù)是準確反映實時降水量的重要依據(jù)，然而在實際使用過程中存在著一定的問題，比如：站點遷址導致其地理位置信息不準確、翻斗式雨量計易受周圍大風影響、雜物堵塞、人為澆水、破壞以及維護不當?shù)纫蛩卦斐捎炅繙y量誤差等，這些都給業(yè)務人員進行數(shù)據(jù)分析帶來了非常大的困擾，無法判斷真實的天氣過程。因此如何準確評估自動氣象站雨量數(shù)據(jù)質量，是一個非常緊迫和亟待解決的問題。

傳統(tǒng)的自動氣象站雨量數(shù)據(jù)質量控制方法主要包括氣候學界限值檢查、區(qū)域界限值檢查、時間一致性檢查和空間一致性檢查等。這些方法是以降水分布滿足時空分布連續(xù)性為前提條件的；而在夏季，我國強對流天氣頻發(fā)，常常會發(fā)生兩類降水系統(tǒng)，分別是局地陣性降水和大面積層狀云中鑲嵌小尺度強對流云團的混合型降水，這樣的降水系統(tǒng)會給現(xiàn)有的自動氣象站雨量數(shù)據(jù)質量控制技術判斷帶來障礙。

隨著雷達遙感技術的發(fā)展以及雷達觀測網(wǎng)的建設，目前國內(nèi)外開始發(fā)展基于雷達定量估測降水的自動氣象站雨量數(shù)據(jù)的質量控制方法。marzen和fuelberg(2005)發(fā)展了一套基于多傳感器降水估計算法修正后的雷達定量估測降水數(shù)據(jù)的雨量數(shù)據(jù)質量控制方法。zhangjian(2016)發(fā)展了一套基于mrms系統(tǒng)中多傳感器的小時雨量實時質量控制系統(tǒng)。叢芳等(2011)利用雷達定量估計降水與自動氣象站站雨量數(shù)據(jù)之間的相似離度關系，進行了將雷達資料應用于自動氣象站雨量數(shù)據(jù)質量控制的試驗性研究。仲凌志等(2013)利用天氣雷達輔助對自動氣象站逐小時降水數(shù)據(jù)進行了實時質量控制及評估。

上述的自動氣象站降水數(shù)據(jù)質量控制方法中，均沒有對雷達數(shù)據(jù)的準確性進行評估，同時雷達定量降水估計中使用的z-r關系系數(shù)也是前人給出的，并不一定適合當前的降水系統(tǒng)。

技術實現(xiàn)要素：

本發(fā)明要解決的技術問題是針對上述現(xiàn)有技術的不足，而提供一種區(qū)域自動氣象站小時雨量數(shù)據(jù)質量控制方法，該區(qū)域自動氣象站小時雨量數(shù)據(jù)質量控制方法同時接入了3類數(shù)據(jù)，分別為天氣雷達數(shù)據(jù)、雨滴譜儀數(shù)據(jù)、自動氣象站數(shù)據(jù)。先利用雨滴譜儀測量的雨滴譜儀數(shù)據(jù)對天氣雷達數(shù)據(jù)進行修正，然后使用修正后的天氣雷達數(shù)據(jù)對自動氣象站數(shù)據(jù)進行自動判別，從而使自動氣象站雨量數(shù)據(jù)更加準確、可靠，提高觀測數(shù)據(jù)質量，為氣象監(jiān)測、預警和預報等業(yè)務的運行奠定良好的基礎。

為解決上述技術問題，本發(fā)明采用的技術方案是：

一種區(qū)域自動氣象站小時雨量數(shù)據(jù)質量控制方法，包括如下步驟。

步驟1，建立判別模型：根據(jù)歷史觀測資料，選擇一段時間統(tǒng)計天氣雷達小時雨量和自動氣象站小時雨量數(shù)據(jù)，建立判別模型，判別模型的公式為：

上式中，f1(x)為可信自動氣象站小時雨量的上邊界，f2(x)為可信自動氣象站小時雨量的下邊界；x表示雷達估測小時雨量；a、k1、t1、k2、t2為方程系數(shù)，根據(jù)經(jīng)驗，a∈[0.5,1]，k1∈[0.5,2]，k2∈[0.5,2]，t1∈[0.1,10]，t2∈[0.1,10]。

步驟2，數(shù)據(jù)讀取及轉換：分別將天氣雷達、雨滴譜儀和區(qū)域自動氣象站的采集數(shù)據(jù)進行讀取，并均轉換成標準格式的數(shù)據(jù)。

步驟3，修正雷達回波強度：利用步驟2轉換后的雨滴譜儀采集數(shù)據(jù)，完成對天氣雷達實測回波強度信息的評估和訂正；評估和訂正的方法如下：

步驟31，通過步驟2轉換后的雨滴譜儀采集數(shù)據(jù)估算出雨滴譜儀回波強度信息；雨滴譜儀回波強度z的計算公式為：

式中，n(di)為單位體積內(nèi)第i級直徑對應的所有降水粒子；di為第i級直徑的粒子。

步驟32，坐標關系映射：將天氣雷達回波數(shù)據(jù)從極坐標映射到笛卡爾坐標。

步驟33，天氣雷達回波強度修正：以雨滴譜儀的空間位置為中心，確定空間匹配窗口；在空間匹配窗口內(nèi)，將步驟32中位于笛卡爾坐標內(nèi)的天氣雷達回波數(shù)據(jù)與步驟31估算出的雨滴譜儀回波強度信息數(shù)據(jù)進行對比，完成空間匹配，計算出天氣雷達回波強度的偏差并進行修正。

步驟4，建立不同降水類型下的天氣雷達回波強度與雨強關系式：利用步驟2中的雨滴譜儀采集數(shù)據(jù)，通過最小二乘法擬合出不同降水類型下的雷達回波強度與雨強關系式。

具體建立方法如下：

步驟41，建立天氣雷達回波強度與雨強的通用z-r關系式，該通用z-r關系式為：

z＝ar^b

式中，z為天氣雷達回波強度；r為雨強；a、b為待定常數(shù)，不同降水類型所對應的a和b值不同。

步驟42，降水類型劃分：先設定一個固定閾值作為劃分降水類型的分割線，然后根據(jù)該固定閾值，將步驟31估算出的雨滴譜儀回波強度z進行劃分，具體劃分方式為：當雨滴譜儀回波強度z大于等于固定閾值時，判斷為對流降水；當雨滴譜儀回波強度z小于固定閾值時，判斷為層狀云降水。

步驟43，計算不同降水類型下的待定常數(shù)a、b；將步驟41建立的通用z-r關系式進行變換，然后采用最小二乘法線性回歸方程得出待定常數(shù)a、b的計算公式，然后根據(jù)待定常數(shù)a、b的計算公式，分別擬合出對流降水和層狀云降水時的常數(shù)a、b值及所對應的z-r關系式；其中，待定常數(shù)a、b的計算公式如下：

式中，xi是雨滴譜儀獲取的雨強數(shù)據(jù)樣本，yi是雨滴譜儀計算出的雷達回波強度數(shù)據(jù)樣本，n為樣本數(shù)量。

步驟5，建立雨量集合：根據(jù)步驟3得到的修正雷達回波強度，代入步驟4建立的不同降水類型下的天氣雷達回波強度與雨強關系式，得到雷達估測小時降水雨量，聯(lián)合區(qū)域自動氣象站小時雨量數(shù)據(jù)構成二維坐標點，得到雨量集合；聯(lián)合雷達覆蓋范圍內(nèi)的區(qū)域自動氣象站小時雨量數(shù)據(jù)，獲取空間上匹配的點對點雨量集合{[x_rg1,y_rr1],[x_rg2,y_rr2],…,[x_rgn,y_rrn]},其中x_rg代表雷達估測的小時降水雨量，y_rr代表區(qū)域自動氣象站的小時雨量，n代表自動氣象站的個數(shù)；

步驟6，自動氣象站小時雨量數(shù)據(jù)判斷：根據(jù)步驟5獲得的雨量集合，分別代入步驟1建立的判別模型中，判斷區(qū)域自動氣象站小時雨量數(shù)據(jù)的質量；當區(qū)域自動氣象站小時雨量數(shù)據(jù)處于上下邊界之間時，即f1(x_rg)＜y_rr＜f2(x_rg)，判斷區(qū)域自動氣象站小時雨量數(shù)據(jù)可信，反之，則判別為異常。

所述步驟2中，雨滴譜儀采集數(shù)據(jù)讀取與轉換后，需進行質量控制；雨滴譜儀采集數(shù)據(jù)質量控制的方法為：

步驟21，建立雨滴譜儀粒子粒徑和下落速度的關系式：雨滴譜儀粒子粒徑和下落速度的關系式為：

v＝-0.002362d⁴+0.07934d³-0.9551d²+4.932d-0.1021

其中d為粒子粒徑，υ為理論下落速度；

步驟22，剔除非真實氣象數(shù)據(jù)：當雨滴譜儀所觀測粒子的下落速度符合如下計算公式時，則將該觀測粒子的下落速度數(shù)據(jù)進行剔除；需要符合的計算公式為：

其中，υ1為觀測粒子的下落速度，υ為步驟21計算出的理論下落速度，a0為設定閾值。

所述步驟2中，設定閾值a0取值為0.6。

所述步驟1中，根據(jù)歷史觀測資料，選擇不少于一個月的天氣雷達小時雨量和自動氣象站小時雨量數(shù)據(jù)，建立判別模型。

所述步驟32，在進行天氣雷達回波數(shù)據(jù)坐標轉換時，需要進行空間分辨率的設置，為保證數(shù)據(jù)不失真，空間分辨率的設置采用天氣雷達距離分辨率的整數(shù)倍，最小分辨率為天氣雷達的一個探測距離單元。

所述步驟42中，在進行降水類型劃分時，所設定的一個固定閾值為35dbz。

本發(fā)明采用上述方法后，能同時接入3類數(shù)據(jù)，分別為天氣雷達數(shù)據(jù)、雨滴譜儀數(shù)據(jù)、自動氣象站數(shù)據(jù)。通過歷史數(shù)據(jù)建立的判別模型及設備間的互補性，實現(xiàn)區(qū)域自動氣象站小時雨量數(shù)據(jù)質量檢驗。首先根據(jù)歷史數(shù)據(jù)建立自動氣象站小時雨量質量檢驗判別模型，其次實時整理雨滴譜儀數(shù)據(jù)、天氣雷達數(shù)據(jù)和自動氣象站雨量數(shù)據(jù)，利用多部雨滴譜儀數(shù)據(jù)，完成天氣雷達回波強度評估和訂正，并擬合出天氣雷達定量估測降水所需的待定參數(shù)，最后建立天氣雷達估測小時雨量與區(qū)域自動氣象站小時雨量集合，通過與預先建立的判別模型進行比較，完成區(qū)域自動氣象站小時雨量數(shù)據(jù)是否異常的判定，從而使自動氣象站雨量數(shù)據(jù)更加準確、可靠，提高觀測數(shù)據(jù)質量，為氣象監(jiān)測、預警和預報等業(yè)務的運行奠定良好的基礎。

附圖說明

圖1顯示了天氣雷達定量估測雨量數(shù)據(jù)與人工核查的自動站雨量數(shù)據(jù)所建立的判別模型。

圖2顯示了本發(fā)明一種區(qū)域自動氣象站小時雨量數(shù)據(jù)質量控制方法的流程示意圖。

具體實施方式

下面結合附圖和具體較佳實施方式對本發(fā)明作進一步詳細的說明。

如圖2所示，一種區(qū)域自動氣象站小時雨量數(shù)據(jù)質量控制方法，包括如下步驟。

步驟1，建立判別模型：根據(jù)歷史觀測資料，選擇一段時間統(tǒng)計天氣雷達小時雨量和自動氣象站小時雨量數(shù)據(jù)，建立判別模型。

判別模型的優(yōu)選建立過程如下所述。

步驟11，歷史數(shù)據(jù)選擇。

在天氣雷達和自動氣象站的歷史觀測數(shù)據(jù)中，選取不少于1個月的天氣雷達估測小時雨量結果rr_old(t)及其覆蓋范圍內(nèi)經(jīng)人工核查的自動氣象站小時雨量數(shù)據(jù)rg_old(t)，其中t為時間采集點數(shù)。

步驟12，建立關統(tǒng)計系圖。

在相同時間尺度的情況下，以雷達估測小時雨量為自變量，自動氣象站小時雨量為因變量，建立如圖1所示的統(tǒng)計關系圖。

步驟13，建立判別模型。

根據(jù)上述統(tǒng)計關系圖，擬合出上下邊界的函數(shù)關系式，建立如下所述的判別模型：

上式中，f1(x)為可信自動氣象站小時雨量的上邊界，f2(x)為可信自動氣象站小時雨量的下邊界；x表示雷達估測小時雨量；a、k1、t1、k2、t2為方程系數(shù)，根據(jù)經(jīng)驗，a∈[0.5,1]，k1∈[0.5,2]，k2∈[0.5,2]，t1∈[0.1,10]，t2∈[0.1,10]；本次實驗中a＝lg5，k1＝1，t1＝1，k2＝1，t2＝5。

另外，上述可信自動氣象站小時雨量是指經(jīng)過人工檢查后的自動氣象站小時雨量數(shù)據(jù)。

步驟2，數(shù)據(jù)讀取及轉換：分別將天氣雷達、雨滴譜儀和區(qū)域自動氣象站的采集數(shù)據(jù)進行讀取，并均轉換成標準格式的數(shù)據(jù)。

在這里，為方便后續(xù)的數(shù)據(jù)處理，首先按照天氣雷達、雨滴譜儀和區(qū)域自動氣象站的數(shù)據(jù)格式完成數(shù)據(jù)讀取和標準化工作，同時需要對雨滴譜儀數(shù)據(jù)進行必要的質量控制，濾除數(shù)據(jù)中異常下落粒子等非真實的氣象數(shù)據(jù)，保留真實有效的采集數(shù)據(jù)結果，使得雨滴譜儀數(shù)據(jù)更為準確、可靠。

雨滴譜儀采集數(shù)據(jù)質量控制的方法為：

步驟21，建立雨滴譜儀粒子粒徑和下落速度的關系式。

雨滴譜儀粒子粒徑和下落速度的關系式為：

v＝-0.002362d⁴+0.07934d³-0.9551d²+4.932d-0.1021

其中d為粒子粒徑，υ為理論下落速度。

步驟22，剔除非真實氣象數(shù)據(jù)。

當雨滴譜儀所觀測粒子的下落速度符合如下計算公式時，則將該觀測粒子的下落速度數(shù)據(jù)進行剔除；需要符合的計算公式為：

其中，υ1為觀測粒子的下落速度，υ為步驟21計算出的理論下落速度，a0為設定閾值，優(yōu)選取值為0.6。

步驟3，修正雷達回波強度：利用步驟2轉換后的雨滴譜儀采集數(shù)據(jù)，完成對天氣雷達實測回波強度信息的評估和訂正。

評估和訂正的方法如下：

步驟31，通過步驟2轉換后的雨滴譜儀采集數(shù)據(jù)估算出雨滴譜儀回波強度信息。

回波強度的定義是單位體積內(nèi)粒子直徑的6次方之和，它同粒子直徑大小密切相關，如下式：

式中，d為粒子直徑；n(d)是直徑為d的粒子數(shù)密度，即降水粒子譜分布。

雨滴譜儀采集數(shù)據(jù)是量化后的結果，因而雨滴譜儀回波強度z的計算公式為：

式中，n(di)為單位體積內(nèi)第i級直徑對應的所有降水粒子；di為第i級直徑的粒子。

步驟32，坐標關系映射：將天氣雷達回波數(shù)據(jù)從極坐標映射到笛卡爾坐標。

天氣雷達的掃描數(shù)據(jù)是以極坐標的形式進行存儲，為獲取空間格點信息，需要進行坐標轉換，設三維網(wǎng)格中任意網(wǎng)格單元的坐標為(ji,wi,hi)，其中ji為緯度，wi為經(jīng)度，hi為高度；天氣雷達天線所在位置坐標為(jr,wr,hr)，其中jr為緯度，wr為經(jīng)度，hr為高度；基于雷達波束傳播和大圓幾何學理論，三維網(wǎng)格任意單元相對于雷達點的極坐標位置(r,az,el)，其中r為斜距，az為方位角，el為仰角，則：

方位az的表達式為：

仰角el的表達式為：

斜距r的表達式為：

r＝sin(s/rm)(rm+hi-hr)/cos(el)

式中，c＝cos(ji)sin(wi-wr)/sin(s/r)，r為地球半徑，rm為等效地球半徑，s＝rcos^-1(sin(jr)sin(ji)+cos(jr)cos(ji)cos(wi-wr))。

在進行天氣雷達回波數(shù)據(jù)坐標轉換時，需要進行空間分辨率的設置，為保證數(shù)據(jù)不失真，空間分辨率的設置采用天氣雷達距離分辨率的整數(shù)倍，最小分辨率為天氣雷達的一個探測距離單元。

步驟33，天氣雷達回波強度修正：以雨滴譜儀的空間位置為中心，確定空間匹配窗口；在空間匹配窗口內(nèi)，將步驟32中位于笛卡爾坐標內(nèi)的天氣雷達回波數(shù)據(jù)與步驟31估算出的雨滴譜儀回波強度信息數(shù)據(jù)進行對比，完成空間匹配，計算出天氣雷達回波強度的偏差并進行修正。

根據(jù)步驟32設置的空間分辨率，來確定本步驟的空間匹配窗口大小m×n，本次具體實施的空間匹配窗口大小采用3×3，假設第i個雨滴譜儀計算得到的回波強度信息為zdi(nj,nw,nh)，其中nj為緯度索引，nw為經(jīng)度索引，nh為高度索引，以雨滴譜儀所在位置為中心，獲取映射到笛卡爾坐標系下的天氣雷達回波強度信息zr_3*3，如下式所示：

式中，zr_3×3代表匹配窗口3×3范圍內(nèi)的天氣雷達回波強度信息，zr代表天氣雷達回波強度。

通過對比集合zr_3×3與zdi(nj,nw,nh)，獲取絕對偏差最小的位置，并記錄下該點與zdi(nj,nw,nh)的差值為δzrdi。按照上述步驟，遍歷所有雨滴譜儀，求取回波偏差值，并計算偏差值的平均值，記為將天氣雷達回波信息zr整體增加偏差值的平均結果zmv，完成對天氣雷達回波強度信息進行修正。

步驟4，建立不同降水類型下的天氣雷達回波強度與雨強關系式。

利用步驟2中的雨滴譜儀采集數(shù)據(jù)，通過最小二乘法擬合出不同降水類型下的雷達回波強度與雨強關系式。

具體建立方法如下：

步驟41，建立天氣雷達回波強度與雨強的通用z-r關系式，該通用z-r關系式為：

z＝ar^b

式中，z為天氣雷達回波強度；r為雨強；a、b為待定常數(shù)，不同降水類型所對應的a和b值不同。

步驟42，降水類型劃分。

本申請中的降水類型分為兩類，一類為對流降水，另一類為層狀云降水。

先設定一個固定閾值作為劃分降水類型的分割線，本具體實施優(yōu)選采用35dbz作為固定閾值。然后根據(jù)該固定閾值，將步驟31估算出的雨滴譜儀回波強度z進行劃分，具體劃分方式為：當雨滴譜儀回波強度z大于等于固定閾值時，也即z≥35dbz時，判斷為對流降水；當雨滴譜儀回波強度z小于固定閾值時，也即z＜35dbz時，判斷為層狀云降水。

步驟43，計算不同降水類型下的待定常數(shù)a、b。

先將步驟41建立的通用z-r關系式進行變換為：lgz＝lg(a)+b×lg(r)。

然后記x＝lg(r)，y＝lg(z)，再采用最小二乘法線性回歸方程得出待定常數(shù)a、b的計算公式，然后根據(jù)待定常數(shù)a、b的計算公式，分別擬合出對流降水和層狀云降水時的常數(shù)a、b值及所對應的z-r關系式。

其中，待定常數(shù)a、b的計算公式如下：

式中，xi是雨滴譜儀獲取的雨強數(shù)據(jù)樣本，yi是雨滴譜儀計算出的雷達回波強度數(shù)據(jù)樣本，n為樣本數(shù)量。

步驟5，建立雨量集合：根據(jù)步驟3得到的修正雷達回波強度，代入步驟4建立的不同降水類型下的天氣雷達回波強度與雨強關系式，得到雷達估測小時降水雨量，聯(lián)合區(qū)域自動氣象站小時雨量數(shù)據(jù)構成二維坐標點，得到雨量集合。聯(lián)合雷達覆蓋范圍內(nèi)的區(qū)域自動氣象站小時雨量數(shù)據(jù)，獲取空間上匹配的點對點雨量集合{[x_rg1,y_rr1],[x_rg2,y_rr2],…,[x_rgn,y_rrn]},其中x_rg代表雷達估測的小時降水雨量，y_rr代表區(qū)域自動氣象站的小時雨量，n代表自動氣象站的個數(shù)。

步驟6，自動氣象站小時雨量數(shù)據(jù)判斷：根據(jù)步驟5獲得的雨量集合，分別代入步驟1建立的判別模型中，判斷區(qū)域自動氣象站小時雨量數(shù)據(jù)的質量；當區(qū)域自動氣象站小時雨量數(shù)據(jù)處于上下邊界之間時，即f1(x_rg)＜y_rr＜f2(x_rg)，判斷區(qū)域自動氣象站小時雨量數(shù)據(jù)可信，反之，則判別為異常。

以上詳細描述了本發(fā)明的優(yōu)選實施方式，但是，本發(fā)明并不限于上述實施方式中的具體細節(jié)，在本發(fā)明的技術構思范圍內(nèi)，可以對本發(fā)明的技術方案進行多種等同變換，這些等同變換均屬于本發(fā)明的保護范圍。