李 諾，韓其飛，馬 英，黃曉東

（1.南京信息工程大學(xué) 地理科學(xué)學(xué)院，江蘇 南京 210041；2.蘭州大學(xué) 草地農(nóng)業(yè)生態(tài)系統(tǒng)國家重點(diǎn)實(shí)驗(yàn)室 草地農(nóng)業(yè)科技學(xué)院，甘肅 蘭州 730000）

0 引言

積雪作為冰凍圈的重要組成成分之一，對氣候變化十分敏感。根據(jù)政府間氣候變化專門委員會（Intergovernmental Panel on Climate Change，IPCC）于2021 年8 月發(fā)布的第六次評估報(bào)告，冰凍圈正呈現(xiàn)加速萎縮狀態(tài)，北半球積雪范圍也在加速減小［1］。積雪因?yàn)槠涮囟ǖ母叻凑章屎偷蛯?dǎo)熱率，對地球能量平衡及水循環(huán)有著深刻的影響［2-3］。青藏高原地處中緯度地區(qū)，作為我國三大穩(wěn)定積雪區(qū)之一，該地區(qū)的積雪影響東亞地區(qū)的大氣環(huán)流及天氣系統(tǒng)，其冰雪融水更是農(nóng)業(yè)灌溉及周邊眾多河流的重要來源［4］。因此準(zhǔn)確獲取青藏高原地區(qū)的積雪信息對研究青藏高原的徑流變化及氣候演變有著重要意義。

遙感技術(shù)相較于傳統(tǒng)的氣象臺站觀測數(shù)據(jù)具有大尺度、獲取信息方便快捷、不受地理?xiàng)l件限制等特征，成為大范圍積雪監(jiān)測的關(guān)鍵手段［5］。搭載在Terra 和Aqua 衛(wèi)星的中分辨率成像光譜儀（Moderate Resolution Imaging Spectroradiometer，MODIS）因其在時間、空間和光譜分辨率的較大優(yōu)勢，在區(qū)域及半球尺度積雪監(jiān)測中應(yīng)用廣泛。MODIS積雪產(chǎn)品第5 版本（V5）提供了二值積雪范圍產(chǎn)品和積雪覆蓋率產(chǎn)品［6］，兩種產(chǎn)品的算法均基于歸一化差值積雪指數(shù)（Normalized Difference Snow Index，NDSI）。二值積雪產(chǎn)品通過設(shè)置NDSI 閾值對每個像元進(jìn)行分類，即一個MODIS 像元如果NDSI 值大于等于0.4 且波段2 反射率大于11%，則該像元被分類為積雪［7］，并引入歸一化植被指數(shù)（Normalized Difference Vegetation Index，NDVI）用于提高M(jìn)ODIS 在林區(qū)的積雪分類精度［8］。評估結(jié)果表明，MODIS 逐日二值積雪產(chǎn)品（V5）在晴空狀態(tài)下的總體分類精度介于85%～99%［9］。積雪覆蓋率（Fractional Snow Cover，F(xiàn)SC）產(chǎn)品通過建立Landsat 獲取的FSC 與NDSI 之間的線性關(guān)系生成［10］。2016年發(fā)布的MODIS 第6版本積雪產(chǎn)品（V6）不再提供積雪二值產(chǎn)品與積雪覆蓋率產(chǎn)品，僅提供網(wǎng)格的NDSI數(shù)據(jù)。V6的產(chǎn)品結(jié)合高程重新定義了地表溫度對高海拔地區(qū)積雪的分類閾值。另外，采用定量圖像恢復(fù)算法重建了Aqua/MODIS 短波紅外波段（波段6），替代了V5 版本使用波段7 計(jì)算NDSI的策略［11］。V6版本較V5版本的積雪產(chǎn)品，質(zhì)量和精度均得到了有效的提升［12］。

云是影響MODIS 逐日積雪遙感產(chǎn)品準(zhǔn)確獲取地表積雪覆蓋范圍的主要因素之一。雖然MODIS逐日積雪產(chǎn)品在晴空狀態(tài)下具有可靠的精度，但由于云的影響，導(dǎo)致該產(chǎn)品在區(qū)域積雪監(jiān)測中受到極大的限制。因此針對MODIS 逐日積雪產(chǎn)品的云下信息恢復(fù)，一直是領(lǐng)域內(nèi)學(xué)者關(guān)注的一個熱點(diǎn)問題。國內(nèi)外許多學(xué)者基于V5 產(chǎn)品針對云下信息恢復(fù)開展了大量的研究。唐志光等［13］利用三次樣條函數(shù)插值法對MODIS 逐日積雪產(chǎn)品MOD10A1 進(jìn)行去云處理，得到青藏高原MODIS 逐日無云的積雪覆蓋率產(chǎn)品，該產(chǎn)品獲取的積雪日數(shù)（Snow-Covered Days，SCD）與地面觀測值得到的SCD 平均一致性為87.0%。侯小剛等［14］通過融合MOD10A1 與MYD10A1 數(shù)據(jù)，發(fā)展了一套適合新疆地區(qū)的逐日無云積雪范圍數(shù)據(jù)集，去云后產(chǎn)品與氣象臺站積雪觀測的一致性為88.1%，總體精度達(dá)到90.6%。黃曉東等［15］利用MODIS 上下午星獲取的逐日積雪產(chǎn)品和被動微波數(shù)據(jù)AMSR-E 雪水當(dāng)量產(chǎn)品，對MODIS 逐日積雪圖像進(jìn)行上下午星合成和鄰近日像元分析，并結(jié)合數(shù)字高程數(shù)據(jù)利用雪線算法進(jìn)行云下信息恢復(fù)，獲取了青藏高原時空連續(xù)的MODIS 積雪逐日產(chǎn)品，總體分類精度達(dá)到90.7%。Gafurov等［16］在阿富汗東北部的Kokcha盆地，綜合了雙星合成、時間濾波、空間濾波等五種方法組成連續(xù)六步去云算法，生成的MODIS 逐日積雪產(chǎn)品，其總體分類精度也達(dá)到90%。Huang 等［17］提出一種基于隱馬爾可夫隨機(jī)場（Hidden Markov Random Field，HMRF）的時空濾波模型，充分利用了立方體內(nèi)的時空信息對MODIS 逐日積雪產(chǎn)品進(jìn)行去云處理，總體分類精度達(dá)88.0%，云量減少到1%。Hall等［18］提出一種基于MOD10C1積雪覆蓋產(chǎn)品的時間濾波去云方法，生成MODIS Cloud-Gap-Filled（CGF）逐日云填補(bǔ)積雪產(chǎn)品。隨著V6 版本積雪產(chǎn)品的發(fā)布，陸續(xù)有學(xué)者基于V6 積雪產(chǎn)品開展了類似研究。邱玉寶等［19］于2021 年11 月在科學(xué)數(shù)據(jù)銀行（Science Data Bank，ScienceDB）更新了基于V6 版本的青藏高原MODIS 逐日無云積雪面積數(shù)據(jù)集，該數(shù)據(jù)集共使用八個連續(xù)步驟進(jìn)行去云。Muhammad 等［20］于2021 年3 月在Data Publisher for Earth & Environmental Science（PANGAEA）發(fā)布了基于MODIS V6產(chǎn)品的亞洲高山區(qū)逐日無云二值積雪數(shù)據(jù)集。郝曉華［21］于2020 年11 月在國家冰川凍土沙漠科學(xué)數(shù)據(jù)中心發(fā)布中國逐日無云500 m 積雪面積產(chǎn)品數(shù)據(jù)集（2000—2020 年），該數(shù)據(jù)集基于MODIS V6 逐日地表反射率產(chǎn)品重建NDSI，通過訓(xùn)練不同土地覆蓋類型下NDSI 識別閾值，然后利用隱馬爾可夫隨機(jī)場時空濾波模型進(jìn)行了去云處理。上述產(chǎn)品自發(fā)布后，在青藏高原還未得到系統(tǒng)驗(yàn)證。

青藏高原作為我國三大積雪分布地區(qū)之一，海拔高且輻射強(qiáng)，導(dǎo)致積雪消融迅速、積雪破碎化嚴(yán)重，使得目前不同積雪產(chǎn)品在該地區(qū)的積雪監(jiān)測精度均不夠理想。因此，從應(yīng)用角度出發(fā)，對MODIS逐日無云積雪產(chǎn)品在青藏高原進(jìn)行系統(tǒng)評估，研究其在積雪監(jiān)測中的精度及影響因素尤為關(guān)鍵。鑒于此，本文利用中高分辨率的Landsat影像數(shù)據(jù)作為參考值，對上述三套MODIS 逐日無云積雪產(chǎn)品進(jìn)行驗(yàn)證，并分析土地覆蓋類型和積雪覆蓋率對MODIS積雪產(chǎn)品精度的影響，為評估青藏高原積雪范圍及其時空動態(tài)變化提供參考依據(jù)。

1 研究區(qū)及數(shù)據(jù)

1.1 研究區(qū)概況

青藏高原位于我國西南部，平均海拔4 000 m左右，被稱為“世界屋脊”“第三極”［22］（圖1）。青藏高原西起帕米爾高原向東至橫斷山脈，自南喜馬拉雅山脈南緣延伸至北昆侖山脈北緣和祁連山脈，總面積約2.57×106km2［23］。青藏高原以冰川、積雪、凍土等形態(tài)儲存了巨大的水資源，我國的黃河、長江、怒江、瀾滄江、雅魯藏布江等均發(fā)源于青藏高原，有“中華水塔”和“亞洲水塔”之稱［24］。青藏高原積雪水儲量關(guān)系著所在區(qū)域及周邊區(qū)域的生活和灌溉用水，影響高原及周邊區(qū)域的植被生長，該地區(qū)的積雪具有重要的水文、生態(tài)及氣候意義［25］。

圖1 青藏高原土地覆蓋類型及Landsat-8影像分布Fig.1 Land cover types in Tibetan Plateau and Landsat-8 images spatial distribution

1.2 MODIS逐日無云積雪產(chǎn)品

（1）中國逐日無云500 m 積雪面積產(chǎn)品數(shù)據(jù)集（MODIS CGF SCE）：來源于國家冰川凍土沙漠科學(xué)數(shù)據(jù)中心（http：//www.ncdc.ac.cn/）。該產(chǎn)品基于MODIS反射率產(chǎn)品MOD09GA/MYD09GA，利用高分辨率Landsat TM 數(shù)據(jù)作為真值，結(jié)合MODIS土地覆蓋分類產(chǎn)品，確定林區(qū)和非林區(qū)下積雪判別的指標(biāo)閾值，利用MODIS 積雪反演算法獲取初級產(chǎn)品，通過隱馬爾科夫隨機(jī)場時空濾波算法進(jìn)行云下信息恢復(fù)。該產(chǎn)品以HDF5 文件格式存儲，每個HDF5文件包含18個數(shù)據(jù)要素，其中包括數(shù)據(jù)值、數(shù)據(jù)起始日期、經(jīng)緯度等。本文采用無云積雪數(shù)據(jù)，為二值影像。其編碼有0、1、2、3、4和255，對應(yīng)的編碼含義為陸地、影像識別積雪、去云插補(bǔ)積雪、雪深插補(bǔ)積雪、水體和填充值［21］。

（2）青藏高原MODIS 逐日無云積雪面積數(shù)據(jù)集（MODIS_Dysno_Cloudfree）：來源于ScienceDB（https：//www.scidb.cn/）。該產(chǎn)品以MODIS V6 積雪產(chǎn)品MOD10A1 和MYD10A1 為基礎(chǔ)，通過上下午星合成、三天合成、“長時間”積雪和陸地法、鄰近四像元法、高程濾波法、修改陰影區(qū)錯誤分類和最大積雪陸地范圍掩膜這7 個連續(xù)步驟，獲得MODIS少云積雪產(chǎn)品，而后將研究區(qū)根據(jù)坡度劃分為“印度平原”“高原北部”“高原腹地”“藏東南山區(qū)”“帕米爾高原”“天山山脈”以及“喜馬拉雅山脈”7 個區(qū)域，針對7 個區(qū)域采用擬合預(yù)期雪線方法去除全部云污染，獲得MODIS 逐日無云積雪產(chǎn)品。本次所用產(chǎn)品為該數(shù)據(jù)集中青藏高原MODIS 逐日無云積雪產(chǎn)品面積數(shù)據(jù)集第二次更新版（2002—2021 年MODIS積雪產(chǎn)品C6.1版，http：//www.csdata.org/p/15/）。該產(chǎn)品為二值產(chǎn)品，影像分類代碼與MOD10A1 一致，其編碼有25、37、40、100、150 和200，對應(yīng)含義為無積雪覆蓋的陸地、湖泊、湖泊不確定、湖冰、湖冰不確定和積雪。該產(chǎn)品利用地面臺站積雪深度觀測數(shù)據(jù)進(jìn)行驗(yàn)證，結(jié)果表明其積雪分類精度為78.4%，去掉雪深≤3 cm 的樣本后，精度達(dá)到89.0%［19］。

（3）亞洲高山區(qū)MODIS 逐日無云積雪數(shù)據(jù)集（M*D10A1GL06）：數(shù)據(jù)來源于PANGAEA（https：//www.pangaea.de/）。該產(chǎn)品以MODIS 積雪產(chǎn)品MOD10A1 和MYD10A1 為基礎(chǔ)，使用8d MOYDGL06*產(chǎn)品（雙星8 天融合產(chǎn)品）分別代替Terra 和Aqua 雪產(chǎn)品中的云像素。然后與全球陸地冰空間測量計(jì)劃（Global Land Ice Measurements from Space，GLIMS）發(fā)布的RGI（Randolph Glacier Inventory）6.0 冰川編目數(shù)據(jù)產(chǎn)品進(jìn)行融合，最終達(dá)到99.99%的去云效果［20］。本產(chǎn)品為二值產(chǎn)品，其編碼含義如表1所示。

表1 M*D10A1GL06積雪產(chǎn)品屬性定義Table 1 Attribute definition of M*D10A1GL06 snow product

1.3 Landsat-8

前期研究大多基于站點(diǎn)雪深觀測數(shù)據(jù)對MODIS積雪產(chǎn)品進(jìn)行精度驗(yàn)證，但由于二者空間尺度的不匹配，導(dǎo)致基于點(diǎn)-面驗(yàn)證的結(jié)果存在較大的不確定性。因此，本研究采用空間分辨率為30 m的Landsat-8 陸地成像儀（Operational Land Imager，OLI）獲取的影像作為驗(yàn)證數(shù)據(jù)。OLI 影像首先要選擇晴空狀態(tài)的數(shù)據(jù)，以避免云雪混淆造成驗(yàn)證數(shù)據(jù)集的不確定性；影像應(yīng)分布在青藏高原積雪分布的主要區(qū)域，且能包含高原主要的土地覆蓋類型；另外，為了避免破碎化積雪對OLI積雪制圖的影響，驗(yàn)證時段選擇積雪分布較為穩(wěn)定的時期。因此，選取時間范圍為2015 年12 月2 日至2016 年1 月23 日共34 景影像作為驗(yàn)證MODIS 逐日積雪產(chǎn)品的參考數(shù)據(jù)。

采用混合像元分析法首先獲取OLI像元尺度的積雪面積比例（Fractional Snow Cover，F(xiàn)SC）。其中端元提取的規(guī)則如表2所示。將30 m的OLI積雪面積比例數(shù)據(jù)聚合成500 m，定義升尺度后的亞像元積雪面積比例閾值0.5，生成二值積雪數(shù)據(jù)。其中喜馬拉雅山中東段的18 景Landsat 積雪數(shù)據(jù)來源于ScienceDB：“2013—2020 年喜馬拉雅山中段和東段Landsat 8 積雪覆蓋范圍數(shù)據(jù)”。該數(shù)據(jù)集采用支持向量機(jī)（Support Vector Machine，SVM）分類方法，選取不同地形、陰影等條件的積雪特征訓(xùn)練樣本進(jìn)行積雪分類，結(jié)合地表水體等輔助數(shù)據(jù)及空間鄰域分析進(jìn)行分類后處理。通過對比Sentinel-2 高分辨率積雪分類數(shù)據(jù)，在900 m×900 m的網(wǎng)格內(nèi)，其積雪覆蓋率相關(guān)系數(shù)在0.95 以上，均方根誤差約0.1%。該數(shù)據(jù)集為二值影像，其編碼有1、2、100、200 和300，對應(yīng)的編碼含義分別為積雪、陸地、河流、冰湖和湖泊［27］。

表2 Landsat-8 OLI端元提取規(guī)則［26］Table 2 Endmember extraction rules for Landsat-8 OLI［26］

1.4 土地覆蓋數(shù)據(jù)

MODIS 土地覆蓋產(chǎn)品（MCD12Q1）用來比較不同土地覆蓋類型下的MODIS 逐日無云積雪產(chǎn)品的精度。采用國際地圈生物圈計(jì)劃（International Geosphere-Biosphere Program，IGBP）分類方案確定的17個土地覆蓋類別，其中包括11個自然植被類別，3個開發(fā)和鑲嵌土地類別，以及3 個非植被土地覆蓋類別［28］。為了避免類型過多造成評價結(jié)果的復(fù)雜性，本文根據(jù)將IGBP 分類方案重新分類為七大類：森林、灌叢、草原、耕地、城市建設(shè)用地，冰川和裸地（圖1）［29］。

2 精度驗(yàn)證方法

使用高分辨率遙感產(chǎn)品驗(yàn)證低分辨率產(chǎn)品的精度是遙感產(chǎn)品驗(yàn)證的常用手段［30-31］。本研究利用高分辨率的Landsat積雪圖對MODIS兩個版本的數(shù)據(jù)進(jìn)行對比驗(yàn)證，采用精度評估指標(biāo)包括總體分類精度（overall accuracy，OA）、低估誤差（underestimated error，UE）、高估誤差（overestimated error，OE），積雪分類精度P和Kappa系數(shù)。其中Kappa一致性檢驗(yàn)是評價兩幅圖像之間一致性的方法［32］。當(dāng)Kappa 系數(shù)為0～0.20 時，表明兩幅圖像的一致性極低；當(dāng)Kappa 系數(shù)為0.20～0.40 時，表明兩幅圖像一致性一般；當(dāng)Kappa 系數(shù)為0.40～0.60 時，表明兩幅圖像一致性中等；當(dāng)Kappa 系數(shù)為0.60～0.80 時，表明兩幅圖像一致性極好；當(dāng)Kappa 系數(shù)為0.80～1.00時，表明兩幅圖像高度一致［33］。

精度評價指標(biāo)計(jì)算公式如下所示：

式中：a為Landsat 和MODIS 均記錄為積雪的樣本數(shù)；b為Landsat記錄有雪而MODIS被分類為非雪的樣本數(shù)，即漏判數(shù)；c為Landsat 和MODIS 均為非雪的樣本數(shù)；d為Landsat記錄為非雪而MODIS被分類為積雪的樣本數(shù)，即為誤判數(shù)；Pe表示偶然一致性，即在偶然機(jī)會下預(yù)測結(jié)果與真實(shí)結(jié)果一致的概率。

3 結(jié)果與討論

3.1 精度評價

驗(yàn)證結(jié)果表明，MODIS CGF SCE、M*D10A1 GL06 和MODIS_Dysno_Cloudfree 產(chǎn)品的總體精度（OA）分別為88.2%、88.2%和89.6%，表明三種積雪產(chǎn)品有著相對較高的精度。其中MODIS CGF SCE低估誤差（UE）最小，僅為9.4%，M*D10A1GL06和MODIS_Dysno_Cloudfree 積雪產(chǎn)品UE 分別為14.7%和18.6%。三套產(chǎn)品的高估誤差（OE）分別為12.2%、8.9%和11.2%。由于MODIS CGF SCE 產(chǎn)品積雪分類精度（P）最高，達(dá)到90.6%，其他兩種產(chǎn)品分別僅為85.3%和81.4%。Kappa 系數(shù)統(tǒng)計(jì)結(jié)果表明三種積雪數(shù)據(jù)與Landsat積雪數(shù)據(jù)的一致性均表現(xiàn)為極好，分別為0.642、0.629 和0.652，MODIS CGF SCE 和MODIS_Dysno_Cloudfree 與Landsat 積雪數(shù)據(jù)一致性最為接近，M*D10A1GL06產(chǎn)品略差。

M*D10A1GL06 和MODIS_Dysno_Cloudfree 產(chǎn)品均是基于MODIS V6 版本生產(chǎn)的去云積雪產(chǎn)品，兩種產(chǎn)品的用戶精度理論上是一致的，造成積雪分類誤差不一致的原因主要是云重分類誤差導(dǎo)致的，去云后M*D10A1GL06 的積雪分類精度優(yōu)于MODIS_Dysno_Cloudfree，說明后者采用的云下信息恢復(fù)算法存在較大的誤分類現(xiàn)象，與站點(diǎn)驗(yàn)證結(jié)論也有一定的差距［20］。MODIS CGF SCE 產(chǎn)品是基于大氣校正的地表反射率數(shù)據(jù)生成的，并且分土地覆蓋類型設(shè)置最優(yōu)NDSI 閾值，在積雪分類時利用雪深等輔助數(shù)據(jù)進(jìn)行插值，因此積雪分類精度最高；而MODIS V6產(chǎn)品是經(jīng)過輻射校正的反射率數(shù)據(jù)計(jì)算NDSI。驗(yàn)證結(jié)果表明，基于大氣校正的地表反射率數(shù)據(jù)在積雪信息提取方面更具優(yōu)勢，積雪分類精度優(yōu)于另外兩種去云產(chǎn)品。

圖2 MODIS逐日無云積雪產(chǎn)品驗(yàn)證結(jié)果Fig.2 Validation results for three kinds of MODIS daily cloud-free snow cover products

3.2 土地覆蓋對MODIS積雪產(chǎn)品精度的影響

表3中統(tǒng)計(jì)了MODIS逐日無云積雪產(chǎn)品在不同土地覆蓋類型下的積雪分類精度。其中草原、裸地和森林的驗(yàn)證樣本（N）分別占比67%（N=2 391 609），20%（N=729 769）和10%（N=370 443），其余驗(yàn)證樣本分布在其他土地覆蓋類型，因樣本占比較小，未統(tǒng)計(jì)精度評價結(jié)果。結(jié)果顯示，三種去云積雪產(chǎn)品的總體分類精度在不同土地覆蓋類型下均大于80%，其中MODIS CGF SCE 在林區(qū)的總體分類精度最優(yōu)，達(dá)到89.6%，MODIS_Dysno_Cloudfree在草原和裸地的總體分類精度最優(yōu)，分別達(dá)到91.6%和84.0%。三種產(chǎn)品在裸地的總體分類精度均相對較差，介于81%～84%之間。MODIS CGF SCE 產(chǎn)品在草原的漏判和誤判誤差相當(dāng)（～10%），在裸地存在嚴(yán)重的高估誤差（20.9%），但在林區(qū)漏判誤差較大（21.5%）。M*D10A1GL06 產(chǎn)品在裸地高估誤差大，在林區(qū)漏判誤差大，分別達(dá)到20.0%和27.2%。MODIS_Dysno_Cloudfree 產(chǎn)品僅在林區(qū)存在較為嚴(yán)重的漏判誤差（21.3%）。MODIS CGF SCE 產(chǎn)品的積雪分類精度在草原和裸地類型下均優(yōu)于另外兩種產(chǎn)品，林區(qū)的積雪識別精度與MODIS_Dys-no_Cloudfree 產(chǎn)品相當(dāng)，分別為78.2%和78.7%，M*D10A1GL06 產(chǎn)品在林區(qū)的積雪識別精度最低，為72.8%。與OLI 積雪數(shù)據(jù)相比，三種產(chǎn)品在草原和裸地類型下，一致性均較高，林區(qū)的一致性較差，但MODIS CGF SCE 產(chǎn)品在林區(qū)與OLI 一致性最高，為0.249。雖然該產(chǎn)品針對不同土地覆蓋類型分別設(shè)置NDSI 閾值，并結(jié)合了林區(qū)積雪指數(shù)（Normalized Difference Forest Snow Index，NDFSI）用于提高林區(qū)積雪分類精度，但驗(yàn)證結(jié)果顯示，相比基于MODIS 標(biāo)準(zhǔn)逐日積雪產(chǎn)品生成的去云產(chǎn)品，對林區(qū)積雪識別精度的改進(jìn)并不明顯，依然存在較高的漏判誤差。

表3 土地覆蓋類型對MODIS積雪產(chǎn)品精度的影響Table 3 Assessment of MODIS cloud-free snow products under different land cover types

3.3 FSC對MODIS積雪產(chǎn)品精度的影響

為了評估積雪覆蓋率（FSC）對MODIS 逐日無云產(chǎn)品積雪分類精度的影響，本研究基于OLI 獲取的FSC，對像元尺度的FSC 進(jìn)行分級，驗(yàn)證不同F(xiàn)SC級別MODIS 逐日無云積雪產(chǎn)品的積雪分類精度。結(jié)果如圖3 所示，三種積雪分類精度隨著積雪覆蓋率的增大在逐漸升高。其中MODIS CGF SCE 逐日無云積雪產(chǎn)品積雪分類精度在NDSI 大于10（原始值/100 代表實(shí)際值）均大于其他兩種產(chǎn)品，Dysno_Cloudfree 產(chǎn)品積雪分類精度相對較差。進(jìn)一步說明基于大氣校正的地表反射率數(shù)據(jù)計(jì)算NDSI，進(jìn)而對積雪信息進(jìn)行提取，相比MODIS 標(biāo)準(zhǔn)逐日積雪產(chǎn)品具有一定的優(yōu)勢。

圖3 不同積雪覆蓋率下的MODIS積雪分類精度Fig.3 Accuracy of snow classification in graded fractional snow cover

4 結(jié)論

近年來，MODIS積雪產(chǎn)品在水文、水資源管理、氣象和氣候變化等眾多領(lǐng)域得到了廣泛的應(yīng)用。在青藏高原地區(qū)，受限于該地區(qū)較強(qiáng)的太陽輻射，積雪消融迅速且破碎化嚴(yán)重，使得MODIS 積雪產(chǎn)品在該地區(qū)的準(zhǔn)確性明顯低于高緯度穩(wěn)定積雪區(qū)。本文選取了青藏高原基于MODIS 數(shù)據(jù)發(fā)展的三套逐日無云積雪產(chǎn)品，其中MODIS CGF SCE 是基于MODIS 地表反射率數(shù)據(jù)生產(chǎn)的，其他兩套產(chǎn)品（M*D10A1GL06 和MODIS_Dysno_Cloudfree）是基于美國雪冰數(shù)據(jù)中心（National Snow and Ice Data Center，NSIDC）發(fā)布的MODIS 逐日積雪產(chǎn)品進(jìn)行去云處理后獲取的。本文利用分辨率較高的Landsat-8 OLI 積雪影像數(shù)據(jù)，對比驗(yàn)證了三套逐日無云積雪產(chǎn)品的精度，并分析了土地覆蓋類型及積雪覆蓋率對產(chǎn)品精度的影響，得出以下結(jié)論：

（1）三套無云積雪產(chǎn)品的總體精度均優(yōu)于80%，MODIS_Dysno_Cloudfree 積雪產(chǎn)品總體分類精度最高，達(dá)到89.6%，Kappa 系數(shù)為0.652。MODIS CGF SCE 低估誤差最小，為9.4%，積雪分類精度也是最高的，達(dá)到90.6%。M*D10A1GL06 高估誤差最低（8.9%）。

（2）土地覆蓋類型對MODIS 積雪產(chǎn)品積雪分類的精度影響較大，尤其在林區(qū)。三種積雪產(chǎn)品在林區(qū)的積雪分類精度均較差，存在較高的漏判誤差。雖然MODIS CGF SCE 產(chǎn)品針對林區(qū)進(jìn)行了NDSI閾值優(yōu)化，但效果甚微。

（3）隨著積雪覆蓋率的增加，三套無云積雪產(chǎn)品的積雪分類精度也隨之增加，其中MODIS CGF SCE 產(chǎn)品在不同NDSI 等級（NDSI＞10）均優(yōu)于其他兩套產(chǎn)品。

綜上，基于MODIS 地表反射率產(chǎn)品生產(chǎn)的逐日無云產(chǎn)品（MODIS CGF SCE），針對不同土地覆蓋類型進(jìn)行NDSI 閾值優(yōu)化，在積雪識別精度方面提升較為明顯，但是林區(qū)積雪識別誤差較大的問題依舊沒有得到有效解決。當(dāng)積雪覆蓋率較低時，積雪識別精度相比未經(jīng)過大氣校正的MODIS 標(biāo)準(zhǔn)積雪產(chǎn)品，精度相對較好，但提升不顯著。因此，在積雪破碎化嚴(yán)重的青藏高原地區(qū)，除了考慮NDSI 閾值的優(yōu)化，應(yīng)進(jìn)一步結(jié)合該地區(qū)的積雪特性，包括環(huán)境特征，進(jìn)一步發(fā)展適合青藏高原的積雪制圖算法。