CN107180158B

CN107180158B - 一种基于温度变化速率的地表温度降尺度方法

Info

Publication number: CN107180158B
Application number: CN201710456288.7A
Authority: CN
Inventors: 陈云浩; 李京; 夏海萍
Original assignee: Beijing Normal University
Current assignee: Beijing Normal University
Priority date: 2017-06-16
Filing date: 2017-06-16
Publication date: 2020-03-10
Anticipated expiration: 2037-06-16
Also published as: CN107180158A

Abstract

本发明公开了一种基于温度变化速率的地表温度降尺度方法，通过静态降尺度和动态降尺度，得到了高时空分辨率的地表温度影像。该方法从温度本身出发，探究地表温度在空间内的变化规律，同时结合了温度年循环模型，将长时间序列的地表温度结合到降尺度模型中，动态地获取了任意时刻的地表温度，使用的降尺度因子是地表温度梯度，该降尺度因子在一年内的变化稳定，相较于NDVI和反照率等降尺度因子，获取更容易，在数据缺失的情况下缺乏时间一致的尺度因子时，可以用天气条件相似的因子代替，并且对结果的影响较小，降尺度方法效果稳定。

Description

一种基于温度变化速率的地表温度降尺度方法

技术领域

本发明涉及一种基于温度变化速率的地表温度降尺度方法。

背景技术

地表温度是估算地表能量通量的重要参数，被广泛应用于各学科领域。这些应用包括土壤水分估计、森林火灾监测、城市热岛效应监测、水文过程研究及气候研究等。

高时间和空间分辨率的地表温度遥感影像在各领域都被广泛需求。但是由于遥感技术的缺陷，同时具有高时间分辨率和高空间分辨率的影像数据难以获取。

航天热红外传感器的出现使得大范围获取地表温度成为可能。然而，受到热红外传感器成像条件的制约，获取高时空分辨率的遥感地表温度仍很难实现。与可见光波段相比，星载热红外传感器的数量很少，其相应数据的空间分辨率也更低。这种低空间分辨率使得热红外影像的混合像元效应更为显著。尽管已有一些卫星发射计划，以获取较高时空分辨率的热红外数据，如小卫星陆地表面热红外成像计划，但是混合像元效应仍无法避免。这使得遥感地表温度降尺度得到了越来越多的关注。

在过去的几十年里，地表温度降尺度的研究已经取得了很大的进展。对地表温度降尺度的方法包括图像融合、统计回归方法、调制像元法和混合法等。这些方法还依赖于除温度以外的其他辅助数据，如NDVI(植被覆盖指数)、反照率等。

发明内容

本发明的目的是解决目前地表温度降尺度算法需要依赖于温度以外的辅助数据，导致算法复杂，稳定性差的技术问题。

为实现以上发明目的，本发明提供一种基于温度变化速率的地表温度降尺度方法，包括如下步骤：

静态降尺度，包括：

a1：根据式T_slope＝arctan(CRLST)＝arctan(ΔT/Δd)，通过输入高分辨率地表温度影像计算高分辨率的地表温度坡度；

式中，T_slope为地表温度坡度，CRLST为地表温度变化速率，ΔT为一幅地表温度影像中两个点之间的温度差异，Δd为一幅地表温度影像中两个点之间的距离；

a2：将高分辨率的地表温度升尺度到低分辨率地表温度作为背景温度；

a3：输入低分辨率的地表温度和高分辨率的地表温度坡度，用移动窗口算法计算降尺度后高分辨率的地表温度。

进一步地，在所述静态降尺度后还包括如下步骤：

动态降尺度，包括：

b1：将温度年循环模型加入到降尺度中，得到

T_high(t)＝T_low(t)+ΔT’(t)＝T_back(t)+ΔT(t)；

其中，T_low(t)＝A sin(2πtf+θ)+B；

式中，T_high(t)是高分辨率的地表温度时间序列，T_low(t)是低分辨率的地表温度时间序列，T_back(t)是一个临时产生的背景温度，ΔT(t)是一幅地表温度影像中两个点之间的温度差异，ΔT’(t)是高、低分辨率的地表温度之间的差值，A是温度的季节振幅，f是频率，θ是相位，B是年平均温度，t是时间；

将地表温度时间序列T_high(t)和T_low(t)输入到非线性回归方程中计算得到温度年循环模型的系数A、θ、B；

b2：根据温度年循环模型的系数A、θ、B计算低分辨率地表温度，作为背景温度；

b3：将高分辨率的温度坡度T_slope用于低分辨率地表温度降尺度。

进一步地，步骤a3中所述移动窗口算法的计算过程如下：

将一个像元分解到9个子像元，则一个像元的地表温度描述为背景温度和温度差异的总和，表示为 T_{high(3i+m，3j+n)}＝T_low(i，j)+ΔT’_{(3i+m，3j+n)}，m，n∈(0，2)；

式中，T_high，T_low和ΔT’分别为高分辨率的地表温度，低分辨率的地表温度以及高、低分辨率的地表温度之间的差值；i和j分别是低分辨率的地表温度的行列号；

降尺度前的地表温度与降尺度后子像元的平均地表温度相等，约束条件如下：

将ΔT’替换为ΔT，则有下式，

T_low(i，j)+ΔT’_{(3i+m，3j+n)}＝T_back(i，j)+ΔT_{(3i+m，3j+n)}；

式中，T_back(i，j)是一个临时产生的背景温度，根据上述约束条件可以计算得到，

将T_back代入式T_{high(3i+m，3j+n)}＝T_low(i，j)+ΔT’_{(3i+m，3j+n)}，则可得降尺度后的高分辨率地表温度T_{high(3i+m，3j+n)}＝T_back(i，j)+ΔT_{(3i+m，3j+n)}，

其中，ΔT_{(3i+m，3j+n)}＝Δd×tan(T_{slope(3i+m，3j+n)})。

进一步地，还对所述高分辨率地表温度 T_{high(3i+m，3j+n)}进行初级平滑处理，以消除格网效应，并最小化温度重新分配产生的误差。

进一步地，所述低分辨率地表温度为3km分辨率的地表温度。

进一步地，所述高分辨率地表温度为1km分辨率的地表温度。

进一步地，步骤b1中所述f取1/365。

进一步地，所述降尺度的精度在白天为2.0K。

进一步地，所述降尺度的精度在夜间为1.0K。

与现有技术相比，本发明的有益效果是：

本算法较为简单，仅从温度本身出发，探究地表温度在空间内的变化规律，同时结合了温度年循环模型，将长时间序列的地表温度结合到降尺度模型中，动态地获取了任意时刻的地表温度，使用的降尺度因子是地表温度梯度，该降尺度因子在一年内的变化稳定，相较于NDVI和反照率等降尺度因子，获取更容易，在数据缺失的情况下缺乏时间一致的尺度因子时，可以用天气条件相似的因子代替，并且对结果的影响较小，降尺度方法效果稳定。

附图说明

图1是静态降尺度算法原理框图；

图2是动态降尺度算法原理框图；

图3是1km分辨率的地表温度遥感数据；

图4是根据图3中的数据计算得出的地表温度的坡度；

图5是3km分辨率的地表温度遥感数据；

图6是对图5中的数据降尺度后的1km分辨率的地表温度；

图7是参考数据和降尺度结果的平均绝对差异；

图8是一年12个月的地表温度数据(3km分辨率)；

图9是一年12个月的地表温度数据(1km分辨率)。

具体实施方式

下面结合附图和具体实施例对本发明作进一步说明。

(一)方法描述

1.1地表温度变化速率

地表温度变化速率的定义类似于坡度在地形中的定义，即温度在特定方向变化的剧烈程度。类似于气象学中的温度梯度，可被表示为 CRLST＝ΔT/Δd. (1)

地表温度的坡度就可由上式计算，

T_slope＝arc tan(CRLST)＝arctan(ΔT/Δd) (2)

式中，T_slope是地表温度的坡度，ΔT是一幅地表温度影像中两点之间的温度差异，Δd代表一幅地表温度影像中两点之间的距离。

1.2地表温度降尺度

一个像元的地表温度可被描述为背景温度和温度差异的总和，表示为 T_high＝T_low+ΔT’； (3)

式中T_high，T_low和ΔT’分别为高分辨率的地表温度，低分辨率的地表温度以及高、低分辨率的地表温度之间的差值。这个式子可被详细描述如下：

T_{high(3i+m，3j+n)}＝T_low(i，j)+ΔT’_{(3i+m，3j+n)}，m，n∈(0，2)； (4) 式中i和j分别是低分辨率的地表温度的行列号；因为此处将像元分解到9 个子像元，所以m和n的范围为0到2。

降尺度前后的温度遵循能量守恒定律，即降尺度前的温度与降尺度后子像元的平均温度相等，约束条件如下：

如ΔT’采用式(2)中的ΔT，则有下式，

T_low(i，j)+ΔT’_{(3i+m，3j+n)}＝T_back(i，j)+ΔT_{(3i+m，3j+n)}； (6)

式中的T_back(i，j)是一个临时产生的背景温度，根据式(5)的约束条件可以计算得到，

将T_back代入式(4)则可得降尺度后温度，

T_{high(3i+m，3j+n)}＝T_back(i，j)+ΔT_{(3i+m，3j+n)}； (8)

ΔT_{(3i+m，3j+n)}＝Δd×tan(T_{slope(3i+m，3j+n)})。 (9)

1.3初级高分辨率LST平滑

对低分辨率的温度数据进行逐像元分解时，会产生明显的格网效应，为了消除格网效应，采用移动窗口对初级高分辨率温度进行平滑。这个平滑的过程不同于以往的平滑方法(如低通滤波)，是对初级高分辨率温度的再次降尺度。

这个过程中的低分辨率像元是初级高分辨率温度在移动窗口内合成的，因此可以在消除格网效应的过程中同时最小化温度重新分配产生的误差。

1.4添加时间序列的温度细节的降尺度

将温度年循环模型(ATCM)加入到降尺度中，得到

T_high(t)＝T_low(t)+ΔT’(t)＝T_back(t)+ΔT(t)； (10)

其中，

T_low(t)＝A sin(2πtf+θ)+B； (11)

式中，T_high(t)是高分辨率的地表温度时间序列，T_low(t)是低分辨率的地表温度时间序列，A是温度的季节振幅，f是频率，在此处取1/365，因为一个行星年是365天，θ是相位，B是年平均温度，t是时间。

(二)方法实现

(1)静态降尺度步骤，参见图1：

步骤一：根据式(2)，将输入的高分辨率地表温度影像用以计算地表温度的坡度(计算结果参见图4)，在该实验中用的高分辨率地表温度影像数据是图3 所示的1km分辨率的地表温度遥感数据；

步骤二：将高分辨率的地表温度升尺度到低分辨率作为背景温度，此处低分辨率的地表温度为图5所示的3km分辨率的地表温度遥感数据；

步骤三：根据1.2中的公式(3)-(9)，输入3km分辨率的地表温度数据和1km 分辨率的地表温度坡度数据，用移动窗口算法计算降尺度后的地表温度，降尺度结果参见图6，参考数据和降尺度结果的平均绝对差异参见图7，此处参考数据是指最初的高分辨率温度影像，本发明是将高分辨率温度升尺度到低分辨率的温度，再对低分辨率温度进行降尺度到高分辨率温度。

(2)动态降尺度步骤，参见图2：

步骤一：将时间序列的温度数据T_high(t)和T_low(t)输入到非线性回归方程中计算温度年循环模型的系数A、θ、B，这个具体计算过程在IDL(交互式数据语言)中实现，计算结果参见表2；

图1和图2中各变量的含义参见表1；

表1.变量说明

表2.计算得到的稳定年循环模型的系数

步骤二：根据表2中的系数，计算低分辨率(如3km分辨率)地表温度数据，作为背景温度；计算结果参见图8；

步骤三：将高分辨率(如1km分辨率)的地表温度坡度用于地表温度降尺度，降尺度结果参见图9。

本降尺度算法方法简单，采用的降尺度方法效果稳定，该方法中采用地表温度坡度作为降尺度因子，相较于NDVI和反照率等降尺度因子，获取更容易。在数据缺失的情况下，可以用类似条件的因子代替，并且对结果的影响较小。

除上述实施方式外，本发明还可以有其他实施方式，凡采用等同替换或等效变换形成的技术方案，均落在本发明的保护范围内。

Claims

1.一种基于温度变化速率的地表温度降尺度方法，其特征在于，包括如下步骤：

静态降尺度，包括：

a1:根据式T_slope＝arctan(CRLST)＝arctan(ΔT/Δd)，通过输入高分辨率地表温度影像计算高分辨率的地表温度坡度；

a2:将高分辨率的地表温度升尺度到低分辨率地表温度作为背景温度；

a3:输入低分辨率的地表温度和高分辨率的地表温度坡度，用移动窗口算法计算降尺度后高分辨率的地表温度；

步骤a3中所述移动窗口算法的计算过程如下：

将一个像元分解到9个子像元，则一个像元的地表温度描述为背景温度和温度差异的总和，表示为

T_{high(3i+m，3j+n)}＝T_low(i，j)+ΔT’_{(3i+m，3j+n)}，m，n∈(0，2)；

将ΔT’替换为ΔT，则有下式，

T_low(i，j)+ΔT’_{(3i+m，3j+n)}＝T_back(i，j)+ΔT_{(3i+m，3j+n)}；

式中，T_back(i,j)是一个临时产生的背景温度，根据上述约束条件可以计算得到，

将T_back代入式T_{high(3i+m,3j+n)}＝T_low(i,j)+ΔT’_(3i+m,3j+n)，则可得降尺度后的高分辨率地表温度T_{high(3i+m,3j+n)}＝T_back(i,j)+ΔT_(3i+m,3j+n)，

其中，ΔT_(3i+m,3j+n)＝Δd×tan(T_{slope(3i+m,3j+n)})。

2.如权利要求1所述的基于温度变化速率的地表温度降尺度方法，其特征在于，在所述静态降尺度后还包括如下步骤：

动态降尺度，包括：

b1:将温度年循环模型加入到降尺度中，得到

T_high(t)＝T_low(t)+ΔT'(t)＝T_back(t)+ΔT(t)；

其中，T_low(t)＝Asin(2πtf+θ)+B；

b2:根据温度年循环模型的系数A、θ、B计算低分辨率地表温度，作为背景温度；

b3:将高分辨率的温度坡度T_slope用于低分辨率地表温度降尺度。

3.如权利要求2所述的基于温度变化速率的地表温度降尺度方法，其特征在于，还对所述高分辨率地表温度T_{high(3i+m,3j+n)}进行初级平滑处理，以消除格网效应，并最小化温度重新分配产生的误差。

4.如权利要求1-3任一项所述的基于温度变化速率的地表温度降尺度方法，其特征在于，所述低分辨率地表温度为3km分辨率的地表温度。

5.如权利要求4所述的基于温度变化速率的地表温度降尺度方法，其特征在于，所述高分辨率地表温度为1km分辨率的地表温度。

6.如权利要求2或3所述的基于温度变化速率的地表温度降尺度方法，其特征在于，步骤b1中所述f取1/365。

7.如权利要求5所述的基于温度变化速率的地表温度降尺度方法，其特征在于，所述降尺度的精度在白天为2.0K。

8.如权利要求5所述的基于温度变化速率的地表温度降尺度方法，其特征在于，所述降尺度的精度在夜间为1.0K。