提示:文章写完后,目录可以自动生成,如何生成可参考右边的帮助文档


文章目录

  • 前言
  • 一、辐射定标的定义
  • 二、DJI M3M/P4M的辐射定标原理(反射率的计算)
  • 1. 相机位置误差校正
  • 2. 相机曝光时间误差校正
  • 3. 黑电平校正
  • 4. 暗角补偿
  • 5. 畸变校正
  • Summary
  • 三、代码实现
  • 参考文献



Github:https://github.com/xhdidi/TIFF/blob/master/Radiometric_calibrations.py

前言

为什么要对航拍数据进行辐射定标

单波段相机获取的是DN值,即像元亮度值,其大小与传感器的辐射分辨率、地物发射频率、大气透过率、散射率有关。无人机载相机作为农业巡航设备,其值大小受太阳角度、双向反射分布函数效应、云的阴影、相机增益及曝光时间等因素影响,不能直接反映作物生长状况,因此需要转化为反射率,即物体表面反射的辐射量/接收的辐射量,才能保证所建立的检测模型在不同设备、不同时间和不同天气条件下的鲁棒性。


一、辐射定标的定义

传感器辐射定标分相对辐射定标绝对辐射定标两种,前者标定和纠正传感器不同探元间的相对响应差异,后者则是利用绝对参照体建立图像DN值(或灰度信息)与实际地物辐亮度的响应关系,将图像DN值转换成辐亮度、反射率等绝对地物信息(具有具体意义,如反射率信息)。

二、DJI M3M/P4M的辐射定标原理(反射率的计算)

图像的像素值代表阳光(入射光)照在地物(target)上的反射值,同时无人机顶端配有光强传感器以获取入射光的信号值。有研究表明,可见光图像的辐射定标方程是非线性的,而多光谱图像的辐射定标方程是线性的[2]。因此,某波段的反射率公式可定义如下,其中如何通过大疆S1学习python 大疆tello python_如何通过大疆S1学习python是调节图像信号与多光谱光强传感器信号之间相互转化的参数,使入射光(多光谱光强传感器信号值)与反射光(图像信号值)的量级(单位)保持统一。
如何通过大疆S1学习python 大疆tello python_反射率_02
由于光强传感器在不同波段存在感度差异,即在相同光源的条件下不同波段的相机和不同波段的多光谱光强传感器会得到不同的信号值,因此需要对相机和感光器各自做感度校准。
在此以NIR波段为基准,所有其他波段的相机都参照 NIR 波段的相机的感度做校准,得到校准参数如何通过大疆S1学习python 大疆tello python_python_03如何通过大疆S1学习python 大疆tello python_python_04,各波段与NIR的调节参数定义如下:
如何通过大疆S1学习python 大疆tello python_如何通过大疆S1学习python_05
则各波段的反射率定义可写作:
如何通过大疆S1学习python 大疆tello python_如何通过大疆S1学习python_06
如何通过大疆S1学习python 大疆tello python_python_03在元数据中保存为:Xmp.drone-dji.SensorGainAdjustment
如何通过大疆S1学习python 大疆tello python_python_08在元数据中保存为:Xmp.drone-dji.Irradiance
如何通过大疆S1学习python 大疆tello python_计算机视觉_09没有给出具体数值,需要参照植被指数计算软件如ENVI等进行手动线性校正,考虑到植被指数的数值本身并无意义,只是作为作物生长状态的特征表达,故此步骤可忽略。
下面对如何通过大疆S1学习python 大疆tello python_计算机视觉_10进行计算(图像校正)。

1. 相机位置误差校正

多个单波段相机拍摄的物理位置存在微小差异,为了能够将不同相机拍摄的图像相匹配(以计算植被指数),需要进行位置误差校正。


如何通过大疆S1学习python 大疆tello python_如何通过大疆S1学习python_11

如何通过大疆S1学习python 大疆tello python_光强_12

在此以NIR相机为基准,TIF文件的元数据保存有该波段相机相对于NIR波段相机的位置偏移量(像素),分别为

Xmp.drone-dji.RelativeOpticalCenterX , Xmp.drone-dji.RelativeOpticalCenterY,直接对原始图像进行padding后裁剪即可。

2. 相机曝光时间误差校正

不同相机的曝光时间存在差异,为了减轻曝光时间差异导致的位置偏移影响,首先对图像进行平滑处理(高斯滤波等),再对图像进行配准(sobel滤波器进行特征提取,Enhanced Correlation Coefficient Maximization做对齐)。

3. 黑电平校正

用归一化的像素值减去黑电平值即可,黑电平值在元数据里保存为Xmp.Camera.BlackCurrent,归一化前值为4096。

4. 暗角补偿

暗角是指画面四角有变暗的现象,从照片来看暗角其实就是指拍摄亮度均匀的场景时,画面四角却出现扇形向外延伸的渐暗区域。

如何通过大疆S1学习python 大疆tello python_计算机视觉_13

暗角补偿公式如下,其中如何通过大疆S1学习python 大疆tello python_python_14为像素坐标,如何通过大疆S1学习python 大疆tello python_python_15为当前像素到暗角补偿中心的距离,如何通过大疆S1学习python 大疆tello python_光强_16为暗角补偿的多项式系数:
如何通过大疆S1学习python 大疆tello python_反射率_17
如何通过大疆S1学习python 大疆tello python_如何通过大疆S1学习python_18
暗角补偿中心坐标在元数据中保存为:Xmp.drone-dji.CalibratedOpticalCenterX, Xmp.drone-dji.CalibratedOpticalCenterY
暗角补偿系数在元数据中保存为:Xmp.drone-dji.VignettingData

5. 畸变校正

相机的径向畸变和切向畸变模型定义如下:
如何通过大疆S1学习python 大疆tello python_光强_19
如何通过大疆S1学习python 大疆tello python_如何通过大疆S1学习python_20

像素坐标系与相机坐标系之间的转换由内参决定(相机参数),相机坐标系到世界坐标系的转换由外参决定(相机位姿)。因此畸变校准只需要计算相机内参。
相机内参定义如下:
如何通过大疆S1学习python 大疆tello python_反射率_21

如何通过大疆S1学习python 大疆tello python_如何通过大疆S1学习python_22在元数据中保存为Xmp.drone-dji.DewarpData,其中如何通过大疆S1学习python 大疆tello python_光强_23为相机参数,如何通过大疆S1学习python 大疆tello python_反射率_24为畸变校准的多项式系数,分别用于校正径向畸变和切向畸变。
采用cv2.undistort函数实现,校正原理在此不详细阐述。

Summary

对以上畸变进行校正,并考虑曝光增益、曝光时间,得到图像信号值定义如下:
如何通过大疆S1学习python 大疆tello python_反射率_25
其中如何通过大疆S1学习python 大疆tello python_计算机视觉_26包括暗角校正和畸变校正,
如何通过大疆S1学习python 大疆tello python_计算机视觉_27在元数据中保存为:Xmp.drone-dji.SensorGain
如何通过大疆S1学习python 大疆tello python_python_28在元数据中保存为:Xmp.drone-dji.ExposureTime

三、代码实现

class TIF_Image():
    def __init__(self, filename):
        self.filename = filename
        self.pyexiv2 = Image(self.filename)
        self.xmp = self.pyexiv2.read_xmp()
        self.exif = self.pyexiv2.read_exif()

        self.W, self.L = int(self.exif['Exif.Image.ImageWidth']), int(self.exif['Exif.Image.ImageLength'])

        t_raster = gdal.Open(filename)
        self.image = t_raster.ReadAsArray()  # <class 'numpy.ndarray'>

    def physical_position_calibration(self):
        '''相机物理位置误差校准'''
        pos_x = int(np.around(eval(self.xmp['Xmp.drone-dji.RelativeOpticalCenterX']), decimals=0, out=None))
        pos_y = int(np.around(eval(self.xmp['Xmp.drone-dji.RelativeOpticalCenterY']), decimals=0, out=None))
        
        # pos_x<0, pos_y>0
        image_pad = np.pad(self.image, ((0, pos_y), (-pos_x, 0)), 'constant', constant_values=(0, 0))
        self.image = image_pad[pos_y:self.L+pos_y, :self.W]

    def exposure_time_calibration(self):
        '''相机曝光时间误差校准(高斯滤波)'''
        self.image = cv2.GaussianBlur(self.image, (3,3), 0)

    def cam_correction(self):
        # 暗角补偿
        center_x = eval(self.xmp['Xmp.drone-dji.CalibratedOpticalCenterX'])
        center_y = eval(self.xmp['Xmp.drone-dji.CalibratedOpticalCenterY'])
        k0, k1, k2, k3, k4, k5 = self.xmp['Xmp.drone-dji.VignettingData'].split(', ')
        k0, k1, k2, k3, k4, k5 = eval(k0), eval(k1), eval(k2), eval(k3), eval(k4), eval(k5)

        max_distance = math.sqrt(
            max((vertex[0] - center_x) ** 2 + (vertex[1] - center_y) ** 2 for vertex in
                [[0, 0], [0, W], [L, 0], [L, W]]))
        for x in range(0, L):
            for y in range(0, W):
                distance = math.sqrt(pow(x - center_x, 2) + pow(y - center_y, 2))
                r = distance / max_distance
                k = k5 * r**6 + k4 * r**5 + k3 * r**4 + k2 * r**3 + k1 * r**2 + k0 * r + 1.0
                self.image[x, y] *= k

        # 畸变校准
        fx, fy, cx, cy, k1, k2, p1, p2, p3 = self.xmp['Xmp.drone-dji.DewarpData'].split(';')[-1].split(',')

        distcoeffs = np.float32([eval(k1), eval(k2), eval(p1), eval(p2), eval(p3)])
        cam_matrix = np.zeros((3, 3))
        cam_matrix[0, 0] = eval(fx)
        cam_matrix[1, 1] = eval(fy)
        cam_matrix[0, 2] = eval(cx) + center_x
        cam_matrix[1, 2] = eval(cy) + center_y
        cam_matrix[2, 2] = 1

        self.image = cv2.undistort(self.image, cam_matrix, distcoeffs)


    def cam(self):
        bit = eval(self.exif['Exif.Image.BitsPerSample'])
        self.image /= pow(2.0, bit)
        # black_current = eval(self.xmp['Xmp.Camera.BlackCurrent']) / 4096
        # self.image -= black_current

        self.cam_correction()

        pcam = eval(self.xmp['Xmp.drone-dji.SensorGainAdjustment'])
        gain = eval(self.xmp['Xmp.drone-dji.SensorGain'])
        etime = eval(self.xmp['Xmp.drone-dji.ExposureTime'])
        denominator = eval(self.xmp['Xmp.drone-dji.Irradiance'])

        self.image *= pcam/(gain*etime/1e6)/denominator

    def calibration(self):
        '''总校准函数'''
        self.physical_position_calibration()
        self.exposure_time_calibration_cv()
        self.cam()

参考文献

[1] P4M图像处理指南 [2] 章凌翔. 基于无人机多光谱影像油菜产量预测及倒伏风险评价[D].西北农林科技大学,2022.DOI:10.27409/d.cnki.gxbnu.2022.000191.