-
遥感在监测积雪方面优势明显, 目前已有一些利用遥感数据获取积雪面积的方法。高空间分辨率的遥感数据, 如TM(Thematic Mapper)、SPOT(Systeme Probatoire d’Observation de la Terre)影像等, 虽然可以得到空间上高精度的积雪面积数据, 但覆盖范围小, 时间分辨率低, 且工作量巨大。中分辨率卫星遥感数据的时间分辨率较高, 但空间分辨率低, 仅仅依靠二值填图难以得到空间分辨率高的积雪面积数据, 因此只有通过中分辨率遥感数据的亚像元雪填图才能得到时间和空间分辨率都高的积雪面积数据。
Metsamaki等基于AVHRR(second series of advanced very high resolution radiometer)数据使用积雪完全覆盖地表反射率和无雪地表反射率间的线性插值, 得到像元内积雪覆盖的百分比值[1, 2]。该算法并非基于可靠物理机理, 所得结果精度不高, 最低为64%。Painter等将多端元光谱混合分析技术引入亚像元雪填图, 利用二流辐射传输模型模拟得到多雪端元, 实测得到多非雪端元。研究表明, 引入包含不同粒径的多种雪端元光谱能够提高混合像元分解雪填图算法的精度[3-5]。Salomonson和Appel基于MODIS(moderate-resolution imaging spectroradiometer)数据利用像元内的积雪面积百分比与归一化雪指数(normalized difference snow index, NDSI)之间的线性关系, 获得亚像元积雪面积数据, 算法简单易于实现[6, 7]。唐志光和王建等将Salomonson的算法引入青藏高原地区的亚像元积雪填图, 获得该区的亚像元积雪面积数据, 但精度不高, 说明Salomonson算法对青藏高原这样的山区来说效果不太好[8]。朱骥等以青藏高原为研究区提出了MODIS数据的多端元自动选取方法, 基于线性混合分解模型实现了MODIS数据的大面积亚像元自动雪填图[9], 该算法简易可靠, 是目前青藏高原地区比较成熟的中分辨率亚像元雪填图算法。亚像元雪填图方面的代表性工作见文献[10-14]。
1998年升空的NOAA-15(即NOAA-K)卫星上的AVHRR/3传感器和1999年发射的EOS/ Terra卫星上的MODIS传感器都不可能提供1998年以前的中分辨率遥感数据。1981年发射的NOAA-7卫星首次使用AVHRR/2传感器, 可以将中分辨率遥感数据延伸到30 a以前, 并提供高时间分辨率数据。由于AVHRR/2传感器缺少AVHRR/3和MODIS中的1.6 μm波段, 目前较可靠的中分辨率亚像元雪填图算法都不能用于实现AVHRR/2数据亚像元雪填图。因此, 需要开展AVHRR/2数据亚像元雪填图研究, 开发相应的算法。
青藏高原地区积雪分布广泛, 对区域水循环和大气变化的影响深刻。但该地区长时段的高精度积雪面积数据目前仍难以获取。地面台站可以得到该地区精确的积雪数据, 但只能观测到局部的一点, 无法开展大面积的监测。而MODIS遥感数据亚像元雪填图无法获得1999年以前的亚像元积雪面积数据。因此为了满足该地区水文、气象以及其他部门的需要, 有必要开展AVHRR/2数据亚像元积雪填图的研究, 以期得到相应亚像元积雪面积的历史数据, 并延长高精度积雪面积数据的时间跨度, 为水文、气象模型提供高精度的时空积雪面积参数。
-
根据AVHRR/2数据, 开展青藏高原亚像元积雪填图研究, 建立AVHRR/2数据亚像元雪填图算法, 关键步骤是确定混合像元分解模型和端元选取方法。
-
本文研究区为青藏高原地区(图 1), 位于25 °N~40 °N, 74°E~104 °E, 是世界平均海拔最高的高原, 平均海拔4 000~5 000 m, 总面积近300万km2, 在中国境内面积约250万km2, 青藏高原地貌景观奇异, 气候最显著特征是高寒, 冰川和常年积雪分布广泛。
-
本研究使用的AVHRR/2数据来自美国地质调查局(United States Geological Survey, USGS)的数据网站, 经过了辐射定标、大气校正和几何配准等预处理, 且是根据10 d内NDVI的最大值得到的AVHRR/2合成数据[15]。AVHRR/2数据总共有5个波段, 第1、2波段为可见光与近红外波段, 第3波段为中红外波段, 第4、5波段为热红外波段, 其空间分辨率约为1 100 m。
本文使用的Landsat/TM数据来自马里兰大学网站, 空间分辨率约为30 m, 在可见光与近红外波段共有6个通道:1、2、3通道分别为可见光的蓝、绿、红通道, 4、5、7通道均为近红外通道(6通道属热红外波段, 本文没有采用)。由于TM数据空间分辨率远小于AVHRR数据, 因此将TM数据的积雪填图数据作为验证AVHRR/2数据雪填图结果的地面参考数据。
-
在现有的混合像元分解模型中, 线性混合模型目前应用最广泛, 可靠性也相对较高[16]。实际应用中同物异谱所带来的光谱多样性和不确定性是混合像元分解算法误差的一个最主要来源。在传统的光谱混合分解中, 一种地物只选取一个端元。不管该端元是否能代表像元中的该地物, 都使用同样的端元, 从而产生分类误差。多端元光谱混合分析考虑同类地物光谱的多样性, 不同像元在进行混合分解计算时选用的同类地物端元可能不同, 可以提高填图算法的精度。
基于此, 本文将多端元线性混合模型作为AVHRR/2数据青藏高原雪填图的算法基础。该模型随机取得各种雪端元和各种非雪端元之间的组合参与计算, 通过最小二乘残差来决定最佳的端元组合, 得到混合像元分解结果。该模型如下:
$${R_i} = f{S_i} + \left( {1 - f} \right){M_i}$$ (1) 式中, f为所解混像元中积雪面积的百分比;Ri为所解混像元第i波段的反射率;Si为第i波段的雪端元反射率;Mi为第i波段的非雪端元反射率。
-
端元选取是混合像元分解过程的关键步骤。为了提高亚像元雪填图的精度和准确性, 必须选取多种类型的雪端元和非雪端元, 建立多种条件下雪和非雪端元的波谱库, 即雪和非雪端元库。
本文主要利用NDVI和AVHRR/2数据某些波段的反射率建立基于知识的回归决策树分析算法, 获取多种雪和非雪端元(图 2)。假定混合像元是由某种特定类型的雪端元和某种非雪端元组成。
1) 雪端元提取规则
根据图 3, 可以发现积雪的第1波段和第2波段反射率非常高, 且第1波段的反射率一般高于第2波段[17]。归一化植被指数是表征地表植被状态的一个有效度量, 其计算公式为:
$$V = ({b_2} - {b_1})/({b_2} + {b_1})$$ (2) 式中, V是NDVI, 即归一化植被指数;b1和b2分别为AVHRR/2数据的第一波段与第二波段的反射率(以下公式中符号的物理意义与此相同)。从式(2) 可以看出, 积雪的NDVI值应当比较小。初步假定积雪端元提取的规则为:
$$\left\{ \begin{array}{l} V \le {m_1}\\ {b_1} > {m_2} \end{array} \right.$$ (3) 式中, m1和m2分别是NDVI和b1的阈值, 都为某一常数。只要确定式(2) 中各阈值的大小, 即可得到提取雪端元的规则。
研究发现, 像元内部积雪的覆盖度越高, NDVI值越小, 可见光波段的反射率则越高[9]。另外, MODIS数据雪端元的NDVI≤-0.04, 第二波段反射率大于0.75[9]。虽然AVHRR/2数据与MODIS数据在光谱上有一定差别, 但这两种传感器的差异不会改变积雪与其他地物的基本光谱特征。
本文以1992年、1993年和1995年中抽出的部分数据作为样本数据, 采取人工方法得到积雪、植被和裸地等不同地物平均光谱曲线, 如图 4所示, 图 4中右边纵坐标只标注NDVI。根据图 4, 可以发现纯雪像元的NDVI值应该略小于0, 而第1通道在80%以上。实验验证后, 确定提取雪端元的规则为:
图 4 AVHRR/2图像中青藏高原地区不同地物的光谱曲线
Figure 4. Spectral Curves of Various Ground Objects in AVHRR/2 Images of Qinghai-Tibet Plateau
$$\left\{ \begin{array}{l} V \le 0.0\\ {b_1} > 0.83 \end{array} \right.$$ (4) 2) 非雪端元提取规则
青藏高原地区的非雪端元只提取植被与裸地两种主要类型, 虽然提取时按植被端元和裸地端元分别提取, 但保存时都放在同一非雪端元库中, 且在库中不再区分植被端元和裸地端元, 视为相同的非雪端元。此外, 非雪端元的提取与雪端元相似, 先采取人工方法得到积雪、植被和裸地等不同地物平均光谱曲线, 如图 4所示(图 4中右边纵坐标只标注NDVI), 然后确定非雪端元提取规则。
(1) 植被端元的提取规则
根据植被的典型光谱特征(图 3), 在AVHRR/2数据中植被的第1波段反射率远小于第2波段, 因此植被的NDVI值相对积雪与裸地而言更高。但在实验中发现单纯使用NDVI提取植被端元光谱容易产生误差, 需要加入其他通道数据作为辅助, 以提高植被端元提取的准确性。此处选择AVHRR/2数据第1通道和第2通道数据参与提取植被端元。最后植被端元的提取模型为:
$$\left\{ \begin{array}{l} V > {m_3}\\ {b_1} < {m_4}\\ {b_2} > {m_5} \end{array} \right.$$ (5) 式中, m3、m4和m5分别是NDVI、b1和b2的阈值常数。从图 4可以看出, 植被的NDVI值一般都在0.4左右, 而第1通道与第2通道的像元值都比较小, 但第2通道相对大些。结合实验结果, 可以确定, 如果
$$\left\{ \begin{array}{l} V > 0.3\\ {b_1} < 0.10\\ {b_2} > 0.11 \end{array} \right.$$ (6) 则可以把植被端元提取出来。
(2) 裸地端元的提取规则
裸地的光谱特征相对来说不太显著, 端元提取的难度也更高。不过, 裸地光谱一般是第1波段反射率与第2波段相差不大, 且第2波段稍稍更高, 其NDVI值较小。采取与植被端元提取相似的措施, 可以得到裸地端元的提取模型和阈值。最后, 裸地的端元提取规则为:
$$\left\{ \begin{array}{l} 0.0 < V < 0.16\\ {b_1} < 0.20\\ {b_2} > 0.30 \end{array} \right.$$ (7) -
1) AVHRR/2原始数据预处理。本文实验中使用USGS处理过的AVHRR/2数据。
2) 水体掩膜。对AVHRR/2图像中的水体在亚像元雪填图前进行掩膜, 使其不参与雪填图计算。
3) 确定混合像元分解模型, 并选取雪和非雪端元。本文使用的线性混合模型参阅文献[16]。
4) 进行AVHRR/2数据亚像元雪填图计算, 并获得青藏高原地区亚像元级积雪面积百分比数据。
5) 验证青藏高原地区AVHRR/2数据亚像元雪填图结果。根据同期的TM数据对该雪填图结果进行验证。
-
根据同期的TM数据对AVHRR/2数据亚像元雪填图结果进行验证。首先, 对TM数据进行几何校正, 再将图像的灰度值转变为反射率值。接着, 根据TM的反射率提取积雪信息, 确定积雪分布面积。本文采取阈值方法将TM图像分为雪和非雪两大类。然后, 将TM的积雪二值图像重采样成与AVHRR/2数据相同空间分辨率的积雪面积百分比模拟图像。最后, 对比TM积雪面积模拟图像与AVHRR/2亚像元积雪面积图像, 以便验证本研究结果的精度与准确性。
选取两景TM数据(TM1、TM2) 作为验证的参考数据。两景数据中心像元的地理位置分别是(33.170 02°N, 94.429 94°E)和(31.740 78°N, 94.040 83°E)。
本研究结果的验证分4级网格尺度进行(表 1), 即1 100×1 100 m、2 200×2 200 m、4 400×4 400 m和8 800×8 800 m。随着网格尺度的变大, 可以消除配准过程带来的误差影响。从表 1可以看出, 当网格为1 100×1 100 m时, 相关系数高于0.83, 均方根误差接近0.1;2 200×2 200 m时, 相关系数高于0.85, 均方根误差低于0.1;4 400×4 400 m时, 相关系数高于0.92, 均方根误差低于0.08;8 800×8 800 m时, 相关系数高于0.95, 均方根误差接近0.05。由此可见, 随着网格尺度增大, 配准误差得到逐步消除。即使在1 100×1 100 m时, 均方根误差也接近0.1。这说明AVHRR/2数据亚像元积雪填图算法是合理可行的。
表 1 本算法的验证结果
Table 1. Validati Results of our Algorithm
TM数据 分辨率/m 相关系数(R) R2 均方根误差(RMSE) TM1 1 100 0.837 2 0.701 0 0.106 1 2 200 0.852 9 0.727 4 0.085 6 4 400 0.923 2 0.852 4 0.062 7 8 800 0.959 6 0.920 8 0.051 2 TM2 1 100 0.846 1 0.715 9 0.117 0 2 200 0.877 9 0.770 8 0.098 5 4 400 0.928 6 0.862 3 0.076 3 8 800 0.962 2 0.925 8 0.056 6 由于10 d合成的AVHRR/2数据与某时刻的TM数据之间存在较大的时间差异, 以及AVHRR/2数据所存在的质量问题, 造成验证结果存在一定的偏差。但总体来说, 本文提出的AVHRR/2数据亚像元雪填图算法是可靠可行的。
-
本文以青藏高原地区为研究区, 针对AVHRR/2的10 d合成数据引入多端元线性光谱混合分析模型, 建立一种亚像元雪填图算法。根据该算法可以得到AVHRR/2数据的亚像元级积雪面积分布百分比数据。
为了得到混合像元分解所需的各种地物端元, 使用了一种直接从AVHRR/2图像中自动选取端元的方法。在选取端元时, 主要采用NDVI、第1波段、第2波段等作为选取端元的指标。雪端元与非雪端元分别选取, 并分别保存于不同的端元数据库中。
基于TM数据, 对该算法的雪填图结果进行了验证。这两种数据的雪填图结果对比所得的均方根误差接近0.1, 对于崎岖地形区来说, 这是比较高的精度。由于所使用的数据为10 d合成的AVHRR/2数据, 所以不能获取每日的AVHRR/2数据亚像元级积雪面积数据。今后, 需加强AVHRR/2数据的前期处理工作, 获得可使用的AVHRR/2逐日数据, 对AVHRR/2逐日数据开展亚像元雪填图, 生产出该区AVHRR/2数据亚像元级积雪面积逐日数据。
Subpixel Snow Mapping Using AVHRR/2 10-Day Compositing Data of Qinghai-Tibet Plateau
-
摘要: 为了满足水文和气象模型对长时段积雪面积数据的需求,基于第二代甚高分辨率辐射计(second series of advanced very high resolution radiometer,AVHRR/2)的10 d合成数据提出了一种青藏高原地区AVHRR/2数据亚像元雪填图算法,将中分辨率遥感数据亚像元级积雪面积数据集延伸至30 a时间跨度。本文算法以多端元线性光谱混合分析模型为基础,以归一化植被指数、第一波段、第二波段等作为选取端元的指标,直接从AVHRR/2图像中自动选取所需雪端元与非雪端元。基于TM数据对该算法的AVHRR/2数据亚像元雪填图结果进行验证,其均方根误差接近0.1,在青藏高原山区具有较高的精度。Abstract: Because there are not AVHRR/3 data and MODIS data before 1998, the time length of their data are not long enough to meet the needs of hydrology and meteorology. However, AVHRR/2 data can extend to 30 years ago. For AVHRR/2 data don't have the data of 1.6 μm waveband, which AVHRR/3 data and MODIS data have, the existing subpixel snow mapping algorithm of moderate resolution data cannot be used to implement the subpixel snow mapping of AVHRR/2 data. Therefore, in order to satisfy the needs of hydrological and meteorological modeling for the long-term data of snow covered area, an algorithm of subpixel snow mapping was established using AVHRR/2 10-day composited data of Qinghai-Tibet Plateau in the work. Moreover, linear spectral mixture model and multiple endmember spectral mixture analysis were introduced into the algorithm. According to the algorithm, the percent data of subpixel snow covered area can be obtained from AVHRR/2 data, and the subpixel snow covered area data can be extended to more than 30 years.For obtaining the required endmembers of ground objects, the method of selecting endmembers was found, in light of which snow endmembers and non-snow endmembers were selected automatically from AVHRR/2 images using the data of NDVI, channel 1 and channel 2 as the indices of selecting the endmembers, and the obtained snow and non-snow endmembers were saved in respective databases. In the work, this algorithm was validated against TM data, and the root-mean-square errors of the validation result are close to 0.1. For such a mountainous region as Qinghai-Tibet Plateau, this algorithm is very reliable and feasible.Because the used AVHRR/2 data are the 10-day composed data, the daily subpixel snow covered area cannot be obtained. In future, it is necessary to implement the pre-processing of AVHRR/2 data, and to obtain usable daily AVHRR/2 data. Then, after the obtained daily AVHRR/2 data were mapped at subpixel scale, the daily subpixel snow area data from AVHRR/2 data will be obtained.
-
Key words:
- snow cover /
- subpixel /
- multiple endmembers /
- AVHRR/2 data /
- Qinghai-Tibet Plateau /
- linear spectral mixture model
-
表 1 本算法的验证结果
Table 1. Validati Results of our Algorithm
TM数据 分辨率/m 相关系数(R) R2 均方根误差(RMSE) TM1 1 100 0.837 2 0.701 0 0.106 1 2 200 0.852 9 0.727 4 0.085 6 4 400 0.923 2 0.852 4 0.062 7 8 800 0.959 6 0.920 8 0.051 2 TM2 1 100 0.846 1 0.715 9 0.117 0 2 200 0.877 9 0.770 8 0.098 5 4 400 0.928 6 0.862 3 0.076 3 8 800 0.962 2 0.925 8 0.056 6 -
[1] Metsmki S, Vepsalainen J, Pulliainen J, et al. Improved Linear Interpolation Method for the Estimation of Snow-Covered Area from Optical Data[J].Remote Sensing of Environment, 2002, 82:64-78 doi: 10.1016/S0034-4257(02)00025-1 [2] Metsmki S, Anttila T S, Markus H J, et al. A Feasible Method for Fractional Snow Cover Mapping in Boreal Zone Based on a Reflectance Model[J].Remote Sensing of Environment, 2005, 95:77-95 doi: 10.1016/j.rse.2004.11.013 [3] Painter T H, Roberts D A, Green R O, et al. The Effect of Grain Size on Spectral Mixture Analysis of Snow-covered Area from AVIRIS[J].Remote Sensing of Environment, 1998, 65:320-332 doi: 10.1016/S0034-4257(98)00041-8 [4] Painter T H, Dozier J, Roberts D A, et al. Retrieval of Subpixel Snow-Covered Area and Grain Size from Imaging Spectrometer Data[J].Remote Sensing of Environment, 2003, 85(1):64-77 doi: 10.1016/S0034-4257(02)00187-6 [5] Painter T H, Rittger K, McKenzie C, et al. Retrieval of Subpixel Snow Covered Area, Grain Size, and Albedo from MODIS[J].Remote Sensing of Environment, 2009, 113:868-879 doi: 10.1016/j.rse.2009.01.001 [6] Salomonson V V, Appel I. Estimating Fractional Snow Cover from MODIS Using the Normalized Difference Snow Index[J].Remote Sensing of Environment, 2004, 89:351-360 doi: 10.1016/j.rse.2003.10.016 [7] Salomonson V V, Appel I. Development of the Aqua MODIS NDSI Fractional Snow Cover Algorithm and Validation Results[J].IEEE Trans Geosci Remot Sen, 2006, 44(7):1747-1756 doi: 10.1109/TGRS.2006.876029 [8] 唐志光, 王建, 彦立利, 等.基于MODIS的青藏高原亚像元积雪覆盖反演[J].干旱区资源与环境, 2013, 27(11):33-38 doi: 10.3969/j.issn.1003-7578.2013.11.006 Tang Zhiguang, Wang Jian, Yan Lili, et al. Estimating Sub-pixel Snow Cover from MODIS in Qinghai-Tibet Plateau[J]. Journal of Arid Land Resources and Environment, 2013, 27(11):33-38 doi: 10.3969/j.issn.1003-7578.2013.11.006 [9] Zhu J, Shi J, Wang Y. Subpixel Snow Mapping of the Qinghai-Tibet Plateau Using MODIS Data[J].International Journal of Applied Earth Observation and Geoinformation, 2012, 18:251-262 doi: 10.1016/j.jag.2012.02.001 [10] 王建.卫星遥感雪盖制图方法对比与分析[J].遥感技术与应用, 1999, 14(4):29-36 doi: 10.11873/j.issn.1004-0323.1999.4.29 Wang Jian. Comparison and Analysis on Methods of Snow Cover Mapping by Using Satellite Remote Sensing Data[J].Remote Sensing Technology and Application, 1999, 14(4):29-36 doi: 10.11873/j.issn.1004-0323.1999.4.29 [11] 延昊, 张国平.混合像元分解法提取积雪盖度[J].应用气象学报, 2004, 15(6):665-671 http://www.cqvip.com/QK/97586X/200406/11507711.html Yan Hao, Zhang Guoping. Unimixing Method Applied to NOAA-AVHRR Data for Snow Cover Estimation[J].Journal of Applied Meteorology, 2004, 15(6):665-671 http://www.cqvip.com/QK/97586X/200406/11507711.html [12] 延昊, 张国平.像元分解法提取积雪边界线[J].山地学报, 2004, 22(1):110-115 Yan Hao, Zhang Guoping. Unimixing Method Applied to Snow Boundary Estimation[J].Mountain Research, 2004, 22(1):110-115 [13] Zhu J, Shi J, Chu H, et al. Approaches to Using End-Members for Sub-pixel Snow Mapping with MODIS Data in Qinghai-Tibet Plateau[C]. The 2010 IEEE International Geoscience and Remote Sensing Symposium (IGARSS 2010), Honolulu, USA, 2010 [14] 陈晓娜, 包安明, 张红利, 等.基于混合像元分解的MODIS积雪面积信息提取及其精度评价[J].资源科学, 2010, 32(9):1761-1768 http://www.oalib.com/paper/1704853 Chen Xiaona, Bao Anming, Zhang Hongli, et al. A Study on Methods and Accuracy Assessment for Extracting Snow Covered Areas from MODIS Images Based on Pixel Unmixing:A Case on the Middle of the Tianshan Mountain[J].Resources Science, 2010, 32(9):1761-1768 http://www.oalib.com/paper/1704853 [15] Eidenshink J C, Faundeen J L. 1-km AVHRR Global Land Dataset:First Stages in Implementation[J].International Journal of Remote Sensing, 1994, 15:3443-3462 doi: 10.1080/01431169408954339 [16] 张洪恩. 青藏高原中分辩率亚像元雪填图算法研究[D]. 北京: 中国科学院遥感应用研究所, 2004 http://cdmd.cnki.com.cn/Article/CDMD-80070-2005121825.htm Zhang Hong'en. Moderate Sub-pixel Snow Mapping Algorithm on Tibetan Plateau[D]. Beijing:Institute of Remote Sensing Applications, CAS, 2004 http://cdmd.cnki.com.cn/Article/CDMD-80070-2005121825.htm [17] 李微, 方圣辉, 佃袁勇.基于光谱分析的云检测算法研究[J].武汉大学学报信息科学版, 2005, 30(5):435-438 doi: 10.3321/j.issn:1671-8860.2005.05.014 Li Wei, Fang Shenghui, Dian Yuanyong.Cloud Detection in MODIS Data Based on Spectrum Analysis[J].Geomatics and Information Science of Wuhan University, 2005, 30(5):435-438 doi: 10.3321/j.issn:1671-8860.2005.05.014 -