本发明属于卫星遥感图像处理与应用的
技术领域：
：。
背景技术：
：：积雪是地球表面最为活跃的自然要素之一，直接影响全球和大陆尺度的地表辐射平衡和能量交换、水文过程和气候变化。传统的积雪覆盖信息的获取方法主要是通过气象站点得到相应的数据，但是地形或物流成本高等因素，在测量方面仍然存在很大的局限性，特别是在发展中国家，这降低了水务机构直接评估积雪特征的能力，并不能大范围的得到准确的雪盖数据。监测雪参数(如积雪面积和雪水当量)是一项具有挑战性的工作。遥感技术为积雪监测的深入勘探和研究提供了极其重要的优势。数十颗卫星被用于监测积雪，包括光学遥感、被动微波和主动微波遥感：积雪具有明显的光谱特征，simica.等对基于光谱特征的modis积雪产品进行精度分析评价，证实modis积雪产品在晴空条件下得到的雪盖信息精确度高。(参见validationofvegetation,modis,andgoes+ssm/isnow-coverproductsovercanadabasedonsurfacesnowdepthobservations)。在光学遥感中，森林遮挡一直是积雪识别的主要影响因素，汪凌霄等分析了红光波段多角度反射率在林区不同积雪覆盖状况下的表现差异，改善林下积雪监测精度(参见天山典型林带积雪的多角度遥感识别)，wang等对比森林地区积雪使用landsat8数据和gf1号数据验证近红外波段(0.85-0.89)比可见光波段更适用于森林积雪遥感，提出ndsfi指数(参见retrievingsnowcoverinforestsofqilianmountainsfromlandsatoptionallandimager)。被动微波积雪产品不受气候条件的影响，主要用于全球范围内的雪深和雪水当量(swe)的研究。研究有yang等基于ssm/i亮温数据对中国地区雪深算法改进(参见improvementofsnowdepthestimationusingssm/ibrightnesstemperatureinchina)。星载主动微波传感器主要为合成孔径雷达sar，常用的有envisat卫星和radarsat-2卫星。在积雪研究中，主要用来雪深和雪水当量的反演，nagler等通过哨兵的光学传感器和sar协同作用，提出一种雪水当量和融雪面积的算法(参见thesentinel-1mission:newopportunitiesforicesheetobservations)到目前为止，国内外研究学者提出了很多积雪提取算法，但仍存在一些明显的缺陷：(1)算法复杂度高；(2)光学传感器受森林冠层遮挡影响，林下积雪覆盖度监测精度低；(3)被动微波数据空间分辨率较低；(4)主动微波遥感在区分土壤和干雪的后向散射信号方面仍然存在问题，没有成熟的雪参数产品，也没有满足一定时空要求的基于星载雷达的全球swe观测卫星。技术实现要素：中国东北地区的下垫面类型可以分为森林区和非森林区。对于非森林区域，利用landsat8oli数据结合snowmap算法得到“真实”的地表积雪信息，使用fy-3b数据计算积雪指数；然后，通过线性回归建立积雪指数和积雪覆盖度的关系；最后，选取精度较高的积雪指数回归方程，得到无云条件下的积雪覆盖产品。针对森林地区，为了解决现有光学卫星数据受到森林冠层遮挡，影响林下积雪监测准确性的问题，本发明采用一种基于计算森林透射率，将林下积雪覆盖度转换为星载多光谱数据的反射率的函数的提取方法，降低森林冠层遮挡的影响，有效提取林下积雪覆盖度信息。本发明采用的技术方案具体步骤如下：步骤一、影像预处理：将极轨卫星fy-3b数据进行辐射定标和几何校正处理，然后根据经纬度信息将图像裁剪选取需要研究的图像区域，根据modis土地分类产品igbp标准将所选取的区域划分成分为森林区域和非森林地区。步骤二、在非森林地区积雪覆盖度测算：将从landsat8获取的数据进行预处理后，使用snowmap算法获取真实的积雪二值图，再经重采样成1km空间分辨率的积雪覆盖率数据；将fy-3b获取的多波段数据进行波段运算，计算归一化差分积雪指数ndsi，计算公式如公式(1)所示：ndsi＝(green-swir)/(green+swir)…………………………(1)其中，green代表fy-3b数据中得到绿色波段的反射率，swir代表fy-3b数据中的短波红外波段的反射率。使用最小二乘法对重采样的积雪覆盖率和归一化差分积雪指数ndsi拟合线性回归方程，建立积雪指数和积雪覆盖度的函数关系，如公式(2)所示，再结合从landsat8获取的数据进行验证：fsc＝a×ndsi+b……………………………………………(2)fsc为积雪覆盖度；a为增益，即拟合函数的斜率；b为偏移量，即拟合函数的截距；将已获得的归一化差分积雪指数ndsi代入公式(2)获得非森林地区积雪覆盖度；步骤三、森林地区积雪覆盖度测算：森林地区积雪覆盖度fsc存在如下关系式：其中，选择研究区域降雪前和降雪后的fy-3b经过预处理的数据获得可见光反射率数据ρλ，obs(fsc)；统计多日非森林地区中有雪像元和无雪像元的卫星可见光反射率数据的平均值用于替代森林下积雪像元的可见光反射率ρλ，snow和森林下无雪像元的可见光反射率ρλ，ground，并设置为固定参数值；其中，代表森林透射率，在非森林地区，透射率接近1，在茂密的森林里，它接近0，通过如下方式获取：定位冬季森林区域中积雪完全覆盖地区(fsc＝1)，获得其对应的卫星可见光反射率ρλ，obs(fsc＝1)，根据公式(4)可以计算森林透射率步骤四、影像后处理：对星载多光谱得到的东北地区积雪覆盖度产品进行均值滤波，结合八邻域内的像元积雪覆盖度信息过滤零散的像元点，如果八邻域内都是无雪像元，将中心像元的值改为无雪像元；如果八邻域都是积雪像元，而中心像元为无雪像元，将中心像元改为八邻域积雪覆盖度的均值。降低由于卫星数据的噪声影响引起的像元影响，得到改善后的东北地区积雪覆盖度产品。在步骤一影像预处理的过程中的预处理步骤包括：(a)辐射定标：利用完整的遥感图像处理平台(envi)读取遥感数据的头文件信息，通过radiometriccalibration工具将数据中的dn值转换为反射率或者亮温数据，也可以通过波段运算，结合数据自身的增益和偏移量的属性信息，获取反射率数据。(b)几何校正：利用空间位置变换关系，采用计算公式和经纬度信息进行的校正。输入经纬度信息后通过envi软件生成glt文件(geographiclookuptable,地理位置查找表文件)。结合查找表文件给原始数据赋予地理信息。(c)根据经纬度信息通过东北地区矢量数据裁剪原始数据，得到感兴趣区域(roi)步骤二、非森林地区的积雪监测过程中：(a)snowmap算法:通过计算landsat8数据中ndsi指数提取积雪范围，设置ndsi的阈值为0.4，大于0.4为积雪像元。结合多光谱的特征去除水体和暗像元干扰(近红外波段>0.11)。得到空间分辨率为30m的积雪二值图。(b)重采样：将空间分辨率为30m的积雪覆盖度重采成1km空间分辨率的积雪覆盖度产品，计算1km×1km范围内的积雪像元个数占总像元的比例，作为积雪覆盖度的“真实”信息。步骤三、森林积雪监测的过程中：(a)数据验证：利用envi平台进行高斯施密特正交变换融合，将辐射校正后的gf-2数据融合成空间分辨率为1m的多光谱遥感数据，观察林下积雪覆盖情况。对比基于森林透射率算法的积雪覆盖度监测结果和传统基于ndsi算法监测结果受森林茂密程度的影响程度。本发明的有益效果：本发明根据积雪指数和森林透射率可以快速、准确的提取东北地区积雪覆盖度，解决了现有积雪覆盖度提取算法复杂度高，忽略多光谱信息，林下积雪受到森林冠层遮挡的影响的问题。此外，本发明通过卫星遥感数据可以大范围连续观测的优势，降低云层干扰，为森林透射率的计算提供大量有效的数据，改善林下积雪监测结果。本发明为中国东北地区积雪覆盖度监测提供有效的算法，为气候变化，水文过程以及雪灾监测提供一定的技术支持。附图说明图1是本发明基于星载多光谱遥感数据东北地区积雪覆盖度监测流程图。图2是本发明实施例1中landsat8数据中snowmap算法监测结果。图3是本发明实施例1中积雪指数和重采样后的积雪覆盖度数据的散点图。图4是本发明实施例1中中国东北地区森林反射率分布情况。图5是本发明实施例1中2月和10月多日融合成的无云遥感数据。图6是本发明实施例1中基于森林透射率积雪监测结果和gf-2数据。图7是本发明实施例1中基于森林透射率积雪监测结果和基于ndsi指数积雪监测结果与森林茂密程度的直方图。图8是本发明实例1中中国东北地区2016年1月4日积雪覆盖度产品。具体实施方式实施例1：如表1.所示，采用2016年1月空间分辨率为1km的多光谱fy-3b数据作为实验数据，根据土地类型将东北地区分为森林地区和非森林地区。对于非森林区域，利用landsat8oli数据结合snowmap算法得到“真实”的地表积雪信息，使用fy-3b数据计算积雪指数；然后，通过线性回归建立积雪指数和积雪覆盖度的关系；最后，选取精度较高的积雪指数回归方程，得到无云条件下的积雪覆盖产品。针对森林地区，利用scamod模型计算森林透射率，将林下积雪覆盖度表示为fy-3b可见光反射率和森林透射率函数，生成更准确的林下积雪覆盖度产品，并且结合gf-2数据进行积雪覆盖度验证。表1.fy-3bvirr光谱特征步骤一：影像预处理(a)辐射定标：利用完整的遥感图像处理平台(envi)读取遥感数据的头文件信息，通过radiometriccalibration工具将数据中的dn值转换为反射率或者亮温数据(b)几何校正：利用空间位置变换关系，采用计算公式和经纬度信息进行的校正。输入经纬度信息后通过envi软件生成glt文件(geographiclookuptable,地理位置查找表文件)。结合查找表文件给原始数据赋予地理信息。步骤二：非森林地区积雪覆盖度监测(a)选取2016年1月5日的一景landsat8数据使用snowmap算法获取真实的积雪二值图，具体流程为通过计算ndsi指数提取积雪范围，如公式(1)所示，green和swir分别为landsat8数据中绿色波段和短波红外波段，设置ndsi的阈值为0.4，结合多光谱的特征去除水体和暗像元干扰(近红外波段b5>0.11)。得到该地区的积雪二值图。经重采样成1km空间分辨率的积雪覆盖率数据，研究区域和积雪监测结果如图2所示。ndsi＝(green-swir)/(green+swir)…………………(1)(b)计算积雪指数ndsi方法通常用于利用卫星数据绘制雪覆盖图,可以很好的鉴别积雪和其他地物。ndsi利用了雪在可见光(绿色)和短波红外(swir)区域的高反射率和低反射率。同时还有其他积雪指数如s3指数，ndfsi指数。计算公式如公式(1)所示ndsi＝(green-swir)/(green+swir)…………………………(1)fy-3b获得的积雪指数与积雪覆盖度散点图如图3所示，随着积雪指数增加，积雪覆盖度变大。大量研究表明，非森林地区积雪覆盖度可以通过基于积雪指数的线性回归方程计算，将样本散点图中积雪覆盖度误差大于15％的样本删除，结合最小二乘法拟合线性回归方程的偏移量和增益，如公式(2)所示fsc＝a×ndsi+b……………………………………(2)fsc为积雪覆盖度，a为增益，b为偏移量，计算求得a＝2.0371，b＝-132.4489。(c)步骤三：森林地区积雪覆盖度监测森林冠层影响林下积雪覆盖度的识别精度，基于ndsi指数得到的积雪覆盖度会有一定程度的低估，对于茂密森林地区，可能会误判成该地区没有积雪覆盖，本文主要使用scamod模型计算森林透射率，改善林下积雪覆盖度监测的方法。scamod模型是基于卫星观测反射率转换成积雪覆盖度的函数的半经验模型。有效的森林透射率是由每个像元的面积和不同地表上的反射率(湿雪，森林冠层和无雪的地面上)作为模型参数，可应用于各种传感器在可见光和近红外范围的数据。卫星观测的森林积雪覆盖区的可见光反射率数据可由二部分可见光反射率信号组成，可以用公式(3)表示，其中森林冠层可见光反射率信号(公式(3)中等式右边第一项)，森林下积雪和无雪混合像元的可见光反射率信号(公式(3)中等式右边第二项)，其中，ρλ，obs(fsc)代表积雪盖度为fsc的森林区域可见光(波长为λ)反射率数据；代表森林透射率，在非森林地区，透射率接近1，在茂密的森林里，它接近0。ρλ，forest代表森林冠层的可见光反射率数据；ρλ，snow代表森林下积雪像元的可见光反射率；ρλ，ground代表森林下无雪像元的可见光反射率。由于森林下积雪像元可见光反射率卫星无法观测，需要采用非森林区域中有雪像元和无雪像元的卫星可见光反射率数据的平均值分别替代ρλ，snow、ρλ，ground。ρλ，forest由降雪前森林区域的卫星可见光反射率数据统计平均获得。经过统计，观测期间ρλ，snow、ρλ，ground和ρλ，forest是固定参数值。定位冬季森林区域中积雪完全覆盖地区(fsc＝1)，获得其对应的卫星可见光反射率ρλ，obs(fsc＝1)，根据公式(4)可以计算森林透射率通过公式(3)和公式(4)，最后可以获得森林积雪覆盖度fsc：其中ρλ，obs(fsc)为卫星观测到的森林区域可见光反射率数据，ρλ,snow＝120、ρλ,ground＝40、ρλ,forest、为已知参数。当大面积观测时，不同地区森林生长状态不同，森林郁闭度和种类都会影响输入参数的准确性。图4为东北地区不同森林反射率分布状况，说明森林分布不均匀，用单一均值不足以代表整个东三省地区森林反射率。因此，本发明通过遥感数据可以连续大范围观测的特点，统计每个像元的森林透射率数据。统计的过程中光学传感器可能受到云的影响，不能得到全部森林地区的反射率参数，因此，使用多天数据融合产品。图5分别为2016年2月1号到10号，2015年10月20号到11月1号融合出东北地区无云fy-3b数据，将其作为scamod模型的主要参数。步骤四：影像后处理gf-2数据和landsat8数据幅宽有限，可能没有云像元干扰监测结果。然而，对中国东北地区积雪覆盖度大范围监测的过程中，存在云干扰的问题，通过考虑空间的邻域信息，将8邻域内的非云像元的积雪覆盖度的均值替代云像元，降低云的干扰，生成更准确的东北地区积雪覆盖度产品，如图8所示。实验结果：在稳定的积雪期，在小范围内认为林下积雪覆盖度是一致的或者值相近的，选取的gf-2数据具有大量的森林像元，且茂密程度不同，将gf-2号数据融合为空间分辨率为1m的影像，从gf-2号数据中人工解译到林下积雪，对比不同方法在该地区林下积雪监测的结果，如图6、图7所示。图6(a)为2016年2月15日gf-2在黑龙江省东部伊春市的假彩色影像，中间区域城镇，四周为森林。图6(b)为研究地区归一化植被指数(ndvi)，(为处理方便，将ndvi指数*100，小于0的区域设为0)，用来表示目标区域中森林的茂密程度。对比图6(c)中modis积雪产品(mod10a1和myd10a1融合)、图6(d)中fy-3b数据基于ndsi指数的回归方程的监测结果和图6(e)中fy-3b数据基于scamod模型对该区域积雪覆盖度监测的结果，图7为基于森林透射率积雪监测结果和基于ndsi指数积雪监测结果与森林茂密程度的直方图。实验结果表明在高空间分辨率gf-2数据观测到部分林下积雪的条件下，modis产品和基于ndsi指数计算的积雪覆盖度结果为无雪或少量雪，积雪覆盖度随着森林茂密程度的增加而减小。相对而言scamod模型监测到了林下积雪，并且监测结果受森林茂密程度的影响较小(如图7所示)。证明了本专利可以为林下积雪监测提供更有效的方法。当前第1页12当前第1页12