CN111563228A - 基于地表入射短波辐射的山地地表反射率地形改正方法 - Google Patents

Abstract

本发明公开了基于地表入射短波辐射的山地地表反射率地形改正方法，包括粗分辨率静止气象卫星数据与高分辨率遥感数据空间匹配；MSG地表入射短波辐射产品空间降尺度；高分辨率遥感数据分层；不同数据层地表反射率数据地形效应改正；不同数据层地表反射率结果综合与最终改正结果生成。本发明通过耦合静止气象卫星数据反演获得粗分辨率地表入射短波辐射数据，通过空间降尺度获得高分辨率地表入射短波辐射数据，然后再对高分辨率卫星地表反射率进行地形辐射改正，通过同步考虑直射辐射和散射辐射，消除了传统方法在直射辐射较低或者缺失地区的过校正或者欠校正问题，为山地地表反射率数据有效开展山地环境监测提供了重要的前期数据处理方法基础。

Description

基于地表入射短波辐射的山地地表反射率地形改正方法

技术领域

本发明涉及基于地表入射短波辐射的山地地表反射率地形改正方法。

背景技术

地表反射率表征地球表面对太阳辐射的反射能力，是揭示地物目标本质的最有用信息，直接关系诸如地表植被指数，叶面积指数，FPAR指数等各种地表生物物理量的估计，对定量遥感研究有着重要意义。

在山地光学遥感影像中，由于地形的存在改变了太阳—地表目标—传感器之间的几何位置关系，使遥感观测的多光谱、多角度、多时空尺度的倾斜地表辐射信息与平坦地表辐射不一致，形成地形影响效应。地形起伏导致不同朝向的坡面接收的太阳辐射不同，致使相同的地物在向阳面与阴面有较大的光谱差异，造成部分像元处于阴影覆盖下，而另一些却处于过度感光状态，导致了“同物异谱”和“同谱异物”现象，直接影响地表参量提取精度。地形校正属于遥感影像辐射和大气校正的范畴，是准确获得地表反射率必不可少的一步，特别是在丘陵、山区等地形起伏区域的遥感监测与应用，消除地形因素的影响是非常必要的。通过地形校正消除地形影响，能够提高地表覆盖分类精度，从而最大程度地提高山区遥感监测的精度。

就地形校正原理而言，其主要是将坡面像元的辐亮度转换为水平像元辐亮度，使两个反射特性相同的地物在同一景遥感图像中具有相同或相近的亮度值，从而减缓或消除地形引起的辐射失真。

随着对遥感机理认知的不断深入，学者们开始尝试从机理角度对太阳辐射的不同分量进行校正，基于辐射传输模型的地形校正方法则综合考虑了太阳指数辐射、大气散射辐射以及周围地形反射辐射，相比此前经验、半经验模型能够取得较好的校正效果。然而，这些模型大多假设地表是各向同性的，即以地表朗伯体特性假设为基础，不符合地表真实情况。为了更加准确的描述山地辐射传输过程，闻建光等结合辐射传输方程，发展出了物理机制更为完善的基于方向反射的山地地表反射率计算模型，并构建了水平地表和坡面地表的反射率转换关系，为地形校正提供了新的方式。然后，其校正结果发现在一些山脊和山谷地区仍然存在一些欠改正的个。

从目前的研究进展来看，虽然目前国内外针对山地光学遥感影像的地形辐射校正问题提出了从简单的经验模型到复杂的大气-BRDF-地形辐射协同校正模型，并在区域尺度得到了很好的示范应用，但也可以看出，当前无论是地形校正模型的求解方法还是在更大尺度的应用研究还远远不够，而且目前没有非常通用的校正模型。

⑴目前现有技术中的C校正模型技术方案为:

C校正模型为基于太阳—地表—传感器的校正模型。C校正是一种基于影像像元和太阳相对入射角余弦值之间线性关系的校正方法，并利用两者之间回归方程的截距和斜率的比值生成新的参数。与其他模型相比，其应用相对简单，且地形校正效果较好。在C校正方法中，仅有一个待确定的变量，即参数C。

C校正模型是目前应用相对较为广泛的地形辐射校正模型。为了修正余弦校正(一个简单的光学函数，校正后像素接受的总辐射与坡面像素接受的总辐射有一个由入射角余弦决定的直接比例关系)中存在的问题，提出将一个半经验系数C应用到余弦校正中，通过研究崎岖山地的影像数据处理发现，每一个波段的反射率ρ_T和光照系数cosγ_i之间存在一个线性相关，即：

ρ_T＝α_λ+b_λcosγ_i,

式中，α_λ、b_λ分别为λ波段的经验系数。

将这个线性相关应用到余弦校正中，得出了C校正模型为

ρ_H＝ρ_T(cosθ_Z+Cλ)/(cosγ_i+C_λ)，

式中，ρ_H为水平地表反射率；θ_Z为太阳天顶角；C_λ＝α_λ/b_λ。

C校正模型在一定程度上避免了余弦校正模型在低光照参数地区的过度校正问题，特别是对于地形背光效果更好。参数C的值也可以通过经验统计回归线的斜率和截距获得，该参数通过提高分母的值在C校正中发挥一个调节作用，以减弱光照像素区域的过度校正问题。

现有的目前C校正模型存在的缺陷为：

①通常将统一的参数C应用于整幅遥感影像，而在一些复杂地形山区，由于坡度、NDVI和地类等变异较大，这种统一参数的C校正方法往往会导致过校正、欠校正等现象。为有效避免整体C校正中通常存在的过校正、欠校正等问题，区分不同的坡度、NDVI区间和不同土地利用/土地覆盖类型分别计算时会花费大量的人力和时间。

②遥感影像的各个波段适应的地形校正策略不同，如LandsatTM影像回归分析表明不同的波段需要采取不同的C校正策略。选择正确的策略需要经过不断的实验，则耗时较久。

③校正参数C无法明确解释天空散射辐射和邻近地形的反射辐射两个辐射分量，对于显著地形效应的陡峭山区，还是难于抑制阴影区域过度校正现象。在中、高山峡谷等复杂地形山区中，因地形起伏在不同坡度、坡向所形成的照度差异对影像光谱信息的影响定量评价和消除则有待进一步研究。

⑵目前现有技术中的SCS+C校正模型技术方案为:

SCS模型在植被覆盖区获得了比余弦校正更好的效果。前面的余弦校正模型和C校正模型都是基于太阳—地表—传感器三者的几何关系来考虑校正方法的，可以统称为STS校正模型。SCS校正模型主要是基于太阳—冠层—传感器三者的几何关系来考虑问题，由于树木的生长是向地性的(垂直于大地水准面)，地形不能控制太阳和树木之间的几何关系，地形影响的只是树木相对地表的位置关系。SCS校正模型是基于冠层的校正模型，与STS模型的着眼点不同，它使光照冠层在从坡面到水平面校正的过程中，光照方向的变化更符合实际情况。假定来自光照冠层的反射辐射因树木的向地生长特性而大大独立于地形，光照冠层的总体反射率与其范围成正比。

在SCS模型中，太阳和冠层之间的几何关系在校正前后保持不变，因此更加符合实际情况，适合于对森林区域的地形辐射校正。然而，SCS校正模型由于散射辐射的影响被忽视了，使得背光区域的坡面存在过度校正问题，于是迫切需要对SCS校正模型进行修正。

SCS校正出现过度校正问题的原因与余弦校正相似、当入射角接近90时，校正系数变得很大而引起过度校正。在C校正中，参数C应用到余弦校正中，仿效天空散射辐射的影响，具有调节过度校正的作用。于是，Scott引入C系数来改进SCS校正模型,提出了SCS+C校正模型，即

ρ_H＝ρ_T(cosθ_Pcosθ_Z+C_λ)/(cosγ_i+C_λ)，

式中，θ_P为坡度角。选择参数C主要是由于它在改进余弦校正模型时有效,同时计算也比较简单。

SCS+C模型具有较优的校正效果，最大程度地改善了影像质量，而且符合实际情况，因而更有发展和应用优势。SCS+C校正在SCS校正框架上引入C系数，减弱了过度校正问题，同是又保留了SCS校正的优点。概括来说，SCS提供了一种更加符合地表实际的地形辐射校正原理，而SCS+C是在这个原理基础上获得的一种具有良好地形辐射校正效果的实践方案。

现有的目前SCS+C校正模型技术存在的缺陷为：

①在SCS+C校正模型中，虽然其它用来说明大气散射辐射和地形散射辐射的参数可能会提供更加精确的结果，但计算量很大且计算时要求输入参数的获取很难实现，与引入的C系数相比没有优势。

②除了小区域的经验统计模型以外，理论上所有的校正模型在校正坡度大于20的陡峭山地影像时都没有很好效果。因为，地形坡度越大的区域，校正方程的分母越小，校正参数值就会越大，从而产生过度校正问题。

综上，目前采用较为广泛的仅考虑太阳直接辐射的经验、半经验校正模型，其适用性在低太阳高度角或者阴影情况下非常低，如何有效提高现有方法的改正效果显得十分重要和急迫，对提高山地遥感应用水平具有重要意义。

发明内容

本发明要解决的技术问题是克服现有的缺陷，提供基于地表入射短波辐射的山地地表反射率地形改正方法，可以有效解决背景技术中的问题。

为了解决上述技术问题，本发明提供了如下的技术方案：

本发明提供基于地表入射短波辐射的山地地表反射率地形改正方法，包括以下步骤：

S1：粗分辨率静止气象卫星数据与高分辨率遥感数据空间匹配；

S2：MSG地表入射短波辐射产品空间降尺度；

S3：高分辨率遥感数据分层；

S4：不同数据层地表反射率数据地形效应改正；

S5：不同数据层地表反射率结果综合与最终改正结果生成。

作为优选方案，为完成不同空间分辨率遥感数据之间的空间尺度转换，步骤S1包括：

通过各像元点自身的经纬度，针对高分影像各像元与粗分影像各像元之间计算距离，并求解出最近像元，即找到了高分像元所对应的粗分像元，进而建立各个像元之间的映射关系，为后面的空间降尺度提供对应关系。

作为优选方案，步骤S2包括：

地表入射短波辐射结构，地表接收的入射短波辐射包括太阳直射辐射、大气散射辐射和邻近地形反射辐射，可以表示为：

R_s＝R_b+R_d+R_r；

地表入射短波辐射空间降尺度，将静止气象卫星数据高时间分辨率反演的地表入射短波辐射降尺度至高分辨率卫星像元尺度。

作为优选方案，所述太阳直射辐射为：

与水平面上类似，入射到起伏地表面上的太阳直接辐射可以表示为：

其中，E₀是太阳常数，由于太阳常数是日地平均距离处给定的值，

当要具体到某一日期时，需要进行订正，当天日地距离订正系数为：

式中，DOY是积日，τ_b是直接辐射大气透过率，θ为太阳直线光线与地表坡面法线的夹角即坡面上的实际太阳入射角，计算公式如下：

cosθ＝cos Z_s cos S+cos Z_s sin S cos(A_s-A)，

其中，Z_s为太阳天顶角，A_s为太阳方位角；S为坡度，A为坡向，实际地形坡度和坡向可利用数字高程模型计算。

作为优选方案，所述天空散射辐射为：

复杂地形下天空散射辐射R_d，是对平坦地形下天空散射辐射R_d,flat，利用天空可视因子SVF进行订正获得：

I_d＝I_d,flat×SVF

I_d,flat＝E₀×dr×τ_d×cos(Z_s)，

其中，τ_b为散射辐射透过率，在晴朗无云条件下天空的散射是一个均质散射，在直接辐射和散射辐射之间存在线性关系：

τ_d＝0.271-0.294×τ_b，

其中，天空可视因子定义为把半球2π空间划分为n等份，目标点上方半球可见部分面积与半球面积之比，在每一片“天空”上沿光线在xy平面上的投影方向追踪，依次计算此方向上每一坡元与起点坡元的高度角，找出最大高度角，记为h_i，每点的SVF可由下式计算：

其中，n取16个方向，搜索半径R可根据实际情况确定。

作为优选方案，所述邻近地形反射辐射为：

对于周围坡面反射产生的附加辐射，可采用Dozier简化的近似计算方法，只考虑地形坡度、天空视角因子以及周围地形的平均反射作用：

Ir＝ρ×C_t×(I_b+I_d)，

其中，C_t为地形结构参数，对于周围地形贡献的反射辐射I_r，C_t包括了反射辐射的各向异性特性，并且也包括坡元与周围可见坡元之间的几何效应，假设下垫面为朗伯体，C_t按下式计算：

C_t＝(1+cosS)/2-SVF，

其中，S是坡度，SVF是天空可视因子，ρ是邻近地形的平均反照率。

作为优选方案，所述各向异性特性为周围坡面有近有远，且坡度、坡向各异；所述坡元为具有一定坡向、坡度的坡面网格单元。

作为优选方案，为了更为准确地获取地表入射短波辐射，步骤S2为：

将静止气象卫星数据高时间分辨率反演的地表入射短波辐射降尺度至高分辨率卫星像元尺度。

作为优选方案，在降尺度过程中，其中一个非常重要的前提即为粗分辨率静止气象卫星数据像元的入射短波辐射估算值是其所对应的高分辨率卫星影像中所有像元的平均值，

因此，粗分辨率入射短波辐射和高分辨率入射短波辐射之间存在如下关系式：

式中，DSR_c为粗分辨率地表入射短波辐射，DSR_h为高分辨率地表入射短波辐射，i代表粗分影像像元所包含的高分影像像元i，对具体每个高分像元i而言，其入射短波辐射DSR_h,i可以表示为如下公式：

式中，在假定已知各高分像元直射和散射辐射透过率的前提下，高分像元DSR可以估算得到。

作为优选方案，考虑到不同地表覆被类型差异引起的地表光谱特征差异及其与地表太阳入射条件响应程度的差异性，

实施步骤S3至S5时，采用归一化植被指数为数据分层的关键指标，按照区段0.25为区间，将待校正影像数据分层4个部分，即<0.25，0.25-0.5，0.5-0.75，和>0.75，然后针对各个数据层，分别采用经验旋转的地形效应改正模型，即如下公式：

SR_C(λ)＝SR(λ)-α(λ)·(DSR-DSR_H)，

式中，SRc为改正后的地表反射率，SR为改正前地表反射率，DSR为地表入射短波辐射，DSRH为平坦地表入射短波辐射，系数α为改正系数；

系数α的估算采用地表反射率与地表入射短波辐射的经验线性拟合获得，即通过线性拟合获取原始影像中地表反射率随着地表入射条件差异带来的地表反射率差异，具体公式如下：

SR(λ)＝α(λ)·DSR+β(λ)，

采用上述改正模型分别应用于每一层数据，获得各层地形效应改正结果，然后集成各层改正结果生成最终的地表反射率地形效应改正结果。

本发明中提供的技术方案，至少具有如下技术效果或者优点：

1、为提高当前山地地表反射率地形效应改正效果，减少目前常用仅考虑太阳直射辐射的地形效应改正方法中在局地低太阳高度角入射条件或者阴影条件下地表反射率过校正或者欠校正问题，本发明提出了基于地表入射短波辐射(同时涵盖太阳直射辐射和散射辐射)的地形效应改正方法。

2、该方法通过耦合静止气象卫星数据反演获得粗分辨率地表入射短波辐射数据，通过空间降尺度获得高分辨率地表入射短波辐射数据，然后再对高分辨率卫星地表反射率进行地形辐射改正。

3、该方法通过同步考虑直射辐射和散射辐射，消除了传统方法在直射辐射较低或者缺失地区的过校正或者欠校正问题，为山地地表反射率数据有效开展山地环境监测提供了重要的前期数据处理方法基础。

附图说明

附图用来提供对本发明的进一步理解，并且构成说明书的一部分，与本发明的实施例一起用于解释本发明，并不构成对本发明的限制。

在附图中：

图1是本发明实施例中基于地表入射短波辐射的山地地表反射率地形改正方法技术流程图。

图2是本发明实施例中复杂地形上地表太阳入射辐射示意图。

图3是本发明目标点天空可视区域示意图。

具体实施方式

为了更好地理解上述技术方案，下面将结合说明书附图以及具体实施方式对上述技术方案进行详细的说明。

实施例：

请参阅图1-3，本实施例提供基于地表入射短波辐射的山地地表反射率地形改正方法，包括以下步骤：

(1)粗分辨率静止气象卫星数据与高分辨率遥感数据空间匹配：

为完成不同空间分辨率遥感数据之间的空间尺度转换，首先需对本发明所采用的的粗分辨率静止气象卫星数据和高分辨率遥感数据之间开展空间匹配，其主要完成过程为通过各像元点自身的经纬度，针对高分影像各像元与粗分影像各像元之间计算距离，并求解出最近像元，即找到了高分像元所对应的粗分像元，进而建立各个像元之间的映射关系，为后面的空间降尺度提供对应关系；

(2)MSG地表入射短波辐射产品空间降尺度：

1)地表入射短波辐射结构组成

地表接收的入射短波辐射(R_s)一般包括太阳直射辐射(R_b)、大气散射辐射(R_d)和邻近地形反射辐射(R_r)，可以表示为:

R_s＝R_b+R_d+R_r (1)

其中，太阳直射辐射为：

与水平面上类似，入射到起伏地表面上的太阳直接辐射可以表示为：

其中E₀是太阳常数，由于太阳常数是日地平均距离处给定的值，当要具体到某一日期时，需要进行订正，当天日地距离订正系数为：

式中DOY是积日；公式(6)中τ_b是直接辐射大气透过率，θ为太阳直线光线与地表坡面法线的夹角即坡面上的实际太阳入射角，计算公式如下：

cosθ＝cos Z_s cos S+cos Z_s sin S cos(A_s-A) (4)

其中Z_s为太阳天顶角；A_s为太阳方位角；S为坡度；A为坡向；实际地形坡度和坡向可以利用数字高程模型计算；

其中，天空散射辐射为：

复杂地形下天空散射辐射R_d，是对平坦地形下天空散射辐射R_d,flat，利用天空可视因子SVF进行订正获得：

I_d＝I_d,flat×SVF

I_d,flat＝E₀×dr×τ_d×cos(Z_s) (5)

其中τ_d为散射辐射透过率，在晴朗无云条件下天空的散射是一个均质散射，在直接辐射和散射辐射之间存在线性关系：

τ_d＝0.271-0.294×τ_b (6)

其中天空可视因子定义为把半球2π空间划分为n等份，目标点上方半球可见部分面积与半球面积之比，具体情况图3所示，在每一片“天空”上沿光线在xy平面上的投影方向追踪，依次计算此方向上每一坡元与起点坡元的高度角，找出最大高度角，记为h_i，每点的SVF可由下式计算：

本研究中n取16个方向，搜索半径R可根据实际情况确定。

其中，邻近地形反射辐射为：

对于周围坡面反射产生的附加辐射，可采用Dozier简化的近似计算方法，只考虑地形坡度、天空视角因子以及周围地形的平均反射作用：

Ir＝ρ×C_t×(I_b+I_d) (8)

其中C_t为地形结构参数。对于周围地形贡献的反射辐射I_r，C_t包括了反射辐射的各向异性特性(周围坡面有近有远，且坡度、坡向各异)，并且也包括坡元(具有一定坡向、坡度的坡面网格单元)与周围可见坡元之间的几何效应(因为对某一坡元的反射辐射有影响的只是它周围的可见坡元)，假设下垫面为朗伯体，C_t按下式计算：

C_t＝(1+cos S)/2-SVF (9)

其中S是坡度，SVF是天空可视因子，ρ是邻近地形的平均反照率；

以上三部分内容简单介绍了地表入射太阳辐射的主要组成及其概念，通过地面观测数据大小可以发现，地表入射太阳辐射中主要由直射辐射和散射辐射，邻近地形反射辐射所占比重很小。因此，在本发明中，邻近地形反射辐射项暂不考虑。

(2)地表入射短波辐射空间降尺度；

目前，高分辨率地表入射短波辐射产品还没有，虽然采用高分辨率卫星数据一定程度上可以反演获得地表短波辐射，但是地表反射率地形效应的影响往往导致地表辐射估算存在较大的不确定性；

为了更为准确地获取地表入射短波辐射，本发明采用空间降尺度的方式，将静止气象卫星数据高时间分辨率反演的地表入射短波辐射降尺度至高分辨率卫星像元尺度；

在降尺度过程中，其中一个非常重要的前提即为粗分辨率静止气象卫星数据像元的入射短波辐射估算值是其所对应的高分辨率卫星影像中所有像元的平均值，因此，粗分辨率入射短波辐射和高分辨率入射短波辐射之间存在如下关系式：

式中，DSR_c为粗分辨率地表入射短波辐射，DSR_h为高分辨率地表入射短波辐射，i代表粗分影像像元所包含的高分影像像元i。对具体每个高分像元i而言，其入射短波辐射DSR_h,i可以表示为如下公式：

式中，粗分像元DSR已知，在假定已知各高分像元直射和散射辐射透过率的前提下，高分像元DSR可以估算得到，为此，本发明采用经验性的直射辐射估算法，如下式所示：

τ_b＝exp(-0.8662·T_LK·m_a·δ_R) (12)

式中，TLK为大气浑浊因子，ma为空气质量比(Air Mass Ratio)，即太阳辐射方向经过的大气质量与太阳天顶直射通过大气质量的比例，可根据以下公式计算得到：

其中

为经过太阳折射纠正后的太阳高度角，可表示为：

另外，公式15中p/p₀可写为高度的函数：

p/p₀＝exp(-z/8434.5) (15)

式中，z为海拔。

与此同时，散射辐射透过率可以采用它与直射辐射之间的经验关系式进行计算，即：

τ_d＝c₁+c₂·τ_b (16)

式中，c1为0.2710，c2为-0.2939。

由于粗分像元地表短波辐射估算未考虑地形影响，因此，基于以上经验关系式，高分辨遥感影像像元的地表入射短波辐射的估算结果仅为高分像元水平条件下的地表入射短波辐射。为进一步获得高分像元不同地形条件下地表入射短波辐射，针对水平估算结果进行地形改正是重要的一步。

依据估算的高分像元水平条件下的地表入射短波辐射，可以反向估算高分像元的直射和散射辐射的透过率，即如下公式所示：

进而获得高分像元水平条件下直射(E_bh)和散射(E_dh)辐射各分量：

最终，结合地形数据(坡度、坡向、天空可视因子和遮蔽因子)，根据以下公式对水平直射和散射进行地形改正，最终得到高分像元真实地形条件下地表入射短波辐射

其中，

(3)地表反射率分层地形效应改正建模：

在获得山地高分像元真实地形环境下地表入射短波辐射的前提下，本发明进而开展地形效应改正建模。此处，考虑到不同地表覆被类型差异引起的地表光谱特征差异及其与地表太阳入射条件响应程度的差异性，

本发明首先对待校正图像进行了数据分层处理。此处，本发明选择采用归一化植被指数(NDVI)为数据分层的关键指标，按照区段0.25为区间，将待校正影像数据分层4个部分，即<0.25，0.25-0.5，0.5-0.75，和>0.75。然后针对各个数据层，分别采用经验旋转的地形效应改正模型，即如下公式：

SR_C(λ)＝SR(λ)-α(λ)·(DSR-DSR_H) (20)

式中，SRc为改正后的地表反射率，SR为改正前地表反射率，DSR为地表入射短波辐射，DSRH为平坦地表入射短波辐射，系数α为改正系数。从公式结构来看，该改正模型的原理即将地表入射短波辐射差异引起的地表反射率变化归一化到同样的入射条件。相比于此前基于太阳直射辐射的改正方法，本改正方法中各像元的地表入射短波辐射都有有效值。因此，公式(20)的改正结果不会出现严重的过校正或者欠校正的情况。其中，找到合适的α值是实现有效改正的关键。

在本发明中，系数α的估算我们采用地表反射率与地表入射短波辐射的经验线性拟合获得，即通过线性拟合获取原始影像中地表反射率随着地表入射条件差异带来的地表反射率差异，具体公式如下：

SR(λ)＝α(λ)·DSR+β(λ) (21)

最终，采用上述改正模型分别应用于每一层数据，获得各层地形效应改正结果，然后集成各层改正结果生成最终的地表反射率地形效应改正结果。

在太阳高度角较小的月份，主要是冬季，地形起伏引起的地形阴影导致研究区的地表反射率呈现非常显著的地形差异，地形效应极为显著。

最后应说明的是：以上所述仅为本发明的优选实施例而已，并不用于限制本发明，尽管参照前述实施例对本发明进行了详细的说明，对于本领域的技术人员来说，其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (10)

1.基于地表入射短波辐射的山地地表反射率地形改正方法，其特征在于：包括以下步骤：

S1：粗分辨率静止气象卫星数据与高分辨率遥感数据空间匹配；

S2：MSG地表入射短波辐射产品空间降尺度；

S3：高分辨率遥感数据分层；

S4：不同数据层地表反射率数据地形效应改正；

S5：不同数据层地表反射率结果综合与最终改正结果生成。

2.根据权利要求1所述的基于地表入射短波辐射的山地地表反射率地形改正方法，其特征在于：为完成不同空间分辨率遥感数据之间的空间尺度转换，步骤S1包括：

通过各像元点自身的经纬度，针对高分影像各像元与粗分影像各像元之间计算距离，并求解出最近像元，即找到了高分像元所对应的粗分像元，进而建立各个像元之间的映射关系，为后面的空间降尺度提供对应关系。

3.根据权利要求2所述的基于地表入射短波辐射的山地地表反射率地形改正方法，其特征在于：步骤S2包括：

地表入射短波辐射结构，地表接收的入射短波辐射包括太阳直射辐射、大气散射辐射和邻近地形反射辐射，可以表示为：

R_s＝R_b+R_d+R_r；

地表入射短波辐射空间降尺度，将静止气象卫星数据高时间分辨率反演的地表入射短波辐射降尺度至高分辨率卫星像元尺度。

4.根据权利要求3所述的基于地表入射短波辐射的山地地表反射率地形改正方法，其特征在于：所述太阳直射辐射为：

与水平面上类似，入射到起伏地表面上的太阳直接辐射可以表示为：

其中，E₀是太阳常数，由于太阳常数是日地平均距离处给定的值，当要具体到某一日期时，需要进行订正，当天日地距离订正系数为：

式中，DOY是积日，τ_b是直接辐射大气透过率，θ为太阳直线光线与地表坡面法线的夹角即坡面上的实际太阳入射角，计算公式如下：

cosθ＝cosZ_scosS+cosZ_ssinScos(A_s-A)，

其中，Z_s为太阳天顶角，A_s为太阳方位角；S为坡度，A为坡向，实际地形坡度和坡向可利用数字高程模型计算。

5.根据权利要求4所述的基于地表入射短波辐射的山地地表反射率地形改正方法，其特征在于：所述天空散射辐射为：

复杂地形下天空散射辐射R_d，是对平坦地形下天空散射辐射R_d,flat，利用天空可视因子SVF进行订正获得：

I_d＝I_d,flat×SVF

I_d,flat＝E₀×dr×τ_d×cos(Z_s)，

其中，τ_b为散射辐射透过率，在晴朗无云条件下天空的散射是一个均质散射，在直接辐射和散射辐射之间存在线性关系：

τ_d＝0.271-0.294×τ_b，

其中，天空可视因子定义为把半球2π空间划分为n等份，目标点上方半球可见部分面积与半球面积之比，在每一片“天空”上沿光线在xy平面上的投影方向追踪，依次计算此方向上每一坡元与起点坡元的高度角，找出最大高度角，记为h_i，每点的SVF可由下式计算：

其中，n取16个方向，搜索半径R可根据实际情况确定。

6.根据权利要求5所述的基于地表入射短波辐射的山地地表反射率地形改正方法，其特征在于：所述邻近地形反射辐射为：

对于周围坡面反射产生的附加辐射，可采用Dozier简化的近似计算方法，只考虑地形坡度、天空视角因子以及周围地形的平均反射作用：

Ir＝ρ×C_t×(I_b+I_d)，

其中，C_t为地形结构参数，对于周围地形贡献的反射辐射I_r，C_t包括了反射辐射的各向异性特性，并且也包括坡元与周围可见坡元之间的几何效应，假设下垫面为朗伯体，C_t按下式计算：

C_t＝(1+cosS)/2-SVF，

其中，S是坡度，SVF是天空可视因子，ρ是邻近地形的平均反照率。

7.根据权利要求5所述的基于地表入射短波辐射的山地地表反射率地形改正方法，其特征在于：所述各向异性特性为周围坡面有近有远，且坡度、坡向各异；所述坡元为具有一定坡向、坡度的坡面网格单元。

8.根据权利要求7所述的基于地表入射短波辐射的山地地表反射率地形改正方法，其特征在于：为了更为准确地获取地表入射短波辐射，步骤S2为：

将静止气象卫星数据高时间分辨率反演的地表入射短波辐射降尺度至高分辨率卫星像元尺度。

9.根据权利要求8所述的基于地表入射短波辐射的山地地表反射率地形改正方法，其特征在于：在降尺度过程中，其中一个非常重要的前提即为粗分辨率静止气象卫星数据像元的入射短波辐射估算值是其所对应的高分辨率卫星影像中所有像元的平均值，

因此，粗分辨率入射短波辐射和高分辨率入射短波辐射之间存在如下关系式：

式中，DSR_c为粗分辨率地表入射短波辐射，DSR_h为高分辨率地表入射短波辐射，i代表粗分影像像元所包含的高分影像像元i，对具体每个高分像元i而言，其入射短波辐射DSR_h,i可以表示为如下公式：

式中，在假定已知各高分像元直射和散射辐射透过率的前提下，高分像元DSR可以估算得到。

10.根据权利要求9所述的基于地表入射短波辐射的山地地表反射率地形改正方法，其特征在于：考虑到不同地表覆被类型差异引起的地表光谱特征差异及其与地表太阳入射条件响应程度的差异性，

实施步骤S3至S5时，采用归一化植被指数为数据分层的关键指标，按照区段0.25为区间，将待校正影像数据分层4个部分，即<0.25，0.25-0.5，0.5-0.75，和>0.75，然后针对各个数据层，分别采用经验旋转的地形效应改正模型，即如下公式：

SR_C(λ)＝SR(λ)-α(λ)·(DSR-DSR_H)，

式中，SRc为改正后的地表反射率，SR为改正前地表反射率，DSR为地表入射短波辐射，DSRH为平坦地表入射短波辐射，系数α为改正系数；

系数α的估算采用地表反射率与地表入射短波辐射的经验线性拟合获得，即通过线性拟合获取原始影像中地表反射率随着地表入射条件差异带来的地表反射率差异，具体公式如下：

SR(λ)＝α(λ)·DSR+β(λ)，

采用上述改正模型分别应用于每一层数据，获得各层地形效应改正结果，然后集成各层改正结果生成最终的地表反射率地形效应改正结果。

Images