中國區域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.