中国区域Landsat地表反射率产品

中国区域Landsat地表反射率产品
作者：彭燕　何国金　张兆明　尹然宇
2020年6月1日

本作品收錄於《中国科学数据》

彭燕, 何国金, 张兆明, 尹然宇. 中国区域Landsat地表反射率产品[J/OL]. 中国科学数据, 2020. (2020-05-25). DOI: 10.11922/csdata.2020.0032.zh.

摘要&关键词

摘要：Landsat数据是记录着人类活动和自然变化最长时间序列的遥感数据集，Landsat系列卫星数据长时序地表反射率产品在森林、水资源、气候变化等领域长时序信息挖掘分析方面具有广泛的用途。本数据集采用基于6S辐射传输模型的大气校正方法，通过辐射定标、模型参数获取、大气校正等数据处理流程，生产了中国1980s–2019年的Landsat系列卫星遥感地表反射率产品，可为遥感用户和科研工作者提供即得即用的高质量地表反射率产品。 最终结果以GeoTiff格式保存，并附带相应的质量文件和元数据文件。

关键词：Landsat 5/7/8；大气校正；6S模型；地表反射率

Abstract & Keywords

Abstract: Landsat satellite remote sensing datasets have been recording the human activities and natural changes since 1970s. And the long-term surface reflectance products of Landsat series data have a wide range of applications like forest monitoring, water resources management, and climate change research. By using the 6S radiative transfer model based atmospheric correction methods，high-quality surface reflectance products over China have been produced from 1980s to 2019, which are Ready To Use (RTU) for researchers and remote sensing application users. The processing steps for such kind of products included radiometric calibration, model input parameter acquisition and atmospheric correction. The products were stored as a format of Geotiff with a corresponding quality assessment files (QA), and metadata files.

Keywords: Landsat 5/7/8; atmospheric correction; 6s model; surface reflectance

数据库（集）基本信息简介

数据库（集）名称	中国区域Landsat地表反射率产品
数据作者	彭燕，何国金，张兆明，尹然宇
数据通信作者	何国金（hegj@radi.ac.cn）
数据时间范围	1980s–2019年
地理区域	中国陆地区域
空间分辨率	30 m
数据量	50 TB
数据格式	*.tiff (GeoTIFF, 16bit interger)
数据服务系统网址	http://databank.casearth.cn http://www.sciencedb.cn/dataSet/handle/984
基金项目	中国科学院A类战略性先导科技专项（XDA19090300）；国家自然科学基金（61731022）。
数据库（集）组成	数据集包括147970景地表反射率产品，以景为单位存放文件夹，每个文件夹中包括每个波段的地表反射率结果，QA文件，缩略图，元数据文件。

Dataset Profile

Title	Landsat surface reflectance products over China
Data corresponding author	He Guojin (hegj@mail.cn)
Data authors	Peng Yan, He Guojin, Zhang Zhaoming, Yin Ranyu
Time range	1980s–2019
Geographical scope	China
Spatial resolution	30 m
Data volume	50 TB
Data format	*.tiff (GeoTIFF, 16bit interger)
Data service system	<http://databank.casearth.cn> <http://www.sciencedb.cn/dataSet/handle/984>
Sources of funding	Strategic Priority Research Program of the Chinese Academy of Sciences (XDA19090300); National Natural Science Foundation of China (61731022).
Database composition	The dataset consists of 147,970 surface reflectance products, which are stored in folders based on imageries. Each folder includes surface reflectance results of each band, QA files, thumbnails and metadata files.

引 言

Landsat系列卫星自1972年发射以来，已获取了大量中空间分辨率的卫星影像，这些影像记录着人类活动和自然变化，成为最长时间序列的星载陆地观测数据集^[1]。地表反射率常用于土地覆盖及变化研究，是很多地表地球物理参数（如，叶面积指数、叶绿素和生物量）反演的关键输入参数。同时，随着对地观测领域进入大数据时代，卫星遥感数据应用用户和研究者不仅希望遥感数据在几何位置上具有一致性，同时也对遥感数据的辐射一致性提出了更高的要求，以更好地进行遥感数据应用和信息挖掘分析。由于卫星遥感数据的大气校正是一项繁琐而专业的工作，因此，研究生产高质量中国区域的Landsat地表反射率产品并公开共享，具有重要的现实意义。

本文提供30 m空间分辨率的中国区域Landsat系列的地表反射率产品，给出了中国区域Landsat数据地表反射率产品生产的技术流程：首先采集覆盖中国区域的Landsat正射影像，然后经过辐射定标，获取6S（Second Simulation of Satellite Signal in the Solar Spectrum）辐射传输模型所需要的几何参数和大气参数并进行大气校正，生成中国区域的Landsat地表反射率产品。

1 数据采集和处理方法

1.1 数据采集方法

本文所生产的中国区域的地表反射率产品主要是在Landsat系列正射影像数据的基础上进行生产的，1980s–2012年为Landsat 5数据，2000–2003年期间为无条带的Landsat 7数据，2013–2019年期间采用Landsat 8数据，这些数据主要由中国遥感卫星地面站接收。对于太阳天顶角大于76°的Landsat数据不进行大气校正处理。

1.2 数据处理方法

中国区域的地表反射率产品是基于6S辐射传输模型进行大气校正而得到的。6S模型是目前比较完善的大气校正模型之一，由Vermote等人在5S模型的基础上改进而来^[2]，适用于0.25–4 μm波长范围内电磁波的大气辐射传输的模拟。假设气溶胶类型为大陆型，并将地表视为平面朗伯体，则其地表反射率可根据表观反射率计算而得到。其公式为^[3]：

\({\rho }_{\mathrm{T}\mathrm{O}\mathrm{A}}={T}_{\mathrm{g}}\left(\mathrm{O}\mathrm{G}\right)×\left[{\rho }_{\mathrm{R}+\mathrm{A}}+{T}_{\mathrm{R}+\mathrm{A}}{T}_{\mathrm{g}}\left({\mathrm{H}}_{2}\mathrm{O}\right)\frac{{\rho }_{\mathrm{s}}}{1-{S}_{\mathrm{R}+\mathrm{A}}{\rho }_{\mathrm{s}}}\right]\) (1)

式中，\({\rho }_{\mathrm{s}}\)为地表反射率；\({\rho }_{\mathrm{T}\mathrm{O}\mathrm{A}}\)为表观反射率，通过辐射定标可计算得到；\(T_g (OG)\))为\({\mathrm{O}}_{2}\)、\({\mathrm{O}}_{3}\)、\(\mathrm{C}{\mathrm{O}}_{2}\)、\(\mathrm{N}{\mathrm{O}}_{2}\)、\(\mathrm{C}{\mathrm{H}}_{4}\)气体总的透过率，\({T}_{\mathrm{g}}\left(\mathrm{O}\mathrm{G}\right)={T}_{\mathrm{g}}\left({\mathrm{O}}_{2}\right)×{T}_{\mathrm{g}}\left({\mathrm{O}}_{3}\right)×{T}_{\mathrm{g}}\left(\mathrm{C}{\mathrm{O}}_{2}\right)×{T}_{\mathrm{g}}\left(\mathrm{N}{\mathrm{O}}_{2}\right)×{T}_{\mathrm{g}}\left(\mathrm{C}{\mathrm{H}}_{4}\right)\)；\({\rho }_{\mathrm{R}+\mathrm{A}}\)为瑞利和气溶胶的反射率；\({T}_{\mathrm{R}+\mathrm{A}}\)为瑞利和气溶胶透过率；\({T}_{\mathrm{g}}\left({\mathrm{H}}_{2}\mathrm{O}\right)\)为水汽的透过率；\({S}_{R+\mathrm{A}}\)为瑞利和气溶胶球面反照率。而\({T}_{\mathrm{g}}\left(\mathrm{O}\mathrm{G}\right)\)，\({\rho }_{\mathrm{R}+\mathrm{A}}\)，\({T}_{\mathrm{R}+\mathrm{A}}\)，\({T}_{\mathrm{g}}\left({\mathrm{H}}_{2}\mathrm{O}\right)\)，\({S}_{\mathrm{R}+\mathrm{A}}\)等大气相关系数均可通过调用6S模型计算得到。而基于6S模型进行大气校正的关键和难点在于大气参数的获取。本数据集是基于公式（1）的计算原理和6S辐射传输模型^[4][3][1]进行生产的，首先需要进行辐射定标，然后获取6S模型所需要的大气参数、DEM以及几何参数进行大气校正，其技术路线图如图1所示。

图片

图1 地表反射率产品生产技术路线图

1.2.1 辐射定标

星上反射率与大气顶层进入卫星传感器的光谱辐射亮度、日地间距离、大气顶层的平均太阳光谱辐照度，以及太阳的天顶角有关。Landsat 5/7的星上反射率计算公式为：

\(\rho =\frac{\mathrm{\pi }{L}_{\mathrm{\lambda }}{d}^{2}}{{E}_{0}\mathit{cos}\left({\theta }_{\mathrm{z}}\right)}\) （2）

式中，ρ为星上反射率，d为日地距离，θz 为太阳天顶角（单位：度，与元数据文件中给出的太阳高度角互为余角），E0 是大气层外相应波长的太阳光谱辐照度，Lλ 为光谱辐射亮度，w/(m2μmsr)，利用定标系数进行线性计算得到。

Landsat 8的星上反射率计算公式为：

\(\rho =\frac{{M}_{\mathrm{\rho }}×DN+{A}_{\rho }}{\mathit{cos}\left({\theta }_{\mathrm{z}}\right)}\) （3）

式中，\({M}_{\mathrm{\rho }}\)和\({A}_{\rho }\)可从元数据文件中获取，分别为\({M}_{\mathrm{\rho }}\)= REFLECTANCE_MULT_BAND_x，\({A}_{\rho }\)= REFLECTANCE_ADD_BAND_x，x为波段号。

1.2.2 模型参数获取

6S模型的输入参数主要包括太阳天顶角、方位角和观测天顶角、观测方位角等观测几何参数以及水汽、臭氧、气溶胶、气压等大气参数和海拔高度信息。在本文中，观测几何参数可根据元数据文件获取，海拔高度信息则是利用0.05°空间分辨率的GCM DEM（Global Climate Model Digital Elevation Model）^[5]得到的。对于Landsat 5/7数据而言，大气校正所需要的水汽数据以及气压数据均来自于NCEP（NOAA National Centers for Environmental Prediction）再分析数据，空间分辨率为2.5°\(×\)2.5°，臭氧数据来自TOMS（Total Ozone Mapping-Spectrometer）数据，空间分辨率为1.25° longitude\(×\)1.00° latitude，气溶胶光学厚度则是利用暗目标法反演得到的^[3]。而对于Landsat 8数据而言，大气校正所需要的水汽数据、臭氧数据、气溶胶光学厚度（AOT）来自于0.05°空间分辨率的MODIS09CMA（MODIS surface reflectance Climate Modeling Grid）数据^[6]。

2 数据样本描述

该数据集是以景为单位存放文件夹，由每个波段的地表反射率产品、质量文件（Quality Assessment, QA）、元数据、缩略图组成。

1）文件夹的命名规则为：卫星-传感器-path-row-成像日期-LSR，如L5-TM-115-026-19840418-LSR。

2）每个波段的地表反射率产品：Landsat 5/7包括波段1，2，3，4，5，7（即分别为蓝、绿、红、近红、中红外1、中红外2），命名规则为：卫星-传感器-path-row-成像日期-LSR-BX（X表示第几个波段）.TIF，如L5-TM-115-026-19840418-LSR-B1.TIF。影像的空间分辨率为30 m，投影坐标系是WGS84 UTM。为了降低存储空间，地表反射率结果由原本的0~1范围内的浮点型均乘以10000变成16位整型，背景填充值为-9999，并进行了“LZW”的无损压缩。

3）质量文件（QA）：包括在原始数据的基础上生成的QA，为PIXEL-QA，主要是对填充值（Fill）、晴空（Clear）、云（Cloud）、云置信度（Cloud Confidence）、云阴影（Cloud Shadow）、冰雪（Snow/Ice）以及水（Water）等信息进行标识，命名规则为：卫星-传感器-path-row-成像日期-PIXEL-QA.TIF；和在地表反射率的基础上生成的QA，主要是对气溶胶相关的信息进行标识，Landsat 5/7为LSR-CLOUD-QA，命名规则为卫星-传感器-path-row-成像日期-LSR-CLOUD-QA.TIF，Landsat 8为LSR-AEROSOL，命名规则为卫星-传感器-path-row-成像日期-LSR-AEROSOL.TIF。Landsat 5/7和Landsat 8的QA的属性表分别如表1–4所示。

表1 Landsat 5/7的PIXEL-QA属性表

位	位值	累计和	属性	像素值
0	1	1	Fill	1
1	2	3	Clear	66,130
2	4	7	Water	68,132
3	8	15	Cloud shadow	72,136
4	16	31	Snow/Ice	80,112,144,176
5	32	63	cloud	96,112,160,176,224
6	64	127	Cloud Confidence00 = none01 = low10 = medium11 =high	Low cloud confidence:66,68,72,80,96,112Medium cloud confidence:130,132,136,144,160,176High cloud confidence:224
7	128	255

表2 Landsat 5/7的LSR-CLOUD-QA的属性表

位	位值	累计和	属性	像素值
0	1	1	Dark Dense Vegetation(DDV)	1,9
1	2	3	CLOUD	2,34
2	4	7	CLOUD Shadow	4,12,20,36,52
3	8	15	Adjacent cloud	8,12,24,40,56
4	16	31	Snow/Ice	16,20,24,48,52,56
5	32	63	water	32,34,36,40,48,52,56

表3 Landsat 8 的PIXEL-QA属性表

位	位值	累计和	属性	像素值
0	1	1	Fill	1
1	2	3	Clear	322,386
2	4	7	Water	324,388,836,900
3	8	15	Cloud shadow	328,392,840,904
4	16	31	Snow/Ice	336, 368, 400, 432, 848, 880, 912, 944
5	32	63	cloud	352, 368, 416, 432, 480, 864, 880, 928, 944, 992
6	64	127	Cloud Confidence00 = none01 = low10 = medium11 =high	Low cloud confidence:322, 324, 328, 336, 352, 368, 832, 836, 840, 848, 864, 880Medium cloud confidence:386, 388, 392, 400, 416, 432, 900, 904, 928, 944High cloud confidence:480, 992
7	128	255
8	256	511	Cirrus Confidence00 = Not set01 = Low10 = medium11 =high	Low cirrus confidence:322, 324, 328, 336, 352, 368, 386, 388, 392, 400, 416, 432, 480High cirrus confidence:832, 836, 840, 848, 864, 880, 900, 904, 912, 928, 944, 992
9	512	1023
10	1024	2047	Terrain Occlusion	>=1024

表4 Landsat 8 的LSR-AEROSOL属性表

位	位值	累计和	属性	像素值
0	1	1	Fill Value	1
1	2	3	Aerosol Valid	1,20,66,74,130,138,194,202
2	4	7	Aerosol Interpolated	4,12,68,76,132,140,196,204
3	8	15	Water	8,12,24,40,56
4	16	31	Climatology level aerosol	2, 4, 8, 10, 12, 16, 24, 28
5	32	63	Low level aerosol	64, 66, 68, 72, 74, 76
6	64	127	Medium level aerosol	128, 130, 132, 136, 138, 140
7	128	255	High level aerosol	192, 194, 196, 200, 202, 204

4）元数据文件：包括原始数据的元数据文件（MTL.txt文件），命名规则为：卫星-传感器-path-row-成像日期-MTL.txt；和地表反射率产品的XML元数据文件，命名规则为：卫星-传感器-path-row-成像日期-LSR.xml。XML元数据文件描述了地表反射率产品的基本信息以及各个波段产品的相关信息。

5）缩略图：包括512像素大小和1024像素大小的缩略图，命名规则分别为：卫星-传感器-path-row-成像日期-LSR-THUMB.JPG和卫星-传感器-path-row-成像日期-LSR-BROWSER.JPG。

图2为2005年中国区域的地表反射率结果示意图，其中图2（a）为大气校正前的DN数值镶嵌结果，图2（b）为地表反射率的镶嵌结果，二者均没有进行任何的色彩调整，从图2中可以看出，进行大气校正后的地表反射率镶嵌结果具有更好的辐射一致性。

图片（a）

图片（b）

图2 中国陆地区域地表反射率产品示意图（a）大气校正前的DN数值镶嵌结果； (b) 地表反射率的镶嵌结果

3 数据质量控制和评估

Feng等人将全球2000年和2005年的Landsat 5/7地表反射率产品分别与相应的MODIS地表反射率产品进行比较分析，结果表明Landsat 5 TM的均方根差（RMSD）约为2.2%–3.5%，Landsat 7 ETM+的均方根差约为1.3%–2.8%^[7]。Peng等人利用2014年6月11日获取的南京地区的实测光谱数据对Landsat 8数据地表反射率产品进行验证，结果显示其地表反射率产品的RMSD约为3%–5%之间^[1]。

4 数据价值

本文推出的1986–2019年中国区域Landsat系列卫星长时序地表反射率产品，采用目前国际上较为成熟的地表反射率算法，精度较高，未来将持续更新，在森林、水资源、气候变化等领域长时序信息挖掘分析方面具有重要的应用价值。

5 数据使用方法和建议

本数据集可通过地球大数据科学工程（CASEarth）Databank在线服务网址（http://databank.casearth.cn）获取数据。用户注册成功并登录系统后，进入平台产品查询界面（图3），产品类型选择“LSR”（图3红框位置所示），然后在所需要的数据下点击下载（如图2中小蓝框位置所示），随后会弹出一个下载对话框（如图2中大蓝框位置所示），该框中列出相应的指数产品结果，QA文件和元数据文件，根据需要下载数据即可。用户可根据行政区、地图选择以及行列号等方式查询所需要的数据。目前共享的产品主要包括1980s–2012年云量小于50%的Landsat 5、2000–2003年云量小于20%的Landsat 7以及2018–2019年的全部Landsat 8的地表反射率产品，后续作者将持续生产我国区域云量大于50%的Landsat 5、云量大于20%的Landsat 7、2013–2017年和2019年以后的Landsat 8地表反射率产品，以提供更好的、持续的数据共享服务。若平台系统中暂时缺少用户所需的数据或者有与本数据相关的其他数据需求，可通过咨询本文作者进行申请。

图片

图3 CASEarth Databank在线服务地表反射率产品查询界面示意图（底图来源：高德地图）

致 谢

衷心感谢王桂周和龙腾飞在产品生产时遇到的存储与效率问题给予的建设性意见。

参考文献

1. ^ ^{1.0 1.1 1.2} PENG Y, HE G J, ZHANG Z M, et al. Study on atmospheric correction approach of Landsat-8 imageries based on 6S model and look-up table[J]. Journal of Applied Remote Sensing, 2016, 10(4). DOI: 10.1117/1.JRS.10.045006.
2. VERMOTE E F, EL SALEOUS N, JUSTICE C O, et al. Atmospheric correction of visible to middle-infrared EOS-MODIS data over land surfaces: background, operational algorithm and validation. J. Geophys. Res, 1997,102(D14): 17131–17141.
3. ^ ^{3.0 3.1 3.2} JEFFREY G M, VERMOTE E F, SALEOUS N E, et al. A Landsat Surface Reflectance Dataset for North America, 1990-2000[J]. IEEE Geoscience and Remote Sensing Letters, 2006, 3(1): 68-72.
4. VERMOTE E, JUSTICE C, CLAVERIE M, et al. Preliminary analysis of the performance of the Landsat 8/OLI land surface reflectance product[J]. Remote Sensing of Environment, 2016,185:46-56.
5. ETOPO5. Data announcement 88-MGG-02, digital relief of the surface of the earth[DB]. Boulder, Colorado: NOAA, National Geophysical Data Center, 1988.
6. GAO B C, AND KAUFMAN Y J. Water vapor retrievals using moderate resolution imaging spectroradiometer (MODIS) near-infrared channels[J]. Journal of Geophysical Research: Atmospheres (1984–2012),2003,108(D13).
7. FENG M, HUANG C Q, CHANNAN S, et al. Quality assessment of Landsat surface reflectance products using MODIS data[J]. Computers & Geosciences, 2012,38(1): 9-22.

数据引用格式

彭燕, 何国金, 张兆明, 尹然宇. 中国区域Landsat地表反射率产品[DB/OL]. Science Data Bank, 2020. (2020-04-29). DOI: 10.11922/sciencedb.984.