免責聲明

首先聲明本人不是氣候方向，對于氣候模式完全不了解，對于訂正技術(shù)完全不了解。

之所以寫這篇文章純粹是有個需求用到，在自己的初步研究下做出來了。

對于結(jié)果不做任何保證。

數(shù)據(jù)介紹

54511站每天一個數(shù)據(jù)，包括多種要素，如氣溫、降水等

1. 不同氣候模式不同情境的模擬，時間為2015-2050

2. 1979-2014實況數(shù)據(jù)

3. 1979-2014氣候模式模擬的歷史數(shù)據(jù)

Quantile Mapping

使用該方法對氣候模式預測數(shù)據(jù)進行訂正，參考這篇文章，尤其要感謝作者不僅實現(xiàn)了該方法，還進行了開源，在github上搜索pycat即可找到，用戶手冊作為重要參考。

其實來說我不是很了解算法實現(xiàn)的細節(jié)，經(jīng)過咨詢相關領域?qū)＜?，大體上是使用實況的概率分布去訂正氣候模式預報的結(jié)果。

程序編寫

開源的pycat不能直接使用，其在用戶手冊中提供的示例以及數(shù)據(jù)均具有一定的特殊性。因為時間關系，至今我尚未搞明白在netcdf文件下寫入投影信息時如何做到的。不過經(jīng)過修改源碼不影響對數(shù)據(jù)進行訂正。

1. 程序使用iris來讀取netcdf數(shù)據(jù)，iris是英國氣象局開源氣象處理軟件，我早期使用過，但后來發(fā)現(xiàn)了xarray就不再使用了...

2. 在pycat中包含了插值過程、匹配組合歷史和氣候模式模擬歷史數(shù)據(jù)等過程。

由于本人的數(shù)據(jù)是站點數(shù)據(jù)，因此需要先寫成netcdf文件，在寫入的過程中需要注意幾個重要的點：

1. 首先保證所有數(shù)據(jù)沒有缺失數(shù)據(jù)，如果有請自動剔除或使用pandas的interpolate函數(shù)進行處理；

2. netcdf文件中coord要增加一些屬性，如axis，其中在x軸要命名為大寫X。

    hist_nc.x.attrs.update({
     'axis':'X',
     'units':'degrees_east',
     'standard_name':'longitude',
     'long_name':'x of longitude'
    })
    hist_nc.y.attrs.update({
     'axis':'Y',
     'units':'degrees_north',
     'standard_name':'latitude',
     'long_name':'x of latitude'
    })
    hist_nc.time.attrs.update({
     'axis':'T',
     'standard_name':'time',
    #         'units':'hours since 1970-01-01 00:00:00',
    #         'calendar':'standard'
    })

稍作解釋，為什么對于time變量沒有寫入untis和calendar，因為xarray在寫入時會自動匹配當前的時間數(shù)據(jù)格式。

3. 由于pycat的作者每次僅對一個變量進行訂正，且考慮變量的概率分布不同還可能采用不同的方法，因此在寫入數(shù)據(jù)時要增加變量的各類屬性。

if iname=='tt':
    standard_name = 'air_temperature'
elif iname == 'mm':
    standard_name = 'precipitation_amount'
else:
    standard_name = 'air_temperature'</pre>

tt為溫度、mm為降水，因為我這里還訂正了其他變量所以其他變量的standard_name均設置為與溫度一致。standard_name會影響使用的訂正算法，我們來看pycat的具體代碼如下：

對應變量的處理函數(shù)

所對應的具體方法如下：

• air_temperature
  :meth:absolute_sdm using normal distribution
• precipitation_amount, surface_downwelling_shortwave_flux_in_air
  :meth:relative_sdm using gamma distribution

4. 變量的單位要按照實際情況進行設置。例如溫度是celsius，降水meters，其他自行查閱。

修改源碼

剛才討論了，pycat示例中涉及投影情況，因此氣候模擬的數(shù)據(jù)包含了3個點的數(shù)據(jù)，在進行訂正前進行插值，而我這里僅僅包含一個點，且不會在netcdf中寫入投影信息讓iris直接能夠獲取，因此需要修改涉及投影這一塊。首先在上一步寫入netcdf中增加x和y的coord，分別設置為116.47和39.8。

修改pycat使之單獨對我們這個需求應用等經(jīng)緯度投影，對pycat的io模塊代碼進行調(diào)整如下：

投影修改

這樣如果在netcdf文件中的投影信息時None，則會自動設置為等經(jīng)緯度投影。

修改提方法中插值部分，如果我們的數(shù)據(jù)僅僅包含1個點，則不需要插值處理。

插值處理

當且僅當模式數(shù)據(jù)點不是1個點時進行插值處理。

運行訂正

from pycat.io import Dataset
from pycat.esd import ScaledDistributionMapping
def esd_func(obsdir, moddir, scedir, work_dir):
 #work-dir即結(jié)果輸出目錄，pycat默認在臨時文件夾
    obs = Dataset(os.path.dirname(obsdir), os.path.basename(obsdir))
    mod = Dataset(os.path.dirname(moddir), os.path.basename(moddir))
    sce = Dataset(os.path.dirname(scedir), os.path.basename(scedir))
    qm = ScaledDistributionMapping(obs, mod, sce, work_dir=work_dir)
    qm.correct()

經(jīng)過上述步驟相信我們得到了不同變量不同模式不同情景的各類文件，只需要將文件路徑傳遞給這個函數(shù)即可。

再次聲明，本人既不了解訂正過程，更沒有對訂正結(jié)果沒有做過驗證。

以上僅僅記錄本人一上午的工作記錄。