北京市平原区2015年—2019年0.8 m地表反射率数据集

北京市平原区2015年—2019年0.8 m地表反射率数据集 Annual 0.8 m surface reflectance data set of Beijing plain area from 2015 to 2019 1 first-author 申茜，1981年生，女，副研究员，研究方向为水环境遥感、高光谱遥感。E-mail： shenqian@aircas.ac.cn 申茜，1981年生，女，副研究员，研究方向为水环境遥感、高光谱遥感。E-mail： shenqian@aircas.ac.cn shenqian@aircas.ac.cn 1 1 1 2 1 2 1 3 中国科学院空天信息创新研究院 数字地球重点实验室, 北京 100094 中国科学院空天信息创新研究院 数字地球重点实验室, 北京 100094 Key Laboratory of Digital Earth Science, Aerospace Information Research Institute, Chinese Academy of Sciences, Beijing 100094, China Key Laboratory of Digital Earth Science, Aerospace Information Research Institute, Chinese Academy of Sciences, Beijing 100094, China 中国遥感卫星地面接收站, 北京 100094 中国遥感卫星地面接收站, 北京 100094 China Remote Sensing Satellite Ground Station, Beijing 100094, China China Remote Sensing Satellite Ground Station, Beijing 100094, China 中国科学院大学 资源与环境学院, 北京 100049 中国科学院大学 资源与环境学院, 北京 100049 College of Resources and Environment, University of Chinese Academy of Sciences, Beijing 100049 ,China College of Resources and Environment, University of Chinese Academy of Sciences, Beijing 100049 ,China 地表反射率产品作为最重要的定量遥感产品，是很多参量化遥感产品的基础数据源，可以被广泛应用于林业、农业、水资源、生态环境、城市环境等典型应用领域。对于米级高分辨率的遥感影像，国内外不提供反射率影像产品。目前主流的国产高分数据源大部分都是蓝、绿、红、近红4个波段的多光谱数据，缺少短波红外波段，难以满足陆地区域浓密植被算法或者清洁水体的短波红外信号近似看作零的假设条件，要完成精确大气校正、生产大范围的GF-2影像拼接的地表反射率产品，是一个挑战。为了更好的推广高分地表反射率的应用，本研究针对北京市平原区，运用全色和多光谱融合、几何精校正和相对辐射一致化算法，以Sentinel-2为参考影像，生成了一套几何精校正后的0.8 m地表反射率影像集。该数据集时间范围为2015年—2019年，每年一期，包含年度产品的覆盖情况及分布矢量，共计184景地表反射率影像，数据量总共为1.63 TB。 As the most important quantitative remote sensing product, surface reflectance products are the basic data source for many parametric remote sensing products and can be widely used in typical applications, such as forestry, agriculture, water resources, ecological environment, and urban environment. However, for meter-level high-resolution remote sensing images, reflectance products are unavailable at home and abroad. Most of the current domestic satellite products use high-resolution multispectral data in the four bands of blue, green, red, and near infrared as the main data source, and accurate atmospheric correction is difficult to achieve mainly due to the lack of short-wave infrared data. A set of spatiotemporal continuous meter-level high-resolution image reflectance products should also be provided for researchers to use. In this study, GF-2 Level-1A standardized products were used as input. The processing steps included panchromatic and multispectral fusion, geometric precision correction and mosaic, and atmospheric correction. A set of annual surface reflectance images with a spatial resolution of 0.8 m covering Beijing plain area from 2015 to 2019 was generated. First, a geometric deviation of subpixel level exists among the four multispectral bands of GF-2 Level-1A products. The Pixel Knife software was used to register each band of GF-2 images accurately, then the images were fused with panchromatic images. Second, with Sentinel-2 as reference images (10 m resolution) and supplementary SRTM DEM data (30 m resolution), we used a regional network adjustment method to achieve geometric precision correction. The geometric accuracy is within 20 m. This method ensures the absolute geometric positioning accuracy and the relative geometric accuracy among images. Third, we used a relative radiation uniformity method to complete atmospheric correction. Regarding Sentinel-2 reflectance images as reference images (10 m resolution), we automatically searched the pseudo invariant points between the reference images and the images to be corrected and built a regression equation band by band. A total of 184 scenes of surface reflectance images with a total data volume of 1.63 TB are acquired. The data set is issued annually, including the coverage and distribution vector of the annual products. In the current data set, for mountain-shaded and relatively clean water bodies, such as Miyun Reservoir water bodies, the retrieval reflectance is often zero or negative. Therefore, this data set is suitable for underlying surface applications, except mountainous areas and relatively clean water bodies. It can prevent the shading phenomenon after the fusion of panchromatic and multispectral images and ensure that the multiscene geometric precision correction images have good geometric consistency at the joints. The validation results show that the atmospheric correction effect of water, road, and vegetation canopy is good. On the premise that the reference image is ready, the processing method proposed in this article can process high-resolution images and output surface reflectance products rapidly, massively, and automatically. Currently, 70% success rate is achieved for atmospheric correction, which should be further improved. 0.8 m地表反射率 GF-2 高分辨率影像 大气校正 北京市平原区 0.8 m surface reflectance data GF-2 high resolution image atmospheric correction plain area of Beijing 1　引言 1 地表反射率（以下简称为反射率）产品作为最重要的定量遥感产品，是很多参量化遥感产品的基础数据源（ 梁顺林 等，2020 梁顺林 等，2020 ， 2016 2016 ；刘良云和张肖，2017； 曹小杰 等，2019 曹小杰 等，2019 ）。可以被广泛应用于林业、农业、水资源、生态环境、城市环境等典型应用领域。目前国际上MODIS、Landsat系列和Sentinel-2等卫星都公开发布了L2A级反射率产品（ Vermote等，2016 Vermote等，2016 ； Richter，1998 Richter，1998 ； Richter和Schläpfer，2015 Richter和Schläpfer，2015 ； Doxani等，2018 Doxani等，2018 ），用户可以直接下载使用。对于米级高分辨率的遥感影像，国内外不提供反射率影像产品。一方面，目前5 m以内的高分数据应用主要集中在土地利用分类、目标解译识别等对反射率要求不高的领域，如果遇到如黑臭水体识别、道路材质识别等对反射率精度要求较高的应用，就需要研究者自己做大量的图像处理工作，无形中升高了高分影像应用的门槛。另一方面，商业公司以发射多颗米级卫星组网实现地球表面全覆盖，但发射的卫星数量越多，成本会越高，如果影像不好用，高成本就难以收回。如何利用好现有的海量多光谱高分影像并降低应用门槛，是所有高分遥感卫星的一个痛点。因此，提供一套时空连续的米级高分辨率影像反射率产品供研究者使用，具有重要的意义。 国产卫星数据反射率产品提供面临主要挑战有2个：一方面，尽管基于图像自身也能够实现气溶胶信息反演，但往往需要图像具备短波红外波段，以满足陆地区域浓密植被算法或者清洁水体的短波红外信号近似看作零的假设条件，这些主流的国产高分数据源基本都是蓝、绿、红、近红4个波段的多光谱数据，缺少短波红外波段，要做到定量化处理，难度较大；另一方面，传感器误差、标准产品的处理误差、定标及真实性检验场的误差（ 吴小丹 等，2014 吴小丹 等，2014 ）、大气校正误差等这一系列误差积累都导致反射率产品的误差较大，需要深入逐一对每个步骤的误差进行严格控制。 本研究主要以GF-2 L1级标准化产品为输入，完成全色和多光谱融合、几何精校正和大气校正等处理，每年度包括多幅GF-2影像并全覆盖北京市平原区，生产出2015年—2019年的年度北京平原区反射率产品。 2　数据集元数据简介 1 北京市平原区2015年—2019年年度0.8 m反射率数据集的名称、作者、地理区域、数据年代、时间分辨率、空间分辨率、数据集组成、数据共享政策等信息见 表1 表1 。 数据集关键词 详细信息 作者信息 申茜，中国科学院空天信息创新研究院，shenqian@aircas.ac.cn 数据量 1.63 TB 数据格式 .shp， .tiff， .xlsx 数据集组成 数据单位为无量纲，除以10000即得到反射率；数据信息表.xlsx包含了所有影像的传感器、中心经度、中心纬度、日期、产品序列号、景号、参考影像、系数选择；名称为“201*镶嵌线+缩略图”的文件夹，包含了2015年—2019年年度镶嵌线矢量和影像缩略图 数据计算环境 C、Python、IDL、ArcGIS 数据共享政策 数据包括影像数据、镶嵌线矢量数据和通过《遥感学报》发表的数据论文。其共享政策如下：（1）“数据”以最便利的方式通过互联网系统免费向全社会开放，用户需要与作者签署书面协议后，获得下载地址，免费浏览、免费下载；（2）最终用户使用“数据”需要按照引用格式在参考文献或适当的位置标注数据来源；（3）不得用于商业目的 3　数据源和处理算法 1 3.1　数据源 2 本文的数据源主要为2015年—2019年的0.8 m的GF-2 L1A级产品共181景，2 m的GF-1 L1A级产品3景作为补充，共计184景。每年度保证覆盖北京平原区一次，由于受云和耀斑影响导致个别区域没有数据覆盖，采用了分辨率2 m的GF-1 L1A级产品补充完整。数据以夏、春、秋季影像优先，个别区域采纳了无冰的冬季影像为补充。 3.2　技术路线 2 中国资源卫星中心分发的GF-2 L1A级产品经过全色和多光谱融合、几何精校正、大气校正3个处理步骤，获得0.8 m 4个波段的反射率产品，技术路线如 图1 图1 所示。 3.3　处理算法 2 3.3.1　全色和多光谱融合 3 中国资源卫星中心分发的GF-2 L1A级产品的4个多光谱波段之间存在亚像元级别的几何偏移。商业软件自带的全色和多光谱融合算法，例如NNDiffuse（ Sun等，2014 Sun等，2014 ）、Gram-Schmidt（ Laben和Brower，2000 Laben和Brower，2000 ）、主成分变换（ Chavez等，1991 Chavez等，1991 ）等算法全都是假定各波段影像是精确配准的，但如果影像不满足该假定就会导致输出的影像呈现出一定程度的晕染现象。对于窄路、细河这种辐射值会发生锐利变化的地物，融合前后的辐射值变化非常大。 这里采用了像素刻刀Pixel-Knife（ 陈甫，2021 陈甫，2021 ）软件，首先通过梯度追踪的精细配准，修正了4个波段之间细微的错位，保证配准误差小于0.1个像元，然后再实现复原融合处理，融合后的影像不仅能够精确地校正波段间的偏离，还能保持高度的波谱保真度。对于窄路、细河这种辐射值会发生锐利变化的地物，融合前后的辐射值几乎不变。采用像素刻刀处理一景GF-2 L1A级影像大约需要8 min（CPU i7 8700及以上，内存32 G及以上），可以大批量处理高分系列多光谱影像（如GF-1/1B/1C/1D PMS、GF-6 PMS等）。 因此，如果采用融合影像做定量化研究，或对辐射信息保真度要求较高，建议先对GF-2多光谱影像各个波段精确配准，然后再与全色影像做融合处理。 3.3.2　几何精校正 3 全色和多光谱融合后GF-2影像空间分辨率约为0.8 m，实现其绝对几何位置的高精度校正需要依赖更高精度的亚米级地面控制点，因此直接采用地面控制点的几何精校正方法是不适用的。事实上，通常非测绘级的实际应用对影像的绝对几何精度要求并不高，关键是要保证整个研究区内、不同时相影像数据的几何一致性。 普通的商业软件如ENVI软件，几何精校正模块可以实现逐一处理单景影像，无法采用区域网平差的方法对多幅影像同时处理。中国遥感卫星地面站正图软件、PixelGrid、GEOWAY CIPS、像素工厂、PCI-GXL等可以采用区域网平差方法同时输出多景正射影像，可同时保证影像的绝对几何定位精度和影像间的相对几何精度。 这里采用了中国遥感卫星地面站正图软件，一方面采用全球公开的Sentinel-2（空间分辨率为10 m）几何影像产品作为参考影像（ European Space Agency，2015 European Space Agency，2015 ），通过影像匹配的方法自动获取控制点（ Long等，2016 Long等，2016 ）；另一方面也在不同影像的重叠区内自动提取连接点，然后综合利用控制点和连接点，在SRTM DEM数据的辅助下（空间分辨率为30 m）进行区域网平差（ Long等，2016 Long等，2016 ；龙腾飞，2016）。区域网平差以GF-2影像自带的RPC参数建立初始几何成像模型，以像方仿射变换作为修正模型（ Grodecki和Dial，2003 Grodecki和Dial，2003 ），采用Huber函数加权的稀疏Levenberg-Marquardt算法（ Kanzow等，2004 Kanzow等，2004 ）进行最优化求解，同时克服粗差点对几何成像模型的影响。最后利用修正后的RPC模型对影像进行几何校正。该方法同时保证了影像的绝对几何定位精度和影像间的相对几何精度。在准备好参考影像的前提下，采用中国遥感卫星地面站正图软件处理GF-2 L1A级影像单景平均需要约10 min，可以大批量处理高分系列多光谱影像，精度基本在20 m以内。但是对于山区、云覆盖区、个别几何偏差本身较大的影像，会存在稍大的偏差。本文逐一检查结果，对1景偏差大于30 m的影像，重新做了单景影像的几何精校正。 3.3.3　大气校正 3 影像成像瞬间通常难以同步测量大气参数，而直接基于实测大气参数进行影像大气校正难以做到大面积和业务化应用。尽管基于图像自身也能够实现气溶胶信息反演，但往往需要有短波红外波段设置，以满足陆地区域浓密植被算法或者清洁水体的短波红外信号近似看作零的假设条件。由于国产高分系列卫星的多光谱卫星数据，往往是蓝、绿、红、近红外4个多光谱波段，缺少短波红外波段，因此基于图像自身的大气校正算法也较难实现。 基于高分系列多光谱影像自身实现大气校正，目前可以采用ENVI提供的FLAASH模块。FLAASH需要以辐亮度数据为输入，因此需要代入高分传感器定标系数将DN值转换为辐亮度；FLAASH是基于MODTRAN5辐射传输模型，人工估算平均能见度或采用暗像元方法估算能见度、辅以其他大气参数实现大气校正。由于定标系数偏差、能见度估算误差、其他参数偏差等原因，导致FLAASH反演的反射率精度有时无法保证。此外，FLAASH不能批量处理影像。 本文采用了一种相对辐射一致化方法（ Shen等，2019 Shen等，2019 ），实现高分系列多光谱影像的大气校正。第1步，仍延续选择与高分影像时空对应的Sentinel-2（空间分辨率为10 m）L2A级地表反射率影像作为参考影像，由于几何精校正参考影像就是采用该影像，再将待处理的高分影像重采样到10 m，保证了待处理影像和参考影像的几何配准精度；第2步，采用多元变化检测算法（ Canty等，2004 Canty等，2004 ），在待处理影像和参考影像中自动搜索出伪不变地物点，即在多时相影像中反射率不随季节和生物周期的变化而变化的地物；第3步，利用这些伪不变地物像元，逐波段构建出Sentinel-2影像地表反射率与高分影像DN值的正交线性回归方程，进而对待处理高分影像进行逐波段地相对辐射一致化处理。采用该算法可以自动化大批量处理高分系列多光谱影像，处理GF-2 L1A级影像单景平均约需要15 min（CPU i7 8700及以上，内存32 G及以上）。 通过观察Sentinel-2影像上北京市平原区范围内的植被生长期和凋谢期，大致划分了4个时段：冬季11月15日—次年4月15日；植被生长期4月15日—5月15日；夏季5月15日—10月15日；植被凋谢期10月15日—11月15日。在4个时段内选择与待处理影像的时相尽可能接近的参考影像，完成大气校正。其中输入184景高分影像，经批量处理后，人工检查发现132景大气校正效果较好，处理后影像的地物与参考影像地物反射率波形基本一致，数值也基本一致。对于52景效果较差的影像，应用以下策略：（1）重新选取感兴趣区域；（2）替换待校正影像或参考影像；（3）如果上下相邻图像正好是同一轨道同一时间拍摄，可以借用邻近图像的回归方程系数完成大气校正；共46景输出了较好的结果；剩余6景影像除水体以外的地物反射率与参考影像地物反射率波形基本一致，水体反射率的其他波段数值基本一致，但水体反射率的红波段、近红外波段往往为零。这6景影像主要覆盖了密云水库以及怀柔区较为清澈的水体，推测是水体较为清澈，高分系列卫星传感器在近红波段响应程度不敏感所导致。此外，所有影像中，部分山区阴影的存在反射率零值情况，在应用时需要注意。 4　算法评价 1 这里分别从全色和多光谱融合、几何精校正、大气校正3个步骤分别介绍评价方法、数据和评价结果。 4.1　全色和多光谱融合评价 2 针对同一景全色和同一景多光谱影像，分别采用了ENVI的NND融合算法和Pixel-Knife复原融合算法，对比了全色和多光谱融合结果。评价数据采用了2016-08-27顺义区（景号：2754656，裁剪区域中心经纬度：40.175789 N，116.5505167 E）的GF-2全色和多光谱影像。 结果如 图2 图2 所示， 图2（a） 图2（a） 中大面积均匀地物效果较好，但是红线圈出的道路两侧、河流边缘有明显的晕染现象，类似于“彩虹”。 图2（b） 图2（b） 是Pixel-Knife的结果，道路两侧、河流等边缘清晰、锐利。这两种融合结果的差异，主要源于GF-2原始影像各波段之间存在细微的几何偏差。商业软件自带的全色和多光谱融合算法都是假定各波段影像的几何位置是精确配准的，但如果假定不满足就会导致输出的影像呈现出晕染现象。同时，如果采用融合影像做定量化研究，仍然建议先对GF-2多光谱影像各个波段精确配准，然后再与全色影像做融合处理。 1,2; 4.2　几何精校正评价 2 首先，针对有重叠区域的两景融合后影像进行评价，分别采用ENVI软件几何精校正模块逐一处理单景影像和采用本文的区域网平差方法同时处理两景影像，对比了几何精校正后结果。其中，数据采用了2018-08-31（景号5469411）和2018-09-05（景号5487784）获取的两景GF-2影像。采用卷帘工具对比结果，如 图3 图3 所示， 图3（a） 图3（a） 和 图3（b） 图3（b） 是左右部分分别为这两景影像， 图3（c） 图3（c） 和 图3（d） 图3（d） 中上下部分分别为这两景影像。可以看到ENVI软件处理结果， 图3（a） 图3（a） 和 图3（c） 图3（c） 中两景几何精校正影像在房屋处存在明显的错位；而本文方法处理结果， 图3（b） 图3（b） 和 图3（d） 图3（d） 中两景几何精校正影像在接边处具有较好的几何一致性。 1,2;3,4; 其次，针对有重叠区域的多景影像进行评价，数据采用2016年北京市平原区37景GF-2影像，先自动匹配得到控制点和连接点，随后通过区域网平差修正影像的RPC模型并对37景影像进行几何精校正。在研究区内均匀选取的30个检查点验证几何校正影像的绝对几何精度，在影像重叠区内均匀选取40个连接检查点验证几何校正影像的接边精度。几何精校正采用的精度评价指标是均方根误差（RMSE），公式如下 <math display="block" id="M1"><mi mathvariant="normal">R</mi><mi mathvariant="normal">M</mi><mi mathvariant="normal">S</mi><mi mathvariant="normal">E</mi><mo>=</mo><mroot><mrow><mfrac><mrow><mn mathvariant="normal">1</mn></mrow><mrow><mi>n</mi></mrow></mfrac><mstyle displaystyle="true"><munderover><mo largeop="true">∑</mo><mrow><mi>i</mi><mo>=</mo><mn mathvariant="normal">1</mn></mrow><mrow><mi>n</mi></mrow></munderover></mstyle><mrow><mo stretchy="false">(</mo><msub><mrow><mover accent="true"><mi>y</mi><mo>̂</mo></mover></mrow><mrow><mi>i</mi></mrow></msub><mo>-</mo><msub><mrow><mi>y</mi></mrow><mrow><mi>i</mi></mrow></msub><msup><mrow><mo stretchy="false">)</mo></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msup><mo>+</mo><mo stretchy="false">(</mo><msub><mrow><mover accent="true"><mi>x</mi><mo>̂</mo></mover></mrow><mrow><mi>i</mi></mrow></msub><mo>-</mo><msub><mrow><mi>x</mi></mrow><mrow><mi>i</mi></mrow></msub><msup><mrow><mo stretchy="false">)</mo></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msup></mrow></mrow><mrow/></mroot></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=22213686&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=22213681&type= 55.88000107 9.56733322 式中， n 表示检查点的个数， <math id="M2"><msub><mrow><mover accent="true"><mi>x</mi><mo>̂</mo></mover></mrow><mrow><mi>i</mi></mrow></msub><mo>,</mo><msub><mrow><mover accent="true"><mi>y</mi><mo>̂</mo></mover></mrow><mrow><mi>i</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=22213692&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=22213689&type= 6.68866634 4.57200003 表示检查点在校正后影像中对应的地面投影坐标， <math id="M3"><msub><mrow><mi>x</mi></mrow><mrow><mi>i</mi></mrow></msub><mo>,</mo><msub><mrow><mi>y</mi></mrow><mrow><mi>i</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=22213697&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=22213694&type= 6.68866634 4.57200003 表示检查点在参考影像中对应的地面投影坐标。 结果显示，几何精校正影像的绝对几何精度（均方根误差）约为18.4 m，接边精度（均方根误差）约为3.8 m。 4.3　大气校正评价 2 首先对比了8种典型地物（水体、植被、人工草坪、红色塑胶跑道、蓝色屋顶、灰色屋顶、裸土、道路）分别在校正后GF-2影像和参考影像中的反射率（无量纲），来评价大气校正算法的有效性。采用的GF-2影像是2018-09-05（景号5487783），Sentinel-2参考影像是2017-10-30。验证样点所在的GF-2影像局部放大图，如 图4（a） 图4（a） 所示；基于验证样点的相对辐射一致化结果定量评价，如 图4（b） 图4（b） 所示。 图4（b） 图4（b） 横轴是参考影像Sentinel-2中8种地物的反射率，纵轴是大气校正后GF-2影像中8种地物的反射率，这些点基本都分布在1∶1线附近，表明数值十分接近，校正效果较好。 1,2; 其次，对比了准同步实验中地物实测等效反射率，和对应点在大气校正后高分影像中的反射率，这里选择了水体、道路、植被冠层3种典型地物，实测数据与影像相差在3 d之内。 采用的实测水体数据包括：2016-09-19和2016-10-09沈阳市城区水体，2018-09-05北京市朝阳区4条城市河流水体，2019-09-18北京市北海公园水体和2019-10-25北京市紫竹院公园水体，共51个同步的水体采样点。对应的同步影像分别为2016-09-19 GF-2 PMS1（景号2827917），2016-10-09 GF-2 PMS2（景号2886011），2018-09-05 GF-2 PMS2（景号5487916），2019-09-18 GF-6 PMS（景号260947）和2019-10-25 GF-6 PMS（景号273296）；前2景沈阳影像采用了2019-09-11 Sentinel-2影像、2018-09-05同步试验采用了2019-10-05 Sentinel-2影像、最后2景影像采用了2019-10-25的Sentinel-2影像作为参考影像，完成了本文的全色和多光谱融合、几何精校正、大气校正处理步骤。由于水体一般采用遥感反射率，单位sr -1 ，所以这里评价水体采用的是遥感反射率，一般可以用地表反射率图像除以π近似获得。实测的水体遥感反射率采用波段等效计算到GF-2的4个波段，用于验证对应点在大气校正后GF-2影像中的遥感反射率产品。 植被冠层实测数据包括：2016-12-07、2018-12-21和2019-04-04北京市小汤山农业示范园冬小麦，共3个准同步的植被冠层采样点。对应的影像分别为2016-12-09 GF-2 PMS2（景号3101468）影像、2018-12-24 GF-1B PMS（景号212295）影像和2019-04-04 GF-1D PMS（景号260863）影像，采用了2017-03-09 Sentinel-2影像作为参考影像。评价植被冠层和道路仍然采用地表反射率，无单位。 道路实测数据包括：2019-07-05北京市邓庄南路附近10个道路采样点。对应的影像是2019-07-03 GF-2 PMS1（景号6612235）影像，采用了2017-05-28 Sentinel-2影像作为参考影像。 结果如 图5 图5 所示，横轴是地物实测等效后的反射率，纵轴是校正后高分影像中地物的反射率，这些点基本都分布在1∶1线附近，校正效果较好。 1,2,3;4; 从 图5 图5 中可以看出，对于水体，大气校正后GF-2影像的遥感反射率与地面实测等效计算的遥感反射率在蓝绿红3个波段都较为接近，在近红波段图像反演值略高于实测值。这可能是由于水体本身在近红波段的信号较低，传感器在该波段的响应程度不敏感的原因，进而造成较大误差。若仅对水体蓝绿红3个波段数据进行分析， R 2 为0.89，MAPE为24.08%，RMSE为0.0039；对于植被冠层，4个波段 R 2 为0.78，除近红波段偏差略大，其他3个波段的反演值和实测值较为接近；对于道路，包括沥青、水泥路、土路、碎石路4种材质， R 2 为0.88，MAPE为20.57%，RMSE为0.016，4个波段的反演值和实测值也都十分接近。通过对比地面实测反射率与大气校正后图像对应点的反射率，表明了3种典型地物的大气校正的效果都较好。 最后，对比了典型地物在多时序高分影像中的反射率，这里选择了水体、植被、屋顶3种典型地物。采用的影像数据为2018-01—2020-08的北京市同一区域的高分影像（包括GF-1 PMS、GF-1B/C/D PMS、GF-2 PMS和GF-6 PMS），时间范围为1—3月、4—6月、7—9月和10—12月每3月中选择一景，选择典型地物感兴趣区的均值，用于对比典型地物在不同时期高分影像中的反射率。同样，水体采用遥感反射率评价，单位sr -1 ；植被和屋顶选择反射率（无量纲）评价。 以平均绝对百分比误差（MAPE）作为大气校正的评价指标， <math display="block" id="M4"><mi mathvariant="normal">M</mi><mi mathvariant="normal">A</mi><mi mathvariant="normal">P</mi><mi mathvariant="normal">E</mi><mo>=</mo><mfrac><mrow><mn mathvariant="normal">100</mn><mi mathvariant="normal">%</mi></mrow><mrow><mi>n</mi></mrow></mfrac><mstyle displaystyle="true"><munderover><mo largeop="true">∑</mo><mrow><mi>i</mi><mo>=</mo><mn mathvariant="normal">1</mn></mrow><mrow><mi>n</mi></mrow></munderover></mstyle><mrow><mfenced open="|" close="|" separators="|"><mrow><mfrac><mrow><msub><mrow><mover accent="true"><mi>y</mi><mo>̂</mo></mover></mrow><mrow><mi>i</mi></mrow></msub><mo>-</mo><msub><mrow><mi>y</mi></mrow><mrow><mi>i</mi></mrow></msub></mrow><mrow><msub><mrow><mi>y</mi></mrow><mrow><mi>i</mi></mrow></msub></mrow></mfrac></mrow></mfenced></mrow></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=22213748&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=22213746&type= 40.72467041 10.24466610 式中， n 表示总数， <math id="M5"><msub><mrow><mover accent="true"><mi>y</mi><mo>̂</mo></mover></mrow><mrow><mi>i</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=22213751&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=22213750&type= 2.37066650 4.57200003 表示校正后影像对应的反射率或遥感反射率， <math id="M6"><msub><mrow><mi>y</mi></mrow><mrow><mi>i</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=22213754&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=22213752&type= 2.37066650 4.57200003 表示实测等效后的反射率或遥感反射率。 结果如 图6 图6 所示。从 图6 图6 中可以看出，3种典型地物多期大气校正后高分的反射率曲线基本一致，数值相差不大，且符合各地物的反射率特征，表明3种典型地物大气校正的效果较好。对于水体，符合典型内陆水体的反射率特征，即绿光波段最高，蓝光和红光波段次之，近红外波段最低。对于植被，经过对比发现2018年1月、4月初和11月、2019年2月和11月、2020年3月符合非生长季植被冠层的反射率特征，即蓝光最低、红光略低于绿光或与绿光持平、近红外最高。2018年9月、2019年5月和8月、2020年5月和8月符合生长季的植被反射率特征，即蓝光波段最低、红光波段低于绿光波段、近红外波段最高；对于屋顶，具有蓝光波段高、绿光和红光波段次之、近红外波段最高的特征，符合蓝顶建筑物的反射率特征。 1,2;3,4;5; 相对辐射一致化方法和参考影像共同完成大气校正仍然存在如下问题：（1）在一些地物变化比较大的区域，不能够有效获取覆盖整个辐射动态范围的伪不变点集合，容易造成回归结果出现偏差；（2）容易受到植被物候的影响，需要尽可能选择和待校正影像时相接近的参考影像；（3）容易受到薄云影响，使得自动化选取伪不变地物存在困难；（4）算法假设影像中所有地物具有统一的线性变换方程，不能够适应大场景内大气条件多变的情况。 5　结语 1 地表反射率产品作为最重要的定量遥感共性产品，是很多参量化遥感产品的基础数据源，可以被广泛应用于林业、农业、水资源、生态环境、城市环境等典型应用领域。该数据集是0.8 m空间分辨率2015年—2019年每年覆盖一遍北京市平原区范围的GF-2地表反射率影像产品。当前的数据集，对于山体阴影、较为清洁的水体如密云水库水体，反演的反射率往往是零值或负值。因此本数据集适合做除山区和较为清洁水体以外的下垫面应用。 本数据集考虑了高分系列中多光谱影像各个波段的精确配准，然后再与全色影像做复原融合处理，避免全色和多光谱融合后影像呈现出晕染现象；几何精校正考虑了区域网平差，保证了多景几何精校正影像在接边处也具有较好的几何一致性；大气校正基于参考图像和相对辐射一致化算法，便于开展业务应用。经验证，8种典型地物分别在校正后的GF-2影像和参考影像中的反射率都十分接近；采用与实测数据相差3 d之内的影像以及多时序高分影像评价大气校正效果，发现水体、道路、植被冠层3种典型地物的大气校正的效果都较好；水体、植被冠层、屋顶3种典型地物在多期大气校正后高分影像中的反射率折线基本一致，且符合各地物的反射率特征，可以满足实际应用需求。 在实际业务应用中，准备好参考影像的前提下，本文的处理方法可以快速、大批量、自动化地处理高分影像，计算输出地表反射率产品。虽然当前仍无法保证大气校正具有100%的成功率，但是约占总数70%的影像可以保证待处理影像与参考影像地物反射率的波形和数值基本一致。剩余30%影像需要人工有针对性的进行替换和尝试后才能输出满意结果，这主要是由于现实情况分布着复杂的下垫面，有时给伪不变地物寻找带来挑战。下一步计划是对于输出的结果图像采用自动化检查替换人工检查，从而自动筛选出合格的结果。 参考文献（References） Canty M J , Nielsen A A and Schmidt M . 2004 . Automatic radiometric normalization of multitemporal satellite imagery . Remote Sensing of Environment , 91 ( 3/4 ): 441 - 451 [ DOI: 10.1016/j.rse.2003.10.024 http://dx.doi.org/10.1016/j.rse.2003.10.024 ] Automatic radiometric normalization of multitemporal satellite imagery Cao X J , Jiang H , Zhang Z M , He G J , Wang M M , Long T F , Wang G Z and Zhang X M . 2019 . Research progress of Landsat series satellite surface reflectance products . Geomatics and Spatial Information Technology , 42 ( 3 ): 93 - 96 Research progress of Landsat series satellite surface reflectance products 曹小杰 , 江洪 , 张兆明 , 何国金 , 王猛猛 , 龙腾飞 , 王桂周 , 张晓美 . 2019 . Landsat系列卫星地表反射率产品研究进展 . 测绘与空间地理信息 , 42 ( 3 ): 93- 96 [ DOI: 10.3969/j.issn.1672-5867.2019.03.025 http://dx.doi.org/10.3969/j.issn.1672-5867.2019.03.025 ] Landsat系列卫星地表反射率产品研究进展 Jr Chavez P S , Sides S C and Anderson J A . 1991 . Comparison of three different methods to merge multiresolution and multispectral data: Landsat TM and SPOT panchromatic . Photogrammetric Engineering and Remote Sensing , 57 ( 3 ): 295 - 303 Comparison of three different methods to merge multiresolution and multispectral data: Landsat TM and SPOT panchromatic Chen F . 2021 . Pixel-Knife [EB/OL]. [ 2020-06-07 ]. https://www.zybuluo.com/novachen/note/426294 ( https://www.zybuluo.com/novachen/note/426294( Pixel-Knife 陈甫 . 2021 . 像素刻刀高分卫星处理软件 [EB/OL]. [ 2020-06-07 ]. https://www.zybuluo.com/novachen/note/426294 https://www.zybuluo.com/novachen/note/426294 像素刻刀高分卫星处理软件 Doxani G , Vermote E , Roger J C , Gascon F , Adriaensen S , Frantz D , Hagolle O , Hollstein A , Kirches G , Li F Q , Louis J , Mangin A , Pahlevan N , Pflug B and Vanhellemont Q . 2018 . Atmospheric correction inter-comparison exercise . Remote Sensing , 10 ( 2 ): 352 [ DOI: 10.3390/re10020352 http://dx.doi.org/10.3390/re10020352 ] Atmospheric correction inter-comparison exercise European Space Agency . 2015 . Sentinel-2 User Handbook, Issue 1 , Revision 2 . Sentinel-2 User Handbook, Issue 1 Grodecki J and Dial G . 2003 . Block adjustment of high-resolution satellite images described by rational polynomials . Photogrammetric Engineering and Remote Sensing , 69 ( 1 ): 59 - 68 [ DOI: 10.14358/PERS.69.1.59 http://dx.doi.org/10.14358/PERS.69.1.59 ] Block adjustment of high-resolution satellite images described by rational polynomials Kanzow C , Yamashita N and Fukushima M . 2004 . Levenberg-Marquardt methods with strong local convergence properties for solving nonlinear equations with convex constraints . Journal of Computational and Applied Mathematics , 172 ( 2 ): 375 - 397 [ DOI: 10.1016/j.cam.2004.02.013 http://dx.doi.org/10.1016/j.cam.2004.02.013 ] Levenberg-Marquardt methods with strong local convergence properties for solving nonlinear equations with convex constraints Laben C A and Brower B V . 2000 . Process for enhancing the spatial resolution of multispectral imagery using pan-sharpening. U.S ., No . 6011875 Process for enhancing the spatial resolution of multispectral imagery using pan-sharpening. U.S Liang S L , Bai R , Chen X N , Cheng J , Fan W J , He T , Jia K , Jiang B , Jiang L M , Jiao Z T , Liu Y B , Ni W J , Qiu F , Song L L , Sun L , Tang B H , Wen J G , Wu G P , Xie D H , Yao Y J , Yuan W P , Zhang Y G , Zhang Y Z , Zhang Y T , Zhang X T , Zhao T J and Zhao X . 2020 . Review of China’s land surface quantitative remote sensing development in 2019 . Journal of Remote Sensing , 24 ( 6 ): 618 - 671 Review of China’s land surface quantitative remote sensing development in 2019 梁顺林 , 白瑞 , 陈晓娜 , 程洁 , 范闻捷 , 何涛 , 贾坤 , 江波 , 蒋玲梅 , 焦子锑 , 刘元波 , 倪文俭 , 邱凤 , 宋柳霖 , 孙林 , 唐伯惠 , 闻建光 , 吴桂平 , 谢东辉 , 姚云军 , 袁文平 , 张永光 , 张玉珍 , 张云腾 , 张晓通 , 赵天杰 , 赵祥 . 2020 . 2019年中国陆表定量遥感发展综述 . 遥感学报 , 24 ( 6 ): 618- 671 [ DOI: 10.11834/jrs.20209476 http://dx.doi.org/10.11834/jrs.20209476 ] 2019年中国陆表定量遥感发展综述 Liang S L , Cheng J , Jia K , Jiang B , Liu Q , Liu S H , Xiao Z Q , Xie X H , Yao Y J , Yuan W P , Zhang X T and Zhao X . 2016 . Recent progress in land surface quantitative remote sensing . Journal of Remote Sensing , 20 ( 5 ): 875 - 898 Recent progress in land surface quantitative remote sensing 梁顺林 , 程洁 , 贾坤 , 江波 , 刘强 , 刘素红 , 肖志强 , 谢先红 , 姚云军 , 袁文平 , 张晓通 , 赵祥 . 2016 . 陆表定量遥感反演方法的发展新动态 . 遥感学报 , 20 ( 5 ): 875- 898 [ DOI: 10.11834/jrs.20166258 http://dx.doi.org/10.11834/jrs.20166258 ] 陆表定量遥感反演方法的发展新动态 Liu L Y and Zhang X . 2017-10-20 . Landsat- 8 land surface reflectance products of China in 2015 . Landsat- The Institute of remote sensing and digital earth , Chinese Academy of Sciences 刘良云 , 张肖 . 2017-10-20 . 2015年全国陆地区域Landsat-8地表反射率产品 . 中国科学院遥感与数字地球研究所) [DOI: 10.12237/casearth.5c19a5690600cf2a3c557bd3] 2015年全国陆地区域Landsat-8地表反射率产品 Long T F . 2016 . Efficient Positioning Technology for Instant Satellite Image Service. Doctoral dissertation, University of Chinese Academy of Sciences . Beijing : University of Chinese Academy of Sciences : 83 (龙腾飞 . 2016. 面向即时卫星影像服务的快速几何定位技术. 北京 : 中国科学院大学 : 83 ) (龙腾飞 Long T F , Jiao W L , He G J and Zhang Z M . 2016 . A fast and reliable matching method for automated georeferencing of remotely-sensed imagery . Remote Sensing , 8 ( 1 ): 56 [ DOI: 10.3390/rs8010056 http://dx.doi.org/10.3390/rs8010056 ] A fast and reliable matching method for automated georeferencing of remotely-sensed imagery Richter R . 1998 . Correction of satellite imagery over mountainous terrain . Applied Optics , 37 ( 18 ): 4004 - 4015 [ DOI: 10.1364/AO.37.004004 http://dx.doi.org/10.1364/AO.37.004004 ] Correction of satellite imagery over mountainous terrain Richter R and Schläpfer D . 2015 . Atmospheric/Topographic Correction for Satellite Imagery . DLR Report DLR-I B565 - 02 / 15 . Wessling, Germany: German Aerospace Center (DLR) Atmospheric/Topographic Correction for Satellite Imagery Shen Q , Yao Y , Li J S , Zhang F F , Wang S L , Wu Y H , Ye H P and Zhang B . 2019 . A CIE color purity algorithm to detect black and odorous water in urban rivers using high-resolution multispectral remote sensing images . IEEE Transactions on Geoscience and Remote Sensing , 57 ( 9 ): 6577 - 6590 [ DOI: 10.1109/TGRS.2019.2907283 http://dx.doi.org/10.1109/TGRS.2019.2907283 ] A CIE color purity algorithm to detect black and odorous water in urban rivers using high-resolution multispectral remote sensing images Sun W H , Chen B and Messinger D . 2014 . Nearest-neighbor diffusion-based pan-sharpening algorithm for spectral images . Optical Engineering , 53 ( 1 ): 013107 [ DOI: 10.1117/1.OE.53.1.013107 http://dx.doi.org/10.1117/1.OE.53.1.013107 ] Nearest-neighbor diffusion-based pan-sharpening algorithm for spectral images Vermote E , Justice C , Claverie M and Franch B . 2016 . Preliminary analysis of the performance of the Landsat 8/OLI land surface reflectance product . Remote Sensing of Environment , 185 : 46 - 56 [ DOI: 10.1016/j.rse.2016.04.008 http://dx.doi.org/10.1016/j.rse.2016.04.008 ] Preliminary analysis of the performance of the Landsat 8/OLI land surface reflectance product Wu X D , Xiao Q , Wen J G , Liu Q , Peng J J and Li X W . 2014 . Advances in uncertainty analysis for the validation of remote sensing products: take leaf area index for example . Journal of Remote Sensing , 18 ( 5 ): 1011 - 1023 Advances in uncertainty analysis for the validation of remote sensing products: take leaf area index for example 吴小丹 , 肖青 , 闻建光 , 刘强 , 彭菁菁 , 李小文 . 2014 . 遥感数据产品真实性检验不确定性分析研究进展 . 遥感学报 , 18 ( 5 ): 1011- 1023 [ DOI: 10.11834/jrs.20143332 http://dx.doi.org/10.11834/jrs.20143332 ] 遥感数据产品真实性检验不确定性分析研究进展 10.11834/jrs.20210407 1007-4619（2021）11-2303-10 申茜，姚月，李利伟，龙腾飞，陈甫，张兵.2021.北京市平原区2015年—2019年0.8 m地表反射率数据集.遥感学报，25（11）： 2303-2312 申茜，姚月，李利伟，龙腾飞，陈甫，张兵.2021.北京市平原区2015年—2019年0.8 m地表反射率数据集.遥感学报，25（11）： 2303-2312 Shen Q， Yao Y， Li L W， Long T F， Chen F and Zhang B. 2021. Annual 0.8 m surface reflectance data set of Beijing plain area from 2015 to 2019. National Remote Sensing Bulletin， 25（11）：2303-2312 Shen Q， Yao Y， Li L W， Long T F， Chen F and Zhang B. 2021. Annual 0.8 m surface reflectance data set of Beijing plain area from 2015 to 2019. National Remote Sensing Bulletin， 25（11）：2303-2312 国家自然科学基金（编号：42030111，41571361）；中国科学院先导专项（编号：XDA23040102） Supported by National Natural Science Foundation of China (No. 42030111, 41571361); Strategic Priority Research Program Project of the Chinese Academy of Sciences (No. XDA23040102) 2020-09-20 2021-11-07 2021-11-26 2021-07-05 遥感学报 11 25