利用无人机搭载的多光谱相机直接测定地表反射率的新方法

引用本文

孙华生, 张远, 史云飞, 赵敏. 利用无人机搭载的多光谱相机直接测定地表反射率的新方法[J]. 光谱学与光谱分析, 2022,42(5): 1581-1587.
SUN Hua-sheng, ZHANG Yuan, SHI Yun-fei, ZHAO Min. A New Method for Direct Measurement of Land Surface Reflectance With UAV-Based Multispectral Cameras[J]. Spectroscopy and Spectral Analysis, 2022,42(5): 1581-1587.
Doi:10.3964/j.issn.1000-0593(2022)05-1581-07 复制到剪切板

Permissions

《光谱学与光谱分析》期刊社所有

利用无人机搭载的多光谱相机直接测定地表反射率的新方法

孙华生¹, 张远^2,^*, 史云飞¹, 赵敏¹

1.临沂大学资源环境学院(山东省水土保持与环境保育研究所), 山东临沂 276000

2.华东师范大学地理科学学院, 崇明生态研究院, 地理信息科学教育部重点实验室, 上海 200241

*通讯作者 e-mail: yzhang@geo.ecnu.edu.cn

作者简介: 孙华生, 1980年生, 临沂大学资源环境学院副教授 e-mail: sunhuasheng@126.com

收稿日期: 2021-04-07 修回日期: 2021-06-10

基金: 上海市科技兴农项目(沪农科创字(2018)第2-1号), 国家自然科学基金项目(41571410, 41977067)资助

摘要

目前, 无人机获取的多光谱遥感数据已被广泛应用于农业、林业、环境等领域的定量监测中。然而, 现有的将多光谱遥感数据转换为地表反射率的方法, 仍然存在一定的缺陷, 如需要依赖地面参考板、无法适应光照条件变化、得到的结果不准确等, 从而影响了多光谱遥感数据定量化应用的效果。为了解决该问题, 提出了一种可以利用无人机搭载的多光谱相机, 直接对地表反射率进行测量的新方法。该方法具有非常强的适应能力, 即使在环境光照强度变化的条件下, 仍然能够得到准确的地表反射率。其中, 如何利用倾斜状态下的光强传感器获取准确的太阳辐照度, 是需要解决的关键问题。对此, 提出了一种利用两个或者更多朝向不同方向的光强传感器, 实现太阳直射和散射辐照度测量的新方法。利用此方法即可将相机记录的数字量化(DN)值直接转换为地表反射率。为了验证本方法的实际效果, 设计了具体的实验验证方案, 对不同日期不同光照条件下获取的无人机遥感数据进行验证。实际测试结果表明: 利用该方法, 得到黑、灰、白三张参考板的反射率在5个多光谱(蓝、绿、红、红边和近红外)波段中最大的平均绝对误差为3.34%, 其对应的标准差为2.11%; 三张参考板在所有波段中最大的平均绝对误差为2.94%, 其对应的标准差为1.84%。由此可见, 在光照强度变化的条件下, 利用该方法实现地表反射率的准确测量是可行性的。该方法极大地简化了无人机遥感数据转换为地表反射率的过程。对多光谱无人机的设计, 以及无人机遥感数据的定量化应用, 都具有重要的参考价值。

关键词: 无人机遥感; 多光谱相机; 地表反射率; 辐射校正; 太阳辐照度

中图分类号:TP73 文献标志码:A

A New Method for Direct Measurement of Land Surface Reflectance With UAV-Based Multispectral Cameras

SUN Hua-sheng¹, ZHANG Yuan^2,^*, SHI Yun-fei¹, ZHAO Min¹

1. Shandong Provincial Key Laboratory of Soil and Water Conservation and Environmental Protection, School of Resources and Environment, Linyi University, Linyi 276000, China

2. Key Laboratory of Geographic Information Science (Ministry of Education), School of Geographic Sciences & Institute of Eco-Chongming, East China Normal University, Shanghai 200241, China

*Corresponding author

Abstract

At present, the multispectral remote sensing data obtained by UAV (Unmanned Aerial Vehicle), have been widely used in the quantitative monitoring of agriculture, forest, environment and other fields. However, the existing methods of converting the multispectral remote sensing data into land surface reflectance still have some defects, e.g. relying on reference boards, unable to adapt to varying solar illumination conditions, or the obtained results are inaccurate, etc., so they may affect the quantitative application effect of the remote sensing data. In order to solve this problem, a novel methodology for direct measurement of land surface reflectance with the multispectral camera is proposed in this study. The method has strong adaptability, and it can get accurate land surface reflectance even under the condition of variable illumination. The key problem that needs to be solved is how to use irradiance sensors to obtain accurate solar irradiance in a tilted state. A new method was proposed to separate the direct and scattering irradiance with two or more irradiance sensors oriented towards different directions to achieve this objective. Therefore, the digital numbers (DNs) recorded by the multispectral bands can be converted into accurate reflectance with the measured irradiance results. UAV-based remote sensing images obtained on the different dates under different illumination conditions were used to validate the actual effect through the validation scheme designed in this study. The experimental results showed that the maximum MAE (mean absolute error) and standard deviation of five multispectral bands (i.e. blue, green, red, red edge and near infrared) for all reference panels is 3.34% and 2.11%, respectively; and the maximum MAE and standard deviation of the three reference panels (i.e. black, gray and white) is 2.94% and 1.84%, respectively. Therefore, the proposed method can obtain accurate reflectance results, even under varying solar illumination conditions. It greatly simplifies the reflectance conversion process of UAV-based remote sensing images. The study results are significant to the UAV-based remote sensing system's design and the quantitative application of multispectral remote sensing data.

Keyword: UAV-based remote sensing; Multispectral camera; Land surface reflectance; Radiometric calibration; Solar irradiance

文章图片

引言

无人机遥感具有操作简便、数据获取成本低、时空分辨率高、可在云下实现数据采集等诸多优势, 弥补了卫星或大型航空遥感系统的缺点, 所以在小范围、高时间和高空间分辨率遥感监测中, 发挥了极为重要的作用。目前, 无人机遥感技术已被广泛应用于农业、林业、环境领域的定量监测^[1]。反射率是地物较为稳定的物理属性, 也是定量遥感模型的重要输入参数之一。在遥感图像定量化应用过程中, 需要将相机记录的数字量化(digital number, DN)值转换为地表反射率。因此, 得到准确的地表反射率是实现遥感数据定量化应用的基础。

在环境光照稳定的条件下, 通常只要在无人机起飞前或航拍完成后, 拍摄两个或更多具有不同反射率的参考板, 通过对每个波段记录的DN值进行线性拉伸, 即可得到每个波段的地表反射率^{[2, 3, 4]}。然而, 实际的状况却是异常复杂的^{[5, 6, 7]}, 这是因为: (1)受云层遮挡和太阳高度角变化的影响, 数据获取时的光照条件可能会随时发生变化; (2)为了得到高质量的图像, 相机曝光时间、光圈大小和ISO值, 需要根据拍摄时的光照条件自动调整; (3)在很多情况下, 在地面上铺设参考板是很难或无法实现的(如森林、水域或沼泽地中), 而且要想在每幅图像中都拍摄到一组参考板也是不现实的。

为了摆脱对参考板的依赖, 一些研究者采用经验线法得到近似的地表反射率^[2]; 还有的建议在无人机上安装一个光强传感器, 来同步记录相机拍摄时每个波段对应的光照强度, 并将该测量结果用于计算地表反射率^{[3, 4, 5, 7, 8, 9]}, 所以获取可靠的地面接收的太阳辐照度是得到准确地表反射率的前提。如果光强传感器不是安装在一个精密的增稳平台上, 受飞行器姿态变化的影响, 导致传感器的测量结果与到达地面的太阳辐射强度并不一致。对此, 一些研究者采用区块校正法来降低传感器的倾斜效应^{[8, 9, 10]}, 但其得到的只是一种近似的结果。还有采用余弦函数对直接测量的结果进行倾斜校正, 但实验测试表明, 余弦校正的结果并不理想^{[3, 6]}。这是因为: 太阳辐射不仅包括直射, 还包括散射, 在不同天气状况下, 太阳辐射的组成有着很大的差异, 而传感器倾斜对太阳直射和散射的影响是完全不同的。为了获取更加准确的太阳辐照度测量结果, 可将光强传感器安装在一个自动调平装置上。例如, Markelin等利用一个光学平衡系统, 可以对光强传感器进行15° 以内的倾斜补偿^[4], 但对多旋翼无人机来说, 在飞行状态下, 其倾斜角度有时甚至会超过30° 。如果只是将光强传感器搭载在一个简易的调平装置上, 并不能完全消除飞行器倾斜造成的影响; 而如果单独为光强传感器安装一个精密的调平装置(如一个三轴增稳云台), 将会大大增加飞行器的负担和制造成本。在实际应用中, 在没有精密调平装置的条件下, 为了使光强传感器的测量结果更加可靠, 常用的做法包括: (1)要求航拍时太阳高度角较大, 并且航线与太阳方位角垂直, 以降低传感器倾斜对太阳辐照度测量结果的影响, 但这种方法并不能消除传感器倾斜的影响; (2)采用悬停方式进行拍摄, 但这种方式只有在没有风的情况下才能使传感器保持水平, 否则, 并不能完全消除传感器倾斜的影响, 而且其作业效率极低。因此, 目前利用无人机光强传感器对太阳辐照度的测量, 仍然存在一定的缺陷, 或者只能在一定的限定条件下使用。

如何在倾斜状态下得到准确的太阳辐照度, 进而得到准确的地表反射率, 是无人机遥感领域多年来一直未能很好解决的技术难题。本研究对该问题进行了深入探讨, 并提出了一种太阳辐照度准确测量的新方法。该方法可以在任意光照条件下, 直接将多光谱图像记录的DN值转换为准确的地表反射率。

1 实验部分

1.1 无人机平台和传感器

利用大疆精灵4 RTK无人机作为遥感数据获取平台。该无人机搭载了一台一体式的多光谱相机, 共有6个1/2.9 英寸的CMOS传感器(详见图1), 包含1个可见光波段的RGB传感器和5个多光谱传感器。 5个多光谱波段参数的描述, 如表1所示。此外, 在该无人机的顶部还集成了一个光强传感器, 可同步检测5个波段的太阳辐照度, 并将其记录在图像的XMP(Extensible Metadata Platform)元数据中^[11]。

	Figure Option View Download New Window
	图1 大疆精灵4 RTK多光谱无人机及其成像系统Fig.1 DJI Phantom 4 RTK multispectral UAV and its imaging system

表1 大疆精灵4 RTK多光谱无人机的传感器参数 Table 1 Sensor parameters of DJI Phantom 4 RTK multispectral UAV

1.2 多光谱图像转换为地表反射率的方法

地表反射率的计算方法如式(1)所示

$R_{λ} = E_{r} (λ) / E_{g} (λ)$ (1)

式(1)中, R_λ为波长为λ 的地表反射率, E_r(λ )为地物对波长为λ 的电磁波在单位波长单位面积的反射辐射功率(W· m^-2· nm^-1), E_g(λ )为到达地面的入射辐射辐照度(W· m^-2· nm^-1)。

对E_r(λ )的定量描述如式(2)所示

$E_{r} (λ) = \int_{0}^{2} \int_{0}^{π} L (λ) \cos θ \sin θ d θ d φ = πL (λ)$ (2)

式(2)中, L为每个波段每个像元对应的地物反射辐照度(W· m^-2· Sr^-1· nm^-1), θ 为太阳辐射的天顶角, φ 为太阳辐射的方位角。

在计算反射率时, 需要首先将DN值转换为对应像元的辐射亮度, 然后再将辐射亮度进一步转换为地表反射率。由于无人机多光谱相机拍摄的原始图像记录的DN值会受曝光时间、 ISO值、暗角效应等因素的影响, 所以需要事先对原始图像记录的DN值进行修正处理, 以得到统一标准下的结果。具体转换方法如式(3)所示

$DN' = (DN - D N_{black}) / N / G_{sensor} / T_{e} \times A \times V$ (3)

式(3)中, DN'为转换后的每个波段的数字量化值, DN为原始图像记录的每个波段的数字量化值, DN_black为黑电平的数字量化值, N为图像可表示的最大灰度值, G_sensor为传感器增益值, T_e为曝光时间, A为传感器敏感度修正系数, V为相机镜头的暗角效应修正值。以上参数都是在实验室中测定的, 其中DN_black, G_sensor, T_e和A可以从XMP元数据中查找到, 而暗角效应修正值V可按照式(4)计算得到

$V (x, y) = (1 + k_{1} r + k_{2} r^{2} + k_{3} r^{3} + k_{4} r^{4} + k_{5} r^{5} + k_{6} r^{6})$ (4)

其中, r= $\sqrt[]{(x - x_{0})^{2} + (y - y_{0})^{2}}$ , x和y分别为每个像元在图像上的列数和行数, x₀和y₀分别为图像中心像元的列数和行数, 而k₁, k₂, k₃, k₄, k₅, k₆也可以从XMP中查找到。

修正后的DN'(λ )与每个波段的辐射亮度L(λ )之间满足简单的线性转换关系, 其表达式如式(5)所示

$L (λ) = G_{λ} \times DN' (λ) + B_{λ}$ (5)

式(5)中, G_λ为增益值, B_λ为偏移值, DN'为修正后的数字量化值。

因此, 根据式(1)、式(2)和式(5), 地表反射率R_λ可表示为式(6)的形式

$\begin{array}{l} R_{λ} = E_{r} (λ) / E_{g} (λ) = π L (λ) / E_{g} (λ) = \\ π [G_{λ} \times DN' (λ) + B_{λ}] / E_{g} (λ) \end{array}$ (6)

如果能够准确地得到图像拍摄时到达地面的入射辐射辐照度E_g(λ ), 通过在一张图像中拍摄到两张或者更多具有不同反射率的参考板(其反射率是已知的), 即可得到每个波段的定标参数G_λ和B_λ, 其计算方法如式(7)所示

$DN' (λ) \times G_{λ} + B_{λ} = R_{λ} E_{g} (λ) / π$ (7)

由此可见, 利用倾斜的光强传感器获取准确的E_g(λ ), 是需要解决的关键问题, 而这也是长期以来一直未能很好解决的技术难题。在本研究中, 提出了一种实现太阳辐照度准确测量的新方法, 对其具体描述如下:

首先, 需要对光强传感器接收的太阳辐射进行定量化描述。其接收的太阳辐射包括: 太阳直射、天空散射和少量的地表反射。对太阳直射来说, 其辐照度测量结果受传感器倾斜的影响, 如式(8)所示

$E'_{dir} = E_{dir} cosz$ (8)

式(8)中, E'_dir为倾斜的光强传感器平面与直射入射方向有一定倾斜角度时的辐照度, E_dir为光强传感器平面与直射入射方向垂直时的辐照度, z为直射辐射在光强传感器平面上的入射角, z可以利用XMP中记录的信息计算出来[详见附注1中式(A4)的描述]。

光强传感器倾斜时接收的天空散射辐照度E'_sca受倾斜的影响, 如式(9)所示

$E'_{sca} = E_{sca} co s^{2} (\frac{s}{2})$ (9)

式(9)中, E'_sca为倾斜平面上接收的天空散射辐照度, E_sca为水平面上接收的天空散射辐照度, s为传感器倾斜平面的坡度, 而s也可以利用XMP中记录的信息计算出来[详见附注1中式(A3)的描述]。

倾斜平面接收的地面反射辐射E'_ref的定量描述如式(10)所示

$E'_{ref} = E_{ref} si n^{2} (s / 2) = {\bar{R}}_{g} E_{g} si n^{2} (s / 2)$ (10)

式(10)中, E_ref为地面反射的辐照度; ${\bar{R}}_{g}$ 为地面的平均反射率, 在可见光和近红外波谱区间, 通常在一般地面上可近似地令 ${\bar{R}}_{g}$ =0.2, 而在雪地上可近似地令 ${\bar{R}}_{g}$ =0.7; E_g为水平地面接收的太阳辐照度, E_g=E_dircosθ +E_dif, θ 为太阳直射在水平地面上的入射角, 而θ 也可以利用XMP中记录的信息计算出来[详见附注1中式(A1)的描述]。

根据以上分析可知, 在光强传感器平面倾斜时, 直接测量的每个波段的太阳辐照度E_m的表达式, 如式(11)所示

$\begin{aligned}E_{m}=& E_{\text {dir }}^{\prime}+E_{\text {sca }}^{\prime}+E_{\text {ref }}^{\prime}=E_{\mathrm{dir}} \cos z+E_{\mathrm{sca}} \cos ^{2}(s / 2)+\\& R_{\mathrm{g}}\left(E_{\mathrm{dir}} \cos \theta+E_{\mathrm{sca}}\right) \sin ^{2}(s / 2)=E_{\mathrm{dir}}[\cos z\\&\left.R_{\mathrm{g}} \cos \theta \sin ^{2}(s / 2)\right]+E_{\mathrm{sca}}\left[\cos ^{2}(s / 2)+R_{\mathrm{g}} \sin ^{2}(s / 2)\right]\end{aligned}$(11)

根据式(11)可知, 如果能够通过某种方法获取朝向两个或者更多方向的太阳辐照度, 就可以直接求出直射E_dir和散射E_sca。需要注意的是: 如果z≥ 90° , 即cos(z)≤ 0, 就会导致直射辐射无法到达光强传感器的正面, 从而使传感器正面处于阴影中。此时, 需要令cos(z)=0。

在得到直射E_dir和散射E_sca之后, 地面接收的太阳辐照度E_g可表示为式(12)的形式

$E_{g} = E_{dir} cosθ + E_{sca}$ (12)

在得出地面接收的太阳辐照度E_g后, 即可根据式(6)得到各个波段每个像素对应的地表反射率R_λ。

实际上, 如果能够得到直射辐射占总太阳辐射的比例 $p (p = E_{dir} / E_{sun}, 而 E_{sun} = E_{dir} + E_{sca}),$ 那么就能将光强传感器测量值$E_m$, 直接转换为地面接收的太阳辐照度$E_g$ 。把$E_{dir}=p E_{sun}$, $E_{sca}=(1-p) E_{sun}$ 代入式(11), 即可得到$E_m$ 和$E_sun$

之间的关系式, 其结果如式(13)所示

$\begin{array}{l} E_{m} = [p \cos z + (1 - p) co s^{2} (s / 2) + \\ {\bar{R}}_{g} (p \cos θ + 1 - p) si n^{2} (s / 2)] E_{s un} \end{array}$ (13)

因此, 地面接收的太阳辐照度E_g可以进一步转换为式(14)的形式

$\begin{array}{l} E_{g} = E_{dir} cosθ + E_{sca} = (p \cos θ + 1 - p) E_{sun} = \\ \frac{p \cos θ + 1 - p}{p \cos z + (1 - p) co s^{2} (s / 2) + {\bar{R}}_{g} (p \cos θ + 1 - p) si n^{2} (s / 2)} E_{m} \end{array}$ (14)

1.3 实验测试

1.3.1 实验区

实验区位于山东省临沂市, 其地理位置为35° 9'53.5″N— 35° 10'02.2″N, 118° 16'00.5″E— 118° 16'08.5″E。该实验区的面积约250 m× 200 m。航线规划是通过DJI GS Pro应用软件实现的, 设定的航高为35 m, 航线重叠度为67%, 旁向重叠度为55%, 共生成14条航线, 每条航线拍摄31张图像。飞行速度约为3.8 m· s^-1, 每次航拍约持续16 min。从2020年7月开始, 一直持续到11月底, 共获取了4个时期的航拍图像。

1.3.2 相机定标

选取黑、灰、白三张参考板(尺寸为18 cm× 25 cm)作为靶标, 以实现无人机多光谱图像的反射率定标和反射率测量结果的验证。各个参考板在不同波段的反射率通过PSR-1100高光谱仪测量得到。该高光谱仪的波谱覆盖范围为320~1 100 nm, 在600 nm处的光谱分辨率为3 nm。黑、灰、白参考板及其高光谱测试结果, 如图2所示。黑、灰、白三张参考板在不同波段的平均反射率, 如表2所示。

	Figure Option View Download New Window
	图2 黑、灰、白参考板及其高光谱测试结果Fig.2 The black, gray and white reference panels and their hyperspectral test results

表2 利用高光谱仪测量的参考板的平均反射率 Table 2 Average reflectances of the reference panels measured by a hyper spectrometer

在得到无人机拍摄的参考板的多光谱图像之后, 根据式(7)构建R_λ和DN'(λ )之间的线性拉伸关系, 即可得到每个波段的G_λ和B_λ的结果。对一台多光谱成像系统来说, 一旦得到其定标参数, 该定标参数通常是较为稳定的, 所以在实际应用中就不必再依靠参考板, 即可得到准确的地表反射率。

1.3.3 验证方案设计

如果飞行器自身携带两个或者更多的朝向不同方向的光强传感器, 即可根据式(11)实现直射E_dir和散射E_sca的求解。然而, 由于本研究采用的无人机只有一个光强传感器, 而且是固定在无人机顶部的, 所以无法直接实现以上解算。为了验证本方法的可行性, 设计了如下的实验方案:

首先, 选取一块开阔的实验场地, 确定其周围没有高大的树木或建筑物。然后, 用手举起飞行器并手动调整其姿态, 以让光强传感器平面朝向不同的方向(在操作时可以使传感器平面分别处于水平、前倾、后倾、左倾、右倾的状态), 并拍照以记录不同方向上各波段接收的太阳辐照度, 从而可实现对各个波段的直射E_dir和散射E_sca的求解。为了检验得到的地表反射率的准确度, 在实验场中还布置了一组黑、灰、白参考板。

分别在不同日期不同光照条件下获取航拍数据, 而每次实验时的天气都较为稳定。在这种情况下, 在一段较短的时间内, 可认为太阳辐射的组分是不变的(即使太阳辐照度有小幅度的变化, 但其组分的变化较小)。因此, 可利用本方法, 得到在航拍前或航拍后的直射比例p, 并根据式(14)计算出所有图像在拍照时刻地面接收的太阳辐照度E_g, 从而可根据该结果得到所有图像的地表反射率。最后, 利用参考板来验证计算结果的准确度。

2 结果与讨论

由于每次实验都是在光照条件较为稳定的情况下进行的(但不同日期的光照条件是不同的), 所以航拍时段的直射比例可以在航拍前或航拍后, 利用本方法来确定。航拍时段不同波段直射比例的计算结果, 如表3所示。

表3 航拍时段不同波段的直射比例p Table 3 Direct solar irradiance proportion (p) of each band during the image collection period

在得到直射比例p之后, 通过式(14)即可将光强传感器直接测量的结果E_m, 转换为地面接收的太阳辐照度E_g。不同日期的具体计算结果如图3所示。

Figure Option
View Download New Window

图3 光强传感器直接测量的太阳辐照度E_m及对其进行校正后得到的地面接收的太阳辐照度E_g的结果
(a): 2020-7-20 14:51— 15:07 E_m; (b): 2020-7-20 14:51— 15:07 E_g; (c): 2020-9-23 16:08— 16:24 E_m; (d): 2020-9-23 16:08— 16:24 E_g; (e): 2020-11-13 13:37— 13:53 E_m; (f): 2020-11-13 13:37— 13:53 E_g; (g): 2020-11-29 13:52— 14:09 E_m; (h): 2020-11-29 13:52— 14:09 E_gFig.3 Measured solar irradiance (E_m) by sensor and calculated the irradiance reaching ground (E_g)
(a): 2020-7-20 14:51— 15:07 E_m; (b): 2020-7-20 14:51— 15:07 E_g; (c): 2020-9-23 16:08— 16:24 E_m; (d): 2020-9-23 16:08— 16:24 E_g; (e): 2020-11-13 13:37— 13:53 E_m; (f): 2020-11-13 13:37— 13:53 E_g; (g): 2020-11-29 13:52— 14:09 E_m; (h): 2020-11-29 13:52— 14:09 E_g

由图3可以看出, 光强传感器的直接测量结果受飞行器自身姿态变化(因航向发生了改变)的影响很大, 尤其是当太阳高度角较小(即入射角较大)时, 其影响极为显著。利用本研究提出的方法可以很好地消除测量结果的剧烈波动。为了进一步检验该校正结果的准确度, 利用事先铺设的黑、灰、白参考板, 来检验无人机多光谱图像得到的反射率的准确度, 其结果如表4所示。

表4 利用无人机多光谱图像得到的参考板在不同波段的平均反射率 Table 4 Average reflectances of the reference panels derived from the UAV-based multispectral images

利用平均绝对误差(MAE), 以及标准差(S)两种指标对测量结果的准确度进行评价。计算结果显示, 蓝(B)、绿(G)、红(R)、红边(RE)和近红外(NIR)波段在所有日期所有参考板的MAE的结果分别为3.34%, 1.36%, 1.45%, 1.81%和1.63%, 而各自的标准差S分别为2.11%, 0.85%, 1.21%, 1.88%和0.59%。黑、灰、白三张参考板在所有日期的所有波段中的MAE的结果分别为0.82%, 1.99%和2.94%, 而各自的标准差S分别为0.55%, 1.31%和1.84%。引起地表反射率测量误差的因素主要包括以下几个方面: (1)参考板不是严格的漫反射体, 得到的反射率会因太阳入射方向和相机观测方向的不同而存在一定的误差; (2)参考板的高光谱测量结果, 以及得到的相机定标参数存在一定的误差; (3)光照条件不是绝对稳定的, 因为在航拍时段内太阳高度角和云量发生了一定的变化, 从而导致真实的直射比例p不是绝对稳定的。

尽管得到的地表反射率存在着一定的误差, 但从以上统计结果可以看出, 在不同光照条件下得到的不同波段的地表反射率, 仍然具有较高的准确度。由此可以推断出, 利用本方法得到的地面接收的太阳辐照度E_g的结果是可靠的。

3 结论

提出了一种全新的可以利用多光谱无人机图像实现地表反射率直接测量的方法。与传统的基于光强传感器实现太阳辐照度测量方法不同的是: 本方法不是直接消除光强传感器倾斜的影响, 而是借助传感器的倾斜效应来实现对光照强度的校正。即使在光照强度变化的条件下, 利用该方法仍然能够得到准确的测量结果, 从而可据此进一步得到准确的地表反射率。本研究解决了多年来一直未能得到很好解决关键技术问题, 对多光谱无人机的设计, 以及无人机遥感数据的定量化应用, 都具有重要意义。受硬件条件的限制, 本研究并没有直接利用无人机搭载的多个传感器, 来获取每个拍照时刻的太阳辐照度。在未来的研究中, 将对飞行器进行硬件上的改进, 使其携带多个朝向不同方向的光强传感器, 以实现每个拍照时刻的直射和散射辐照度的准确测量, 并检验其实际测量效果。

附注1: 几个角度的计算方法

(1)太阳直射的天顶角θ 可以通过式(A1)来确定

$\cos (θ) = \sin (Ψ) \sin (δ) + \cos (Ψ) \cos (δ) \cos (H)$ (A1)

式(A1)中, Ψ 为某一点的纬度(弧度); δ 为太阳赤纬, δ =-23.44/180π cos(2π /365(N+10)), N为某日期在一年中的天数; H为时角(以南为参考)。 H可以通过式(A2)计算出来

$H = π / 12 (T_{solar} - 12)$ (A2)

式(A2)中, T_solar为地方太阳时, 即根据太阳的具体位置所确定的时刻, 与地方标准时间T_Lst的关系为T_solar=T_Lst+12/π (Lon_loc-Lon_st)+E, Lon_loc为某一点的经度(弧度), Lon_st为地方标准时间所采用的标准经度(弧度), E为地球绕太阳公转时的进动和转速变化对地方太阳时的修正, E=[9.87sin(2B)-7.53cos(B)-1.5sin(B)]/60, 而B=2π (N-81)/364。

(2)传感器倾斜的坡度s可以利用式(A3)来确定

$\cos (s) = \frac{{v^{T}}_{Z} v_{n}}{| v_{Z} | | v_{n} |}$ (A3)

式(A3)中, v_Z为Z轴的方向向量(v_Z=(0, 0, 1)^T); v_n为斜面的法线向量, v_n=Rv_Z(R为旋转矩阵, 可以利用俯仰角、侧滚角和偏航角计算出来, 它们都记录在XMP文件中)。

(3)太阳直射在倾斜传感器坡面上的入射角z可以利用式(A4)来确定

$\cos (z) = \cos (θ) \cos (s) + \sin (θ) \sin (s) \cos (φ - a)$ (A4)

式(A4)中, θ 为太阳直射的天顶角; s为传感器倾斜的坡度; φ 为太阳方位角, a为传感器倾斜的坡向。而φ 和a可以分别通过式(A5)和式(A6)计算出来。注: φ 和a的参考方向需要与式(A2)中的时角H一致, 即以南为参考。

$φ = \arctan 2 (\sin φ, \cos φ)$ (A5)

式(A5)中, $\sin φ = \cos (δ) \sin (Η) / \sin (θ), \cos φ = \frac{\cos (θ) \sin (Ψ) - \sin (δ)}{\sin (θ) \cos (Ψ)} 。$

$a = \arctan 2 (\sin a, \cos a) + π$ (A6)

式(A6)中, $\sin a = v_{n (x)} / \sqrt[]{v_{n (x)}^{2} + v_{n (y)}^{2}}, \cos a = v_{n (y)} / \sqrt[]{v_{n (x)}^{2} + v_{n (y)}^{2}}$ , $v_{n(x)} $ 和$v_{n(y)} $ 分别为$v_n$ 的第1和第2个元素。注:当传感器处于绝对水平状态时 $, (即 v_{n (x)}^{2} + v_{n (y)}^{2} = 0),$ 则无坡向。

参考文献

文献列表

[1]	TmušiĆ G, Manfreda S, Aasen H, et al. Remote Sens. , 2020, 12(1001): 1. [本文引用:1]
[2]	von Bueren S K, Burkart A, Hueni A, et al. Biogeosciences, 2015, 12: 163. [本文引用:2]
[3]	Hakala T, Markelin L, Honkavaara E, et al. Sensors, 2018, 18(1417): 1. [本文引用:3]
[4]	Markelin L, Suomalainen J, Hakala T, et al. Int. Arch. Photogramm. Remote Sens. Spatial Inf. Sci. , 2018, XLII-1: 283. [本文引用:3]
[5]	Hakala T, Honkavaara E, Saari H, et al. Int. Arch. Photogramm. Remote Sens. Spatial Inf. Sci. , 2013, XL-1/W2: 189. [本文引用:2]
[6]	Honkavaara E, Hakala T, Markelin L, et al. Int. Arch. Photogramm. Remote Sens. Spatial Inf. Sci. , 2014, XL-1: 155. [本文引用:2]
[7]	Burkart A, Hecht V L, Kraska T, et al. Precis. Agric. , 2017, 19(1): 134. [本文引用:2]
[8]	Honkavaara E, Saari H, Kaivosoja J, et al. Remote Sens. , 2013, 5: 5006. [本文引用:2]
[9]	Miyoshi G T, Imai N N, Tommaselli A M G, et al. Int. J. Remote Sens. , 2018, 39(15-16): 4910. [本文引用:2]
[10]	Honkavaara E, Khoramshahi E. Remote Sens. , 2018, 10(2): 256. [本文引用:1]
[11]	DJ-Innovations(大疆创新科技有限公司, DJI). P4 Multispectral Image Processing Guide(P4 Multispectral图像处理指南), 2020. [本文引用:1]

2020

0.0

... 目前, 无人机遥感技术已被广泛应用于农业、林业、环境领域的定量监测^[1] ...

2015

0.0

... 在环境光照稳定的条件下, 通常只要在无人机起飞前或航拍完成后, 拍摄两个或更多具有不同反射率的参考板, 通过对每个波段记录的DN值进行线性拉伸, 即可得到每个波段的地表反射率^[2,3,4] ...

... 为了摆脱对参考板的依赖, 一些研究者采用经验线法得到近似的地表反射率^[2] ...

2018

0.0

... 还有的建议在无人机上安装一个光强传感器, 来同步记录相机拍摄时每个波段对应的光照强度, 并将该测量结果用于计算地表反射率^{[3,4,5, 7,8,9]}, 所以获取可靠的地面接收的太阳辐照度是得到准确地表反射率的前提 ...

... 还有采用余弦函数对直接测量的结果进行倾斜校正, 但实验测试表明, 余弦校正的结果并不理想^{[3, 6]} ...

2018

0.0

... 以内的倾斜补偿^[4], 但对多旋翼无人机来说, 在飞行状态下, 其倾斜角度有时甚至会超过30#cod#x000b0 ...

2013

0.0

... 然而, 实际的状况却是异常复杂的^[5,6,7], 这是因为: (1)受云层遮挡和太阳高度角变化的影响, 数据获取时的光照条件可能会随时发生变化 ...

2014

0.0

... 然而, 实际的状况却是异常复杂的^[5,6,7], 这是因为: (1)受云层遮挡和太阳高度角变化的影响, 数据获取时的光照条件可能会随时发生变化 ...

... 还有采用余弦函数对直接测量的结果进行倾斜校正, 但实验测试表明, 余弦校正的结果并不理想^{[3, 6]} ...

2017

0.0

... 然而, 实际的状况却是异常复杂的^[5,6,7], 这是因为: (1)受云层遮挡和太阳高度角变化的影响, 数据获取时的光照条件可能会随时发生变化 ...

2013

0.0

... 对此, 一些研究者采用区块校正法来降低传感器的倾斜效应^[8,9,10], 但其得到的只是一种近似的结果 ...

2018

0.0

... 对此, 一些研究者采用区块校正法来降低传感器的倾斜效应^[8,9,10], 但其得到的只是一种近似的结果 ...

2018

0.0

... 对此, 一些研究者采用区块校正法来降低传感器的倾斜效应^[8,9,10], 但其得到的只是一种近似的结果 ...

0.0

... 此外, 在该无人机的顶部还集成了一个光强传感器, 可同步检测5个波段的太阳辐照度, 并将其记录在图像的XMP(Extensible Metadata Platform)元数据中^[11] ...