曹 娟，姚曉軍**，靳惠安,2，張調(diào)風(fēng)，高永鵬，張大弘，趙全寧

(1：西北師范大學(xué)地理與環(huán)境科學(xué)學(xué)院，蘭州 730070)(2：甘肅林業(yè)職業(yè)技術(shù)學(xué)院，天水 741020)(3：青海省氣候中心，西寧 810000)(4：云南大學(xué)國際河流與生態(tài)安全研究院，昆明 650091)(5：青海省氣象科學(xué)研究所，西寧 810000)

政府間氣候變化專業(yè)委員會(huì)(IPCC)第五次評(píng)估報(bào)告(AR5)指出，過去130多年全球平均氣溫升高了0.85℃，主要升溫時(shí)間段為1983-2012年的30年間[1]. 作為大氣圈、冰凍圈和水圈的聯(lián)系紐帶，湖泊通過地表水的匯聚和蒸發(fā)以及地下水匯/滲流等過程參與地表水循環(huán)，并對(duì)氣候變化響應(yīng)敏感[2-3]. 青藏高原地區(qū)98個(gè)氣象站點(diǎn)觀測(cè)數(shù)據(jù)表明，1982-2012年平均氣溫升幅高達(dá)1.9℃，是全球平均升溫幅度的2倍[4]. 青藏高原的自然環(huán)境和生態(tài)系統(tǒng)在全球占有特殊地位，該地區(qū)也被視作氣候變化敏感區(qū)[5-6]，如青藏高原部分湖泊湖冰物候特征近期表現(xiàn)為開始凍結(jié)和完全凍結(jié)日期推遲、封凍期減少趨勢(shì)[7-10]. 青海湖是維系青藏高原東北部生態(tài)安全的重要水體，在保護(hù)脆弱的高原湖泊濕地資源、阻擋西部荒漠化向東蔓延和保護(hù)鳥類及裸鯉等生物資源等方面發(fā)揮著不可替代作用. 隨著西部大開發(fā)和絲綢之路經(jīng)濟(jì)帶的建設(shè)，尤其是蘭新高速鐵路的開通，青海湖已成為旅游熱點(diǎn)景區(qū)，而在冬季時(shí)有車輛掉入湖中事故發(fā)生. 因此，開展青海湖湖冰厚度研究不僅對(duì)于認(rèn)識(shí)氣候變暖背景下的湖冰響應(yīng)規(guī)律具有重要的理論價(jià)值，而且對(duì)于制定科學(xué)合理的冬季青海湖管理措施具有現(xiàn)實(shí)意義.

作為湖泊凍結(jié)期的重要物理參數(shù)之一，湖冰厚度不僅反映了水-氣界面能量交換強(qiáng)度和物質(zhì)遷移過程，而且具有重要的生態(tài)價(jià)值和經(jīng)濟(jì)價(jià)值[11]，較其時(shí)間屬性(如封凍期和消融期)難以監(jiān)測(cè)[12]. 傳統(tǒng)的冰厚測(cè)量主要采用接觸式方法，如鉆孔取冰、電熱絲融冰以及壓力傳感器測(cè)厚法等，此類方法雖測(cè)量精度較高，但存在效率低、數(shù)據(jù)量少、損壞冰層等缺點(diǎn)[13]. 隨著遙感技術(shù)快速發(fā)展和測(cè)量儀器設(shè)備的改進(jìn)，一些非接觸式冰厚測(cè)量方法開始出現(xiàn)，如基于冰和水的電導(dǎo)率差別原理采用電磁感應(yīng)儀測(cè)量冰厚，在水下安裝仰視聲吶通過發(fā)射聲脈沖根據(jù)回波時(shí)延差確定冰厚，該類方法僅能用于定點(diǎn)測(cè)量，難以反映湖面冰厚整體分布情況[14-15]；衛(wèi)星數(shù)據(jù)(如CryoSat-2衛(wèi)星搭載的Ku波段SIRAL)雖可用于監(jiān)測(cè)大范圍湖泊冰厚，但對(duì)薄冰反演效果較差且難以剔除雪對(duì)反照率的影響[16]. 從熱力學(xué)角度來看，湖冰是開闊水域在大氣和水的相互作用下通過內(nèi)能和熱能轉(zhuǎn)化達(dá)到臨界狀態(tài)的產(chǎn)物，因此其厚度亦可借助熱力學(xué)模型來模擬，如度日法湖冰生長模型[17]、CLIMo模型[18]、Zubov模型[19]、Mylake模型[20]等. 目前國際上對(duì)于湖冰厚度的研究主要集中在長時(shí)間序列模擬結(jié)果分析及與實(shí)測(cè)數(shù)據(jù)的對(duì)比，即模擬精度驗(yàn)證，或單點(diǎn)測(cè)量短時(shí)間序列的冰厚特征分析. 如Dibike等[20]基于NARR數(shù)據(jù)采用Mylake模型模擬了北美地區(qū)(40°～75°N)的湖泊冰厚，結(jié)果表明最大冰厚普遍分布在40 cm左右且隨著湖深的增加冰厚呈現(xiàn)減小趨勢(shì)，其中位于加拿大東西部及較高海拔地區(qū)的湖泊冰厚更大. Duguay等[18]以Barrow(北極)、Poker Flat Research Range(亞北極)、Churchill(高北森林)3個(gè)代表性地區(qū)的湖泊為研究對(duì)象，基于研究區(qū)附近的氣象站點(diǎn)數(shù)據(jù)采用CLIMO模型進(jìn)行冰厚模擬，發(fā)現(xiàn)隨雪深變化模擬的最大冰層厚度介于159～227 cm，整個(gè)冬季平均絕對(duì)誤差為2 cm. 就青海湖而言，目前系統(tǒng)性的冰厚觀測(cè)仍很少，僅在湖區(qū)東南側(cè)的下設(shè)水文站有冰厚觀測(cè)記錄，因人工觀測(cè)點(diǎn)的單一性、湖冰形成的時(shí)空異質(zhì)性，以及人工觀測(cè)頻次較低和觀測(cè)位置不定，對(duì)于青海湖湖冰厚度時(shí)空分布及演化過程仍知之甚少. 本文基于實(shí)測(cè)鉆孔冰厚數(shù)據(jù)、飛航式測(cè)冰雷達(dá)數(shù)據(jù)和較高時(shí)空分辨率的Terra MODIS遙感影像，采用度日法湖冰生長模型模擬青海湖冰厚變化，開展冰厚精度評(píng)價(jià)，并分析青海湖2000-2019年冰厚時(shí)空變化規(guī)律，從而為深入研究青海湖湖冰生消機(jī)理提供參考.

1 研究區(qū)概況

青海湖地處青藏高原東北部(36°32′～37°15′N，99°36′～100°46′E)，跨青海省剛察、海晏和共和3縣，東西長約109 km，南北寬約40 km(圖1). 根據(jù)青海省第一次水利普查公布的2013年測(cè)量數(shù)據(jù)，當(dāng)水位在3193.5 m時(shí)，青海湖水域面積為4294 km2，平均水深18.3 m，最大水深26.6 m，蓄水量785.2×108m3，湖水呈弱堿性，pH值為9.23，含鹽量為14.13 g/L[21]，是中國最大的內(nèi)陸湖和咸水湖. 青海湖是新構(gòu)造運(yùn)動(dòng)導(dǎo)致周圍山地上升隆起，外泄通道堵塞形成的構(gòu)造斷陷湖，湖中現(xiàn)有海心山和三塊石兩個(gè)湖心島，湖區(qū)東面自北向南依次分布著尕海、沙島湖、海晏灣和耳海4個(gè)子湖. 青海湖入湖徑流較大的河流主要有布哈河、沙柳河以及哈爾蓋河等，其中布哈河徑流量最大，約占入湖徑流總量的60%[22]. 湖區(qū)處于我國東部季風(fēng)區(qū)、西北干旱區(qū)和西南高寒區(qū)的交匯地帶，屬典型的高原半干旱高寒氣候[23]，年平均氣溫介于-1.0～1.5℃，并由東南向西北逐漸降低，年降水量為350～400 mm，其中夏季降水達(dá)60%以上，年蒸發(fā)量可達(dá)1300～2000 mm. 因湖水中含有無機(jī)鹽類，青海湖凍結(jié)溫度比0℃稍低，通常于每年12月中旬至1月初開始凍結(jié)，12月底或次年1月初完全被冰層覆蓋，3月下旬冰層逐漸消融，至4月初或中旬完全消融. 1958-1983年沙陀寺水文站觀測(cè)結(jié)果表明，青海湖冰厚一般為0.5 m左右，最厚達(dá)0.7 m，且由岸邊向湖心方向逐漸變薄[24]. 基于遙感數(shù)據(jù)的監(jiān)測(cè)結(jié)果表明2000-2016年青海湖封凍期和完全封凍期平均為(88±20)和(77±20)d,湖冰存在期和湖冰消融期持續(xù)時(shí)間平均為(108±18)和(10±6)d[10].

圖1 青海湖及湖冰厚度測(cè)區(qū)(點(diǎn))分布

2 數(shù)據(jù)與方法

2.1 湖冰鉆孔數(shù)據(jù)與飛航式雷達(dá)測(cè)厚數(shù)據(jù)

2019年1月13日-3月24日期間，在青海湖開展湖冰鉆孔測(cè)厚和飛航式測(cè)冰雷達(dá)測(cè)厚試驗(yàn)，測(cè)厚樣點(diǎn)分布及測(cè)飛區(qū)域如圖1所示. 出于工作人員安全考慮和受無人機(jī)電池電量限制，測(cè)厚樣點(diǎn)及測(cè)飛區(qū)域繞湖區(qū)周邊布設(shè). 為保障鉆孔數(shù)據(jù)與飛航式測(cè)冰雷達(dá)測(cè)厚數(shù)據(jù)時(shí)空一致性，兩項(xiàng)試驗(yàn)同時(shí)同地開展，其中鉆孔取冰厚度采用L型數(shù)字式量冰尺(精度為0.01 mm)測(cè)量(圖2a)；飛航式測(cè)冰雷達(dá)測(cè)厚采用IGPR-30雷達(dá)系統(tǒng)，以無人機(jī)為載具平臺(tái)，并配備激光測(cè)距、GPS差分定位和高清攝像頭等模塊，設(shè)備雷達(dá)中心頻率為400 MHz，采樣時(shí)間間隔為2 ps，探測(cè)精度可達(dá)mm級(jí)，冰層探測(cè)厚度>6 m(圖2b). 該設(shè)備的冰厚測(cè)量原理是雷達(dá)向下發(fā)射高頻電磁波，當(dāng)?shù)竭_(dá)空氣-冰界面和冰-水界面時(shí)發(fā)生反射，反射的電磁波信號(hào)被天線接收，假設(shè)雷達(dá)波發(fā)射至各界面和返回被天線接收所有時(shí)間相同，則冰厚可根據(jù)雙程走時(shí)、冰層介電常數(shù)及電磁波在空氣中的傳播速度計(jì)算得到. 2019年1月13日-14日，在湖區(qū)西北側(cè)及鳥島附近各布設(shè)一個(gè)200 m×200 m觀測(cè)場(chǎng)，并按50 m等間距在冰面布設(shè)標(biāo)記點(diǎn)，首先采用飛航式測(cè)冰雷達(dá)對(duì)各標(biāo)記點(diǎn)冰厚進(jìn)行測(cè)量，然后鉆孔取冰測(cè)量冰厚并用于校正飛航式測(cè)冰雷達(dá)模擬冰厚的介電常數(shù). 2019年1月30日-31日，在仙女灣、哈爾蓋河入湖口、尕海以及湖區(qū)東部進(jìn)行鉆孔取冰測(cè)厚(5個(gè)測(cè)點(diǎn))和雷達(dá)測(cè)厚(6個(gè)測(cè)區(qū))試驗(yàn). 2019年3月1日-3日，在泉吉、鳥島、布哈河口以及湖區(qū)西南部泉灣、峽灣、黑馬河、三江源祭海臺(tái)、二郎劍等地共測(cè)點(diǎn)10個(gè)，雷達(dá)測(cè)厚區(qū)共8個(gè). 2019年3月22日-24日，在二郎劍、黑馬河、布哈河口、鳥島等地共測(cè)飛7個(gè)區(qū)域. 為更好地了解不同方向湖冰厚度分布情況，飛航式測(cè)冰雷達(dá)測(cè)量時(shí)隨機(jī)采用多種航線(圖2c～f)，并在每個(gè)測(cè)區(qū)各測(cè)飛2次冰厚. 以上湖冰厚度實(shí)測(cè)數(shù)據(jù)用于擬合度日法湖冰生長模型參數(shù)及模擬結(jié)果誤差分析.

圖2 青海湖冰厚測(cè)量照片與無人機(jī)航飛路線

圖3 ERA5氣候數(shù)據(jù)與觀測(cè)值的溫度對(duì)比

2.2 ERA5氣候再分析資料

采用2000-2019年ERA5 Climate Reanalysis數(shù)據(jù)集中青海湖地區(qū)距地面2 m處氣溫格點(diǎn)數(shù)據(jù)作為研究區(qū)氣溫背景資料，參考青海湖歷年開始凍結(jié)和湖冰完全消融日期[10]，截取該數(shù)據(jù)集中的每年10月1日至翌年4月30日數(shù)據(jù). ERA5 Climate Reanalysis數(shù)據(jù)集是歐洲中期天氣預(yù)報(bào)中心(ECMWF)于2017年發(fā)布的全球氣候再分析資料第五代產(chǎn)品，提供了對(duì)全球近期氣候的綜合數(shù)字描述，時(shí)間分辨率為12 h，空間分辨率為0.125°×0.125°，該數(shù)據(jù)從歐洲中期天氣預(yù)報(bào)中心網(wǎng)站(https://www.ecmwf.int)下載獲得. 因ERA5 Climate Reanalysis數(shù)據(jù)集為NC格式且文件較多，利用Python語言編寫腳本程序，將文件格式批量轉(zhuǎn)換為GRID格式并計(jì)算逐格點(diǎn)日均溫?cái)?shù)據(jù).

距青海湖最近的國家基準(zhǔn)站僅有剛察氣象站(圖1)，通過將2000-2019年每年10月1日至翌年4月30日站點(diǎn)實(shí)測(cè)氣溫?cái)?shù)據(jù)與該站點(diǎn)所在ERA5氣候再分析數(shù)據(jù)格點(diǎn)對(duì)應(yīng)分析，發(fā)現(xiàn)二者趨勢(shì)較為一致(圖3). 分析結(jié)果表明實(shí)測(cè)數(shù)據(jù)與ERA5氣候數(shù)據(jù)相關(guān)系數(shù)為0.88(P<0.05)，均方根誤差(RMSE)為5.17℃，兩組數(shù)據(jù)間具有較好的相關(guān)性且相對(duì)誤差較小. 相較于用該氣象站點(diǎn)數(shù)據(jù)來表征青海湖整個(gè)湖面氣溫?cái)?shù)據(jù)，ERA5氣候再分析數(shù)據(jù)分辨率較高，能夠更客觀地反映青海湖湖面氣溫變化.

2.3 光學(xué)衛(wèi)星遙感影像

為了解青海湖湖面凍結(jié)狀況及湖冰空間分布情況，從美國國家航空航天局陸面分布式數(shù)據(jù)中心(NASALPDA)網(wǎng)站(https://data.giss.nasa.gov)共下載2888景MODIS MOD09GQ數(shù)據(jù)產(chǎn)品，該數(shù)據(jù)時(shí)間分辨率為1 d，空間分辨率為250 m. 此外，從美國地質(zhì)勘探局(USGS)網(wǎng)站(http://glovis.usgs.gov)共下載19景2001-2019年青海湖封凍期且云量覆蓋較小的Landsat TM/ETM+/OLI遙感影像，用于人工目視解譯提取青海湖湖冰范圍，以評(píng)估基于MODIS MOD09GQ數(shù)據(jù)提取的湖冰面積占比精度. 鑒于MODIS MOD09GQ數(shù)據(jù)量較大且需要進(jìn)行影像拼接，首先利用MRT工具對(duì)其進(jìn)行批量拼接和幾何校正，然后通過編寫的Python腳本程序提取湖冰并計(jì)算其面積占比.

2.4 度日法湖冰生長模型

湖冰可被視作水-氣界面的絕緣膜，湖冰生消改變大氣與湖水之間的熱量、動(dòng)量和物質(zhì)交換[25]，如湖冰在生長過程中由于熱能損失可改變湖面周圍空氣的熱量. 早在19世紀(jì)，Stefan[26]就認(rèn)為冰的生長過程是一種特殊的熱傳導(dǎo)過程，冰底結(jié)冰釋放的熱量在線性溫度梯度情況下通過冰傳導(dǎo)出去. Lepp?ranta[25]在其研究中指出，Stefan提出的理想狀態(tài)下度日法湖冰生長模型的應(yīng)用前提是不考慮熱慣性和內(nèi)部熱源，冰層底部水中沒有熱通量且冰層頂部的溫度已知，在模擬過程中只考慮溫度變化對(duì)冰厚的影響，則冰體每日的生長速率為：

ρi·L·dH/dt=ki(Tf-T0)/H

(1)

式中，H為冰體厚度；ρi為冰體密度；L為冰體融化潛熱；ki為冰體導(dǎo)熱系數(shù)；T0為冰體上表面溫度；Tf為冰體下表面溫度. 若設(shè)t=0時(shí)H=H0為初始條件，則式(1)的解析式為：

(2)

(3)

(4)

若時(shí)間間隔t為1 d，S則被定義為水體處于凍結(jié)期的負(fù)積溫(℃·d)，冰厚H計(jì)算可簡(jiǎn)化為：

(5)

3 結(jié)果與分析

3.1 基于實(shí)測(cè)數(shù)據(jù)的青海湖冰厚時(shí)空變化特征

表1為2019年1-3月青海湖實(shí)測(cè)冰厚數(shù)據(jù)，測(cè)區(qū)編號(hào)如圖1所示. 由表1可知，在觀測(cè)時(shí)段內(nèi)，各采樣點(diǎn)的冰厚總體呈增長趨勢(shì)，但各月份冰厚增長速率不同. 1月30日測(cè)量數(shù)據(jù)顯示湖區(qū)平均冰厚為40.2 cm，3月2日平均冰厚為43.5 cm，3月22日為51.4 cm，整個(gè)3月平均增長速率為0.30 cm/d，較2月份(0.12 cm/d)增長迅速. 但結(jié)果受該時(shí)間段內(nèi)所測(cè)區(qū)域冰厚影響較大(可能所測(cè)地區(qū)均為冰厚較大或較小區(qū)域). 空間尺度上，實(shí)測(cè)冰厚大體呈現(xiàn)北厚南薄、東厚西薄特點(diǎn). 就各測(cè)區(qū)平均增長速率而言，測(cè)區(qū)1(2019-01-30)、3(2019-01-31)、7(2019-03-01)僅有同一時(shí)間段內(nèi)采樣記錄(未計(jì)算平均增長速率)，平均冰厚分別為42.7、46.4、49.2 cm. 有兩段時(shí)間采樣記錄的測(cè)區(qū)有9個(gè)，可分為3個(gè)主要變化區(qū)：(1)冰厚快速增長區(qū)包括測(cè)區(qū)6、13、9，冰厚平均增長速率均大于0.30 cm/d，分別為0.46、0.34、0.32 cm/d. (2)冰厚平緩增長區(qū)平均增長速率均分布在0.25 cm/d左右，包括測(cè)區(qū)11(0.28 cm/d)、10(0.25 cm/d)、8(0.21 cm/d). (3)冰厚緩慢增長區(qū)包括測(cè)區(qū)12、14，冰厚平均增長速率僅分別為0.15、0.11 cm/d. 測(cè)區(qū)2、4、5、15將3個(gè)采樣時(shí)間均包括在內(nèi)，其對(duì)應(yīng)平均冰厚增長速率依次減小，分別為0.30、0.24、0.12、0.08 cm/d.

表1 2019年1-3月青海湖實(shí)測(cè)冰厚數(shù)據(jù)

“航線長度”表示飛航式測(cè)冰雷達(dá)測(cè)飛長度，“-”表示鉆孔測(cè)點(diǎn)無長度.

3.2 2019年青海湖冰厚模擬結(jié)果及精度評(píng)估

圖4為2018-2019年冬季青海湖逐日氣溫與冰厚模擬值變化. 由該圖可知，2018年11月8日青海湖開始凍結(jié)，當(dāng)日湖面平均氣溫為-6.44℃，2019年3月29日冰厚達(dá)到最大值時(shí)氣溫為-1.21℃. 連續(xù)負(fù)值氣溫持續(xù)天數(shù)為141 d，累計(jì)負(fù)積溫-1537.64℃，最大冰厚為51.05 cm. 將湖區(qū)相同位置的冰厚實(shí)測(cè)值與模擬值對(duì)應(yīng)分析，兩者之間的相關(guān)系數(shù)R2為0.92(P<0.05). 在模擬期內(nèi)，青海湖冰厚平均增長速率為0.34 cm/d，負(fù)積溫絕對(duì)值每平均累加10℃冰厚增加0.25 cm，其中2018年11月8日-11月21日及2018年12月1日-12月11日的冰厚日平均增長速率分別為0.71、0.59 cm/d，均大于0.50 cm/d；2019年2月21日-3月29日的平均增長速率為0.14 cm/d；其余時(shí)間段內(nèi)湖冰厚度日平均增長速率介于0.20～0.50 cm/d之間.

圖4 2018-2019年冬季青海湖氣溫變化與冰厚模擬值

為檢驗(yàn)度日法湖冰生長模型模擬的青海湖冰厚精度，對(duì)模擬結(jié)果與校正后的雷達(dá)測(cè)冰厚度數(shù)據(jù)進(jìn)行對(duì)比(圖5)，結(jié)果表明平均誤差為±2 cm，但湖區(qū)不同部位冰厚模擬精度有所差異，其中誤差較大區(qū)域主要是布哈河和哈爾蓋河入湖口處，其次是湖區(qū)南側(cè)和西南側(cè). 如1月13日鳥島附近冰厚誤差可達(dá)10.45 cm，這可能與湖冰形成初期厚度較薄以及布哈河淡水補(bǔ)給有關(guān)，布哈河水流匯入湖中不僅通過熱量交換使湖水溫度有所上升，而且可通過動(dòng)力過程(如漩渦、湍流)使冰層底部不斷被侵蝕，從而導(dǎo)致模擬數(shù)值大于湖冰實(shí)際厚度. 隨著湖冰厚度增加，模擬精度有所改善，如鳥島附近同一測(cè)區(qū)在1月30日和3月1日誤差分別為0.42 cm和1.35 cm；湖區(qū)其他部位冰厚誤差大多介于1～3 cm之間. 圖5的另一個(gè)顯著特征是，度日法湖冰生長模型模擬的冰厚值在3月中旬前普遍高估，而之后則存在低估現(xiàn)象，如哈爾蓋河入湖口處冰厚誤差在1月31日為3.28 cm，3月24日則為-2.23 cm，這表明湖冰生消不僅與氣溫相關(guān)，而且受湖面風(fēng)速、太陽輻射、水深、湖流等要素影響[10]. 整體而言，在缺乏多要素觀測(cè)條件下，采用度日法湖冰生長模型模擬的湖冰厚度雖與實(shí)測(cè)冰厚有所差異，但其精度亦在可接受范圍內(nèi)，這為基于氣溫觀測(cè)數(shù)據(jù)預(yù)測(cè)青海湖湖冰厚度提供了一種捷徑.

圖5 2019年青海湖冰厚模擬及誤差

3.3 2000-2019年青海湖冰厚時(shí)空變化特征

3.3.1 青海湖湖冰厚度時(shí)間變化特征 圖6為2000-2019年青海湖年際(當(dāng)年11月至翌年3月)及各月平均冰厚變化圖. 在年際尺度上，近20年青海湖平均冰厚呈波動(dòng)變化趨勢(shì)，平均冰厚介于32～37 cm. 根據(jù)平均冰厚大小，可將此時(shí)段劃分為3個(gè)階段：(1)2000-2007年冰厚穩(wěn)定期，此時(shí)段青海湖平均冰厚為34.65 cm，最大變幅為2.26 cm；(2)2008-2016年冰厚劇烈波動(dòng)期，平均冰厚由2008年的32.78 cm快速增長至2010年的36.43 cm；之后湖冰厚度穩(wěn)定處于高值，2010-2014年青海湖平均厚度為36.18 cm，最大變幅僅1.38 cm，其中2013年湖冰平均厚度(37.05 cm)達(dá)到近20年最大值；2015-2016年湖冰平均厚度出現(xiàn)驟降，兩年平均冰厚僅32.87 cm，其中2016年為研究時(shí)段平均冰厚最小(32.08 cm)的年份；(3)2017-2018年冰厚恢復(fù)期，此時(shí)段青海湖湖冰平均厚度為34.79 cm，略高于第一階段.

由圖6亦可知，青海湖在11月份平均冰厚較小(10.61 cm)，12月份和1月份平均冰厚快速增加，可分別 達(dá)24.45和37.33 cm，增厚速率分別為0.45和0.41 cm/d；之后，2月份(0.29 cm/d)和3月份(0.14 cm/d)冰厚變化速率較小，這兩個(gè)月的平均冰厚分別為45.55和50.03 cm. 這表明青海湖在凍結(jié)前期(12月-翌年1月)冰厚增長迅速，而中期(2-3月)冰厚增長變得緩慢. 湖冰厚度在各月的差異直接導(dǎo)致其平均冰厚在年際上的不同，如2010-2014年平均冰厚最大主要與這一時(shí)段1-3月冰厚普遍較大有關(guān).

圖6 2000-2019年青海湖封凍期年際及各月平均冰厚變化

3.3.2 青海湖湖冰厚度空間變化特征 基于2000-2019年ERA5 Climate Reanalysis氣溫?cái)?shù)據(jù)，采用度日法湖冰生長模型計(jì)算同時(shí)段青海湖湖冰厚度可知，近20年間青海湖多年平均冰厚在空間上分布不均，總體上介于23～41 cm之間. 受ERA5 Climate Reanalysis氣溫?cái)?shù)據(jù)空間分辨率影響，模擬的湖冰平均厚度與實(shí)際冰厚必然存在差異，如后者在空間上應(yīng)為連續(xù)變化，但其仍能反映出青海湖湖冰厚度在空間上的顯著差異性，即整體呈現(xiàn)北厚南薄、東厚西薄的特點(diǎn)(圖7). 青海湖冰厚空間分布格局與其凍結(jié)消融空間模式密切相關(guān). 已有研究表明，青海湖首先在鄰近海晏灣的東側(cè)邊緣開始凍結(jié)，之后湖區(qū)東北部、西北部開始封凍，湖冰由湖岸逐漸擴(kuò)張至湖心[10]. 相應(yīng)地，凍結(jié)越早和凍結(jié)期越長的區(qū)域其湖冰越厚. 此外，湖冰厚度還與湖泊水深、入湖河流徑流及湖下熱流分布有關(guān)，如湖區(qū)西側(cè)雖凍結(jié)較早，但布哈河河水匯入使該區(qū)域冰厚較?。欢^(qū)西南側(cè)(三江源祭海臺(tái)周邊)湖冰最薄的原因主要是該區(qū)域分布有溫泉，高溫水流涌入使湖水溫度高于其他湖區(qū)，并不斷融化湖冰，使湖冰厚度難以增大.

圖7 2000-2019年青海湖平均冰厚空間分布及變化趨勢(shì)

青海湖各區(qū)域2000-2019年平均冰厚的變化幅度不同. 以湖泊質(zhì)心為中心將青海湖劃分為4個(gè)區(qū)域，從左上方開始按順時(shí)針方向分別稱為Ⅰ、Ⅱ、Ⅲ、Ⅳ區(qū)(圖7). 就多年平均變化速率而言，Ⅱ區(qū)的平均變化速率最快(-0.046 cm/a)，Ⅳ區(qū)變化速率最慢(-0.026 cm/a)，Ⅰ區(qū)與Ⅲ區(qū)平均變化速率分別為-0.031、-0.041 cm/a. 年際變化速率最大的區(qū)域位于湖區(qū)東南角(-0.051 cm/a)，最小的區(qū)域則位于西南角(-0.019 cm/a)，各區(qū)的變化趨勢(shì)均為負(fù)值，表明2000-2019年青海湖各區(qū)域平均冰厚呈減小趨勢(shì)，整體而言湖區(qū)西部變化幅度小于東部. 這是由于青海湖湖區(qū)東側(cè)沒有受到淡水補(bǔ)給以及盛行西北風(fēng)的擾動(dòng)，其冰厚的變化更能體現(xiàn)出隨氣候變化的較大變化幅度，而湖區(qū)西側(cè)冰層受到此類干擾因素較多，在一定程度上減緩了氣候變化所導(dǎo)致的冰厚變化，從而表現(xiàn)出更為穩(wěn)定的變化幅度.

3.4 青海湖冰厚變化與湖冰物候特征的關(guān)系

3.4.1 2000-2019年青海湖湖冰物候特征 圖8為2000-2019年每年11月至翌年4月青海湖逐日湖冰面積占比及日均溫變化曲線，顯然湖泊凍結(jié)過程較為緩慢，而湖冰消融過程則較為迅速. 具體而言，每年11月份青海湖湖冰面積占比介于1%～5%且波動(dòng)劇烈，這主要與湖冰初生過程中冰層較薄，新生湖冰或碎冰極易受氣溫驟升或大風(fēng)天氣影響有關(guān)，導(dǎo)致湖冰在形成初期出現(xiàn)反復(fù)凍結(jié)消融現(xiàn)象[11]. 至12月初，當(dāng)湖冰面積占比達(dá)到5%時(shí)，湖冰趨于穩(wěn)定并隨著溫度降低湖冰范圍逐漸擴(kuò)大直至完全凍結(jié)(12月下旬或1月上旬). 3月中下旬，隨著氣溫上升湖冰開裂消融，至4月上旬或下旬湖冰完全消融.

圖8 2000-2019年青海湖湖冰面積占比及日均溫變化

前人在研究青海湖湖冰物候特征時(shí)以湖冰面積占比達(dá)到10%作為湖泊開始凍結(jié)和完全消融日期判別依據(jù)[10]，此閾值最早由Reed等[28]在研究北美阿拉斯加西南地區(qū)湖冰物候特征時(shí)所提出，其認(rèn)為當(dāng)湖冰面積占比小于10%時(shí)湖冰會(huì)出現(xiàn)反復(fù)消融凍結(jié)現(xiàn)象. 由圖8可知，當(dāng)青海湖的湖冰面積占比接近5%時(shí)，湖冰凍結(jié)或消融發(fā)展趨勢(shì)將不會(huì)改變，因此將湖冰面積占比連續(xù)3 d達(dá)到5%定義為湖泊開始凍結(jié)日期和湖冰完全消融日期，由開始凍結(jié)日期至湖冰面積占比首次達(dá)到95%的時(shí)間間隔定義為凍結(jié)時(shí)長，由湖冰面積占比最后達(dá)到95%日期至完全消融日期定義為消融時(shí)長，湖冰面積占比介于95%～100%的時(shí)長則為完全封凍時(shí)長. 統(tǒng)計(jì)表明，近20年青海湖平均開始凍結(jié)和完全凍結(jié)日期分別為12月9日和1月6日，平均開始消融和完全消融日期分別為3月23日和4月10日. 但是，在個(gè)別年份青海湖凍結(jié)和消融日期存在較大差異，如2016-2017年開始凍結(jié)和完全凍結(jié)均最晚，分別為2016年12月25日和2017年1月23日，這應(yīng)與該年11-12月氣溫偏高有關(guān)，根據(jù)距青海湖最近的剛察氣象站觀測(cè)資料，2016年11-12月平均氣溫為-7.77℃，較2000-2019年同期氣溫偏高2.63℃. 開始凍結(jié)和完全凍結(jié)最早日期分別為2012年11月30日和2005年12月25日，開始消融和完全消融最早日期均出現(xiàn)在2016年，分別為3月15日和3月22日，最晚日期分別為2012年3月31日和2011年4月16日. 盡管2013-2017年青海湖開始凍結(jié)時(shí)間較晚(多為12月中旬)、開始消融和完全消融時(shí)間較早(分別出現(xiàn)在3月中旬和3月末)，但整體來看并未發(fā)生顯著變化，4個(gè)時(shí)間節(jié)點(diǎn)的變幅均介于4～6 d.

3.4.2 冰厚變化與湖冰物候特征的關(guān)系 從凍結(jié)時(shí)長、完全凍結(jié)時(shí)長和消融時(shí)長變化(圖9)來看，青海湖凍結(jié)時(shí)長明顯大于消融時(shí)長，其平均時(shí)長分別為29和11 d，完全凍結(jié)時(shí)長平均為75 d，最大和最小時(shí)長分別為94 d(2005-2006年)和50 d(2016-2017年). 除個(gè)別年份外，青海湖凍結(jié)時(shí)長和消融時(shí)長雖有所波動(dòng)，但整體變化不大；完全凍結(jié)時(shí)長在2014-2017年有明顯的減少趨勢(shì)，且2004-2006年也較小，與之相對(duì)應(yīng)的剛察氣象站1-3月平均氣溫也較高(普遍高于-7.80℃)，較2000-2019年1-3月平均氣溫高0.66℃. 與采用相同數(shù)據(jù)和方法而閾值為10%的研究結(jié)果[10]對(duì)比，兩種結(jié)果表現(xiàn)出一定差異性，具體表現(xiàn)為平均開始凍結(jié)日期有所提前(相差7 d)，完全凍結(jié)日期和開始消融日期基本一致(相差1 d)，完全消融日期推后(相差6 d). 平均凍結(jié)時(shí)長較10%閾值的分析結(jié)果長14 d，完全凍結(jié)時(shí)長短2 d，消融時(shí)長僅長1 d.

湖冰物候(封凍期長短)是影響湖泊冰厚變化的重要因素之一. 根據(jù)分析青海湖各年凍結(jié)時(shí)長、完全凍結(jié)時(shí)長、消融時(shí)長、封凍期(湖泊開始凍結(jié)至完全消融[29])和湖冰年平均厚度可知(圖9)，湖冰厚度與4個(gè)時(shí)期相關(guān)系數(shù)R2分別為-0.05、0.13、0.06、0.23(P值均小于0.05)，其中冰厚與完全凍結(jié)時(shí)長及封凍期的相關(guān)性較大，這是因?yàn)楹嬖诘臅r(shí)期越長，意味著湖泊水面被冰層覆蓋的時(shí)間越長，除冰-氣間的熱量交換外，湖冰表面的高反射率使冰體吸收的太陽輻射減少，并進(jìn)一步加劇近地面氣溫降低和冰-氣熱量交換增強(qiáng)，最終導(dǎo)致湖冰凍結(jié)深度更大. 表明冰厚垂直方向的增厚主要集中于完全封凍期，而凍結(jié)期與消融期則表現(xiàn)為湖冰在湖面水平方向的凍結(jié)擴(kuò)張或消融破裂，故而與厚度變化相關(guān)性較低.

圖9 2000-2019年青海湖湖冰物候與平均冰厚的關(guān)系

4 討論

當(dāng)湖泊表面被冰覆蓋時(shí)，湖冰厚度增長主要取決于熱力條件和水力條件，其中下冰面-湖水界面的熱交換主要是湖冰底部與水流之間的湍流熱交換，冰體內(nèi)部通過導(dǎo)熱完成熱量傳遞，冰體表面則與大氣發(fā)生熱交換，并受氣溫、風(fēng)速、短波和長波輻射的吸收反射及云量等條件影響. 研究表明，冰層垂向熱傳導(dǎo)通量是驅(qū)動(dòng)冰層生長的主要能量源，隨著冰厚的增加冰內(nèi)垂向熱傳導(dǎo)通量逐漸下降. 當(dāng)冰厚為12～30 cm時(shí)，湖冰上冰層垂向熱傳導(dǎo)通量隨氣溫的變化存在明顯高頻波動(dòng)特性，而冰厚大于30 cm的中下冰層波動(dòng)特性明顯衰減[30]. 青海湖每年11月至翌年1月氣溫溫度梯度較大，冰層較薄，大氣與湖水之間的熱量傳輸劇烈，冰內(nèi)垂向熱傳導(dǎo)通量波動(dòng)明顯，這導(dǎo)致湖冰的厚度增長速率更快；而2-3月溫度梯度變化較小，當(dāng)冰厚大于35 cm時(shí)，冰內(nèi)垂向熱傳導(dǎo)通量波動(dòng)明顯衰減導(dǎo)致冰層生長速率減緩.

湖泊水域范圍通常具有波動(dòng)性，如在豐水期水位上升、面積增大，而在枯水期則呈相反趨勢(shì). 受降水和冰川融水增加以及蒸發(fā)減少等影響，青藏高原大多數(shù)地區(qū)湖泊面積近期呈增加趨勢(shì)[9, 31]. 已有研究表明，1973-2018年青海湖水位和面積總體呈現(xiàn)先下降(減少)后上升(增加)趨勢(shì)，并大致以2004年為界，2004-2018年湖泊水位上升2.07 m，面積增加19.59 km2，2018年的水域面積基本與20世紀(jì)70年代末期的面積持平[32]. 根據(jù)2019年夏季考察，青海湖面積仍在持續(xù)增加，鳥島附近道路、停車場(chǎng)及周邊草場(chǎng)均被淹沒. 通常，湖泊水深越大，其相應(yīng)的熱容也越大，并直接影響湖冰的生消周期[9]. 通過對(duì)比本研究時(shí)段內(nèi)各年青海湖冰厚與水位變化趨勢(shì)發(fā)現(xiàn)，兩者之間并不存在相關(guān)性，如2000-2004年青海湖水位持續(xù)下降，2004年之后水位快速上升，但其冰厚在對(duì)應(yīng)時(shí)段內(nèi)并沒有發(fā)生類似的變化趨勢(shì)，這可能與青海湖鹽度較小及平均水深較淺有關(guān)，并導(dǎo)致冰下與湖水熱量交換強(qiáng)度微弱.

隨著青海湖面積變化，近年來湖岸線亦發(fā)生動(dòng)態(tài)變化. 已有研究表明，近45年來岸線發(fā)生變化最大的區(qū)域主要是青海湖東岸的沙島、西岸的鳥島和鐵布卡灣及北岸沙柳河入口區(qū)域，變化時(shí)間點(diǎn)與湖泊面積退縮及擴(kuò)張時(shí)間一致[32]. 以2004年為界，前期湖岸線向湖內(nèi)推進(jìn)，尤其導(dǎo)致沙島湖與青海湖主體脫離，2004年之后湖岸線持續(xù)向外后退，其中2017-2018年是相鄰年間岸線變化最大的年份[31]. 青海湖凍結(jié)初期較薄的湖冰會(huì)被吹至湖岸邊，沿著湖岸形成岸冰，而這主要與湖岸線形狀和局地盛行風(fēng)向有關(guān)[28]. 湖岸線較為復(fù)雜的區(qū)域(如湖區(qū)東部海晏灣、尕海、沙島及湖區(qū)西部的鐵卜卡灣)最先出現(xiàn)岸冰，且最早形成連續(xù)的固定冰，上述區(qū)域亦為青海湖空間平均冰厚較厚的區(qū)域. 雖然研究時(shí)段內(nèi)岸線變化與平均冰厚變化并沒有相同的變化趨勢(shì)，但岸線變化最大的年份(2017-2018年)平均冰厚較前一年增加2.46 cm，亦是近20年來冰厚變化較大的年份(除2010-2011年較前一年增加2.93 cm).

湖冰厚度變化不僅影響溫室氣體排放強(qiáng)度[33]，而且通過控制湖面反射率和湖-氣間物質(zhì)與能量交換影響湖區(qū)熱量和能量收支平衡[3,34-35]，并通過改變湖水溫度分層控制湖泊生物地球化學(xué)的水生生物食物鏈進(jìn)程[36]. 同時(shí)對(duì)于交通運(yùn)輸，漁業(yè)冬捕、冬季旅游和地方文化認(rèn)同具有重要意義[37-38]. 近10余年青海湖冬季旅游持續(xù)增熱，如地方部門舉辦冰面徒步穿越青海湖活動(dòng)，個(gè)別游客也會(huì)私自下湖游玩. 盡管青海湖國家級(jí)自然保護(hù)區(qū)管理局已出臺(tái)相關(guān)政策規(guī)范冬季冰上旅游安全，但每年仍有車輛陷入湖中事故發(fā)生，這一方面與青海湖范圍廣闊難以有效監(jiān)管有關(guān)，另一方面也與游客對(duì)青海湖冰厚認(rèn)識(shí)不足有關(guān). 上節(jié)分析表明青海湖通常在每年1月中下旬平均冰厚可達(dá)35 cm左右，其承載能力也會(huì)得到極大增強(qiáng)，因此從1月下旬至3月上旬開展冰面活動(dòng)較為安全. 但需要注意的是，因青海湖冰厚在空間上的異質(zhì)性，在一些湖區(qū)應(yīng)加強(qiáng)游客管理，尤其是在湖區(qū)西南側(cè).

5 結(jié)論

本文基于ERA5 Climate Reanalysis氣溫?cái)?shù)據(jù)集、MODIS MOD09GQ數(shù)據(jù)產(chǎn)品和2019年湖冰鉆孔測(cè)厚數(shù)據(jù)及雷達(dá)測(cè)厚數(shù)據(jù)，應(yīng)用度日法湖冰生長模型計(jì)算青海湖2019年湖冰厚度，并對(duì)2000-2019年青海湖冰厚時(shí)空變化進(jìn)行分析，得到以下結(jié)論：

1)分析2019年1-3月的實(shí)測(cè)數(shù)據(jù)可知，1月30日湖區(qū)平均冰厚為40.2 cm，3月2日平均冰厚測(cè)量值為43.5 cm，3月22日為51.4 cm，3月平均增長速率為0.30 cm/d，較2月份(0.12 cm/d)增長迅速. 實(shí)測(cè)冰厚空間整體呈現(xiàn)北厚南薄，東厚西薄的特點(diǎn). 就冰厚平均增長速率而言，測(cè)區(qū)1、3、7僅有同一時(shí)間段內(nèi)觀測(cè)記錄，有兩段時(shí)間內(nèi)觀測(cè)記錄的9個(gè)測(cè)區(qū)，主要有冰厚快速增長區(qū)(測(cè)區(qū)6、13、9)、冰厚平緩增長區(qū)(測(cè)區(qū)11、10、8)、冰厚較慢增長區(qū)(測(cè)區(qū)12、14). 測(cè)區(qū)2、4、5、15包含3段采樣時(shí)間，其對(duì)應(yīng)平均冰厚增長速率依次減小.

2)2018年11月-2019年3月湖冰厚度平均增長速率為0.34 cm/d，與實(shí)測(cè)數(shù)據(jù)相比模擬冰厚平均誤差為±2 cm，但在湖區(qū)不同區(qū)域模擬精度有所差異，其中誤差較大區(qū)域?yàn)椴脊涌?、哈爾蓋河入湖處，以及湖區(qū)西南部以及湖區(qū)南部. 同時(shí)，冰厚模擬數(shù)值存在3月中旬前高估而之后低估現(xiàn)象.

3)2000-2019年青海湖年平均冰厚為32～37 cm，年際變化呈波動(dòng)趨勢(shì)，其中2000-2007年為冰厚穩(wěn)定期，2008-2016年為冰厚劇烈波動(dòng)期，2017-2018年為冰厚恢復(fù)期. 青海湖凍結(jié)前期冰厚增長迅速，12月與1月湖冰增長速率分別高達(dá)0.45和0.41 cm/d；2月后冰厚增長速率放緩，2月和3月的增長速率分別為0.29和0.14 cm/d，這與冰層垂向熱傳導(dǎo)通量在不同深度冰層內(nèi)隨氣溫變化的波動(dòng)頻率密切相關(guān). 青海湖湖冰平均厚度在空間上存在顯著差異，整體表現(xiàn)為北厚南薄、東厚西薄的空間格局，且湖區(qū)西部冰厚變化幅度小于東部.

4)近20年青海湖平均開始凍結(jié)和完全凍結(jié)日期分別為12月9日和1月6日，平均開始消融和完全消融日期分別為3月23日和4月10日，平均凍結(jié)時(shí)長和消融時(shí)長分別為29和11 d，完全封凍時(shí)長平均為75 d. 湖冰平均厚度變化與完全封凍時(shí)長以及封凍期相關(guān)性較高，而與凍結(jié)時(shí)長及消融時(shí)長相關(guān)性較小，表明冰層增厚主要集中于完全封凍期，凍結(jié)期與消融期則為湖冰在湖面的凍結(jié)擴(kuò)張或消融破裂.

湖冰演化是一個(gè)極為復(fù)雜的過程，盡管本文模擬了2000-2019年青海湖冰厚空間分布特征及變化趨勢(shì)，但因度日法湖冰生長模型僅考慮氣溫對(duì)冰厚的影響，并未考慮冰體內(nèi)部導(dǎo)熱過程及冰下的熱力與水力作用，不能模擬消融期的冰厚變化及未能考慮積雪存在對(duì)冰厚模擬的影響，這使得模擬結(jié)果與真實(shí)湖冰厚度存在一定差異，在未來模擬青海湖湖冰厚度時(shí)應(yīng)加以考慮. 同時(shí)，受試驗(yàn)觀測(cè)區(qū)域限制，對(duì)青海湖中間大片區(qū)域的冰厚模擬結(jié)果仍有待驗(yàn)證.

致謝：在開展湖冰測(cè)厚試驗(yàn)期間，青海湖國家級(jí)自然保護(hù)區(qū)管理局和剛察縣氣象局工作人員、大連中?？萍及l(fā)展有限公司陳潔博士及云南大學(xué)祁苗苗博士給予了大量幫助，在此致以謝忱.