中国区域大气-地形-BRDF 全要素地表反射率产品

张兆明，唐朝，何国金*，彭燕，龙腾飞，冷宛春

1. 中国科学院空天信息创新研究院，北京 100094

2. 海南省地球观测重点实验室，海南三亚 572029

3. 中国科学院大学电子电气与通信工程学院，北京 100049

引 言

对地观测领域已经步入大数据时代，每天接收的卫星遥感数据呈爆炸式增长。通常，卫星遥感数据以原始影像数据产品（Digital Number 值，即DN 值）的形式 提供给用户。DN 值没有明确的物理意义，不同时间、不同地点、不同传感器的DN 值不具可比性，遥感数据用户难以直接使用。

为充分挖掘海量卫星数据的应用价值，迫切需要在原始影像数据的基础上进行辐射定量化处理，研发用户能直接使用的“即得即用”（Ready To Use，RTU)产品[1]。地表反射率（Land Surface Reflectance，LSR）是最基本的遥感RTU 产品之一，是开展遥感土地利用分类、长时序地表覆盖变化检测等应用的基础。

为了从原始卫星影像数据得到精确的地表反射率，需要综合考虑传感器辐射定标、大气散射和吸收、地形起伏以及地表二向反射（Bidirectional Reflectance Distribution Function，BRDF）的影响。目前已有的地表反射率产品通常只考虑了传感器辐射定标和大气散射及吸收的影响[2]。在地形起伏地区，地表接收的辐射能量和传感器所接收的信号要受到地形起伏的影响。表现在遥感图像上，阴坡上的像元接收到的入射辐照度较弱、亮度值较低；与此相反，阳坡上的像元接收到的入射辐照度较强、亮度值较高。这样，处在阳坡和阴坡的同类地物的像元亮度值并不相同，而不同地物却可能具有相同或相近的亮度值。这种光谱信息的失真，严重影响了山区遥感图像的信息提取精度。地形校正的目的是减弱或消除地形起伏对传感器接收信号的影响，恢复地物在水平地表条件下的真实反射率。而由于BRDF 效应的存在，在大气和地形条件相同的情况下，同一卫星传感器不同时相的遥感影像由于太阳天顶角和太阳方位角的变化，以及同一时间不同卫星传感器影像之间由于观测天顶角和观测方位角的不同，这些因素都会导致同一地物的地表反射率出现差异，为了得到更精确的地表反射率RTU 产品，需要研发大气-地形-BRDF 耦合校正技术，从而生成大气-地形-BRDF 全要素地表反射率产品。

1 数据采集和处理方法

1.1 数据采集方法

选择中国遥感卫星地面站接收的Landsat 5 卫星数据，覆盖中国全境，数据获取时间为2009 年，尽可能选择无云、高质量的影像数据，所有数据均做过高精度的正射校正处理。

1.2 数据处理方法

在对卫星原始影像各波段数据做高精度辐射定标的基础上，针对卫星传感器参数特点，利用6S辐射传输模型逐像元计算的方式模拟大气散射和吸收影响。6S 模型需要的输入参数利用NCEP（National Centers for Environmental Prediction，美国国家环境预报中心）再分析资料或者与卫星过境时间同步的MODIS（Moderate Resolution Imaging Spectroradiometer，中分辨率成像光谱仪）大气参数产品（大气水蒸汽含量、气溶胶光学厚度、臭氧含量等），具体处理过程参见文献[3]。地形校正利用改进C 地形校正算法，该算法具有较好的地形影响去除效果[4]。BRDF 校正利用C 因子法[5]，将地表反射率调整到同一“太阳-地表-传感器”观测几何下，从而减小太阳角度及卫星观测角度变化对地表反射率反演的影响。通过上述处理，最终实现传感器辐射定标、大气、地形和BRDF 全要素耦合校正，得到地表反射率。

假设太阳天顶角为z，太阳方位角为φ，倾斜坡面的坡度和坡向分别为S和A，坡面太阳入射角为i，改进C 地形校正模型计算过程如下[4]：

其中，T为地形校正系数，h为天空因子，h0为cosi等于0 时的天空因子，c为回归系数，c=b/m, b 和m 可以通过在图像上自动选取样点，对地形校正前地表反射率和cosi之间的关系进行回归分析得到。

在实际计算时，利用NDVI（Normalized Difference Vegetation Index，归一化差值植被指数）区分植被和非植被地表，以NDVI=0.4 为阈值，把影像分为植被区和非植被区，同时将坡度分级，分别统计c系数。

BRDF 校正模型表达式如下：

ρ为对应于波长λ、观测角度Ω（观测天顶角和方位角）和太阳照射角度Ω΄（太阳天顶角和方位角）的光谱反射率。Kvol(Ω, Ω΄)和Kgeo(Ω, Ω΄)分别为体散射核（Ross Thick kernel）和几何光学散射核（Li Sparse Reciprocal Kernel），这两项只与太阳-传感器观测几何(Ω, Ω΄)有关，是太阳天顶角、观测天顶角及相对方位角的三角函数。fiso为各向同性散射系数，等于太阳天顶照射、传感器天顶观测时的地表反射率；fvol和fgeo为权重系数，分别表示体散射和几何光学散射所占的比例，后三个参数和波长相关。具体校正过程参见文献[5]。

2 数据样本描述

数据产品按照影像行列号（path/row）和获取时间存放在不同的文件夹中，文件夹中共包含9 个文件，其中数据文件为TIFF 格式，分波段存储，共6 个文件，文件名末尾为B1-B5 以及B7。数据质量文件为TIFF 格式，以QA（Quality Attribute）标识。元数据文件为txt 格式。数据文件的空间分辨率为30 m，为了降低存储空间，将地表反射率反演结果（浮点型）乘以10000 变成16 位整型，背景填充值为-9999。例如文件夹名为L5-TM-119-036-20090113-LSR-BRDF-TC，表示轨道号为119-036，获取时间为2009 年1 月13 日的Landsat 5 全要素地表反射率产品，L5-TM-119-036-20090113-LSRBRDF-TC-B1.TIF 表示第一波段的数据文件，L5-TM-119-036-20090113-LSR-CLOUD-QA.TIF 表示云质量文件，L5-TM-119-036-20090113-PIXEL-QA.TIF 表示像元质量文件，L5-TM-119-036-20090113-MTL.txt 表示元数据文件。

3 数据质量控制和评估

本节以具体例子通过对比分析来展示大气-地形-BRDF 全要素地表反射率校正的效果。

通过图1 对比可以看出，经过地形校正，地形起伏对地表反射率的影响得到很大程度的抑制，山体阴阳坡同种地物的地表反射率近似相等，地表反射率图像变得更加“平坦”和真实。

图1 地形校正效果对比

为了检验BRDF 校正的效果，选取地表反射率不随季节变化的地表（大片均质的裸地等，简称不变地表），选择冬季（1 月）和夏季（7 月）的卫星影像，由于冬季和夏季太阳高度角和卫星观测角的变化，这些不变地表的地表反射率随季节出现一定的波动（图2a），经过BRDF 校正，减小了太阳高度角及卫星观测角变化对地表反射率反演的影响，冬季和夏季影像上这些不变地表的地表反射率变得更为接近。

图2 BRDF 校正效果对比

4 数据价值

中国全境大气-地形-BRDF 全要素地表反射率产品可为地表覆盖变化调查、遥感地表参数反演、全球变化研究、生态环境监测等遥感应用提供高精度的基础数据产品。

数据作者责任分工

张兆明（1980—），男，河南省郑州市人，博士，正高级工程师，研究方向为遥感数据智能处理与信息挖掘。主要承担工作：研究思路与研究方案设计、算法研究、论文撰写。

唐朝（1996—），男，浙江省湖州市人，硕士研究生，研究方向为遥感数据智能处理与信息挖掘。主要承担工作：算法研究、数据处理。

何国金（1968—），男，福建省武平县人，博士，研究员，研究方向为遥感数据智能处理与信息挖掘。主要承担工作：总体指导，研究思路与研究方案设计。

彭燕（1988—），女，湖南省郴州市人，硕士，工程师，研究方向为遥感数据智能处理与信息挖掘。主要承担工作：数据挑选、数据处理、程序编写、产品生产。

龙腾飞（1986—），男，湖北省武汉市人，博士，副研究员，研究方向为遥感数据智能处理与信息挖掘。主要承担工作：数据处理、程序编写。

冷宛春（1994—），女，辽宁省锦州市人，硕士研究生，研究方向为遥感数据智能处理与信息挖掘。主要承担工作：数据处理、精度评价。