CN114092835A - 基于不同时空分辨率的归一化植被指数数据时空融合方法 - Google Patents

基于不同时空分辨率的归一化植被指数数据时空融合方法 Download PDF

Info

Publication number
CN114092835A
CN114092835A CN202210079071.XA CN202210079071A CN114092835A CN 114092835 A CN114092835 A CN 114092835A CN 202210079071 A CN202210079071 A CN 202210079071A CN 114092835 A CN114092835 A CN 114092835A
Authority
CN
China
Prior art keywords
ndvi
time
remote sensing
sensing image
pixel
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202210079071.XA
Other languages
English (en)
Other versions
CN114092835B (zh
Inventor
车向红
孙擎
刘纪平
王勇
徐胜华
罗安
杜凯旋
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Chinese Academy of Surveying and Mapping
Original Assignee
Chinese Academy of Surveying and Mapping
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Chinese Academy of Surveying and Mapping filed Critical Chinese Academy of Surveying and Mapping
Priority to CN202210079071.XA priority Critical patent/CN114092835B/zh
Publication of CN114092835A publication Critical patent/CN114092835A/zh
Application granted granted Critical
Publication of CN114092835B publication Critical patent/CN114092835B/zh
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/25—Fusion techniques

Landscapes

  • Engineering & Computer Science (AREA)
  • Data Mining & Analysis (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Artificial Intelligence (AREA)
  • Evolutionary Biology (AREA)
  • Evolutionary Computation (AREA)
  • Physics & Mathematics (AREA)
  • General Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Image Processing (AREA)

Abstract

一种基于不同时空分辨率的归一化植被指数数据时空融合方法，将不同时空分辨率的归一化植被指数数据进行时空融合，利用基于相邻NDVI观测的噪点二次过滤和高时间分辨率的NDVI中值的噪点过滤；结合线性插值和空间过滤，对不同时空分辨率的NDVI数据进行自动化融合生成高时间高空间分辨率NDVI数据；并利用空间过滤算法消除MODIS像元边界效应，最终得到目标区域的高空间分辨率高时间分辨率NDVI数据。本发明减少后续融合处理的不确定性，不用考虑输入数据是否晴空无云；通过结合与目标像元NDVI值的接近程度和与目标像元的欧式距离，对融合的目标像元NDVI值进行校正，提高融合影像的连续性。

Description

基于不同时空分辨率的归一化植被指数数据时空融合方法
技术领域
本申请涉及一种归一化植被指数数据的获取方法，具体的，涉及利用不同时空分辨率的归一化植被指数数进行时空融合后获得高时空分辨率的归一化植被指数数据的方法。
背景技术
植被指数是由卫星不同波段数据组合而成，能简单、有效地反映植物生长状况，与植被覆盖度之间存在线性关系，是标识植被生长状况和覆盖度的重要指标。当前已经发展了数十种植被指数，其中，归一化植被指数NDVI（Normalized Vegetation Index）可以消除与仪器定标、太阳角、地形、云、阴影和大气条件有关的大多数影响，并且增加了对植被的响应能力，在植被监测的遥感应用最为广泛。目前常用的NDVI 数据有EOS/MODIS NDVI数据，该数据时间分辨率高（逐日），能反映植被的动态变化，然而空间分辨率较粗（250m-1000m），使得其获取的像元信息存在大量的混合像元现象，这对于植被类型复杂的中小尺度区域，无法准确反映植被生长状况和覆盖度。Landsat系列卫星数据能够提供高空间分辨率（30m）的NDVI数据，然而其8-16天的重访周期很难监测地表植被的动态变化。为了得到高空间分辨率和高时间分辨率遥感数据，已提出了多种多源遥感数据融合方法，常用的主要分为两类：基于线性光谱混合模型和自适应模型。
（1）基于线性光谱混合模型融合方法：根据线性光谱混合模型，低空间分辨率遥感数据的像元值可以当成中高空间分辨率遥感数据像元值的线性组合。假设同一类别中高空间分辨率遥感数据像元值一样，利用从中高空间分辨率遥感数据获得的各类别丰度矩阵，采用最小二乘法，可以从低空间分辨率像元反射率推算出中高空间分辨率像元值。然而上述基于线性光谱混合模型的方法，解算的中高空间分辨率像元值都只是类别平均值，而不是真正的像元值。计算丰度矩阵时同一类型反射率一样的假设对于空间异质性较高的地表植被具有局限性。
（2）自适应融合方法：为避免解算线性光谱混合模型，同时考虑像元反射率的空间可变性，Gao等提出了一种自适应遥感图像时空融合方法（Spatial and TemporalAdaptive Reflectance Fusion Model，STARFM），忽略了不同传感器光谱响应函数的差异，假设地表覆盖类型和系统误差不随时间发生变化，将MODIS光谱融合到Landsat影像中，STARFM方法对空间异质性区域模拟效果有待提高。为此，在STARFM方法基础上，一系列的改进方法相继提取，如考虑了土地覆盖变化的STAARCH（Spatial Temporal AdaptiveAlgorithm for mapping Reflectance Change）、需较高数据质量Landsat影像的ESTARFM（Enhanced Spatial and Temporal Adaptive Reflectance Fusion Model）、FSDAF（Flexible Spatiotemporal DAta Fusion）等。上述数据融合方法有效地提高了影像的时空分辨率，但也存在局限，对数据融合的输入要求较高：需要影像质量较好且离目标时刻较近的MODIS-Landsat影像对，但由于条带坏数据或云/云阴影等无效值影响，这些影像很难用于融合，难以实现自动且高效地融合来进行大区域长时间序列的地表植被监测。
综上，现有的遥感数据融合方法一方面在地表植被复杂的区域融合效果不理想。另一方面，融合方法对于输入的数据要求较为严格，均是针对晴空无云的遥感数据，对于云覆盖区域仍存在瓶颈。
因此，如何提供一种时空融合方法，以获取高时间分辨率和高空间分辨率的NDVI数据，并以解决云/云阴影污染区域难以融合等问题，并减少对于输入NDVI数据源的质量好坏的依赖性，成为现有技术亟需解决的技术。
发明内容
本发明的目的在于提出一种基于不同时空分辨率的归一化植被指数数据时空融合方法，用于缓解现有技术无法对地表植被NDVI进行长时间序列精细化监测的技术问题。
为达此目的，本发明采用以下技术方案：
一种基于不同时空分辨率的归一化植被指数数据时空融合方法，包括如下步骤：
不同时空分辨率的归一化植被指数数据获取步骤S110：
分别获取目标区域的高空间分辨率低时间分辨率的第一遥感图像和低空间分辨率高时间分辨率第二遥感图像，通过地表反射率数据分别计算第一和第二遥感图像的归一化植被指数数据NDVI，生成第一遥感图像和第二遥感图像的归一化植被指数数据NDVI时间序列；
对归一化植被指数数据进行优化处理步骤S120：
分别对第一遥感图像和第二遥感图像的归一化植被指数数据NDVI时间序列进行噪点去除，得到稳定的第一遥感图像和第二遥感图像的归一化植被指数数据NDVI时间序列；
初步的高时空分辨率归一化植被指数数据生成步骤S130：
基于第一和第二遥感图像在时间上的线性对应，根据第一遥感图像的低时间分辨率性，得到低时间分辨率的每个像元的NDVI差值时间序列，对上述NDVI差值时间序列进行线性差值得到高时间分辨率的NDVI差值时间序列，而后与每个像元在高时间分辨率点上的第二遥感图像的归一化植被指数数据NDVI时间序列求和得到初步的高空间分辨率高时间分辨率NDVI时间序列数据；
相似像元集计算获取步骤S140：
选择离融合时间点
Figure 655267DEST_PATH_IMAGE001
最近的第一遥感图像的NDVI数值，通过计算目标窗口内每个像元与目标像元的归一化植被指数插值，来寻找目标像元的相似像元集；
最终的高时空分辨率的归一化植被指数数据计算步骤S150：
计算相似像元与目标像元的欧式距离，并将所述欧式距离转换为相似像元权重值，将相似像元在初步的高空间分辨率高时间分辨率归一化植被指数上对应的值与所述权重值用加权求和得到目标像元最终的NDVI值，逐像元遍历研究区为每一个插值NDVI像元进行空间过滤，生成最终的高空间分辨率高时间分辨率NDVI数值。
可选的，对于目标区域的高空间分辨率低时间分辨率的第一遥感图像，计算并生成第一遥感图像的归一化植被指数数据NDVI时间序列具体为：
对目标区域的高空间分辨率低时间分辨率的第一遥感图像进行大气校正，生成地表反射率数据；同时采用FMASK云检测算法，对所述第一遥感图像数据中的云/云阴影、水体、冰雪的像元进行检测过滤；利用公式1计算出每个过滤后像元的归一化植被指数数据NDVI，形成第一遥感图像的归一化植被指数数据NDVI时间序列，
Figure 90928DEST_PATH_IMAGE002
（1）
其中，
Figure 252787DEST_PATH_IMAGE003
为近红外波段反射率，
Figure 672267DEST_PATH_IMAGE004
为红光波段反射率；
生成第二遥感图像的归一化植被指数数据NDVI时间序列具体为：获取目标区域的低空间分辨率高时间分辨率第二遥感图像，当第二遥感图像为地表反射率数据时，直接采用公式（1），计算得到第二遥感图像的归一化植被指数数据NDVI时间序列。
可选的，在步骤S120中：
对于第一遥感图像的归一化植被指数数据NDVI时间序列，将冬季中NDVI<0.1的序列点进行过滤，然后将同时满足公式2和公式3的像元判定为污染像元，进行过滤，得到优化后的第一遥感图像的归一化植被指数数据NDVI时间序列，
Figure 777277DEST_PATH_IMAGE005
(2)
Figure 333023DEST_PATH_IMAGE006
(3)
其中，NDVI i+j 表示以NDVI i 为中心的四个相邻时间点 NDVI 值，j = -2、-1、1 和2，AVERAGE是计算相邻四个NDVI值的平均值，STD是对应的标准差。
可选的，在步骤S120中，
对于第二遥感图像的归一化植被指数数据NDVI时间序列，过滤冬季NDVI<0.1的时间点，利用中值法进行过滤，以每日目标时间点为中心，将9天的窗口内的第二遥感图像的归一化植被指数数据NDVI数据求中值，减少时间序列噪点，将过滤后的第二遥感图像NDVI转投影并重采样为第一遥感图像的投影和空间分辨率，生成优化的第二遥感图像的归一化植被指数数据NDVI时间序列。
可选的，在步骤S130中，
寻找一年中k个时间对应的MODIS和Landsat NDVI影像对，影像对相应位置的关系如公式4所示，利用公式5将影像对对应的NDVI逐像元求差，生成低时间分辨率的每个像元NDVI差值时间序列，含k个值，
Figure 744413DEST_PATH_IMAGE007
（4）
其中，(x, y)为第一遥感图像的 (NDVIL)和第二遥感图像(NDVIM)相对应的像元位置，T j (j=1,2,…k)为第一遥感图像和第二遥感图像对应的影像获取时间；
Figure 634877DEST_PATH_IMAGE008
是(x, y)位置的已知时间点的低时间分辨率的第一遥感图像与第二遥感图像的NDVI差值，
Figure 258756DEST_PATH_IMAGE009
（5）
将已知时间点的低时间分辨率的NDVI差值时间序列
Figure 200168DEST_PATH_IMAGE008
进行线性插值，获取高时间分辨率的按日的NDVI差值时间序列
Figure 2033DEST_PATH_IMAGE008
，每个像元上任意一个时间点上的高空间分辨率高时间分辨率归一化植被指数
Figure 130526DEST_PATH_IMAGE010
则为该时间点的第二遥感图像归一化植被指数NDVI值与该时间点的高时间分辨率NDVI差值
Figure 558096DEST_PATH_IMAGE008
的和，利用公式6得到初步的高空间分辨率高时间分辨率归一化植被指数（NDVI）数据
Figure 603281DEST_PATH_IMAGE010
，
Figure 356474DEST_PATH_IMAGE011
(6)
其中，
Figure 972263DEST_PATH_IMAGE012
为按日的任意时间点
Figure 203524DEST_PATH_IMAGE001
的第二遥感图像NDVI值，
Figure 805013DEST_PATH_IMAGE008
为插值形成的按日的任意时间点
Figure 729107DEST_PATH_IMAGE001
的高时间分辨率的NDVI差值时间序列。
可选的，在步骤S140中，
选择离融合时间点
Figure 347039DEST_PATH_IMAGE001
最近的第一遥感图像的NDVI数值，对应时间点为
Figure 850832DEST_PATH_IMAGE013
，假设
Figure 90184DEST_PATH_IMAGE001
到
Figure 670332DEST_PATH_IMAGE013
时间段地表植被变化较小，使用最近时间
Figure 260713DEST_PATH_IMAGE013
的第一遥感图像的NDVI来寻找目标像元的相似像元集，每一个目标像元对应的相似像元选择依据：离目标像元距离越近，和目标像元的NDVI值越接近，则与目标像元越接近；首先在以目标像元为中心，利用公式7在m×m的目标窗口内计算每个像元与目标像元的NDVI差值，将NDVI差值从小到大排序，选择前t个差值对应的像元作为目标像元的相似像元集，t为目标像元对应的相似像元的个数，如果在最近时间
Figure 833777DEST_PATH_IMAGE013
的第一遥感图像 NDVI数据上找不到相似像元，就依次循环寻找下一个离融合时间点
Figure 442482DEST_PATH_IMAGE001
最近的时间点，在其第一遥感图像 NDVI数据上按上述步骤查找相似像元集，直到找到相似像元集终止循环，
Figure 708378DEST_PATH_IMAGE014
(7)
其中，n为目标窗口内像元位置标识，n =1, 2, …m×m，
Figure 786056DEST_PATH_IMAGE015
目标像元
Figure 928191DEST_PATH_IMAGE016
在
Figure 142135DEST_PATH_IMAGE013
时刻的第一遥感图像 NDVI值，
Figure 578932DEST_PATH_IMAGE017
为目标窗口内的每个像元
Figure 393174DEST_PATH_IMAGE018
在
Figure 308040DEST_PATH_IMAGE013
时刻的第一遥感图像的 NDVI,
Figure 642069DEST_PATH_IMAGE019
为目标窗口内计算每个像元与目标像元的NDVI差值。
可选的，在步骤S150中，
采用公式8计算相似像元到目标像元的欧氏距离，然后采用公式9将欧式距离转为相似像元权重值，采用公式10将相似像元在初步的高空间分辨率高时间分辨率归一化植被指数
Figure 249768DEST_PATH_IMAGE020
上对应的值采与相似像元权重值加权求和得到目标像元最终的NDVI值，逐像元遍历研究区为每一个插值NDVI像元进行空间过滤，生成最终的高空间分辨率高时间分辨率NDVI数值，
Figure 787191DEST_PATH_IMAGE021
(8)
Figure 505748DEST_PATH_IMAGE022
(9)
Figure 959863DEST_PATH_IMAGE023
(10)
其中，ω取值为融合数据的空间分辨率，例如ω=30，(x, y) 为目标像元的位置坐标，
Figure 987731DEST_PATH_IMAGE024
为相似像元的位置坐标，
Figure 527297DEST_PATH_IMAGE025
为过滤后融合NDVI数值。
可选的，第一遥感图像为Landsat序列遥感图像，第二遥感图像为逐日的MODIS序列遥感图像。
本发明进一步公开了一种存储介质，用于存储计算机可执行指令，其特征在于：
所述计算机可执行指令在被处理器执行时执行上述的基于不同时空分辨率的归一化植被指数数据时空融合方法。
本发明具有如下的优点：1、在用FMASK对高空间分辨率低时间分辨率LandsatNDVI数据过滤后，设计了基于相邻NDVI值的噪点二次过滤和MODIS NDVI 中值的噪点过滤，使用自动化过滤方法，保证输入数据的准确性，减少后续融合处理的不确定性。
2、本发明将线性插值与空间过滤结合，不用考虑诸如Landsat第一遥感图像和诸如MODIS的第二遥感图像的输入数据是否晴空无云，对第一遥感图像和第二遥感图像进行自动化融合生成高时间高空间分辨率NDVI数据。
3、对于空间过滤算法，通过结合与目标像元NDVI值的接近程度和与目标像元的欧式距离，获取每一个融合像元的相似像元，对融合的目标像元NDVI值进行校正，减弱因融合粗空间分辨率影像引入的像元边界效应，提高融合影像的连续性。
附图说明
图1是根据本发明具体实施例的基于不同时空分辨率的归一化植被指数数据时空融合方法的流程图；
图2是根据本发明具体实施例的归一化植被指数数据时空融合方法的逻辑流程图；
图3是根据本发明具体实施例的归一化植被指数数据时空融合方法的结果图；
图4是根据本发明具体实施例的归一化植被指数数据时空融合方法的结果图中子区域1的详细展示图；
图5是根据本发明具体实施例的归一化植被指数数据时空融合方法的结果图中子区域2的详细展示图。
具体实施方式
下面结合附图和实施例对本发明作进一步的详细说明。可以理解的是，此处所描述的具体实施例仅仅用于解释本发明，而非对本发明的限定。另外还需要说明的是，为了便于描述，附图中仅示出了与本发明相关的部分而非全部结构。
本发明主要在于：将不同时空分辨率的归一化植被指数数据进行时空融合，利用设计了基于相邻NDVI观测的噪点二次过滤和高时间分辨率的NDVI中值的噪点过滤；结合线性插值和空间过滤，对不同时空分辨率的NDVI数据进行自动化融合生成高时间高空间分辨率NDVI数据；并利用空间过滤算法消除MODIS像元边界效应，最终得到目标区域的高空间分辨率高时间分辨率NDVI数据。
具体的，参见图1，示出了根据本发明具体实施例的基于不同时空分辨率的归一化植被指数数据时空融合方法的流程图，图2示出了归一化植被指数数据时空融合方法的逻辑流程图，该时空融合方法包括如下步骤：
不同时空分辨率的归一化植被指数数据获取步骤S110：
分别获取目标区域的高空间分辨率低时间分辨率的第一遥感图像和低空间分辨率高时间分辨率第二遥感图像，通过地表反射率数据分别计算第一和第二遥感图像的归一化植被指数数据NDVI，生成第一遥感图像和第二遥感图像的归一化植被指数数据NDVI时间序列。
在本发明中，对于高空间分辨率低时间分辨率的第一遥感图像，以Landsat序列遥感图像为例，对于低空间分辨率高时间分辨率第二遥感图像，以逐日的MODIS序列遥感图像为例。但本发明不限于此，也可以是具有上述不同的时间和空间分辨率的序列的其它遥感图像，均可以利用本发明的方法进行计算。
具体的，在步骤S110中，对于目标区域的高空间分辨率低时间分辨率的第一遥感图像，计算并生成第一遥感图像的归一化植被指数数据NDVI时间序列具体为：
目标区域的高空间分辨率低时间分辨率的第一遥感图像，利用LEDAPS（LandsatEcosystem Disturbance Adaptive Processing System）或LaSRC（Landsat 8 SurfaceReflectance Code）算法对第一遥感图像数据，例如Landsat数据，进行大气校正，生成地表反射率数据；同时采用FMASK云检测算法，对所述第一遥感图像数据中的云/云阴影、水体、冰雪的像元进行检测过滤；利用公式1计算出每个过滤后像元的归一化植被指数数据NDVI，形成第一遥感图像的归一化植被指数数据NDVI时间序列。
Figure 783966DEST_PATH_IMAGE026
（1）
其中，
Figure 92587DEST_PATH_IMAGE027
为近红外波段反射率，
Figure 789892DEST_PATH_IMAGE004
为红光波段反射率。
生成第二遥感图像的归一化植被指数数据NDVI时间序列具体为：
获取目标区域的低空间分辨率高时间分辨率第二遥感图像，例如逐日的MODIS序列遥感图像，该第二遥感图像为地表反射率数据，直接采用公式（1），计算得到第二遥感图像的归一化植被指数数据NDVI时间序列。
对归一化植被指数数据进行优化处理步骤S120：
分别对第一遥感图像和第二遥感图像的归一化植被指数数据NDVI时间序列进行噪点去除，得到稳定的第一遥感图像和第二遥感图像的归一化植被指数数据NDVI时间序列。
具体的，对于第一遥感图像的归一化植被指数数据NDVI时间序列，将冬季中NDVI<0.1的序列点进行过滤，然后将同时满足公式2和公式3的像元判定为污染像元，进行过滤，得到优化后的第一遥感图像的归一化植被指数数据NDVI时间序列，
Figure 285595DEST_PATH_IMAGE005
(2)
Figure 611534DEST_PATH_IMAGE006
(3)
其中，NDVI i+j 表示以NDVI i 为中心的四个相邻时间点 NDVI 值，j = -2、-1、1 和2，AVERAGE是计算相邻四个NDVI值的平均值，STD是对应的标准差。
对于第二遥感图像的归一化植被指数数据NDVI时间序列，过滤冬季NDVI<0.1的时间点，因为第二遥感图像的空间分辨率高，其NDVI时间序列，例如MODIS NDVI为逐日的高时间分辨率数据，时间序列密度大，只以简单的中值法进行过滤，以每日目标时间点为中心，将9天的窗口内的MODIS NDVI数据求中值，减少时间序列噪点，为了与Landsat NDVI在像元尺度进行融合，将过滤后的第二遥感图像MODISNDVI转投影并重采样为第一遥感图像Landsat数据的投影和空间分辨率，生成优化的第二遥感图像的归一化植被指数数据NDVI时间序列，例如MODIS NDVI时间序列。
初步的高时空分辨率归一化植被指数数据生成步骤S130：
设计时空融合算法将诸如 Landsat 和 MODIS NDVI 数据的第一遥感图像和第二遥感图像的归一化植被指数数据NDVI时间序列融合为高空间分辨率高时间分辨率NDVI时间序列，融合算法包括时间序列线性插值和空间过滤两步，本步骤主要涉及时间序列线性插值。
基于第一和第二遥感图像在时间上的线性对应，根据第一遥感图像的低时间分辨率性，得到低时间分辨率的每个像元的NDVI差值时间序列，对上述NDVI差值时间序列进行线性差值得到高时间分辨率的NDVI差值时间序列，而后与每个像元在高时间分辨率点上的第二遥感图像的归一化植被指数数据NDVI时间序列求和得到初步的高空间分辨率高时间分辨率NDVI时间序列数据。
具体的，寻找一年中k个时间对应的MODIS和Landsat NDVI影像对，影像对相应位置的关系如公式4所示，利用公式5将影像对对应的NDVI逐像元求差，生成低时间分辨率的每个像元NDVI差值时间序列，含k个值，此时的k值应当与第一遥感图像，即Landsat的时间分辨率有关。
Figure 23930DEST_PATH_IMAGE007
（4）
其中，(x, y)为第一遥感图像Landsat的 (NDVIL)和第二遥感图像MODIS (NDVIM)相对应的像元位置，T j (j=1,2,…k)为第一遥感图像Landsat和第二遥感图像MODIS对应的影像获取时间；
Figure 144333DEST_PATH_IMAGE008
是(x, y)位置的已知时间点的低时间分辨率的第一遥感图像与第二遥感图像的NDVI差值，主要是由于观测误差、太阳与观测角变化、噪点或传感器系统误差引起的。
Figure 127332DEST_PATH_IMAGE028
（5）
将已知时间点的低时间分辨率的NDVI差值时间序列
Figure 742115DEST_PATH_IMAGE008
进行线性插值，获取高时间分辨率的NDVI差值时间序列
Figure 759750DEST_PATH_IMAGE008
，例如获取逐日的NDVI差值时间序列
Figure 316633DEST_PATH_IMAGE008
，每个像元上任意一个时间点上的高空间分辨率高时间分辨率归一化植被指数（NDVI）
Figure 786929DEST_PATH_IMAGE010
则为该时间点的第二遥感图像归一化植被指数NDVI值与该时间点的高时间分辨率NDVI差值
Figure 438359DEST_PATH_IMAGE008
的和，具体如公式6，形成初步的高空间分辨率高时间分辨率归一化植被指数（NDVI）数据
Figure 576079DEST_PATH_IMAGE010
，
Figure 303864DEST_PATH_IMAGE029
(6)
其中，
Figure 261455DEST_PATH_IMAGE030
为按日的任意时间点
Figure 220971DEST_PATH_IMAGE031
的第二遥感图像NDVI值，
Figure 213198DEST_PATH_IMAGE008
为插值形成的按日的任意时间点
Figure 846304DEST_PATH_IMAGE031
的高时间分辨率的NDVI差值时间序列。
接下来，将基于高空间分辨率高时间分辨率NDVI数据
Figure 806039DEST_PATH_IMAGE010
，设计空间过滤算法。第一遥感图像和第二遥感图像的 NDVI 图像空间分辨率的不一致意味着插值预测主要是与 MODIS 像素范围对应，这将形成块状效应，示例性的，图2中的（e）部分。为了处理插值融合数据
Figure 815583DEST_PATH_IMAGE010
中的块状效应，设计空间过滤算法，寻找需空间过滤的目标像元的相似像元集来消除
Figure 396737DEST_PATH_IMAGE010
块状效应。
因此，本发明继续包括：
相似像元集计算获取步骤S140：
在本步骤中，融合时间点
Figure 466324DEST_PATH_IMAGE031
生成的高空间分辨率高时间分辨率
Figure 149241DEST_PATH_IMAGE010
影像由于融合了第二遥感图像MODIS的 NDVI，已经存在块状效应，因此选择为目标像元选择相似像元集不可基于
Figure 962476DEST_PATH_IMAGE010
的初步融合结果选择，而应当选择离融合时间点
Figure 929295DEST_PATH_IMAGE031
最近的第一遥感图像Landsat的 NDVI数值，
因此，本步骤选择离融合时间点
Figure 887893DEST_PATH_IMAGE031
最近的第一遥感图像的NDVI数值，通过计算目标窗口内每个像元与目标像元的归一化植被指数插值，来寻找目标像元的相似像元集。
具体为：选择离融合时间点
Figure 572952DEST_PATH_IMAGE031
最近的第一遥感图像Landsat的NDVI图像，对应时间点为
Figure 924299DEST_PATH_IMAGE013
，假设
Figure 745624DEST_PATH_IMAGE031
到
Figure 373658DEST_PATH_IMAGE013
时间段地表植被变化较小，使用最近时间
Figure 280434DEST_PATH_IMAGE013
的第一遥感图像LandsatNDVI来寻找目标像元的相似像元集，每一个目标像元对应的相似像元选择依据：离目标像元距离越近，和目标像元NDVI值越接近，则与目标像元越接近；首先在以目标像元为中心，利用公式7在m×m的目标窗口内计算每个像元与目标像元的NDVI差值，将NDVI差值从小到大排序，选择前t个差值对应的像元作为目标像元的相似像元集，t为目标像元对应的相似像元的个数，如果在最近时间
Figure 169893DEST_PATH_IMAGE013
的第一遥感图像 Landsat NDVI数据上找不到相似像元，就依次循环寻找下一个离融合时间点
Figure 94993DEST_PATH_IMAGE031
最近的时间点，在其第一遥感图像 Landsat NDVI数据上按上述步骤查找相似像元集，直到找到相似像元集终止循环。
Figure 677284DEST_PATH_IMAGE032
(7)
其中，n为目标窗口内像元位置标识，n =1, 2, …m×m，
Figure 71356DEST_PATH_IMAGE033
为目标像元(x, y)在
Figure 498926DEST_PATH_IMAGE013
时刻的第一遥感图像 NDVI值，
Figure 45576DEST_PATH_IMAGE034
为目标窗口内的每个像元
Figure 533189DEST_PATH_IMAGE035
在
Figure 148979DEST_PATH_IMAGE013
时刻的第一遥感图像 NDVI,
Figure 895087DEST_PATH_IMAGE036
为目标窗口内计算每个像元与目标像元的NDVI差值。
最终的高时空分辨率的归一化植被指数数据计算步骤S150：
计算相似像元与目标像元的欧式距离，并将所述欧式距离转换为相似像元权重值，将相似像元在初步的高空间分辨率高时间分辨率归一化植被指数
Figure 545511DEST_PATH_IMAGE037
上对应的值与所述权重值用加权求和得到目标像元最终的NDVI值，逐像元遍历研究区为每一个插值NDVI像元进行空间过滤，生成最终的高空间分辨率高时间分辨率NDVI数值。
采用公式8计算相似像元到目标像元的欧氏距离，然后采用公式9将欧式距离转为相似像元权重值，采用公式10将相似像元在初步的高空间分辨率高时间分辨率归一化植被指数
Figure 204025DEST_PATH_IMAGE037
上对应的值采与权重值加权求和得到目标像元最终的NDVI值，逐像元遍历研究区为每一个插值NDVI像元进行空间过滤，生成最终的高空间分辨率高时间分辨率NDVI数值，
Figure 572690DEST_PATH_IMAGE021
(8)
Figure 349585DEST_PATH_IMAGE022
(9)
Figure 854516DEST_PATH_IMAGE023
(10)
其中，ω取值为融合数据的空间分辨率，例如ω=30，(x, y) 为目标像元的位置坐标，
Figure 683931DEST_PATH_IMAGE024
为相似像元的位置坐标，
Figure 523580DEST_PATH_IMAGE038
为过滤后融合NDVI数值。
实施例：
本实施时空融合方法具体包括如下的步骤，其中附图3示出了归一化植被指数数据时空融合方法的结果图，图3中的a部分为，b-c部分为获取目标区域的高空间分辨率低时间分辨率Landsat序列遥感图像，对该遥感图像进行大气校正、云检测预处理，分别计算Landsat影像的NDVI，生成植被NDVI时间序列，图3中的d部分为获取目标区域的低空间分辨率高时间分辨率MODIS序列遥感图像，计算得到MODIS影像的NDVI，生成植被NDVI时间序列。图3中的e部分为对Landsat和MODIS NDVI时间序列进行优化，生成稳定的Landsat和MODISNDVI时间序列；基于Landsat与MODIS NDVI差异在时间序列上是线性的，通过线性插值获取高时间分辨率的NDVI差值，与MODIS NDVI求和，生成初步的高空间分辨率高时间分辨率NDVI数据。图3中的f部分为，通过空间算法消除MODIS像元边界效应，最终得到的目标区域的高空间分辨率高时间分辨率NDVI数据图。
图4，图5的顺序与图2部分的顺序相同，只是对应图3中的区域1和区域2部分的详细放大展示图。
步骤1：筛选2018年1月-2018年12月所有Landsat-8大气顶层反射率数据，该数据空间分辨率为30m，时间分辨率为16天回访周期。该研究区处于p034r032图幅上。
步骤2：基于所选的Landsat-8大气顶层反射率数据，采用LaSRC开源包将其处理为地表反射率数据，然后基于FMASK开源包，对所述Landsat-8时间序列数据进行云/云阴影提取，将受云、云阴影、水体、冰雪污染的像元过滤。Landsat-8云检测结果信息具体如下：
1 = 陆地；
2 = 水体；
3 = 云/云阴影；
4 = 冰雪；
步骤3：筛选2018年1月-2018年12月所有MODIS NBAR数据（MCD43A4），该数据为经过BRDF（Bidirectional Reflectance Distribution Function)）校正的地表反射率数据，时间分辨率为逐天，所在图幅为H09V04。
步骤4：基于过滤的Landsat-8地表反射率数据和MCD43A4数据，利用公式1计算Landsat和MODIS NDVI植被指数序列，计算公式如下：
Figure 362223DEST_PATH_IMAGE039
（1）
其中，NDVI为归一化植被指数，
Figure 987240DEST_PATH_IMAGE040
为近红外波段反射率，
Figure 253136DEST_PATH_IMAGE041
为红光波段反射率。
步骤5：尽管FMASK算法结果可以过滤Landsat云/云阴影、水体和冰雪等污染像元，但仍然存在未检测到污染像元。因为污染像元NDVI值大多比较小，鉴于此，首先将冬季（1-3月和11-12月）NDVI<0.1的序列点进行过滤，然后计算每个目标时间点周围四个时间点NDVI值的平均值和标准差，如果目标观测满足相邻4个观测值的最小值（即满足公式2），且小于上述均值减去标准差（即满足公式3），则该目标时间点对应的NDVI为为污染像元，对其进行过滤，生成最终的Landsat NDVI时间序列。
Figure 815967DEST_PATH_IMAGE005
（2）
Figure 192721DEST_PATH_IMAGE006
（3）
其中，NDVI i+j 表示以NDVI i 为中心的四个相邻时间点 NDVI 值（j = -2、-1、1 和2）。 AVERAGE是计算相邻四个NDVI值的平均值，STD是对应的标准差。
步骤6：将MODIS NDVI冬季（1-3月和11-12月）NDVI<0.1的序列点进行过滤，由于逐日的高时间分辨率，MODIS NDVI以简单的中值法进行过滤，以每日目标时间点为中心，将9天的窗口内的MODIS NDVI数据求中值，减少时间序列噪点。为了与Landsat NDVI在像元尺度进行融合，将过滤后的MODIS影像转投影并重采样为Landsat数据的投影和空间分辨率，生成最终的MODIS NDVI时间序列。
步骤7：基于步骤5和步骤6生成的Landsat与MODIS NDVI时间序列，寻找一年中k个时间对应的MODIS和Landsat NDVI影像对，影像对相应位置的关系如公式4所示，利用公式5将影像对对应的NDVI逐像元求差，生成每个像元NDVI差值时间序列，含k个值。
Figure 672244DEST_PATH_IMAGE042
（4）
其中，(x, y)为Landsat ((NDVIL))和MODIS (NDVIM)相对应的像元位置，T j (j=1,2,…k)为低时间分辨率下的Landsat和MODIS影像获取时间；
Figure 109042DEST_PATH_IMAGE043
是(x, y)位置Landsat与MODIS NDVI差值，主要是由于观测误差、太阳与观测角变化、噪点或传感器系统误差引起的。
Figure 923283DEST_PATH_IMAGE044
（5）
步骤8：将NDVI差值时间序列
Figure 838149DEST_PATH_IMAGE045
进行线性插值，获取高时间分辨率的逐日的NDVI差值时间序列
Figure 437758DEST_PATH_IMAGE045
，i=1,2,…365。
步骤9；每个像元上高时间分辨率的逐日的任意一个时间点
Figure 793260DEST_PATH_IMAGE031
上的高空间分辨率高时间分辨率NDVI
Figure 579950DEST_PATH_IMAGE046
则为该时间点的MODIS NDVI值与该时间点的NDVI差值求和(具体公式6)，形成初步的高空间分辨率高时间分辨率NDVI数据
Figure 298507DEST_PATH_IMAGE046
。
Figure 752622DEST_PATH_IMAGE047
(6)
其中，
Figure 780490DEST_PATH_IMAGE048
为高时间分辨率的逐日的任意时间点
Figure 54477DEST_PATH_IMAGE031
的MODIS NDVI值，
Figure 311146DEST_PATH_IMAGE049
为插值形成的任意时间点
Figure 370500DEST_PATH_IMAGE031
的Landsat与MODIS NDVI差值。
步骤10：选择离融合时间点
Figure 320001DEST_PATH_IMAGE031
最近的Landsat NDVI影像，对应时间点为
Figure 81284DEST_PATH_IMAGE013
，假设
Figure 141644DEST_PATH_IMAGE031
到
Figure 819618DEST_PATH_IMAGE013
时间段地表植被变化较小，基于最近时间
Figure 674442DEST_PATH_IMAGE013
的 Landsat NDVI来寻找目标像元的相似像元集。每一个目标像元对应的相似像元选择依据：离目标像元距离越近，和目标像元NDVI值越接近，则与目标像元越接近。
步骤11: 以目标像元为中心，利用公式7在m×m的目标窗口内计算每个像元与目标像元的NDVI差值，将NDVI差值从小到大排序，选择前t个差值对应的像元作为目标像元的相似像元集，本实施中设定m设定为5，t设定为10。如果在最近时间
Figure 923021DEST_PATH_IMAGE013
的 Landsat NDVI数据上找不到相似像元，就返回步骤10继续寻找下一个离融合时间点
Figure 493065DEST_PATH_IMAGE031
最近的时间点，在其Landsat NDVI影像上按上述步骤查找相似像元集，直到找到相似像元集进入下一个步骤。
Figure 776279DEST_PATH_IMAGE050
(7)
其中，本实施中设定m=5，n为目标窗口内像元位置标识，n =1, 2, …m×m，
Figure 425173DEST_PATH_IMAGE051
为目标像元(x, y)在
Figure 301993DEST_PATH_IMAGE013
时刻的Landsat NDVI值，
Figure 704155DEST_PATH_IMAGE052
为目标窗口内的每个像元
Figure 828494DEST_PATH_IMAGE018
在
Figure 25120DEST_PATH_IMAGE013
时刻的Landsat NDVI,
Figure 982712DEST_PATH_IMAGE019
为目标窗口内计算每个像元与目标像元的NDVI差值
步骤12：基于公式8计算相似像元到目标像元的欧氏距离，然后采用公式9将欧式距离转为相似像元的权重值。
Figure 703412DEST_PATH_IMAGE021
(8)
Figure 430059DEST_PATH_IMAGE022
(9)
其中，距离Dj需要限制在一个合适的范围内，才能对权重Wj产生合理的影响。因此，公式8中使用了1和
Figure 63166DEST_PATH_IMAGE053
，并且Dj的范围相应地从1到
Figure 258786DEST_PATH_IMAGE054
。算法中ω取值为融合数据的空间分辨率，即ω=30。t为目标像元对应的相似像元的个数，(x, y) 为目标像元的位置坐标，
Figure 268331DEST_PATH_IMAGE055
为相似像元的位置坐标。
步骤13：利用公式10把相似像元在
Figure 115064DEST_PATH_IMAGE056
影像上对应的值采用加权求和即为目标像元最终的NDVI值。
Figure 919072DEST_PATH_IMAGE057
(10)
步骤14：逐像元遍历研究区，重复步骤5-步骤13，为每一NDVI像元进行时间线性插值和空间过滤，生成最终的高空间分辨率高时间分辨率NDVI影像，提取结果如图2的f部分所示。
本发明进一步公开了一种存储介质，用于存储计算机可执行指令，所述计算机可执行指令在被处理器执行时执行上述的基于不同时空分辨率的归一化植被指数数据时空融合方法。
本发明具有如下的优点：1、在用FMASK对高空间分辨率低时间分辨率LandsatNDVI数据过滤后，设计了基于相邻NDVI值的噪点二次过滤和MODIS NDVI 中值的噪点过滤，使用自动化过滤方法，保证输入数据的准确性，减少后续融合处理的不确定性。
2、本发明将线性插值与空间过滤结合，不用考虑诸如Landsat第一遥感图像和诸如MODIS的第二遥感图像的输入数据是否晴空无云，对第一遥感图像和第二遥感图像进行自动化融合生成高时间高空间分辨率NDVI数据。
3、对于空间过滤算法，通过结合与目标像元NDVI值的接近程度和与目标像元的欧式距离，获取每一个融合像元的相似像元，对融合的目标像元NDVI值进行校正，减弱因融合粗空间分辨率影像引入的像元边界效应，提高融合影像的连续性。
显然，本领域技术人员应该明白，上述的本发明的各单元或各步骤可以用通用的计算装置来实现，它们可以集中在单个计算装置上,可选地，他们可以用计算机装置可执行的程序代码来实现，从而可以将它们存储在存储装置中由计算装置来执行，或者将它们分别制作成各个集成电路模块，或者将它们中的多个模块或步骤制作成单个集成电路模块来实现。这样，本发明不限制于任何特定的硬件和软件的结合。
以上内容是结合具体的优选实施方式对本发明所作的进一步详细说明，不能认定本发明的具体实施方式仅限于此，对于本发明所属技术领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干简单的推演或替换，都应当视为属于本发明由所提交的权利要求书确定保护范围。

Claims (9)

1.一种基于不同时空分辨率的归一化植被指数数据时空融合方法，其特征在于，包括如下步骤：
不同时空分辨率的归一化植被指数数据获取步骤S110：
分别获取目标区域的高空间分辨率低时间分辨率的第一遥感图像和低空间分辨率高时间分辨率第二遥感图像，通过地表反射率数据分别计算第一和第二遥感图像的归一化植被指数数据NDVI，生成第一遥感图像和第二遥感图像的归一化植被指数数据NDVI时间序列；
对归一化植被指数数据进行优化处理步骤S120：
分别对第一遥感图像和第二遥感图像的归一化植被指数数据NDVI时间序列进行噪点去除，得到稳定的第一遥感图像和第二遥感图像的归一化植被指数数据NDVI时间序列；
初步的高时空分辨率归一化植被指数数据生成步骤S130：
基于第一和第二遥感图像在时间上的线性对应，根据第一遥感图像的低时间分辨率性，得到低时间分辨率的每个像元的NDVI差值时间序列，对上述NDVI差值时间序列进行线性差值得到高时间分辨率的NDVI差值时间序列，而后与每个像元在高时间分辨率点上的第二遥感图像的归一化植被指数数据NDVI时间序列求和得到初步的高空间分辨率高时间分辨率NDVI时间序列数据；
相似像元集计算获取步骤S140：
选择离融合时间点
Figure DEST_PATH_IMAGE001
最近的第一遥感图像的NDVI数值，通过计算目标窗口内每个像元与目标像元的归一化植被指数插值，来寻找目标像元的相似像元集；
最终的高时空分辨率的归一化植被指数数据计算步骤S150：
计算相似像元与目标像元的欧式距离，并将所述欧式距离转换为相似像元权重值，将相似像元在初步的高空间分辨率高时间分辨率归一化植被指数上对应的值与所述权重值用加权求和得到目标像元最终的NDVI值，逐像元遍历研究区为每一个插值NDVI像元进行空间过滤，生成最终的高空间分辨率高时间分辨率NDVI数值。
2.根据权利要求1所述的归一化植被指数数据时空融合方法，其特征在于：
在步骤S110中，对于目标区域的高空间分辨率低时间分辨率的第一遥感图像，计算并生成第一遥感图像的归一化植被指数数据NDVI时间序列具体为：
对目标区域的高空间分辨率低时间分辨率的第一遥感图像进行大气校正，生成地表反射率数据；同时采用FMASK云检测算法，对所述第一遥感图像数据中的云/云阴影、水体、冰雪的像元进行检测过滤；利用公式1计算出每个过滤后像元的归一化植被指数数据NDVI，形成第一遥感图像的归一化植被指数数据NDVI时间序列，
Figure 20078DEST_PATH_IMAGE002
其中，
Figure DEST_PATH_IMAGE003
为近红外波段反射率，
Figure 410084DEST_PATH_IMAGE004
为红光波段反射率；
生成第二遥感图像的归一化植被指数数据NDVI时间序列具体为：获取目标区域的低空间分辨率高时间分辨率第二遥感图像，当第二遥感图像为地表反射率数据时，直接采用公式（1），计算得到第二遥感图像的归一化植被指数数据NDVI时间序列。
3.根据权利要求2所述的归一化植被指数数据时空融合方法，其特征在于：
在步骤S120中：
对于第一遥感图像的归一化植被指数数据NDVI时间序列，将冬季中NDVI<0.1的序列点进行过滤，然后将同时满足公式2和公式3的像元判定为污染像元，进行过滤，得到优化后的第一遥感图像的归一化植被指数数据NDVI时间序列，
Figure DEST_PATH_IMAGE005
其中，NDVI i+j 表示以NDVI i 为中心的四个相邻时间点 NDVI 值，j = -2、-1、1 和 2，AVERAGE是计算相邻四个NDVI值的平均值，STD是对应的标准差。
4.根据权利要求3所述的归一化植被指数数据时空融合方法，其特征在于：
在步骤S120中，
对于第二遥感图像的归一化植被指数数据NDVI时间序列，过滤冬季NDVI<0.1的时间点，利用中值法进行过滤，以每日目标时间点为中心，将9天的窗口内的第二遥感图像的归一化植被指数数据NDVI数据求中值，减少时间序列噪点，将过滤后的第二遥感图像NDVI转投影并重采样为第一遥感图像的投影和空间分辨率，生成优化的第二遥感图像的归一化植被指数数据NDVI时间序列。
5.根据权利要求4所述的归一化植被指数数据时空融合方法，其特征在于：
在步骤S130中，
寻找一年中k个时间对应的MODIS和Landsat NDVI影像对，影像对相应位置的关系如公式4所示，利用公式5将影像对对应的NDVI逐像元求差，生成低时间分辨率的每个像元NDVI差值时间序列，含k个值，
Figure 982010DEST_PATH_IMAGE006
其中，(x, y)为第一遥感图像的
Figure DEST_PATH_IMAGE007
和第二遥感图像
Figure 220225DEST_PATH_IMAGE008
相对应的像元位置，T j (j=1,2,…k)为第一遥感图像和第二遥感图像对应的影像获取时间；
Figure DEST_PATH_IMAGE009
是(x, y)位置的已知时间点的低时间分辨率的第一遥感图像与第二遥感图像的NDVI差值，
Figure 725156DEST_PATH_IMAGE010
将已知时间点的低时间分辨率的NDVI差值时间序列
Figure DEST_PATH_IMAGE011
进行线性插值，获取高时间分辨率的按日的NDVI差值时间序列
Figure 23413DEST_PATH_IMAGE011
，每个像元上任意一个时间点上的高空间分辨率高时间分辨率归一化植被指数
Figure 551477DEST_PATH_IMAGE012
则为该时间点的第二遥感图像归一化植被指数NDVI值与该时间点的高时间分辨率NDVI差值
Figure 718016DEST_PATH_IMAGE011
的和，利用公式6得到初步的高空间分辨率高时间分辨率归一化植被指数（NDVI）数据
Figure 15137DEST_PATH_IMAGE012
，
Figure DEST_PATH_IMAGE013
其中，
Figure 15454DEST_PATH_IMAGE014
为按日的任意时间点
Figure 155448DEST_PATH_IMAGE001
的第二遥感图像NDVI值，
Figure 461097DEST_PATH_IMAGE011
为插值形成的按日的任意时间点
Figure 2937DEST_PATH_IMAGE001
的高时间分辨率的NDVI差值时间序列。
6.根据权利要求5所述的归一化植被指数数据时空融合方法，其特征在于：
在步骤S140中，
选择离融合时间点
Figure 377417DEST_PATH_IMAGE001
最近的第一遥感图像的NDVI数值，对应时间点为
Figure DEST_PATH_IMAGE015
，假设
Figure 411232DEST_PATH_IMAGE001
到
Figure 653995DEST_PATH_IMAGE015
时间段地表植被变化较小，使用最近时间
Figure 925707DEST_PATH_IMAGE015
的第一遥感图像的NDVI来寻找目标像元的相似像元集，每一个目标像元对应的相似像元选择依据：离目标像元距离越近，和目标像元的NDVI值越接近，则与目标像元越接近；首先在以目标像元为中心，利用公式7在m×m的目标窗口内计算每个像元与目标像元的NDVI差值，将NDVI差值从小到大排序，选择前t个差值对应的像元作为目标像元的相似像元集，t为目标像元对应的相似像元的个数，如果在最近时间
Figure 595723DEST_PATH_IMAGE015
的第一遥感图像 NDVI数据上找不到相似像元，就依次循环寻找下一个离融合时间点
Figure 585676DEST_PATH_IMAGE001
最近的时间点，在其第一遥感图像 NDVI数据上按上述步骤查找相似像元集，直到找到相似像元集终止循环，
Figure 632129DEST_PATH_IMAGE016
其中，n为目标窗口内像元位置标识，n =1, 2, …m×m，
Figure DEST_PATH_IMAGE017
为目标像元
Figure 289507DEST_PATH_IMAGE018
在
Figure 740211DEST_PATH_IMAGE015
时刻的第一遥感图像 NDVI值，
Figure DEST_PATH_IMAGE019
为目标窗口内的每个像元
Figure 748618DEST_PATH_IMAGE020
在
Figure 333183DEST_PATH_IMAGE015
时刻的第一遥感图像的 NDVI,
Figure DEST_PATH_IMAGE021
为目标窗口内计算每个像元与目标像元的NDVI差值。
7.根据权利要求6所述的归一化植被指数数据时空融合方法，其特征在于：
在步骤S150中，
采用公式8计算相似像元到目标像元的欧氏距离，然后采用公式9将欧式距离转为相似像元权重值，采用公式10将相似像元在初步的高空间分辨率高时间分辨率归一化植被指数
Figure 842137DEST_PATH_IMAGE022
上对应的值采与相似像元权重值加权求和得到目标像元最终的NDVI值，逐像元遍历研究区为每一个插值NDVI像元进行空间过滤，生成最终的高空间分辨率高时间分辨率NDVI数值，
Figure DEST_PATH_IMAGE023
(8)
Figure 260480DEST_PATH_IMAGE024
(9)
Figure 349659DEST_PATH_IMAGE025
(10)
其中，ω取值为融合数据的空间分辨率，例如ω=30，(x, y) 为目标像元的位置坐标，
Figure 347702DEST_PATH_IMAGE026
为相似像元的位置坐标，
Figure 573147DEST_PATH_IMAGE027
为过滤后融合NDVI数值。
8.根据权利要求1-7中任意一项所述的归一化植被指数数据时空融合方法，其特征在于：
第一遥感图像为Landsat序列遥感图像，第二遥感图像为逐日的MODIS序列遥感图像。
9.一种存储介质，用于存储计算机可执行指令，其特征在于：
所述计算机可执行指令在被处理器执行时执行权利要求1-8中任意一项所述的基于不同时空分辨率的归一化植被指数数据时空融合方法。
CN202210079071.XA 2022-01-24 2022-01-24 基于不同时空分辨率的归一化植被指数数据时空融合方法 Expired - Fee Related CN114092835B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210079071.XA CN114092835B (zh) 2022-01-24 2022-01-24 基于不同时空分辨率的归一化植被指数数据时空融合方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210079071.XA CN114092835B (zh) 2022-01-24 2022-01-24 基于不同时空分辨率的归一化植被指数数据时空融合方法

Publications (2)

Publication Number Publication Date
CN114092835A true CN114092835A (zh) 2022-02-25
CN114092835B CN114092835B (zh) 2022-04-19

Family

ID=80309214

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210079071.XA Expired - Fee Related CN114092835B (zh) 2022-01-24 2022-01-24 基于不同时空分辨率的归一化植被指数数据时空融合方法

Country Status (1)

Country Link
CN (1) CN114092835B (zh)

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115359369A (zh) * 2022-10-19 2022-11-18 中国科学院、水利部成都山地灾害与环境研究所 一种基于时相自适应的山地卫星影像融合方法和系统
CN116071644A (zh) * 2022-12-20 2023-05-05 中化现代农业有限公司 逐日叶面积指数数据反演方法、装置、设备及存储介质
CN116187848A (zh) * 2023-02-23 2023-05-30 安徽星太宇科技有限公司 基于数据分析的遥感载荷综合效能评估系统
CN116975784A (zh) * 2023-09-19 2023-10-31 四川省水利科学研究院 一种高时空分辨率mpdi数据集构建方法、系统及存储介质

Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108613933A (zh) * 2018-06-13 2018-10-02 中南林业科技大学 基于多源遥感数据融合的林地干旱时空动态监测方法
CN111369483A (zh) * 2020-03-05 2020-07-03 北京师范大学 一种融合多源遥感数据生成高时空分辨率遥感数据的方法
CN111666896A (zh) * 2020-06-09 2020-09-15 中国科学院地理科学与资源研究所 一种基于线性融合模型的遥感影像时空融合方法
CN111915694A (zh) * 2020-07-01 2020-11-10 河海大学 一种顾及时空特征的云覆盖像元地表温度重建方法
CN112052720A (zh) * 2020-07-16 2020-12-08 贵州师范学院 一种基于直方图聚类的高时空归一化植被指数ndvi的融合模型
CN113469145A (zh) * 2021-09-01 2021-10-01 中国测绘科学研究院 一种基于高时空分辨率遥感数据的植被物候提取方法
CN113609899A (zh) * 2021-06-24 2021-11-05 河海大学 基于遥感时序分析的退耕地信息定位显示系统

Patent Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108613933A (zh) * 2018-06-13 2018-10-02 中南林业科技大学 基于多源遥感数据融合的林地干旱时空动态监测方法
CN111369483A (zh) * 2020-03-05 2020-07-03 北京师范大学 一种融合多源遥感数据生成高时空分辨率遥感数据的方法
CN111666896A (zh) * 2020-06-09 2020-09-15 中国科学院地理科学与资源研究所 一种基于线性融合模型的遥感影像时空融合方法
CN111915694A (zh) * 2020-07-01 2020-11-10 河海大学 一种顾及时空特征的云覆盖像元地表温度重建方法
CN112052720A (zh) * 2020-07-16 2020-12-08 贵州师范学院 一种基于直方图聚类的高时空归一化植被指数ndvi的融合模型
CN113609899A (zh) * 2021-06-24 2021-11-05 河海大学 基于遥感时序分析的退耕地信息定位显示系统
CN113469145A (zh) * 2021-09-01 2021-10-01 中国测绘科学研究院 一种基于高时空分辨率遥感数据的植被物候提取方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
XUEGANG XING,ET AL.: "An Effective High Spatiotemporal Resolution NDVI Fusion Model Based on Histogram Clustering", 《REMOTE SENSING》 *

Cited By (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115359369A (zh) * 2022-10-19 2022-11-18 中国科学院、水利部成都山地灾害与环境研究所 一种基于时相自适应的山地卫星影像融合方法和系统
CN116071644A (zh) * 2022-12-20 2023-05-05 中化现代农业有限公司 逐日叶面积指数数据反演方法、装置、设备及存储介质
CN116071644B (zh) * 2022-12-20 2023-08-08 中化现代农业有限公司 逐日叶面积指数数据反演方法、装置、设备及存储介质
CN116187848A (zh) * 2023-02-23 2023-05-30 安徽星太宇科技有限公司 基于数据分析的遥感载荷综合效能评估系统
CN116187848B (zh) * 2023-02-23 2024-01-30 安徽星太宇科技有限公司 基于数据分析的遥感载荷综合效能评估系统
CN116975784A (zh) * 2023-09-19 2023-10-31 四川省水利科学研究院 一种高时空分辨率mpdi数据集构建方法、系统及存储介质
CN116975784B (zh) * 2023-09-19 2023-12-29 四川省水利科学研究院 一种高时空分辨率mpdi数据集构建方法、系统及存储介质

Also Published As

Publication number Publication date
CN114092835B (zh) 2022-04-19

Similar Documents

Publication Publication Date Title
CN114092835B (zh) 基于不同时空分辨率的归一化植被指数数据时空融合方法
Henriksen et al. Extracting accurate and precise topography from LROC narrow angle camera stereo observations
Gruszczyński et al. Comparison of low-altitude UAV photogrammetry with terrestrial laser scanning as data-source methods for terrain covered in low vegetation
CN109840553B (zh) 耕地作物类型的提取方法及系统、存储介质、电子设备
KR100678395B1 (ko) 위성영상과 항공영상을 이용한 지리정보 데이터의 실시간위치보정 시스템 및 그 방법
CN102800074B (zh) 基于轮廓波变换的sar图像变化检测差异图生成方法
CN111832518B (zh) 基于时空融合的tsa遥感影像土地利用方法
CN107917881B (zh) 基于冠层内辐射传输机理的光学遥感影像地形校正方法
CN116245757B (zh) 多模态数据的多场景通用性遥感影像云修复方法和系统
CN113744249B (zh) 一种海洋生态环境损害调查方法
CN112364289B (zh) 一种通过数据融合提取水体信息的方法
CN111666896A (zh) 一种基于线性融合模型的遥感影像时空融合方法
Chen et al. A stepwise framework for interpolating land surface temperature under cloudy conditions based on the solar-cloud-satellite geometry
CN113628098A (zh) 一种多时相无云卫星影像重构方法
Qin et al. A coarse elevation map-based registration method for super-resolution of three-line scanner images
CN115761321A (zh) 面向光学遥感影像水体缺失数据重建方法及系统
CN114882169A (zh) 一种基于三维数据的电网工程大数据智能分析系统及方法
Liu et al. Automatic building height estimation with shadow correction over heterogeneous compact cities using stereo Gaofen-7 data at sub-meter resolution
CN114494039A (zh) 一种水下高光谱推扫图像几何校正的方法
CN110991567A (zh) 一种由点到面的遥感瞬时地表温度数据检测方法
Sarabandi et al. Building inventory compilation for disaster management: Application of remote sensing and statistical modeling
CN115205682B (zh) 一种ndvi最大值遥感数据产品无缝生产处理方法
Oshio et al. Generating DTM from DSM Using a Conditional GAN in Built-Up Areas
CN115100537B (zh) 一种基于遥感影像的潮汐能资源评估方法
CN117036987B (zh) 一种基于小波域交叉配对的遥感影像时空融合方法及系统

Legal Events

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
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20220419

CF01 Termination of patent right due to non-payment of annual fee