CN104360040A - 一种基于starfm融合技术的遥感墒情监测方法 - Google Patents

Abstract

本发明公开了一种基于STARFM融合技术的遥感墒情监测方法，该方法主要包括：步骤一、根据TM、ETM数据对HJ星CCD数据进行预处理；步骤二、利用大气辐射传输模型MODTRAN4实现HJ星的逐像元大气校正；步骤三、对MODIS数据进行预处理；步骤四、搜索相似像素，计算HJ星和MODIS数据的相似像素权重、转化系数；步骤五、基于STARFM融合算法进行预测。本发明通过STARFM融合方法综合利用HJ星和MODIS数据进行融合，获得时空较高分辨率的干旱监测结果，克服了单一卫星MODIS数据空间分辨率低和HJ星不能实现逐日观测的难题，使得反演的时间分辨率提高到1天，墒情监测的空间分辨率提高到30米，在农业、水利、环保与生态监测领域具有重要应用价值。

Description

一种基于STARFM融合技术的遥感墒情监测方法

技术领域

本发明涉及一种基于STARFM融合技术的遥感墒情监测方法，具体涉及一种综合利用HJ星和MODIS数据获得较高时空分辨率的干旱监测方法，属于环境监测技术领域。

背景技术

干旱是主要的自然灾害之一，由于其影响范围大、持续时间长，并且后续影响久等特点使得干旱造成的经济损失远远超过了其他气象灾害，因此对干旱状况进行监测，及时了解旱情信息，掌握干旱动向，为相关部门提供科学合理的依据，对减少因干旱造成的损失具有十分重要的意义。

传统的干旱监测方法主要有土壤湿度计法、称重法、电阻法等，虽然具有精度较高的优点，但是由于监测花费的人力物力大、监测范围有限、监测速度较慢，难以满足大范围的实时监测。随着科学技术的进步，遥感技术因其宽广的覆盖范围、较高的时空分辨率以及获取数据快捷等优点，在干旱的监测研究工作中作为大范围内信息获取的有效手段得到了广泛的应用，较传统的干旱监测方法有不可比拟的潜在优势。

目前应用较为广泛的旱情遥感监测评估方法主要有：基于土壤热惯量的旱情监测方法、基于蒸散模型的旱情监测方法和基于植被指数的旱情监测方法。然而上述利用单一卫星遥感技术的监测干旱方法虽然取得了一定成功，但由于单一卫星本身的缺陷，如分辨率低，周期长等，无法完全满足干旱监测的需要。而基于多卫星遥感的监测利用不同卫星各自的特点，具有其优势：可以有目的对各变量参数采取遥感融合模型，提高反演的时空精度，可以增大观测范围、增强观测质量、提高观测频率，弥补了单个卫星监测造成区域旱情监测不准的缺陷。因此，基于多卫星融合进行干旱监测，将能很大程度上提高干旱监测的准确度。但是，目前基于多卫星融合进行干旱监测的方案并不多见，反演的时间精度低，墒情监测的空间分辨率精度不高，干旱监测结果并不理想。

发明内容

本发明所要解决的技术问题是针对前述背景技术中的缺陷，综合HJ星和MODIS数据，提出一种耦合两者数据的方法进行干旱监测，能够实现较高时空分辨率的干旱监测结果，提高干旱监测水平。

本发明为解决上述技术问题采用以下技术方案：

一种基于STARFM融合技术的遥感墒情监测方法，通过STARFM融合方法综合利用HJ星和MODIS数据进行融合，得到干旱监测结果，具体包括以下步骤：

步骤1、根据TM、ETM数据对HJ星CCD数据进行预处理；

步骤2、利用大气辐射传输模型MODTRAN4实现HJ星的逐像元大气校正；

步骤3、MODIS数据预处理；

步骤4、搜索相似像素，计算HJ星和MODIS数据的相似像素权重、转化系数；

步骤5、基于STARFM融合算法进行预测：基于STARFM融合方法获取同一时间的较高分辨率的HJ星CCD数据与较低分辨率MODIS数据，通过计算影像间空间分布的差异，结合另一时间较低分辨率的MODIS数据进行相应时间较高分辨率的HJ星CCD数据的预测。

进一步的，本发明的基于STARFM融合技术的遥感墒情监测方法，所述步骤1具体包括步骤：

(1.1)对TM、ETM、HJ星CCD数据进行镶嵌处理并裁剪；

(1.2)根据裁剪出的TM、ETM底图作为基准图像校正HJ星CCD数据，控制校正误差；

(1.3)对HJ星CCD数据进行辐射定标，并计算表观反射率。

进一步的，本发明的基于STARFM融合技术的遥感墒情监测方法，所述步骤2具体包括：利用MODTRAN 4模型模拟辐射传输过程及对应大气参数下的表观辐亮度，建立大气校正查找表，对HJ数据进行逐像元校正，其中辐射传输方程为：

$L_m=L_0+\frac{\mathrm{\ρ}}{1-\mathrm{s\ρ}}\·\frac{{\mathrm{TF}}_d}{\mathrm{\π}}$

其中L_m表示接收的表观辐亮度，L₀表示路径辐射既程辐射，T表示地表到传感器的透过率，s为大气底层球面反照率，ρ为像元反射率，F_d为到达地表的下行辐射通量密度。

进一步详细介绍步骤2：

(2.1)根据输入的大气参数假设三个像元反射率ρ，MODTRAN运行三次得到三个模拟的星上幅亮度L_m；根据得到的结果ρ和L_m，利用辐射传输方程公式建立三个以L₀，s及TF_d/π为变量的三元一次方程组，求解方程组得出一组与输入大气参数对应的L₀，s及TF_d/π。

(2.2)按照上述步骤输入不同大气参数建立通用查找表。

(2.3)读取HJ星的CCD数据，若是DN值将其图像转化为表观幅亮度图像；读取气溶胶光学厚度(AOT)数据的每个像元，根据每个像元的气溶胶光学厚度值，找到与该大气参数对应的三个参数：L₀，s及TF_d/π。

(2.4)逐像元校正CCD或HSI数据。

进一步的，本发明的基于STARFM融合技术的遥感墒情监测方法，其特征在于，所述步骤3中MODIS数据的预处理包括几何校正、图像裁剪、大气校正。

进一步的，本发明的基于STARFM融合技术的遥感墒情监测方法，所述步骤4中搜索相似像素包括搜索HJ卫星和MODIS数据成像时间相同或相近的图像，计算HJ星和MODIS数据的相似像素权重函数如下：

$W_{\mathrm{ijk}}=(1/C_{\mathrm{ijk}})/{\mathrm{\Σ}}_{i=1}^W{\mathrm{\Σ}}_{j=1}^W{\mathrm{\Σ}}_{k=1}^n(1/C_{\mathrm{ijk}})$

其中，W_ijk为权重系数，C_ijk是由窗口中心点的预测像元与窗口中其它像元的光谱距离、时间距离、空间距离计算的结果，W是窗口的大小，n代表输入已知图像对的数目。

进一步的，本发明的基于STARFM融合技术的遥感墒情监测方法，所述步骤5基于STARFM融合算法进行预测，其计算公式如下：

$L(x_{W/2},y_{W/2},T_2)={\mathrm{\Σ}}_{i=1}^W{\mathrm{\Σ}}_{j=1}^W{\mathrm{\Σ}}_{k=1}^n{W}_{\mathrm{ijk}}(M(x_i,y_i,T_2)+L(x_i,y_i,T_1)-M(x_i,y_i,T_1))$

其中，L(x_W/2,y_W/2,T₂)是预测的T₂时刻的高分辨像元值，W是窗口的大小，n代表输入已知图像对的数目，W_ijk为权重系数，(x_W/2,y_W/2)是窗口中间的像元，M(x_i,y_i,T₂)是窗口位置(x_i,y_i)处在T₂时刻的低分辨率像元值，L(x_i,y_i,T₁)和M(x_i,y_i,T₁)则分别是高分辨率蒸散发ET和低分辨率蒸散发ET在T₁时刻的相应像元值。

本发明采用以上技术方案与现有技术相比，具有以下技术效果：

本发明利用STARFM融合方法通过每天的MODIS数据得到每天较高分辨率的HJ星数据，解决了单一卫星监测得到的MODIS数据空间分辨率较低以及HJ星不能每天观测的缺陷，将HJ星与MODIS相结合，得到时空较高分辨率的合成影像，使得反演的时间精度提高到一天，墒情监测的空间分辨率精度提高到30m，从而可以得到较高时空分辨率的干旱监测结果。

附图说明

图1是本发明的一种基于STARFM融合技术的遥感墒情监测方法的流程图。

具体实施方式

下面结合附图对本发明的技术方案做进一步的详细说明：

本技术领域技术人员可以理解的是，除非另外定义，这里使用的所有术语(包括技术术语和科学术语)具有与本发明所属领域中的普通技术人员的一般理解相同的意义。还应该理解的是，诸如通用字典中定义的那些术语应该被理解为具有与现有技术的上下文中的意义一致的意义，并且除非像这里一样定义，不会用理想化或过于正式的含义来解释。

首先，介绍一下HJ星(环境减灾卫星)的基本概念：HJ星是国家2008年9月发射的用于环境与灾害监测预报的小卫星星座。主要有HJ-1-A星和HJ-1-B星两颗星座。HJ-1-A卫星和HJ-1-B卫星的轨道完全相同，相位相差180度。HJ-1-A星搭载了2台CCD相机和1台超光谱成像仪(HSI)，HJ-1-B星搭载2台CCD相机和1台红外相机(IRS)。两卫星上装载的两台CCD相机以星下点对称放置联合完成对地观测幅宽度700公里的推扫成像，具有30米的地面象元分辨率，两台CCD相机组网后重访周期仅为2天。

此外，A卫星装载有一台超光谱成像仪，完成宽度为50公里对地观测、100米的地面象元分辨率和110～128个光谱谱段的推扫成像。装载在B卫星的一台红外相机，完成720公里的对地幅宽、150米/300米的地面像元分辨率、近短中长4个光谱谱段的成像，如表1。

表1环境减灾A、B星有效载荷指标一览表

接着，介绍MODIS的基本概念：MODIS是美国国家航空航天局发射的搭载在Terra和Aqua上的主要的探测器，也是实行全球免费接收并无偿使用的传感器。它具有36个光谱波段，包括可见光～热红外(0.405～14.385μm)，不同通道根据其光谱的反应，可以用于获得不同的目标信息。数据空间分辨率有250m、500m、1000m三个尺度，扫描观测宽度达2330km，每1～2天对地球表面观测一次，提供全球较高辐射度分辨率的图像数据。

TM、ETM简介：

TM是一种遥感器，搭载在美国陆地卫星Landsat系列卫星上。TM影像是指美国陆地卫星4～5号专题制图仪(thematic mapper)所获取的多波段扫描影像。有7个波段，其波谱范围：TM-1为0.45～0.52微米，TM－2为0.52～0.60微米，TM－3为0.63～0.69微米，以上为可见光波段；TM-4为0.76～0.90微米，为近红外波段；TM-5为1.55～1.75微米，TM-7为2.08～2.35微米，为中红外波段；TM-6为10.40～12.50微米，为热红外波段。影像空间分辨率除热红外波段为120米外，其余均为30米，像幅185×185公里。每波段像元数达61662个(TM-6为15422个)。一景TM影像总信息量为230兆字节，约相当于MSS影像的7倍。因TM影像具较高空间分辨率、波谱分辨率、极为丰富的信息量和较高定位精度，成为20世纪80年代中后期得到世界各国广泛应用的重要的地球资源与环境遥感数据源。能满足有关农、林、水、土、地质、地理、测绘、区域规划、环境监测等专题分析和编制1∶10万或更大比例尺专题图，修测中小比例尺地图的要求。

星上携带专题制图仪ETM，ETM具有8个波段，其中1-5波段和7波段是多光谱波段，空间分辨率是30米，第六波段是热红外波段，空间分辨率是120米，第8波段为全色波段，分辨率为15米。景宽185公里，景面积为34225平方公里。此星一直在正常运转，2000年遥感卫星地面站开始接收数据。

本实施例采用2011年3月25日、4月11日的MODIS与HJ星数据，以及2011年5月17日的MODIS数据，来得到2011年5月17日的融合结果，如图1所示，具体实现步骤为：

步骤1、根据TM、ETM数据对HJ星CCD数据进行预处理；

步骤2、利用大气辐射传输模型MODTRAN4实现HJ星的逐像元大气校正，其计算流程如下：

(1)根据输入的大气参数假设三个地表反射率(ρ＝0，0.5，0.8)，MODTRAN运行三次得到三个模拟的星上幅亮度L_m；根据得到的结果ρ和L_m，利用辐射传输方程公式建立三个以L₀，s及TF_d/π为变量的三元一次方程组，求解方程组得出一组与输入大气参数对应的L₀，s及TF_d/π。

(2)按照上述步骤输入不同大气参数建立通用查找表。

(3)读取HJ星的CCD数据，若是DN值将其图像转化为表观幅亮度图像。读取气溶胶光学厚度(AOT)数据的每个像元，根据每个像元的气溶胶光学厚度值，找到与该大气参数对应的三个参数：L₀，s及TF_d/π。

(4)逐像元校正CCD或HSI数据。

步骤3、MODIS数据的预处理，包括：几何校正、图像裁剪、大气校正等。

步骤4、搜索相似像素，计算HJ星和MODIS数据的计算权重函数如下：

$W_{\mathrm{ijk}}=(1/C_{\mathrm{ijk}})/{\mathrm{\Σ}}_{i=1}^W{\mathrm{\Σ}}_{j=1}^W{\mathrm{\Σ}}_{k=1}^n(1/C_{\mathrm{ijk}})$

其中，C_ijk是由窗口中心点的预测像元与窗口中其他的像元(包括多时相数据)的光谱距离、时间距离、空间距离计算的结果，W是窗口的大小，n代表输入已知图像对的数目。

步骤5、基于STARFM融合方法获取同一时间(T₁)的较高分辨率的HJ星CCD数据与较低分辨率MODIS数据，通过计算影像间空间分布的差异，结合另一时间(T₂)较低分辨率的MODIS数据进行相应时间较高分辨率的HJ星CCD的预测，其计算公式如下：

$L(x_{W/2},y_{W/2},T_2)={\mathrm{\Σ}}_{i=1}^W{\mathrm{\Σ}}_{j=1}^W{\mathrm{\Σ}}_{k=1}^n{W}_{\mathrm{ijk}}(M(x_i,y_i,T_2)+L(x_i,y_i,T_1)-M(x_i,y_i,T_1))$

其中，L(x_W/2,y_W/2,T₂)是预测的T₂时刻的高分辨像元值，W是窗口的大小，W_ijk为权重系数，(x_W/2,y_W/2)是窗口中间的像元，M(x_i,y_i,T₂)是窗口位置(x_i,y_i)处在T₂时刻的像元值，L(x_i,y_i,T₁)和M(x_i,y_i,T₁)则是高分辨率ET和低分辨率ET在T₁时刻的相应像元值。

本技术领域技术人员可以理解的是，本发明中已经讨论过的各种操作、方法、流程中的步骤、措施、方案可以被交替、更改、组合或删除。进一步地，具有本发明中已经讨论过的各种操作、方法、流程中的其他步骤、措施、方案也可以被交替、更改、重排、分解、组合或删除。进一步地，现有技术中的具有与本发明中公开的各种操作、方法、流程中的步骤、措施、方案也可以被交替、更改、重排、分解、组合或删除。

以上只是对本发明的优选实施方式进行了描述。对该技术领域的普通技术人员来说，根据以上实施方式可以很容易地联想到其它的优点和变形。因此，本发明并不局限于上述实施方式，其仅仅作为例子对本发明的一种形态进行详细、示范性的说明。在不背离本发明宗旨的范围内，本领域普通技术人员在本发明技术的方案范围内进行的通常变化和替换，都应包含在本发明的保护范围之内。

Claims (6)

1.一种基于STARFM融合技术的遥感墒情监测方法，其特征在于，通过STARFM融合方法综合利用HJ星和MODIS数据进行融合，得到干旱监测结果，具体包括以下步骤：

步骤1、根据TM、ETM数据对HJ星CCD数据进行预处理；

步骤2、利用大气辐射传输模型MODTRAN4实现HJ星的逐像元大气校正；

步骤3、MODIS数据预处理；

步骤4、搜索相似像素，计算HJ星和MODIS数据的相似像素权重、转化系数；

步骤5、基于STARFM融合算法进行预测：基于STARFM融合方法获取同一时间的较高分辨率的HJ星CCD数据与较低分辨率MODIS数据，通过计算影像间空间分布的差异，结合另一时间较低分辨率的MODIS数据进行相应时间较高分辨率的HJ星CCD数据的预测。

2.根据权利要求1所述的基于STARFM融合技术的遥感墒情监测方法，其特征在于，所述步骤1具体包括步骤：

(1.1)对TM、ETM、HJ星CCD数据进行镶嵌处理并裁剪；

(1.2)根据裁剪出的TM、ETM底图作为基准图像校正HJ星CCD数据，控制校正误差；

(1.3)对HJ星CCD数据进行辐射定标，并计算表观反射率。

3.根据权利要求1所述的基于STARFM融合技术的遥感墒情监测方法，其特征在于，所述步骤2具体包括：利用MODTRAN 4模型模拟辐射传输过程及对应大气参数下的表观辐亮度，建立大气校正查找表，对HJ数据进行逐像元校正，其中辐射传输方程为：

$L_m=L_0+\frac{\mathrm{\ρ}}{1-\mathrm{s\ρ}}\·\frac{{\mathrm{TF}}_d}{\mathrm{\π}}$

其中L_m表示接收的表观辐亮度，L₀表示路径辐射既程辐射，T表示地表到传感器的透过率，s为大气底层球面反照率，ρ为像元反射率，F_d为到达地表的下行辐射通量密度。

4.根据权利要求1所述的基于STARFM融合技术的遥感墒情监测方法，其特征在于，所述步骤3中MODIS数据的预处理包括几何校正、图像裁剪、大气校正。

5.根据权利1所述基于STARFM融合技术的遥感墒情监测方法，其特征在于，所述步骤4中搜索相似像素包括搜索HJ卫星和MODIS数据成像时间相同或相近的图像，计算HJ星和MODIS数据的相似像素权重函数如下：

$W_{\mathrm{ijk}}=(1/C_{\mathrm{ijk}})/{\mathrm{\Σ}}_{i=1}^W{\mathrm{\Σ}}_{j=1}^W{\mathrm{\Σ}}_{k=1}^n(1/C_{\mathrm{ijk}})$

其中，W_ijk为权重系数，C_ijk是由窗口中心点的预测像元与窗口中其它像元的光谱距离、时间距离、空间距离计算的结果，W是窗口的大小，n代表输入已知图像对的数目。

6.根据权利1所述基于STARFM融合技术的遥感墒情监测方法，其特征在于，所述步骤5基于STARFM融合算法进行预测，其计算公式如下：

$L(x_{W/2},y_{W/2},T_2)={\mathrm{\Σ}}_{i=1}^W{\mathrm{\Σ}}_{j=1}^W{\mathrm{\Σ}}_{k=1}^n{W}_{\mathrm{ijk}}(M(x_i,y_i,T_2)+L(x_i,y_i,T_1)-M(x_i,y_i,T_1))$

其中，L(x_W/2,y_W/2,T₂)是预测的T₂时刻的高分辨像元值，W是窗口的大小，n代表输入已知图像对的数目，W_ijk为权重系数，(x_W/2,y_W/2)是窗口中间的像元，M(x_i,y_i,T₂)是窗口位置(x_i,y_i)处在T₂时刻的低分辨率像元值，L(x_i,y_i,T₁)和M(x_i,y_i,T₁)则分别是高分辨率蒸散发ET和低分辨率蒸散发ET在T₁时刻的相应像元值。