




技术领域technical field
本发明涉及卫星传感器数据处理与应用技术领域,尤其涉及一种用于多卫星传感器的积雪覆盖度混合像元分解方法。The invention relates to the technical field of satellite sensor data processing and application, in particular to a snow coverage mixed pixel decomposition method for multi-satellite sensors.
背景技术Background technique
积雪有助于地球的辐射能量平衡,并作为广泛的蓄水层,影响着各种气候和水文过程。常规的地面站监测不能准确地获取大范围的监测结果,地面站的监测范围是否具有代表性直接影响最终的监测精度。随着卫星技术的发展,光学遥感中传感器丰富的波段信息可以提供准确的积雪覆盖信息,同时实现大范围的观测。Snow cover contributes to the Earth's radiative energy balance and acts as an extensive aquifer that affects various climatic and hydrological processes. Conventional ground station monitoring cannot accurately obtain large-scale monitoring results, and whether the monitoring range of ground stations is representative directly affects the final monitoring accuracy. With the development of satellite technology, the rich band information of sensors in optical remote sensing can provide accurate snow cover information and realize large-scale observation at the same time.
发明内容Contents of the invention
(一)要解决的技术问题(1) Technical problems to be solved
鉴于现有技术的上述缺点、不足,本发明提供一种用于多卫星传感器的积雪覆盖度混合像元分解方法。In view of the above-mentioned shortcomings and deficiencies of the prior art, the present invention provides a snow coverage mixed pixel decomposition method for multi-satellite sensors.
(二)技术方案(2) Technical solution
为了达到上述目的,本发明采用的主要技术方案包括:In order to achieve the above object, the main technical solutions adopted in the present invention include:
本发明实施例提供一种用于多卫星传感器的积雪覆盖度混合像元分解方法,包括:An embodiment of the present invention provides a snow coverage mixed pixel decomposition method for multi-satellite sensors, including:
S1、基于预先设定的多个遥感卫星传感器,获取相应的影像数据;S1. Obtain corresponding image data based on multiple preset remote sensing satellite sensors;
所述影像数据包括:多光谱地表反射率数据、成像几何以及掩膜数据;The image data includes: multi-spectral surface reflectance data, imaging geometry and mask data;
S2、基于所述多光谱反射率数据和预先设定的影像端元提取规则进行多种地类的端元提取,获取地类的端元集合;S2. Based on the multi-spectral reflectance data and the preset image endmember extraction rules, perform endmember extraction of various land types, and obtain endmember sets of land types;
所述地类的端元集合包括多个地类的端元;The endmember set of the land category includes endmembers of a plurality of land categories;
当提取失败则从预先获取的辅助数据中调取地物的地面测量光谱曲线,并根据预先获取的遥感卫星传感器的光谱响应函数,计算匹配传感器波段光谱响应的端元,获取地面测量端元库;When the extraction fails, the ground measurement spectral curve of the ground object is retrieved from the pre-acquired auxiliary data, and according to the pre-acquired spectral response function of the remote sensing satellite sensor, the endmember matching the sensor band spectral response is calculated to obtain the ground measurement endmember library ;
所述预先获取辅助数据包括积雪、植被、土壤光谱曲线,以及光谱响应数据;The pre-acquired auxiliary data includes snow cover, vegetation, soil spectral curves, and spectral response data;
S3、基于所述地类的端元集合,获取典型地类端元库;S3. Based on the endmember set of the land type, obtain a typical endmember library of the land type;
S4、基于所述典型地类端元库和地面测量端元库,进行多端元光谱混合分析获取积雪覆盖度;S4. Based on the typical land-type endmember library and the ground measurement endmember library, perform multi-endmember spectral mixing analysis to obtain snow coverage;
S5、计算卫星像元积雪覆盖度产品。S5. Calculate the snow coverage product of the satellite pixel.
优选的,preferred,
所述影像数据包括:MODIS影像数据、Landsat TM/ETM+/OLI影像数据、Sentinel-2MSI影像数据、FY-4A AGRI影像数据以及Himawari-8 AHI影像数据。The image data includes: MODIS image data, Landsat TM/ETM+/OLI image data, Sentinel-2MSI image data, FY-4A AGRI image data and Himawari-8 AHI image data.
优选的,preferred,
所述步骤S2中预先设定的影像端元提取规则为:The preset image endmember extraction rule in the step S2 is:
所述MODIS影像数据、Sentinel-2/MSI影像数据、Landsat TM/ETM+/OLI影像数据中的积雪端元提取规则为:NDVI<-0.035且NDSI>0.75且R0.55>0.7;植被端元提取规则为:NDVI>0.7且NDSI<-0.4;土壤和岩石端元的提取规则为:0<NDVI<0.15且NDSI<-0.4;水体端元提取规则为:NDWI>0.2且R0.86<0.2;The snow cover endmember extraction rules in the MODIS image data, Sentinel-2/MSI image data, and Landsat TM/ETM+/OLI image data are: NDVI<-0.035 and NDSI>0.75 and R0.55 >0.7; vegetation endmember extraction The rules are: NDVI>0.7 and NDSI<-0.4; the extraction rules of soil and rock endmembers are: 0<NDVI<0.15 and NDSI<-0.4; the extraction rules of water endmembers are: NDWI>0.2 and R0.86 <0.2;
NDVI是指归一化植被指数;NDSI是指归一化积雪指数;R0.55指的是0.55微米波段的地表反射率;NDVI refers to the normalized difference vegetation index; NDSI refers to the normalized difference snow cover index; R0.55 refers to the surface reflectance in the 0.55 micron band;
Himawari-8/AHI影像数据和FY-4A/AGRI影像数据中的积雪端元提取规则为:NDVI<-0.03且NDSI>0.74且R0.55>0.5;植被端元提取规则为:NDVI>0.65且NDSI<-0.4;土壤和岩石端元的提取规则为:0<NDVI<0.15且NDSI<-0.4。The snow endmember extraction rules in Himawari-8/AHI image data and FY-4A/AGRI image data are: NDVI<-0.03 and NDSI>0.74 and R0.55 >0.5; the vegetation endmember extraction rules are: NDVI>0.65 and NDSI<-0.4; the extraction rules of soil and rock endmembers are: 0<NDVI<0.15 and NDSI<-0.4.
优选的,所述步骤S3包括:Preferably, said step S3 includes:
基于所述地类的端元集合,对所述地类的端元集合中的影像端元进行筛选处理,获取典型地类端元库。Based on the endmember set of the land type, the image endmembers in the endmember set of the land type are screened to obtain a typical endmember library of the land type.
优选的,所述步骤S3具体包括:Preferably, the step S3 specifically includes:
S31、采用公式(1)计算每一条端元光谱的矢量长度;S31, using formula (1) to calculate the vector length of each endmember spectrum;
其中,||r||为端元光谱的矢量长度;rk为第k波段的反射率;N为波段数目;Among them, ||r|| is the vector length of the endmember spectrum; rk is the reflectivity of the kth band; N is the number of bands;
S32、将积雪、植被、土壤地物每一类端元内的所有矢量长度进行升序排列,获取积雪、植被、土壤地物每一类端元相应的序列;S32. Arranging all vector lengths in each type of endmembers of snow cover, vegetation, and soil features in ascending order, and obtaining the corresponding sequence of each type of endmembers of snow cover, vegetation, and soil features;
S33、针对所述积雪、植被、土壤地物每一类端元相应的序列,将每一类内端元按照等间隔分割为n个子集;S33. For the sequences corresponding to each type of endmembers of the snow cover, vegetation, and soil features, divide each type of endmembers into n subsets at equal intervals;
其中,所述间隔为||Ri||;Wherein, the interval is ||Ri ||;
n为预先设定值;||R||max像元最大反射率;||R||min像元最小反射率;n is a preset value; ||R||max is the maximum reflectance of a pixel; ||R||min is the minimum reflectance of a pixel;
S34、对每一个子集内的所有光谱取中值或均值,作为所述子集的典型端元;S34. Taking the median or mean value of all spectra in each subset as a typical end member of the subset;
S35、基于所述每一子集的典型端元,获取典型地类端元库;S35. Based on the typical endmembers of each subset, obtain a typical endmember-like library;
所述典型地类端元库库中包括:积雪类地物的n个典型端元、植被类地物的n个典型端元以及土壤类地物的n个典型端元。The typical endmember library includes: n typical endmembers of snow-like features, n typical endmembers of vegetation-like features, and n typical endmembers of soil-like features.
优选的,所述步骤S4包括:Preferably, said step S4 includes:
S41、针对基于所述典型地类端元库和地面测量端元库中,将NDSI<0或1.6μm波段的反射率大于0.3的像元标记为无雪,对NDSI>0.7的像元标记为纯雪;S41. Based on the typical ground-like endmember library and the ground measurement endmember library, mark the pixels with NDSI<0 or 1.6 μm band reflectivity greater than 0.3 as no snow, and mark the pixels with NDSI>0.7 as pure snow;
S42、对提取积雪端元成功的像元标记为纯雪,对提取非雪的像元标记为无雪;S42. Mark the pixels successfully extracted from snow cover as pure snow, and mark the pixels extracted from non-snow as no snow;
S43、根据端元库的积雪、土壤和植被光谱,以及预先设定的端元面积比例约束,建立混合像元光谱数据库;S43. Establish a mixed pixel spectral database according to the snow cover, soil and vegetation spectra of the endmember library and the preset endmember area ratio constraints;
S44、根据所述混合像元光谱数据库,仅针对混合像元进行逐个处理,确定积雪-土壤二分模型与积雪-植被二分模型;S44. According to the mixed pixel spectral database, only the mixed pixel is processed one by one, and the snow-soil dichotomous model and the snow-vegetation dichotomous model are determined;
所述混合像元为标记为纯雪和无雪的像元外的像元;The mixed pixel is a pixel marked as pure snow and a pixel other than snow-free;
在像元的NDVI超过0.3时仅使用积雪-植被模型,否则按顺序使用积雪-土壤模型与积雪-植被模型;Only use the snow-vegetation model when the NDVI of the pixel exceeds 0.3, otherwise use the snow-soil model and the snow-vegetation model in sequence;
S45、计算混合像元光谱与不同积雪覆盖度对应的混合像元光谱数据的差异,记录通过残差检验的RMSE最小的结果以及对应的积雪覆盖度;S45. Calculate the difference between the mixed pixel spectrum and the mixed pixel spectral data corresponding to different snow coverage degrees, and record the result of the smallest RMSE passing the residual error test and the corresponding snow coverage degree;
其中,RMSE是均方根误差(root mean square error);Among them, RMSE is the root mean square error (root mean square error);
若均不符合残差检验且当前像元NDVI不超过0.3时,记录通过残差检验的RMSE,以积雪-植被模型重复步骤S45,得到新的积雪覆盖度解算的RMSE,对比二者选择较小者对应的积雪覆盖度作为最终结果。If none of them meet the residual test and the NDVI of the current pixel does not exceed 0.3, record the RMSE that passed the residual test, and repeat step S45 with the snow cover-vegetation model to obtain the new RMSE calculated by the snow cover degree, and compare the two The snow coverage corresponding to the smaller one is selected as the final result.
优选的,所述步骤S43中预先设定的端元面积比例约束为:Preferably, the preset end-member area ratio constraint in step S43 is:
约束条件为:The constraints are:
Fi≥0;Fi ≥ 0;
其中,Rλ是波长λ处混合像元的反射率;Riλ和Fi分别为第i个端元的反射率和丰度;ελ为拟合残差;M为端元数量;Among them, Rλ is the reflectivity of the mixed pixel at wavelength λ; Riλ and Fi are the reflectivity and abundance of the i-th endmember, respectively; ελ is the fitting residual; M is the number of endmembers;
在解算时,满足预先设定的第一条件、第二条件、第三条件、第四条件、第五条件;When solving, satisfy the preset first condition, second condition, third condition, fourth condition and fifth condition;
所述第一条件为:端元的面积比例符合“非负”及“和为一”约束,但允许1%的误差;The first condition is: the area ratio of the end members conforms to the constraints of "non-negative" and "sum to one", but an error of 1% is allowed;
所述第二条件为:在像元的NDVI超过0.3时使用积雪-植被两种地类混合模型,否则按顺序使用积雪-土壤模型与积雪-植被模型;The second condition is: when the NDVI of the pixel exceeds 0.3, use the snow cover-vegetation two land type mixed models, otherwise use the snow cover-soil model and the snow cover-vegetation model in sequence;
所述第三条件为:残差检验时,所有波段的残差的RMSE在不得超过0.025;The third condition is: when the residual error is checked, the RMSE of the residual error of all bands shall not exceed 0.025;
所述第四条件为:残差检验时,一半波段数量的连续波段的残差不超过0.025;The fourth condition is: when the residual error is checked, the residual error of the continuous bands with half the number of bands does not exceed 0.025;
所述第五条件为:RMSE最小的结果作为各端元的面积比例。The fifth condition is: the minimum RMSE is used as the area ratio of each end member.
优选的,所述步骤S5包括:Preferably, said step S5 includes:
S51、在掩膜图层基础上对积雪覆盖度反演结果图层进行掩膜;S51. Masking the snow coverage inversion result layer on the basis of the mask layer;
S52、选定针对积雪覆盖度值域分布范围相应的颜色表,采用无符号整型、LZW压缩方式的Geotiff文件,数据文件同时保存影像原始投影信息,该数据文件为计算的积雪覆盖度产品。S52, select the color table corresponding to the distribution range of the snow cover degree, adopt the Geotiff file of unsigned integer type and LZW compression mode, the data file saves the original projection information of the image at the same time, and the data file is the calculated snow cover degree product.
(三)有益效果(3) Beneficial effects
本发明的有益效果是:本发明的一种用于多卫星传感器的积雪覆盖度混合像元分解方法,适用于常用遥感光学影像包括MODIS、Sentinel-2/MSI、Landsat TM/ETM+/OLI影像,以及静止气象卫星Himawari-8/AHI、FY-4A/AGRI影像,可扩展性好、为积雪覆盖遥感制图以及多源积雪覆盖结果融合、对比分析等研究提供了实用工具。The beneficial effects of the present invention are: a method for decomposing mixed pixels of snow coverage for multi-satellite sensors of the present invention is applicable to commonly used remote sensing optical images including MODIS, Sentinel-2/MSI, Landsat TM/ETM+/OLI images , and Himawari-8/AHI, FY-4A/AGRI images of geostationary meteorological satellites, which have good scalability and provide practical tools for remote sensing mapping of snow cover, fusion of multi-source snow cover results, comparative analysis and other research.
基于多端元光谱混合分析理论和影像端元库自动构建技术,算法具有物理基础,精度高。Based on the theory of multi-endmember spectral mixing analysis and the automatic construction technology of image endmember library, the algorithm has a physical basis and high precision.
利用基于影像的端元自动提取和选择技术,以及积雪覆盖度解算的查找表技术,计算效率高,可根据本方法实现积雪覆盖的遥感实时监测。Using the image-based automatic extraction and selection technology of endmembers and the look-up table technology for snow coverage calculation, the calculation efficiency is high, and the remote sensing real-time monitoring of snow coverage can be realized according to this method.
附图说明Description of drawings
图1为本发明的一种用于多卫星传感器的积雪覆盖度混合像元分解方法流程图;Fig. 1 is a kind of flow chart of the method for decomposing snow cover mixed pixel for multi-satellite sensor of the present invention;
图2为本发明实施例中用于多卫星传感器的积雪覆盖度混合像元分解方法示意图;FIG. 2 is a schematic diagram of a snow coverage mixed pixel decomposition method for multi-satellite sensors in an embodiment of the present invention;
图3为本发明实施例中MODIS影像积雪类型自动端元选择(a)与端元提取(b)结果示例图;Fig. 3 is an example diagram of the results of MODIS image snow type automatic endmember selection (a) and endmember extraction (b) in the embodiment of the present invention;
图4为本发明实施例中MODIS影像土壤类型自动端元选择(c)与端元提取(d)结果示例图;Fig. 4 is an example diagram of the results of MODIS image soil type automatic endmember selection (c) and endmember extraction (d) in the embodiment of the present invention;
图5为本发明实施例中MODIS影像植被类型自动端元选择(e)与端元提取(f)结果示例图。Fig. 5 is an example diagram of the results of automatic endmember selection (e) and endmember extraction (f) of MODIS image vegetation types in the embodiment of the present invention.
具体实施方式Detailed ways
为了更好的解释本发明,以便于理解,下面结合附图,通过具体实施方式,对本发明作详细描述。In order to better explain the present invention and facilitate understanding, the present invention will be described in detail below through specific embodiments in conjunction with the accompanying drawings.
为了更好的理解上述技术方案,下面将参照附图更详细地描述本发明的示例性实施例。虽然附图中显示了本发明的示例性实施例,然而应当理解,可以以各种形式实现本发明而不应被这里阐述的实施例所限制。相反,提供这些实施例是为了能够更清楚、透彻地理解本发明,并且能够将本发明的范围完整的传达给本领域的技术人员。In order to better understand the above technical solutions, the following will describe exemplary embodiments of the present invention in more detail with reference to the accompanying drawings. Although exemplary embodiments of the present invention are shown in the drawings, it should be understood that the invention may be embodied in various forms and should not be limited to the embodiments set forth herein. Rather, these embodiments are provided so that the present invention can be more clearly and thoroughly understood, and the scope of the present invention can be fully conveyed to those skilled in the art.
参见图1,本实施例提供一种用于多卫星传感器的积雪覆盖度混合像元分解方法,包括:Referring to Fig. 1, the present embodiment provides a snow coverage mixed pixel decomposition method for multi-satellite sensors, including:
S1、基于预先设定的多个遥感卫星传感器,获取相应的影像数据;S1. Obtain corresponding image data based on multiple preset remote sensing satellite sensors;
所述影像数据包括:多光谱地表反射率数据、成像几何以及掩膜数据;The image data includes: multi-spectral surface reflectance data, imaging geometry and mask data;
S2、基于所述多光谱反射率数据和预先设定的影像端元提取规则进行多种地类的端元提取,获取地类的端元集合;S2. Based on the multi-spectral reflectance data and the preset image endmember extraction rules, perform endmember extraction of various land types, and obtain endmember sets of land types;
所述地类的端元集合包括多个地类的端元;The endmember set of the land category includes endmembers of a plurality of land categories;
当提取失败则从预先获取的辅助数据中调取地物的地面测量光谱曲线,并根据预先获取的遥感卫星传感器的光谱响应函数,计算匹配传感器波段光谱响应的端元,获取地面测量端元库;When the extraction fails, the ground measurement spectral curve of the ground object is retrieved from the pre-acquired auxiliary data, and according to the pre-acquired spectral response function of the remote sensing satellite sensor, the endmember matching the sensor band spectral response is calculated to obtain the ground measurement endmember library ;
所述预先获取辅助数据包括积雪、植被、土壤光谱曲线,以及光谱响应数据;The pre-acquired auxiliary data includes snow cover, vegetation, soil spectral curves, and spectral response data;
本实施例中辅助数据包括从约翰斯霍普金斯大学光谱库提取的积雪、植被、土壤光谱,以及从各传感器影像产品官方获取的光谱响应数据。The auxiliary data in this embodiment include snow cover, vegetation, and soil spectra extracted from the spectral library of Johns Hopkins University, and spectral response data officially obtained from image products of various sensors.
S3、基于所述地类的端元集合,获取典型地类端元库;S3. Based on the endmember set of the land type, obtain a typical endmember library of the land type;
S4、基于所述典型地类端元库和地面测量端元库,进行多端元光谱混合分析获取积雪覆盖度;S4. Based on the typical land-type endmember library and the ground measurement endmember library, perform multi-endmember spectral mixing analysis to obtain snow coverage;
S5、计算卫星像元积雪覆盖度产品并保存到文件。S5. Calculate the satellite pixel snow coverage product and save it to a file.
优选的,preferred,
所述影像数据包括:MODIS影像数据、Landsat TM/ETM+/OLI影像数据、Sentinel-2MSI影像数据、FY-4A AGRI影像数据以及Himawari-8 AHI影像数据。The image data includes: MODIS image data, Landsat TM/ETM+/OLI image data, Sentinel-2MSI image data, FY-4A AGRI image data and Himawari-8 AHI image data.
本实施例中MODIS影像数据(NASA提供的卫星像元地表反射率数据(MOD09GA与MYD09GA),其中包含成像几何和掩膜数据)、Landsat TM/ETM+/OLI影像数据(USGS提供的地表反射率数据)、Sentinel-2 MSI影像数据(欧空局提供的L2A地表反射率产品)、FY-4AAGRI影像数据(中国气象局提供的AGRI中国区多通道影像、云掩膜及成像几何产品)以及Himawari-8 AHI影像数据(日本气象厅生产的AHI多通道影像包括成像几何、云掩膜及成像几何产品)。In this embodiment, MODIS image data (satellite pixel surface albedo data (MOD09GA and MYD09GA) provided by NASA, which includes imaging geometry and mask data), Landsat TM/ETM+/OLI image data (surface albedo data provided by USGS ), Sentinel-2 MSI image data (L2A surface reflectance product provided by ESA), FY-4AAGRI image data (multi-channel image, cloud mask and imaging geometry product of AGRI China area provided by China Meteorological Administration), and Himawari- 8 AHI image data (the AHI multi-channel image produced by the Japan Meteorological Agency includes imaging geometry, cloud mask and imaging geometry products).
优选的,preferred,
所述步骤S2中预先设定的影像端元提取规则为:The preset image endmember extraction rule in the step S2 is:
所述MODIS影像数据、Sentinel-2/MSI影像数据、Landsat TM/ETM+/OLI影像数据中的积雪端元提取规则为:NDVI<-0.035且NDSI>0.75且R0.55>0.7;植被端元提取规则为:NDVI>0.7且NDSI<-0.4;土壤和岩石端元的提取规则为:0<NDVI<0.15且NDSI<-0.4;水体端元提取规则为:NDWI>0.2且R0.86<0.2。The snow cover endmember extraction rules in the MODIS image data, Sentinel-2/MSI image data, and Landsat TM/ETM+/OLI image data are: NDVI<-0.035 and NDSI>0.75 and R0.55 >0.7; vegetation endmember extraction The rules are: NDVI>0.7 and NDSI<-0.4; the extraction rules of soil and rock endmembers are: 0<NDVI<0.15 and NDSI<-0.4; the extraction rules of water endmembers are: NDWI>0.2 and R0.86 <0.2.
Himawari-8/AHI影像数据和FY-4A/AGRI影像数据中的积雪端元提取规则为:NDVI<-0.03且NDSI>0.74且R0.55>0.5;植被端元提取规则为:NDVI>0.65且NDSI<-0.4;土壤和岩石端元的提取规则为:0<NDVI<0.15且NDSI<-0.4。The snow endmember extraction rules in Himawari-8/AHI image data and FY-4A/AGRI image data are: NDVI<-0.03 and NDSI>0.74 and R0.55 >0.5; the vegetation endmember extraction rules are: NDVI>0.65 and NDSI<-0.4; the extraction rules of soil and rock endmembers are: 0<NDVI<0.15 and NDSI<-0.4.
本实施例中优选的,所述步骤S3包括:基于所述地类的端元集合,对所述地类的端元集合中的影像端元进行筛选处理,获取典型地类端元库(也就是图2中的影像自动提取与选择端元库)。Preferably in this embodiment, the step S3 includes: based on the endmember set of the land type, screening the image endmembers in the endmember set of the land type to obtain a typical endmember library (also It is the image automatic extraction and selection endmember library in Fig. 2).
本实施例中优选的,所述步骤S3具体包括:Preferably in this embodiment, the step S3 specifically includes:
S31、采用公式(1)计算每一条端元光谱的矢量长度;S31, using formula (1) to calculate the vector length of each endmember spectrum;
其中,||r||为端元光谱的矢量长度;rk为第k波段的反射率;N为波段数目。Among them, ||r|| is the vector length of the endmember spectrum; rk is the reflectivity of the kth band; N is the number of bands.
S32、将积雪、植被、土壤地物每一类端元内的所有矢量长度进行升序排列,获取积雪、植被、土壤地物每一类端元相应的序列。S32. Arrange in ascending order all vector lengths in each type of endmembers of snow cover, vegetation, and soil features, and obtain a sequence corresponding to each type of endmembers of snow cover, vegetation, and soil features.
S33、针对所述积雪、植被、土壤地物每一类端元相应的序列,将每一类内端元按照等间隔分割为n个子集;S33. For the sequences corresponding to each type of endmembers of the snow cover, vegetation, and soil features, divide each type of endmembers into n subsets at equal intervals;
其中,所述间隔为||Ri||;Wherein, the interval is ||Ri ||;
n为预先设定值;||R||max像元最大反射率;||R||min像元最小反射率。n is a preset value; ||R||max is the maximum reflectance of a pixel; ||R||min is the minimum reflectance of a pixel.
图3、图4、图5显示了一景MODIS影像中积雪、土壤以及植被类型的自动端元选择与端元提取结果,从影像中直接提取的植被端元数量巨大,经过端元选择,可以提取基本代表了类内光谱差异的少数代表性较好的端元。Figure 3, Figure 4, and Figure 5 show the results of automatic endmember selection and endmember extraction of snow cover, soil, and vegetation types in a MODIS image. The number of vegetation endmembers directly extracted from the image is huge. After endmember selection, A small number of well-represented endmembers that basically represent the spectral differences within the class can be extracted.
S34、对每一个子集内的所有光谱取中值或均值,作为所述子集的典型端元。S34. Take a median or mean value of all spectra in each subset, and use it as a typical end member of the subset.
S35、基于所述每一子集的典型端元,获取典型地类端元库。S35. Based on the typical endmembers of each subset, obtain a typical endmember library.
所述典型地类端元库库中包括:积雪类地物的n个典型端元、植被类地物的n个典型端元以及土壤类地物的n个典型端元。The typical endmember library includes: n typical endmembers of snow-like features, n typical endmembers of vegetation-like features, and n typical endmembers of soil-like features.
本实施例中优选的,所述步骤S4包括:Preferably in this embodiment, the step S4 includes:
S41、针对基于所述典型地类端元库和地面测量端元库中,将NDSI<0或1.6μm波段的反射率大于0.3的像元标记为无雪,对NDSI>0.7的像元标记为纯雪。S41. Based on the typical ground-like endmember library and the ground measurement endmember library, mark the pixels with NDSI<0 or 1.6 μm band reflectivity greater than 0.3 as no snow, and mark the pixels with NDSI>0.7 as pure snow.
S42、对提取积雪端元成功的像元标记为纯雪,对提取非雪的像元标记为无雪。S42. Mark the pixels whose snow cover end members are successfully extracted as pure snow, and mark the pixels whose non-snow is extracted as no snow.
S43、根据端元库的积雪、土壤和植被光谱,以及预先设定的端元面积比例约束,建立混合像元光谱数据库。S43. Establish a mixed pixel spectral database according to the snow cover, soil and vegetation spectra of the endmember library and the preset endmember area ratio constraints.
S44、根据所述混合像元光谱数据库,仅针对混合像元进行逐个处理,确定积雪-土壤二分模型与积雪-植被二分模型。S44. According to the mixed pixel spectral database, only the mixed pixel is processed one by one to determine a snow-soil dichotomous model and a snow-vegetation dichotomous model.
所述混合像元为标记为纯雪和无雪的像元外的像元。The mixed cells are cells other than those marked as pure snow and no snow.
在像元的NDVI超过0.3时仅使用积雪-植被模型,否则按顺序使用积雪-土壤模型与积雪-植被模型。When the NDVI of the pixel exceeds 0.3, only the snow cover-vegetation model is used, otherwise the snow cover-soil model and the snow cover-vegetation model are used in sequence.
S45、计算混合像元光谱与不同积雪覆盖度对应的混合像元光谱数据的差异,记录通过残差检验的RMSE最小的结果以及对应的积雪覆盖度;S45. Calculate the difference between the mixed pixel spectrum and the mixed pixel spectral data corresponding to different snow coverage degrees, and record the result of the smallest RMSE passing the residual error test and the corresponding snow coverage degree;
若均不符合残差检验且当前像元NDVI不超过0.3时,记录通过残差检验的RMSE,以积雪-植被模型重复步骤S45,得到新的积雪覆盖度解算的RMSE,对比二者选择较小者对应的积雪覆盖度作为最终结果。If none of them meet the residual test and the NDVI of the current pixel does not exceed 0.3, record the RMSE that passed the residual test, and repeat step S45 with the snow cover-vegetation model to obtain the new RMSE calculated by the snow cover degree, and compare the two The snow coverage corresponding to the smaller one is selected as the final result.
在本实施的实际应用中多端元光谱混合分析具体步骤如下:In the actual application of this implementation, the specific steps of multi-terminal spectral mixing analysis are as follows:
步骤1,对NDSI<0或1.6μm波段的反射率大于0.3的像元标记为无雪,对NDSI>0.7的像元标记为纯雪。In
步骤2,对提取积雪端元成功的像元标记为纯雪,对提取林地、草地等非雪的像元标记为无雪。Step 2: Mark the pixels that successfully extract the snow endmembers as pure snow, and mark the pixels that extract non-snow such as woodland and grassland as no snow.
步骤3,构建查找表(一幅影像构建一个查找表,根据像元混合光谱,在查找表中找到与之差异最小的“混合光谱”,满足其条件计算的积雪覆盖度就是估计的目标结果)。以1%为积雪覆盖度的步长,考虑1%的估算误差,根据端元库的积雪、土壤和植被光谱,以及前述端元面积比例约束,建立积雪-土壤二分模型与积雪-植被二分模型下的不同积雪覆盖度对应的混合像元光谱数据库。Step 3, build a lookup table (construct a lookup table for one image, find the "mixed spectrum" with the smallest difference in the lookup table according to the mixed spectrum of the pixel, and the calculated snow coverage that satisfies its conditions is the estimated target result ). Taking 1% as the step size of snow coverage and considering the estimation error of 1%, according to the snow cover, soil and vegetation spectra of the end-member pool, and the aforementioned end-member area ratio constraints, a snow-soil dichotomy model and snow cover model were established. - Mixed pixel spectral database corresponding to different snow cover degrees under the vegetation dichotomy model.
步骤4,开始计算,仅针对混合像元(即标记为纯雪和无雪的像元外)的像元进行逐个处理,首先确定混合像元模型。在像元的NDVI超过0.3时仅使用积雪-植被模型,否则按顺序使用积雪-土壤模型与积雪-植被模型。Step 4, start the calculation, and only process the mixed pixels (i.e., except the pixels marked as pure snow and snow-free) one by one, and firstly determine the mixed pixel model. When the NDVI of the pixel exceeds 0.3, only the snow-vegetation model is used, otherwise the snow-soil model and the snow-vegetation model are used in sequence.
步骤5,计算混合像元光谱与不同积雪覆盖度对应的混合像元光谱数据的差异,记录通过残差检验的RMSE最小的结果以及对应的积雪覆盖度。
步骤6,若均不符合残差检验且当前像元NDVI不超过0.3时,记录步骤5的RMSE,以积雪-植被模型重复步骤5,得到新的积雪覆盖度解算的RMSE,对比二者选择较小者对应的积雪覆盖度作为最终结果。Step 6. If none of them meet the residual test and the NDVI of the current pixel does not exceed 0.3, record the RMSE of
优选的,所述步骤S43中预先设定的端元面积比例约束为:Preferably, the preset end-member area ratio constraint in step S43 is:
约束条件为:The constraints are:
Fi≥0;Fi ≥ 0;
其中,Rλ是波长λ处混合像元的反射率;Riλ和Fi分别为第i个端元的反射率和丰度;ελ为拟合残差;M为端元数量。Among them, Rλ is the reflectivity of the mixed pixel at wavelength λ; Riλ and Fi are the reflectivity and abundance of the i-th end member, respectively; ελ is the fitting residual; M is the number of end members.
在解算时,满足预先设定的第一条件、第二条件、第三条件、第四条件、第五条件。During the calculation, the preset first condition, second condition, third condition, fourth condition and fifth condition are satisfied.
所述第一条件为:端元的面积比例符合“非负”及“和为一”约束,但允许1%的误差。The first condition is: the area ratio of the end members meets the constraints of "non-negative" and "sum-to-unity", but an error of 1% is allowed.
所述第二条件为:在像元的NDVI超过0.3时使用积雪-植被两种地类混合模型,否则按顺序使用积雪-土壤模型与积雪-植被模型。The second condition is: when the NDVI of the pixel exceeds 0.3, the mixed snow-vegetation model is used; otherwise, the snow-soil model and the snow-vegetation model are used in sequence.
所述第三条件为:残差检验时,所有波段的残差的RMS在不得超过0.025。The third condition is: when the residual is checked, the RMS of the residual of all bands shall not exceed 0.025.
所述第四条件为:残差检验时,一半波段数量的连续波段的残差不超过0.025。The fourth condition is: when checking the residuals, the residuals of half the number of continuous bands do not exceed 0.025.
所述第五条件为:RMSE最小的结果作为各端元的面积比例。The fifth condition is: the minimum RMSE is used as the area ratio of each end member.
优选的,所述步骤S5包括:Preferably, said step S5 includes:
S51、在掩膜图层基础上对积雪覆盖度反演结果图层进行掩膜。S51. Masking the snow cover degree inversion result layer based on the mask layer.
S52、选定针对积雪覆盖度值域分布范围相应的颜色表(积雪覆盖度为蓝-青-黄渐变色),采用无符号整型、LZW压缩方式的Geotiff文件,数据文件同时保存影像原始投影信息,该数据文件为计算的积雪覆盖度产品。S52, select the color table corresponding to the distribution range of the snow cover degree value range (the snow cover degree is blue-green-yellow gradient color), adopt the Geotiff file of unsigned integer type and LZW compression mode, and save the image at the same time in the data file Raw projection information, this data file is the product of calculated snow cover.
本实施例中的输出文件中保留积雪覆盖估算误差和成像几何图层,为积雪覆盖度的质量评估提供参考。The snow cover estimation error and the imaging geometry layer are retained in the output file of this embodiment, which provides a reference for the quality assessment of the snow cover degree.
本实施例中用于多卫星传感器的积雪覆盖度混合像元分解方法,根据不同卫星光学传感器光谱波段设置确定端元提取规则;根据多光谱影像自动提取积雪、土壤、植被等端元,建立典型端元库,并根据端元光谱矢量选择典型端元;根据典型端元以及最小二乘法对混合像元解算,获取最优解,实现积雪覆盖度反演。发明实施例提供的多卫星传感器的积雪覆盖度混合像元分解方法,适用于多种卫星影像包括MODIS、Sentinel-2/MSI、Landsat TM/ETM+/OLI影像,以及静止气象卫星Himawari-8/AHI、FY-4A/AGRI影像,通过多端元光谱混合分析理论和基于影像的端元自动提取和选择技术,实现复杂地表条件下混合像元积雪覆盖度的自动解算,其结果可用于山区水文模型、水资源管理、数值天气预报、陆面模式以及气候变化研究。In this embodiment, the snow coverage mixed pixel decomposition method used for multi-satellite sensors determines the endmember extraction rules according to the spectral band settings of different satellite optical sensors; automatically extracts snow cover, soil, vegetation and other endmembers based on multispectral images, Establish a typical endmember library, and select typical endmembers according to the endmember spectral vector; calculate the mixed pixels according to the typical endmembers and the least square method, obtain the optimal solution, and realize the inversion of snow coverage. The snow coverage mixed pixel decomposition method of the multi-satellite sensor provided by the embodiment of the invention is applicable to various satellite images including MODIS, Sentinel-2/MSI, Landsat TM/ETM+/OLI images, and stationary meteorological satellite Himawari-8/ AHI, FY-4A/AGRI images, through multi-endmember spectral hybrid analysis theory and image-based endmember automatic extraction and selection technology, realize automatic calculation of mixed pixel snow coverage under complex surface conditions, and the results can be used in mountainous areas Hydrological models, water resource management, numerical weather prediction, land surface models, and climate change research.
本实施例的一种用于多卫星传感器的积雪覆盖度混合像元分解方法,适用于常用遥感光学影像包括MODIS、Sentinel-2/MSI、Landsat TM/ETM+/OLI影像,以及静止气象卫星Himawari-8/AHI、FY-4A/AGRI影像,可扩展性好、为积雪覆盖遥感制图以及多源积雪覆盖结果融合、对比分析等研究提供了实用工具。基于多端元光谱混合分析理论和影像端元库自动构建技术,算法具有物理基础,精度高。利用基于影像的端元自动提取和选择技术,以及积雪覆盖度解算的查找表技术,计算效率高,可根据本方法实现积雪覆盖的遥感实时监测。A snow coverage mixed pixel decomposition method for multi-satellite sensors in this embodiment is applicable to commonly used remote sensing optical images including MODIS, Sentinel-2/MSI, Landsat TM/ETM+/OLI images, and stationary weather satellite Himawari -8/AHI, FY-4A/AGRI images have good scalability and provide practical tools for snow cover remote sensing mapping, multi-source snow cover result fusion, comparative analysis and other research. Based on the theory of multi-endmember spectral mixing analysis and the automatic construction technology of image endmember library, the algorithm has a physical basis and high precision. Using the image-based automatic extraction and selection technology of endmembers and the look-up table technology for snow coverage calculation, the calculation efficiency is high, and the remote sensing real-time monitoring of snow coverage can be realized according to this method.
在本发明的描述中,需要理解的是,术语“第一”、“第二”仅用于描述目的,而不能理解为指示或暗示相对重要性或者隐含指明所指示的技术特征的数量。由此,限定有“第一”、“第二”的特征可以明示或者隐含地包括一个或者更多个该特征。在本发明的描述中,“多个”的含义是两个或两个以上,除非另有明确具体的限定。In the description of the present invention, it should be understood that the terms "first" and "second" are used for description purposes only, and cannot be interpreted as indicating or implying relative importance or implicitly indicating the quantity of indicated technical features. Thus, a feature defined as "first" and "second" may explicitly or implicitly include one or more of these features. In the description of the present invention, "plurality" means two or more, unless otherwise specifically defined.
在本发明中,除非另有明确的规定和限定,术语“安装”、“相连”、“连接”、“固定”等术语应做广义理解,例如,可以是固定连接,也可以是可拆卸连接,或成一体;可以是机械连接,也可以是电连接;可以是直接相连,也可以通过中间媒介间接相连;可以是两个元件内部的连通或两个元件的相互作用关系。对于本领域的普通技术人员而言,可以根据具体情况理解上述术语在本发明中的具体含义。In the present invention, unless otherwise clearly specified and limited, terms such as "installation", "connection", "connection" and "fixation" should be understood in a broad sense, for example, it can be a fixed connection or a detachable connection , or integrated; it can be mechanically connected or electrically connected; it can be directly connected or indirectly connected through an intermediary; it can be the internal communication of two elements or the interaction relationship between two elements. Those of ordinary skill in the art can understand the specific meanings of the above terms in the present invention according to specific situations.
在本发明中,除非另有明确的规定和限定,第一特征在第二特征“上”或“下”,可以是第一和第二特征直接接触,或第一和第二特征通过中间媒介间接接触。而且,第一特征在第二特征“之上”、“上方”和“上面”,可以是第一特征在第二特征正上方或斜上方,或仅仅表示第一特征水平高度高于第二特征。第一特征在第二特征“之下”、“下方”和“下面”,可以是第一特征在第二特征正下方或斜下方,或仅仅表示第一特征水平高度低于第二特征。In the present invention, unless otherwise clearly specified and limited, the first feature is "on" or "under" the second feature, which may be that the first and second features are in direct contact, or the first and second features pass through an intermediary indirect contact. Moreover, the first feature being "above", "above" and "above" the second feature may mean that the first feature is directly above or obliquely above the second feature, or simply means that the level of the first feature is higher than that of the second feature . "Below", "below" and "under" the first feature may mean that the first feature is directly below or obliquely below the second feature, or it just means that the first feature is lower in level than the second feature.
在本说明书的描述中,术语“一个实施例”、“一些实施例”、“实施例”、“示例”、“具体示例”或“一些示例”等的描述,是指结合该实施例或示例描述的具体特征、结构、材料或者特点包含于本发明的至少一个实施例或示例中。在本说明书中,对上述术语的示意性表述不必须针对的是相同的实施例或示例。而且,描述的具体特征、结构、材料或者特点可以在任一个或多个实施例或示例中以合适的方式结合。此外,在不相互矛盾的情况下,本领域的技术人员可以将本说明书中描述的不同实施例或示例以及不同实施例或示例的特征进行结合和组合。In the description of this specification, the description of the terms "one embodiment", "some embodiments", "embodiment", "example", "specific example" or "some examples" refers to the A particular feature, structure, material, or characteristic is described as included in at least one embodiment or example of the invention. In this specification, the schematic representations of the above terms are not necessarily directed to the same embodiment or example. Furthermore, the described specific features, structures, materials or characteristics may be combined in any suitable manner in any one or more embodiments or examples. In addition, those skilled in the art can combine and combine different embodiments or examples and features of different embodiments or examples described in this specification without conflicting with each other.
尽管上面已经示出和描述了本发明的实施例,可以理解的是,上述实施例是示例性的,不能理解为对本发明的限制,本领域的普通技术人员在本发明的范围内可以对上述实施例进行改动、修改、替换和变型。Although the embodiments of the present invention have been shown and described above, it can be understood that the above embodiments are exemplary and should not be construed as limiting the present invention, those skilled in the art can make the above-mentioned The embodiments are subject to alterations, modifications, substitutions and variations.
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202011295856.8ACN112395989B (en) | 2020-11-18 | 2020-11-18 | A Mixed Pixel Decomposition Method for Snow Coverage for Multi-satellite Sensors |
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202011295856.8ACN112395989B (en) | 2020-11-18 | 2020-11-18 | A Mixed Pixel Decomposition Method for Snow Coverage for Multi-satellite Sensors |
| Publication Number | Publication Date |
|---|---|
| CN112395989A CN112395989A (en) | 2021-02-23 |
| CN112395989Btrue CN112395989B (en) | 2023-07-14 |
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202011295856.8AActiveCN112395989B (en) | 2020-11-18 | 2020-11-18 | A Mixed Pixel Decomposition Method for Snow Coverage for Multi-satellite Sensors |
| Country | Link |
|---|---|
| CN (1) | CN112395989B (en) |
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN114926732A (en)* | 2022-04-12 | 2022-08-19 | 北京艾尔思时代科技有限公司 | Multi-sensor fusion crop deep learning identification method and system |
| CN114757923B (en)* | 2022-04-19 | 2025-06-27 | 南京大学 | Method and device for automatic extraction of spectral endmembers |
| CN114780904B (en)* | 2022-06-17 | 2022-09-27 | 中国科学院、水利部成都山地灾害与环境研究所 | An end-member adaptive remote sensing inversion method for mountain vegetation coverage |
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101963664A (en)* | 2010-09-28 | 2011-02-02 | 中国科学院东北地理与农业生态研究所 | Microwave remote sensing pixel element decomposing method based on land and water living beings classifying information |
| CN103808736A (en)* | 2014-02-13 | 2014-05-21 | 吉林大学 | Saline-alkali soil characteristic detection method based on passive microwave mixed pixel decomposition technology |
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US5323317A (en)* | 1991-03-05 | 1994-06-21 | Hampton Terry L | Method and apparatus for determining runoff using remote geographic sensing |
| CN102608592B (en)* | 2012-04-05 | 2013-06-12 | 吉林大学 | Snow passive microwave mixed pixel decomposition method based on classified information of five types of ground features |
| CN102636779B (en)* | 2012-05-07 | 2013-08-21 | 武汉大学 | Extraction method for coverage rate of sub-pixel accumulated snow based on resampling regression analysis |
| CN104142142B (en)* | 2014-07-01 | 2016-08-24 | 北京师范大学 | Whole world vegetation fraction estimation method |
| CN105930817A (en)* | 2016-05-05 | 2016-09-07 | 中国科学院寒区旱区环境与工程研究所 | A road snow disaster monitoring and early warning method based on multi-source remote sensing data |
| CN106372434B (en)* | 2016-08-31 | 2020-03-06 | 中国科学院遥感与数字地球研究所 | A method and device for estimating instantaneous surface emissivity for passive microwave remote sensing |
| CN106951909A (en)* | 2016-11-16 | 2017-07-14 | 中国科学院遥感与数字地球研究所 | A Snow Cover Recognition Method for GF‑4 Satellite Remote Sensing Images |
| CN109376742A (en)* | 2018-09-19 | 2019-02-22 | 中国科学院东北地理与农业生态研究所 | A snow extraction method and system based on remote sensing images |
| CN110136194B (en)* | 2019-05-21 | 2022-11-11 | 吉林大学 | Snow coverage measuring and calculating method based on satellite-borne multispectral remote sensing data |
| CN111611965B (en)* | 2020-05-29 | 2020-12-22 | 中国水利水电科学研究院 | A method of land surface water extraction based on Sentinel-2 images |
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101963664A (en)* | 2010-09-28 | 2011-02-02 | 中国科学院东北地理与农业生态研究所 | Microwave remote sensing pixel element decomposing method based on land and water living beings classifying information |
| CN103808736A (en)* | 2014-02-13 | 2014-05-21 | 吉林大学 | Saline-alkali soil characteristic detection method based on passive microwave mixed pixel decomposition technology |
| Publication number | Publication date |
|---|---|
| CN112395989A (en) | 2021-02-23 |
| Publication | Publication Date | Title |
|---|---|---|
| CN112395989B (en) | A Mixed Pixel Decomposition Method for Snow Coverage for Multi-satellite Sensors | |
| CN107389036B (en) | A large-scale spatial-scale vegetation coverage calculation method combined with UAV images | |
| CN109993237B (en) | Method and system for rapid water extraction based on high-resolution satellite optical remote sensing data | |
| CN101963664B (en) | Microwave remote sensing pixel element decomposing method based on land and water living beings classifying information | |
| Liu et al. | Analysis of the urban heat island effect in Shijiazhuang, China using satellite and airborne data | |
| CN111046613B (en) | Optimal river channel calculation method based on path tracking and river network extraction method based on multi-temporal remote sensing image | |
| Liu et al. | Large-scale mapping of gully-affected areas: An approach integrating Google Earth images and terrain skeleton information | |
| CN111178169B (en) | Urban surface covering fine classification method and device based on remote sensing image | |
| CN103544477B (en) | Vegetation Coverage Estimation Method Based on Improved Linear Spectral Mixture Model | |
| Liu et al. | Comparative analysis of fractional vegetation cover estimation based on multi-sensor data in a semi-arid sandy area | |
| CN105930772A (en) | City impervious surface extraction method based on fusion of SAR image and optical remote sensing image | |
| CN109919250B (en) | Soil moisture-considered evapotranspiration space-time characteristic fusion method and device | |
| CN111007039B (en) | Automatic extraction method and system for sub-pixel level water body of medium-low resolution remote sensing image | |
| WO2025007773A1 (en) | Aerosol optical depth inversion method | |
| Bakar et al. | Spatial assessment of land surface temperature and land use/land cover in Langkawi Island | |
| CN115343226B (en) | A multi-scale vegetation coverage remote sensing calculation method based on UAV | |
| CN108764326A (en) | Urban impervious surface extracting method based on depth confidence network | |
| CN114813651A (en) | Remote sensing water quality inversion method combining difference learning rate and spectrum geometric characteristics | |
| Feng et al. | A hierarchical extraction method of impervious surface based on NDVI thresholding integrated with multispectral and high-resolution remote sensing imageries | |
| CN107688776A (en) | A kind of urban water-body extracting method | |
| Zhang et al. | GLC_FCS10: a global 10-m land-cover dataset with a fine classification system from Sentinel-1 and Sentinel-2 time-series data in Google Earth Engine | |
| CN119229279A (en) | Refined inversion method and system for forest and grass vegetation coverage based on multi-scale remote sensing | |
| Noshadi et al. | Prediction of surface soil color using ETM+ satellite images and artificial neural network approach | |
| He et al. | Linear spectral mixture analysis of Landsat TM data for monitoring invasive exotic plants in estuarine wetlands | |
| Zhang et al. | Mapping of circular or elliptical vegetation community patches: a comparative use of SPOT-5, ALOS and ZY-3 imagery |
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |