楊雅茹，趙春雷，李弘毅，邵東航，紀(jì)文政

1.中國科學(xué)院西北生態(tài)環(huán)境資源研究院，蘭州 730000

2.中國科學(xué)院大學(xué)，北京 100000

3.河北省氣象科學(xué)研究所，石家莊 050021

4.河北省氣象與生態(tài)環(huán)境重點(diǎn)實(shí)驗(yàn)室，石家莊 050021

引 言

隨著2022年冬奧會(huì)的熱潮，冰雪運(yùn)動(dòng)日益興起[1]，但I(xiàn)PCC（Intergovernmental Panel on Climate Change，聯(lián)合國政府間氣候變化專門委員會(huì)）第五次評(píng)估報(bào)告指出全球氣候變暖的加劇，加快了積雪、冰川的消融[2-3]。在此背景下未來三十年冬奧賽場的自然積雪條件是否足以支持冰雪運(yùn)動(dòng)的開展逐漸成為研究熱點(diǎn)。積雪物候是描述積雪動(dòng)態(tài)變化的重要指標(biāo)，且利用積雪日數(shù)、積雪初日、終日、積雪期長度等變量能夠反映地區(qū)的水文、積雪、氣候狀況[4]。張家口賽區(qū)作為北京2022年冬奧會(huì)冰雪運(yùn)動(dòng)主要競賽場地之一，該地區(qū)未來三十年的積雪物候狀況需引起關(guān)注。

目前已有眾多學(xué)者針對(duì)冰雪運(yùn)動(dòng)場館的積雪模擬開展了一系列研究，主要包括滑雪資源評(píng)估和積雪期預(yù)測(cè)兩個(gè)方面。其中，大量研究側(cè)重于利用模型預(yù)測(cè)未來氣候變化下滑雪場的自然降雪或人工造雪量。FANG等通過SkiSim 2.0模型預(yù)估了在未來氣候變化背景下我國116個(gè)滑雪場的積雪季長度、造雪要求等積雪條件[5]；PONS等利用GRENBLS模型分析了比利牛斯山滑雪場在未來氣候變化不同的排放情景下積雪變化的脆弱性及可持續(xù)性[6]；WOBUS等將UEB模型結(jié)合未來造雪條件，模擬了美國247個(gè)滑雪場的自然積雪，并對(duì)人工造雪條件進(jìn)行了評(píng)估[7]；除此之外，國外許多學(xué)者還重點(diǎn)研究了未來氣候條件下滑雪場積雪資源對(duì)滑雪季長度的影響[8-9]，而國內(nèi)對(duì)于滑雪場的積雪狀況研究大多建立在對(duì)歷史資料數(shù)據(jù)的分析上，或集中于區(qū)域滑雪的氣候條件評(píng)估[10-14]。在以往的研究中，缺乏針對(duì)滑雪場未來積雪分布及物候的預(yù)測(cè)數(shù)據(jù)集，而本數(shù)據(jù)集旨在提供張家口賽區(qū)未來積雪物候數(shù)據(jù)。

本數(shù)據(jù)集以CMIP6（Coupled Model Intercomparison Project Phase 6）氣候模式預(yù)報(bào)數(shù)據(jù)為驅(qū)動(dòng)數(shù)據(jù)，結(jié)合GBEHM（Geomorphology Based Eco-Hydrological Model）寒區(qū)生態(tài)水文模型制備了2021-2050年北京冬奧會(huì)張家口賽區(qū)積雪范圍數(shù)據(jù)和積雪物候數(shù)據(jù)，本數(shù)據(jù)集有望為滑雪場提供人工造雪和儲(chǔ)雪的參考，指導(dǎo)未來冬奧賽場建設(shè)與運(yùn)營、積雪資源管理，以及利于開展積雪物候與氣象變化的相關(guān)研究。

1 數(shù)據(jù)采集和處理方法

1.1 研究區(qū)

2022 年北京冬奧會(huì)張家口賽區(qū)位于張家口市崇禮區(qū)，地理范圍為 40°52′-41°10′N，115°18′-115°32′E。崇禮區(qū)包括兩大賽區(qū)，分別為云頂滑雪公園和國家跳臺(tái)滑雪中心。賽區(qū)位置示意圖見圖1。

該地區(qū)屬于溫帶半濕潤與半干旱氣候，冬季平均氣溫約-12 °C，年降水量約490 mm，降雪充足，年均降雪量為63.5 mm，雪季為11月到次年2月，積雪期可達(dá)159 天，積雪深度最大可達(dá)1500 mm，滑雪期長達(dá)100天左右。賽區(qū)海拔約1160-2175 m，形成山脈縱橫的地貌特征[15]，其獨(dú)特的氣候和地形有利于冰雪運(yùn)動(dòng)開展[16]。

圖1 北京冬奧會(huì)張家口賽區(qū)示意圖Figure 1 The schematic diagram of Zhangjiakou Venue Cluster of the Olympics Winter Games Beijing 2022

1.2 數(shù)據(jù)源

本數(shù)據(jù)集基于 GBEHM模型[17-21]結(jié)合 CMIP6未來氣候預(yù)報(bào)數(shù)據(jù)[22]制備了北京冬奧會(huì)張家口賽區(qū)未來三十年的積雪范圍數(shù)據(jù)，并在此基礎(chǔ)上制作了相應(yīng)地區(qū)和時(shí)間序列下的積雪物候數(shù)據(jù)。

CMIP6氣候預(yù)報(bào)數(shù)據(jù)由世界氣候研究計(jì)劃（World Climate Research Programme，WCRP）耦合模擬工作組（Working Group on Coupled Modeling，WGCM）組織生產(chǎn)[22]，在其官方網(wǎng)站https://esgfnode.llnl.gov/search/cmip6可進(jìn)行參量選擇和下載。GBEHM模型所需的CMIP6氣候預(yù)報(bào)參量為：下行短波輻射通量、下行長波輻射通量、大氣壓、相對(duì)濕度、比濕、降水、氣溫、風(fēng)速（表1）。

表1 模型所用CMIP6數(shù)據(jù)參量信息Table 1 CMIP6 data parameter information used in GBEHM

1.3 積雪物候定義

積雪物候是表示積雪動(dòng)態(tài)變化的重要指標(biāo)。本研究使用積雪范圍、積雪日數(shù)、積雪初日、積雪終日、連續(xù)積雪天數(shù)、積雪期長度、積雪物候變化趨勢(shì)、積雪物候變化趨勢(shì)顯著性等指標(biāo)進(jìn)行刻畫。根據(jù)以往的研究[23]，各變量可定義如下：

（1）積雪范圍（Snow-Covered Area, SCA）是將識(shí)別為雪的柵格數(shù)與像元大小的乘積；

（2）積雪日數(shù)（Snow-Covered Days, SCD）是在積雪期內(nèi)積雪深度達(dá)到或超過1 cm的日數(shù)之和；

（3）積雪初日（Snow Onset Day, SOD）是在一個(gè)水文年內(nèi)（9月1日至次年8月31日）首次連續(xù)5天被記錄為雪的首日對(duì)應(yīng)日期；

（4）積雪終日（Snow End Day, SED）是在一個(gè)水文年內(nèi)（9月1日至次年8月31日）最后連續(xù)5天被記錄為雪的末日對(duì)應(yīng)日期；

（5）連續(xù)積雪天數(shù)（Snow Duration Days, SDD）是記錄每個(gè)柵格有連續(xù)積雪的天數(shù)，本數(shù)據(jù)從中選取最大值作為該柵格的最大連續(xù)積雪天數(shù)；

（6）積雪期長度（Length of Snow Period, SPL）是積雪終日與積雪初日之差，即為：

1.4 數(shù)據(jù)處理方法

1.4.1 數(shù)據(jù)預(yù)處理

CMIP6氣候預(yù)報(bào)參量作為模型驅(qū)動(dòng)數(shù)據(jù)，需降尺度處理成時(shí)間分辨率為逐小時(shí)，空間分辨率1 km的數(shù)據(jù)。本文采用統(tǒng)計(jì)降尺度方法，具體步驟如下：

（1）大尺度預(yù)報(bào)因子的選擇：選擇近地表氣溫、近地表風(fēng)速、下行短波輻射、上行長波輻射、氣壓、降水、比濕、相對(duì)濕度作為大尺度預(yù)報(bào)因子，利用克里金插值方法將所有大尺度預(yù)報(bào)因子重采樣到1 km的空間分辨率。

（2）建立統(tǒng)計(jì)降尺度模型：采用主成分分析方法對(duì)大尺度預(yù)報(bào)因子參量進(jìn)行降維和濾波處理后，利用線性回歸模型對(duì)參量進(jìn)行迭代回歸訓(xùn)練，進(jìn)而建立最優(yōu)的統(tǒng)計(jì)降尺度模型。

（3）基于獨(dú)立站點(diǎn)觀測(cè)的近地表氣溫、近地表風(fēng)速、下行短波輻射、上行長波輻射、氣壓、降水、比濕、相對(duì)濕度對(duì)統(tǒng)計(jì)降尺度模型進(jìn)行性能檢驗(yàn)。

（4）CMIP6降尺度數(shù)據(jù)生成，將通過檢驗(yàn)的各參量數(shù)據(jù)輸入建立的統(tǒng)計(jì)降尺度模型，最終得到高質(zhì)量降尺度預(yù)報(bào)數(shù)據(jù)。

1.4.2 數(shù)據(jù)生產(chǎn)

（1）積雪范圍數(shù)據(jù)

積雪范圍數(shù)據(jù)的生產(chǎn)將GBEHM水文模型[17-21]作為基礎(chǔ)，將CMIP6氣候預(yù)報(bào)參量降尺度結(jié)果作為驅(qū)動(dòng)數(shù)據(jù)，輸入GBEHM模型得到2021-2050年北京冬奧會(huì)張家口賽區(qū)逐日的積雪范圍數(shù)據(jù)產(chǎn)品。GBEHM寒區(qū)生態(tài)水文模型耦合了分布式水文模型GBHM[17]、CoLM[24]和CLM[25]陸面過程模型，并在風(fēng)吹雪部分、參數(shù)標(biāo)定以及融水徑流追蹤方面做出改進(jìn)[21]。改進(jìn)后的GBEHM水文模型充分考慮了積雪的積累和消融過程，能夠較為精準(zhǔn)地模擬積雪過程。模型框架如圖2。

圖2 GBEHM水文模型框架Figure 2 Basic framework of GBEHM hydrological model

積雪范圍數(shù)據(jù)基于模型中積雪聚集-消融過程模塊，其原理為以積雪質(zhì)能平衡來模擬冰川、積雪、凍土的聚集和消融過程。在GBEHM模型中，積雪能量平衡方程[21]表示為：

式中，Cs是積雪的比熱(J/kg K)，Ts為特定雪層溫度(K)，Tf為水的凍結(jié)溫度(K)，z是節(jié)點(diǎn)到達(dá)地表的距離(m)，θi為冰的體積分?jǐn)?shù)(m3/m3)，ρi為冰的密度(m3/m3)，Lil為冰的熔化潛熱(J/kg)，Ks為雪的導(dǎo)熱系數(shù)(W/K m)，IR為輻射通量(W/m2)。Esur為積雪表面能量交換(W/m2)。

積雪質(zhì)量平衡方程為：

式中，θ1是液態(tài)水的體積分?jǐn)?shù)(m3/m3)，ρ1是液態(tài)水的密度(m3/m3)，U1是液態(tài)水通量(kg·m-2·s-1)。

（2）積雪物候數(shù)據(jù)

積雪物候數(shù)據(jù)的生產(chǎn)是在積雪范圍數(shù)據(jù)產(chǎn)品的基礎(chǔ)上，根據(jù)各個(gè)參數(shù)的定義在Python程序中逐一實(shí)現(xiàn)。積雪日數(shù)由各柵格在積雪期的積雪二值數(shù)據(jù)累加所得；在一個(gè)水文年內(nèi)首次連續(xù)5天發(fā)生降雪時(shí)，5天內(nèi)首日對(duì)應(yīng)的日期賦予柵格，即為積雪初日；積雪終日的生產(chǎn)類似，水文年內(nèi)末次連續(xù)5天發(fā)生降雪時(shí)的最后一天的日期賦予柵格，即可得到積雪終日數(shù)據(jù)。在一個(gè)水文年內(nèi)將連續(xù)降雪的日數(shù)中的最大值賦予柵格即為最大連續(xù)積雪天數(shù)數(shù)據(jù)；積雪期長度數(shù)據(jù)的生產(chǎn)只需將日期對(duì)應(yīng)的積雪終日和積雪初日數(shù)據(jù)作差后賦予柵格。

1.4.3 精度檢驗(yàn)

本數(shù)據(jù)集利用混淆矩陣來驗(yàn)證產(chǎn)品的精度，三個(gè)驗(yàn)證指標(biāo)分別為總體精度（Overall Accuracy，OA）、漏分精度（Omission Error，OE）、多測(cè)精度（Commission Error，CE）[23]。混淆矩陣具體如表2：

表2 混淆矩陣Table 2 Confusion matrix

各指標(biāo)的公式為：

式中，S表示雪，N表示非雪，SS指產(chǎn)品識(shí)別為雪的柵格真值也識(shí)別為雪，SN指產(chǎn)品識(shí)別為雪的柵格真值卻為非雪，NN指產(chǎn)品和真值都識(shí)別為非雪，NS指產(chǎn)品識(shí)別為非雪但真值為雪。

2 數(shù)據(jù)樣本描述

本數(shù)據(jù)集共包括5個(gè)部分，分別為：

（1）積雪范圍：共30幅影像，格式為NetCDF，時(shí)間分辨率為逐日，空間分辨率為1 km。包括2021-2050年張家口賽區(qū)的積雪范圍數(shù)據(jù)。命名方式如：zhangjiakou_snow_binary_2050.nc表示張家口賽區(qū)2050年的積雪范圍分布。圖3為張家口賽區(qū)2021年12月1日積雪分布情況。

圖3 2021年12月1日張家口賽區(qū)積雪范圍Figure 3 The snow cover area of Zhangjiakou Venue Cluster on December 1st, 2021

（2）積雪日數(shù)：共29幅影像，格式為NetCDF，時(shí)間分辨率為日，空間分辨率為1 km。包括各水文年張家口賽區(qū)的積雪日數(shù)影像。命名方式如：2049_2050_SCD.nc表示2049年9月1日至2050年8月31日積雪日數(shù)總和。其中：2049-2050表示本數(shù)據(jù)的時(shí)間序列為2049年9月1日至2050年8月31日；SCD表示積雪日數(shù)。圖4為2049年水文年張家口賽區(qū)累計(jì)積雪日數(shù)分布。

圖4 2049年水文年張家口賽區(qū)累計(jì)積雪日數(shù)分布Figure 4 The accumulated snow cover days in Zhangjiakou Venue Cluster in the hydrological year of 2049

（3）積雪初日終日：共58幅影像，格式為NetCDF，時(shí)間分辨率為日，空間分辨率1 km。包括積雪初日和積雪終日兩部分，命名方式如：2049_2050_SOD.nc表示賽區(qū)2049年9月1日至2050年8月31日（一個(gè)水文年）的積雪初日影像。其中：2049_2050的命名與積雪日數(shù)相同；SOD表示積雪初日；SED表示積雪終日。圖5為2049年水文年張家口賽區(qū)積雪初日分布，圖6為2049年水文年張家口賽區(qū)積雪終日分布。

圖5 2049年水文年張家口賽區(qū)積雪初日分布Figure 5 The snow onset day in Zhangjiakou Venue Cluster in the hydrological year of 2049

圖6 2049年水文年張家口賽區(qū)積雪終日分布Figure 6 The snow end Day in Zhangjiakou Venue Cluster in the hydrological year of 2049

（4）最大連續(xù)積雪天數(shù)：共29幅影像，格式為NetCDF。包括2021-2050年各水文年張家口賽區(qū)最大連續(xù)積雪天數(shù)影像。文件命名方式如：2049_2050_MAX_SDD.nc表示2049年9月1日至2050年8月31日連續(xù)積雪天數(shù)最大值分布影像數(shù)據(jù)。其中：2049_2050的命名與積雪日數(shù)相同；MAX_SDD表示最大連續(xù)積雪天數(shù)。圖7為2049年水文年張家口賽區(qū)最大連續(xù)積雪天數(shù)分布。

圖7 2049年水文年張家口賽區(qū)最大連續(xù)積雪天數(shù)分布Figure 7 The maximum snow duration days in Zhangjiakou Venue Cluster in the hydrological year of 2049

（5）積雪期長度：共29幅影像，格式為NetCDF。包括2021-2050年各水文年張家口賽區(qū)積雪期長度影像。命名方式如：2049_2050_SPL.nc表示2049年9月1日至2050年8月31日積雪期長度。其中：2049_2050的命名與積雪日數(shù)相同；SPL表示積雪期長度。圖8為2049年水文年張家口賽區(qū)積雪期長度分布。

圖8 2049年水文年張家口賽區(qū)積雪期長度分布Figure 8 The length of snow period in Zhangjiakou Venue Cluster in the hydrological year of 2049

3 數(shù)據(jù)質(zhì)量控制和評(píng)估

利用SNOMAP算法和Landsat8 OLI影像獲得2021年張家口賽區(qū)積雪范圍二值數(shù)據(jù)作為參考數(shù)據(jù)進(jìn)行質(zhì)量評(píng)估。SNOMAP是Hall[26]等提出的積雪識(shí)別算法，在此基礎(chǔ)上WANG等[27]對(duì)其進(jìn)行了改進(jìn)，改進(jìn)后的算法在林區(qū)和陰影區(qū)積雪提取精度上有了很大提高，總體分類精度可達(dá)97.5%，適用于森林茂密的張家口賽區(qū)。Landsat8 OLI影像的選擇條件為圖像覆蓋張家口賽區(qū)且云量小于10%，最終篩選出滿足條件的六幅影像，對(duì)應(yīng)日期為2021年2月2日、2021年4月7日和23日。

由于驗(yàn)證數(shù)據(jù)與積雪范圍數(shù)據(jù)空間分辨率不一致，為了避免重采樣帶來的誤差，驗(yàn)證點(diǎn)選取賽區(qū)內(nèi)氣象站點(diǎn)（兩個(gè)賽區(qū)內(nèi)共22個(gè)氣象站，平均海拔1200 m）所在柵格。分別提取積雪范圍產(chǎn)品與驗(yàn)證數(shù)據(jù)對(duì)應(yīng)氣象站點(diǎn)的柵格值，3天共有66對(duì)數(shù)據(jù)可用于驗(yàn)證。檢驗(yàn)方法使用1.4.3中的混淆矩陣方法，驗(yàn)證情況如表3所示。

表3 驗(yàn)證情況Table 3 Verification

4 數(shù)據(jù)價(jià)值

張家口賽區(qū)是2022年北京冬奧會(huì)主要競賽場地之一，該地區(qū)積雪物候變化對(duì)未來冰雪賽事和賽區(qū)布置有著重要影響。隨著氣候變暖的加劇，積雪、冰川的融化也隨之加劇[1]，過去二十年來張家口賽區(qū)積雪面積呈現(xiàn)下降趨勢(shì)，滑雪產(chǎn)業(yè)也受到?jīng)_擊，這與冰雪運(yùn)動(dòng)日益興起的趨勢(shì)不相適應(yīng)。本數(shù)據(jù)集則為進(jìn)一步了解未來張家口賽區(qū)積雪物候變化提供了參考，可為冰雪賽事的開展提供人工降雪量的指導(dǎo)，有助于未來滑雪場的建設(shè)、運(yùn)營與維護(hù)。

5 數(shù)據(jù)使用方法和建議

數(shù)據(jù)文件夾中為各積雪物侯數(shù)據(jù)，可根據(jù)需要選取參數(shù)和相應(yīng)年份。可使用Panoply、PyNcView等軟件查看數(shù)據(jù)，本數(shù)據(jù)集也提供將NetCDF轉(zhuǎn)為TIFF的代碼，可以在ArcGIS中進(jìn)行柵格運(yùn)算以及生成地圖，進(jìn)行未來積雪物候趨勢(shì)分析等操作。本數(shù)據(jù)集還可以用于氣候響應(yīng)等相關(guān)研究。