本发明属于地表温度遥感反演的技术领域，特别是涉及一种有效的modis地表温度角度校正方法。
背景技术：
地表温度是表征地表过程变化的重要的特征物理量，是地表—大气能量交换的直接驱动因子，综合了地表与大气的相互作用以及大气和陆地之间能量交换的结果。精确的地表温度估算能够提供地表能量平衡状态的时空变化信息，不仅有助于评估地表能量与水文平衡、热惯量和土壤湿度，而且有助于获取全球表面温度及长期变化情况，因此，地表温度被广泛应用于地表能量平衡、气候变化和资源环境监测等研究领域。而由于地表温度的空间异质性，要获取区域和全球尺度上地表温度的时空分布，常规的地面定点观测难以实现，卫星遥感是非常有效的手段。
目前通过遥感手段获取地表温度时存在的一个主要问题是：通常将地表作为均质体，在反演过程中往往忽略地表温度本身的方向性。对于均一地表，这一假设是合理的；但对于非均匀地表覆盖的混合像元，其内部具有复杂组合结构且存在不同的组分，像元内部各组分在不同观测角度下的比例不同，同时各个组分间的多次散射也使得像元发射率依赖于观测方向，而导致同温像元的热辐射也表现出一定的角度效应。对于大角度观测传感器(如modis、avhrr)来说，不同像元对应不同的观测角度，温度的角度效应将导致观测数据反演的同一景影像不同像元的温度缺乏可比性，从而影响遥感获取的植被和土壤温度的使用。
技术实现要素：
为了解决上述技术问题中的不足之处，本发明提供了一种有效的modis地表温度角度校正方法。
为了解决上述技术问题，本发明采用的技术方案是：一种有效的modis地表温度角度校正方法，包括以下步骤：
步骤1：初步构建输入数据集；
步骤2：选择基于端元信息模型，对未校正的混合像元温度进行分解，得到植被和土壤组分的温度；
步骤3：确定反射率计算模型，进行反射率角度校正；
步骤4：在反射率角度校正的基础上，得到角度校正后的植被覆盖率；
步骤5：利用未校正混合像元温度分解得到的植被和土壤温度，结合角度校正后的植被覆盖率，得到角度校正后混合像元的地表温度。
所述的地表温度角度校正方法，其中，所述步骤1中，所述数据集包括气象数据(大气温度ta、饱和水汽压、入射太阳短波辐射rg、下行长波辐射lsky)和遥感数据(modis反射率数据mod09ga、brdf数据mcd43d、地表温度数据mod11_l2)。
所述的地表温度角度校正方法，其中，所述步骤2中，所述基于端元信息模型，根据土壤和植被含水量的极限值，模型包括四个极限端元：1)干燥裸土端元；2)干燥植被端元；3)湿润裸土端元；4)湿润植被端元。
所述的地表温度角度校正方法，其中，所述步骤2中，所述对未校正的混合像元温度进行分解，得到植被和土壤组分的温度，具体步骤为：
步骤2011：进行土壤-植被系统能量分解，得到土壤和植被组分的能量；
步骤2012：在得到土壤和植被组分能量的基础上，计算基于端元信息模型四个极限端元的地表温度；
步骤2013：利用极限端元地表温度和未校正的植被覆盖率，计算土壤表层含水量为0，植被开始受水分胁迫像元的临界地表温度；
步骤2014：根据未校正混合像元地表温度与临界地表温度的关系，分解得到混合像元中植被和土壤组分的温度。
所述的地表温度角度校正方法，其中，所述步骤2011中，进行植被和土壤能量分解的过程为：土壤热通量g＝gs；净辐射rn分解为土壤净辐射rn,s和植被净辐射rn,v，公式一：
rn＝rn,s+rn,v
所述的地表温度角度校正方法，其中，所述步骤2012中，在得到土壤和植被组分能量的基础上，计算基于端元信息模型四个极限端元的地表温度tsd，tvd，tsw和tvw，公式二：
其中，ρ为空气密度；cp为定压比热；γ为干湿球常数；δ为饱和水汽压差-温度的斜率；vpd为水气压亏缺；ta为近地表气温；rvw和rvd分别为供水充足和干燥时的冠层阻抗；rav和ras分别为植被和土壤上层的空气动力学阻抗。
所述的地表温度角度校正方法，其中，所述步骤2013中，通过端元地表温度tsd，tvw和未校正的植被覆盖率fv定义土壤表层含水量为0，植被开始受水分胁迫像元的临界地表温度(t*)，公式三：
所述的地表温度角度校正方法，其中，所述步骤2014中，所述根据未校正混合像元地表温度(tr)与临界地表温度(t*)的关系，计算得到混合像元中植被和土壤组分的温度的方法：
a)当tr≤t*时，计算对应的植被和土壤组分的温度，公式四：
b)当tr≥t*时，计算对应的植被和土壤组分的温度，公式五：
所述的地表温度角度校正方法，其中，所述未校正的植被覆盖率fv，计算公式为，公式六：
其中，ndvi为归一化植被指数，ndvimin和ndvimax分别为对应裸土的最小ndvi以及对应全植被覆盖的最大ndvi。ndvi表达式为，公式七：
其中，rnir为近红外反射率，rred为红光反射率，均从mod09ga数据中得到。
所述的地表温度角度校正方法，其中，所述步骤3中，所述反射率计算模型为ross-librdf核函数模型，模型中地表反射率通过三种基本散射类型的核函数之和表示，其中三种类型的散射为：(1)各项同性散射；(2)辐射传输体散射；(3)几何光学散射。模型的表达式为，公式八：
其中，r为反射率；θ为太阳天顶角；
所述的地表温度角度校正方法，其中，所述步骤3中，所述反射率角度校正方法为：通过设定反射率计算模型(公式八)核函数中太阳天顶角θ、观测天顶角
所述的地表温度角度校正方法，其中，所述反射率角度校正方法中，确定角度校正后反射率计算模型中的核函数，公式九：
其中，o为光照与阴影面的重叠；ξ、θ'、
所述的地表温度角度校正方法，其中，所述反射率角度校正方法中，确定核函数的系数：波段λ的核函数系数fiso(λ)，fvol(λ)和fgeo(λ)均通过遥感mcd43d数据的dn值乘以比例因子获得，数据质量通过mcd43d31文件控制。据此，可确定近红外和红外波段的核函数系数，得到角度校正后像元的近红外和红光波段的反射率值(rnir,c和rred,c)。
所述的地表温度角度校正方法，其中，所述步骤4中，所述植被覆盖率的角度校正步骤为：
步骤4011：在步骤3校正得到近红外反射率rnir,c和红光反射率值rred,c的基础上，得到角度校正的归一化植被指数ndvic，公式十：
步骤4012：在步骤4011得到角度校正后的植被指数的基础上，得到角度校正后的植被覆盖率fv,c，公式十一：
所述的地表温度角度校正方法，其中，所述步骤5中，所述进行混合像元地表温度角度校正的方法：在步骤2分解未校正混合像元温度得到植被组分温度(tv)和土壤组分温度(ts)的基础上，结合步骤4得到的角度校正后的植被覆盖率(fv,c)，进行混合像元地表温度的角度校正，得到校正后温度(tr,c)，公式十二：
tr,c＝(1-fv,c)ts+fv,ctv
综上，通过未校正混合像元分解得到的植被和土壤组分温度，结合角度校正后的反射率得到的植被覆盖率，得到了角度校正后混合像元的地表温度。
本发明主要发展了有效的modis地表温度角度校正方法，通过基于端元信息模型，对未校正的混合像元温度进行分解，得到植被和土壤组分的温度；通过对modis数据驱动的反射率计算模型进行角度校正，得到校正后的反射率值，进而得到角度校正后的植被覆盖率；结合分解得到的植被和土壤组分温度与角度校正的植被覆盖率，得到角度校正后混合像元的地表温度值。
本发明具有以下有益效果：能够对卫星传感器通过任意观测角度观测的反射率值校正为以目标角度观测的反射率值，进而可实现植被覆盖率的角度校正；结合分解得到的植被和土壤组分温度，可实现混合像元地表温度的角度校正，使得同一景影像不同像元得到的温度具有可比性。本发明提高了遥感反演地表温度的实用价值，有助于提供更准确的地表能量平衡状态的时空变化信息，对地表能量平衡、气候变化和资源环境监测等的研究具有重要意义。
具体实施方式
本发明提出一种有效的modis地表温度角度校正方法，下面对本发明作进一步说明：
本发明的具体实施方案如下：
步骤1：初步构建输入数据集，包括气象数据(大气温度ta、饱和水汽压、入射太阳短波辐射rg、下行长波辐射lsky)和遥感数据(modis反射率数据mod09ga、brdf数据mcd43d、地表温度数据mod11_l2)的构建。
步骤2：选择基于端元信息模型，对未校正的混合像元温度进行分解，得到植被和土壤组分的温度。根据土壤和植被含水量的极限值，基于端元信息模型包括四个极限端元：1)干燥裸土端元；2)干燥植被端元；3)湿润裸土端元；4)湿润植被端元。对未校正的混合像元温度进行分解，得到植被和土壤组分的温度，具体步骤为：
(1)首先分解土壤-植被系统中植被和土壤组分的能量：土壤热通量g＝gs；净辐射分解为土壤净辐射rn,s和植被净辐射rn,v，公式如下：
rn＝rn,s+rn,v
(2)在得到土壤和植被组分能量的基础上，计算基于端元信息模型四个极限端元的地表温度tsd，tvd，tsw和tvw，公式如下：
其中，ρ为空气密度；cp为定压比热；γ为干湿球常数；δ为饱和水汽压差-温度的斜率；vpd为水气压亏缺；ta为近地表气温；rvw和rvd分别为供水充足和干燥时的冠层阻抗；rav和ras分别为植被和土壤上层的空气动力学阻抗。
(3)然后通过端元地表温度tsd，tvw和未校正的植被覆盖率fv计算土壤表层含水量为0，植被开始受水分胁迫像元的临界地表温度(t*)，公式如下：
(4)最后根据未校正混合像元地表温度(tr)与临界地表温度(t*)的关系，分解得到混合像元中植被和土壤组分的温度：
a)当tr≤t*时，计算对应的植被和土壤组分的温度，公式如下：
b)当tr≥t*时，计算对应的植被和土壤组分的温度，公式如下：
其中未经过角度校正的植被覆盖率fv，计算公式如下：
式中，ndvi为归一化植被指数，ndvimin和ndvimax分别为对应裸土的最小ndvi以及对应全植被覆盖的最大ndvi。ndvi表达式如下：
式中，rnir为近红外反射率，rred为红光反射率，均从mod09ga数据中得到。
步骤3：确定反射率计算模型，进行反射率角度校正。选择ross-librdf核函数模型计算反射率。模型中地表反射率通过代表三种基本散射类型的核函数之和得到，其中三种类型的散射为：(1)各项同性散射；(2)辐射传输体散射；(3)几何光学散射。模型的表达式如下：
其中，r为反射率；θ为太阳天顶角；
反射率角度校正方法为：通过设定反射率计算模型核函数中太阳天顶角θ、观测天顶角
(1)首先确定角度校正后反射率计算模型中的核函数：
其中，o为光照与阴影面的叠合；ξ、θ'、
(2)确定核函数的系数：波段λ的核函数系数fiso(λ)，fvol(λ)和fgeo(λ)均通过遥感mcd43d数据的dn值乘以比例因子获得，数据质量通过mcd43d31文件控制。据此，可确定近红外和红外波段的核函数系数，得到角度校正后像元近红外和红光波段的反射率值(rnir,c和rred,c)。
步骤4：在反射率角度校正的基础上，得到角度校正后的植被覆盖率。在校正得到近红外反射率rnir,c和红光反射率值rred,c的基础上，首先计算归一化植被指数ndvic，然后计算植被覆盖率fv,c，公式如下：
步骤5：在步骤2分解未校正混合像元温度得到植被组分温度(tv)和土壤组分温度(ts)的基础上，结合步骤4得到的角度校正后的植被覆盖率(fv,c)，进行混合像元地表温度的角度校正，得到角校正后的温度(tr,c)，公式如下：
tr,c＝(1-fv,c)ts+fv,ctv
综上，通过未校正混合像元分解得到的植被和土壤组分温度，结合角度校正后的反射率得到的植被覆盖率，得到了角度校正后混合像元的地表温度。
本发明通过基于端元信息模型，对未校正的混合像元温度进行分解，得到植被和土壤组分的温度；通过对modis数据驱动的反射率计算模型进行角度校正，得到校正后的反射率值，进而得到校正后的植被覆盖率；结合分解得到的植被和土壤组分的温度与角度校正的植被覆盖率，得到角度校正后混合像元的地表温度值。本发明使得同一景影像不同像元得到的温度值均为角度校正后的观测值，具有可比性，提高了遥感获取地表温度的实用价值。
以上显示和描述了本发明的基本原理和主要特征及本发明的优点。本行业的技术人员应该了解，上述实施方式并非是对本发明的限制，本发明也并不仅限于上述举例，本技术领域的技术人员在本发明的技术方案范围内所做出的变化、改型、添加或替换，也均属于本发明的保护范围。本发明中的温度计算模型不限制于基于端元模型，采用其他涉及植被覆盖率的模型计算地表温度均在本发明保护范围内。