一种结合全局与局域信息的遥感地表温度年变化模拟方法

1.本发明涉及一种结合时域全局与局域信息的遥感地表温度年变化模拟方法，属于遥感地表温度年变化模拟的研究领域。


背景技术：

2.地表温度是衡量地表——大气间能量交换的关键参数，在各类地学应用中扮演着至关重要的角色，例如城市热环境监测、土壤水分和蒸散发估算、森林火灾检测等。热红外遥感是获取大范围地表温度的主要方式，然而由于卫星重访周期和云遮挡等因素的影响，遥感观测的地表温度在时间和空间维度上均存在较多的缺失，这极大地限制了地表温度数据的应用。因此通过时序离散的卫星观测地表温度重建时空连续的地表温度，是热红外遥感研究的热点问题。
3.年内地表温度变化模型作为典型的时域全局算法，可以通过少量离散的观测地表温度重建年内时间连续的地表温度。同时，这类模型作为重要的地表温度重建算法，是连接地表温度年际和日内变化模型的纽带和探索地表温度时间尺度演变模式的关键，已经受到学界的密切关注，例如时空无缝地表温度产品生产，地表温度和城市热岛长时间序列分析、地表温度时空降尺度等。学界对年内地表温度重建开展了一系列的探究，现有的年内地表温度重建模型主要可以分为三类。第一类是统计或傅里叶级数型模型，主要通过纯粹的数学统计方法对地表温度进行时间插值，例如傅里叶分析、统计型时间序列分析。这类模型可以较好地模拟地表温度年变化，但模型较为复杂，泛化能力较低。第二类是物理或准物理模型，这类模型从地表能量平衡方程出发，能体现地表热性质的变化机制，且模型各参数在量化地表性质方面具有重要价值。例如三参数地表温度年内变化模型和五参数地表温度年内变化模型等。与统计或傅里叶级数型模型相比，模型形式更简单，自由参数较少，更容易通过少量数据点构建模型，且模型中的自由参数可量化地标性质。然而，受限于简单的形式，这类模型精度有待提高。第三类是统计型与物理型相结合的半物理模型，这类模型是在物理模型的基础上，旨在通过引入统计方法或地表实测辅助数据等辅助数据(例如空气温度、ndvi等)反映地表温度逐日波动情况。例如考虑了空气温度和ndvi的改进的地表温度年变化模型(atce模型)、考虑了植被物候影响的基于植被的地表温度年变化模型(patc模型)以及考虑了多种辅助参数和谐波数量的混合地表温度年变化模型(atch模型)。相较于上述两类模型，这类模型提高了对地表温度逐日波动的模拟精度，但上述模型都是时域全局模型，并未考虑地表温度在局域尺度的相关性。
4.上述研究表明，通过谐波的拓展和地表温度相关辅助参数的耦合，年内地表温度变化模型的精度可以得到提高。然而，现有模型仍存在以下问题：(1)难以反映地表温度空间相关性。年内地表温度变化模型本质上属于时域全局算法，其仅从时域全局角度考虑地表温度的特点，并未考虑局域时空尺度上地表温度的相关性，这会造成地表温度缺值较少时模型精度下降。(2)难以揭示地表温度与相关辅助参数间的复杂关系。除了太阳辐射外，地表温度还受天气条件、高程等因素的影响，然而它们之间并非简单的线性回归关系，这会
导致相关参数难以较好刻画地表温度的逐日波动情况。因此，需要研究新的能够用时考虑时域全局与局域时空信息的地表温度重建方法，进一步提高模型的精度和适用性。


技术实现要素：

5.本发明的目的在于：克服上述现有年内地表温度变化模型对时空局域信息的忽略，提出一种结合时域全局与局域信息的遥感地表温度年变化模拟方法。
6.地表温度在时域全局与时空局域尺度上都有较强的相关性。地表温度年内变化可以分解为两个分量：其一是由太阳辐射变化引起的年内周期性变化(即全局变化)，其二是由天气条件与地表属性等因素引起的逐日波动(即局域变化)。其中，地表温度的年周期变化主要是通过地表温度的全局信息进行模拟，而地表温度的逐日变化主要通过地表温度的局域信息进行刻画。因此，本发明同时耦合了地表温度的全局与局域信息，充分考虑了两个地表温度年内分量的模拟，提出一种结合全局与局域信息的遥感地表温度年变化模拟方法(atc_gl模型)。该方法可以通过公式(1)对目标区域进行地表温度重建，
7.lst(d)＝lst_global(d) lst_local(d)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(1)
8.式中，lst(d)表示相对春分点第d天的重建地表温度，lst_global(d)和lst_local(d)分别表示相对春分点第d天地表温度的年周期变化和逐日波动。
9.在时域全局尺度，本发明通过五参数地表温度年变化模型(acp5模型)模拟地表温度的年周期变化。由于模型的简洁性以及对太阳辐射周期性的充分考虑，五参数地表温度年变化模型可以较好模拟地表温度年周期变化。地表温度年周期变化的模拟可以表示为：
[0010][0011]
其中mask表示年均地表温度，yast1和yast2分别表示年内地表温度的两次振幅，d表示相对于春分点的天数，θ表示相对春分点的相位移，k1和k2是两个常量，分别为2π/365和4π/365。
[0012]
在时空局域尺度，通过随机森林算法建立地表温度与相关参数的统计回归关系，进而刻画缺失的地表温度逐日波动。由于固定的时空窗口难以同时兼顾模型的精度和适用性，因此本发明以缺失像元为中心，建立动态滑动时空窗口，并在动态滑动时空窗口中实现地表温度逐日波动的重建。
[0013]
本发明考虑了地表温度在时空尺度上的相关性，因此令x和y分别表示空间维度上横轴和纵轴的距离，令d为时间维度上的距离。假设整个研究为a，可以表示为：
[0014][0015]
其中z
x
，zy和zd分别表示研究区域a在空间横轴、空间纵轴和时间维度上覆盖的范围，x，y和d分别表示a在三个维度上的最大值。
[0016]
在研究区域内，任一像元点可以表示为p＝p[x,y,d]，其中x和y分别表示p在研究区域内的横坐标和纵坐标，d表示该天相对于春分点的偏移天数。令某地表温度缺失像元为p0，则p0＝p0[x0,y0,d0]我们以缺失像元p0为中心，建立初始时空窗口。初始时空窗口的空间横轴偏移量为α0，空间纵轴偏移量为β0，时间维度偏移量为γ0。因此初始时空窗口b1可以表示为：
[0017]
b1＝{p[x,y,z]∈a:x
0-a0≤x≤x0 a0,y
0-β0≤y≤y0 β0,d
0-γ0≤d≤d0 γ0}
ꢀꢀꢀ
(4)
[0018]
为探究初始时空窗口的大小，本发明使用交叉验证法，对最佳初始时空窗口的选择进行探究(如附图3所示)通过探究我们的初始空间横轴(α0)和初始空间纵轴的偏移量(β0)均设置为4，初始时间偏移量(γ0)设置为0。
[0019]
此外，初始时空窗口内已知地表温度的数量和分布也会影响模型的精度。因此本发明设置一个阈值ω＝10，当初始时空窗口内的已知地表温度小于ω时，同时扩大时空窗口，令每次在空间横轴、空间纵轴和时间维度上的扩大步调为i。因此，当初始时空窗口内的已知地表温度少于阈值ω时，动态时空窗口b可以表示为：
[0020][0021]
其中，扩大步调每次的增量为1，直至时空窗口内的已知地表温度数量大于阈值ω。
[0022]
当时空窗口b内的已知地表温度数量大于阈值时，本发明通过随机森林算法建立地表温度与相关辅助参数的统计回归关系，进而重建缺失的地表温度逐日波动。鉴于地表温度和空气温度的逐日波动具有较高的相似性，且地表温度的变化还受地表覆盖类型的变化和高程差异的影响，因此本发明选取是空气温度逐日波动、归一化植被指数(ndvi)以及高程(dem)。值得注意的是，当地表温度数据缺失时，该像元该天的地表温度逐日波动缺失。因此重建缺失的地表温度逐日波动可以表示为：
[0023]
lst_local(d)＝f
rf
(sat_local(d),ndvi,dem)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(6)
[0024]
其中，空气温度逐日波动(sat_local(d))可以表示为：
[0025][0026]
式中，sat_local(d)表示相对春分点第d天的空气温度逐日波动，sat(d)表示观测的空气温度，sat_global(d)表示通过五参数地表温度年变化模型模拟的空气温度年内变化，参数mask1、yast3、yast4、θ1通过最小二乘法，利用各个像元年内大于四个的观测地表温度进行拟合。
[0027]
值得注意的是，本发明提出的atc_gl模型是一个地表温度重建框架，可以根据不同的研究需求，引入更多不同组合的辅助参数。这可以在不同场景中同时考虑模型的精度的效率。
[0028]
本发明同时提高了模型精度和适用性，可精确重建地表温度，有望适用于高精度地表温度产品的生产。
附图说明
[0029]
结合附图对本发明作进一步的说明。
[0030]
图1为本发明遥感地表温度年变化模拟方法的流程图。
[0031]
图2为本发明提出的atc_gl模型的重建结果与原始观测数据图。
[0032]
图3为初始空间窗口的选择对模型精度的影响图，w1
–
w8分别表示初始空间横轴、空间纵轴偏移分别为2、3、4、5、6、7、8、9km。
具体实施方式
[0033]
下面结合附图和具体实施例对本发明做进一步说明。
[0034]
如图1所示为本发明实施例结合全局与局域信息的遥感地表温度年变化模拟方法的流程图。
[0035]
本发明方法包括以下步骤：
[0036]
步骤1、通过下式计算地表温度年周期变化
[0037]
lst_global(d)＝acp5(d)
[0038]
＝mask yast1·
sin(k1(d θ)) yast2·
sin(k2(d θ))
[0039]
lst_local(d)分别表示相对春分点第d天地表温度的年周期变化，mask表示年均地表温度，yast1和yast2分别表示年内地表温度的两次振幅，d表示相对于春分点的天数，θ表示相对春分点的相位移，k1和k2是两个常量。本步骤中，k1＝2π/365，k2＝4π/365。
[0040]
步骤2、以缺失像元为中心，建立动态滑动时空窗口，并在动态滑动时空窗口中实现地表温度逐日波动的重建，具体步骤如下：
[0041]
步骤2.1、令x和y分别表示空间维度上横轴和纵轴距离，令d为时间维度上的距离，假设整个研究为a，可以表示为：
[0042][0043]
其中z
x
，zy和zd分别表示研究区域a在空间横轴、空间纵轴和时间维度上覆盖的范围，x，y和d分别表示研究区域a在三个维度上的最大值；在研究区域a内，任一像元点可以表示为p＝p[x,y,d]，其中x和y分别表示p在研究区域内的横坐标和纵坐标，d表示该天相对于春分点的偏移天数；
[0044]
步骤2.2、令某地表温度缺失像元为p0，则p0＝p0[x0,y0,d0]；以缺失像元p0为中心，建立初始时空窗口，初始时空窗口的空间横轴偏移量为α0，空间纵轴偏移量为β0，时间维度偏移量为γ0，初始时空窗口b1表示为：
[0045]
b1＝{p[x,y,z]∈a:x
0-a0≤x≤x0 a0,y
0-β0≤y≤y0 β0,d
0-γ0≤d≤d0 γ0}
[0046]
步骤2.3、设置初始空间横轴偏移量α0、初始空间纵轴的偏移量β0和初始时间偏移量γ0，当初始时空窗口内的已知地表温度数量小于阈值ω，则以步调为i(本实施例中，扩大步调i＝1)扩大时空窗口，直到时空窗口内的已知地表温度数量不小于阈值ω，时空窗口
b由下式：
[0047][0048]
为探究初始时空窗口的大小，本发明使用交叉验证法，对最佳初始时空窗口的选择进行探究(如附图3所示)通过探究我们的初始空间横轴(α0)和初始空间纵轴的偏移量(β0)均设置为4，初始时间偏移量(γ0)设置为0。此外，初始时空窗口内已知地表温度的数量和分布也会影响模型的精度。因此本发明设置一个阈值ω＝10，当初始时空窗口内的已知地表温度小于ω时，同时扩大时空窗口。
[0049]
步骤2.4、当时空窗口b内的已知地表温度数量大于阈值时，则通过随机森林算法建立地表温度与相关辅助参数的统计回归关系，进而重建缺失的地表温度逐日波动。相关辅助参数包括空气温度逐日波动sat_local(d)、归一化植被指数(ndvi)和数字高程模型(dem)。相关辅助参数需重采样至与地表温度同一空间分辨率。
[0050]
鉴于地表温度和空气温度的逐日波动具有较高的相似性，且地表温度的变化还受地表覆盖类型的变化和高程差异的影响，因此本发明选取是空气温度逐日波动、归一化植被指数以及高程。值得注意的是，当地表温度数据缺失时，该像元该天的地表温度逐日波动缺失。因此，重建缺失的地表温度逐日波动表示为：
[0051]
lst_local(d)＝f
rf
(sat_local(d),ndvi,dem)
[0052]
其中，空气温度逐日波动sat_local(d)通过下式计算：
[0053]
sat_local(d)＝sat(d)-sat_global(d)
[0054]
＝sat(d)-{mask1 yast3·
sin(k1(d θ1)) yast4·
sin(k2(d θ1))}
[0055]
式中，sat_local(d)表示相对春分点第d天的空气温度逐日波动，sat(d)表示相对春分点第d天的观测的空气温度，sat_global(d)表示通过五参数地表温度年变化模型模拟的相对春分点第d天的空气温度年内变化，参数mask1、yast3、yast4、θ1通过各个像元年内大于四个的观测地表温度进行最小二乘法拟合获得。
[0056]
步骤4、本发明通过五参数地表温度年变化模型(acp5模型)模拟地表温度的年周期变化。由于模型的简洁性以及对太阳辐射周期性的充分考虑，五参数地表温度年变化模型可以较好模拟地表温度年周期变化。
[0057]
通过下式进行对目标区域的地表温度重建。
[0058]
lst(d)＝lst_global(d) lst_local(d)
[0059]
式中，lst(d)表示相对春分点第d天的重建地表温度，lst_global(d)和lst_local(d)分别表示相对春分点第d天地表温度的年周期变化和逐日波动。
[0060]
通过本发明重建目标区域的地表温度后，随机其中一个像元点，展示本发明重建结果与原始地表温度图(如图2所示)。
[0061]
除上述实施例外，本发明还可以有其他实施方式。凡采用等同替换或等效变换形成的技术方案，均落在本发明要求的保护范围。