乔 凯,黄石生,智喜洋,孙 晅,赵 明
(1.北京跟踪通信技术研究所,北京 100094;2.哈尔滨工业大学 空间光学工程研究中心,黑龙江 哈尔滨 150001;3.大连海事大学 信息科学技术学院,辽宁 大连 116026)
高分辨率遥感卫星能够迅速、准确地获取地物信息,在军民应用领域发挥着重要作用。特别是以KH-12为代表的美国军用遥感卫星,地面分辨率优于0.1 m,以WorldView-3为代表的商业遥感卫星,可提供分辨率优于0.4 m全色图像[1-3]。
但现有卫星在轨图像质量仍不够理想,主要体现为图像模糊、有雾状感、边缘锐度不强,且存在动态范围偏窄、有效量化和有效信息不足的问题。立足于光学遥感全链路,通过图像复原的方法可消除相机光学系统、探测器、成像电子学以及卫星平台环境和大气路径等因素引起的图像模糊、噪声等像质退化。关于模糊、噪声处理算法的研究很多[4-7]。然而关于图像动态范围展宽方面,多是通过直方图调整、多曝光图像融合、多尺度图像融合等地面处理方式,以增强图像清晰度和对比度[8-11]。这些方法仍无法从根本上解决图像有效量化不足的问题,且图像中暗场景和边缘区域的处理性能会受到制约。实际上,目前相机动态范围(灰阶)可达到211以上,理论上能够满足通常成像条件下场景信息的获取需求,但国内在轨卫星(高分二号、吉林一号等)数据统计分析结果显示,有效量化仅达到了6~8位,这主要是由于卫星在轨运行过程中设置的成像参数和实时拍摄场景的动态范围存在不适配性。
因此,依据卫星在轨成像实时拍摄的动态场景进行相机成像参数的设计与动态调整,对进一步提升暗或亮场景细节信息的分辨能力,具有重要意义和实际工程应用价值。本文从卫星在轨成像动态范围的影响机理出发,研究卫星在轨动态场景实时匹配方法,同时利用无人机开展飞行试验验证。
图1 不同卫星图像直方图拉伸前后分布对比Fig.1 Distribution comparison of different satellite image histograms before and after stretching
从低动态范围卫星在轨图像直方图(如图1所示)拉伸前后结果来看,与高动态卫星图像相比,其卫星图像灰度层次不够丰富,灰度分布在相对较窄范围,而且直方图经拉伸后出现灰阶断层。这也是导致图像整体偏暗或高亮度场景区域出现饱和、景物纹理细腻度不够的主要原因。
光学卫星通常采用线阵推扫式TDICCD相机,通过多级TDI-CCD对地面同一场景进行多次曝光,从而获取高信噪比、高质量的图像。为了避免高动态场景区域图像灰度的过饱和,卫星在轨运行中一般采用固定的积分级数与增益[12-13]。设置增益的作用是使芯片输出电压与AD转换器电压相适应,从而提高输出量化位数。假定模数转换器的最大电压为2 V,曝光时间固定时,模拟电压为0.5~1.1 V,通过增益放大2倍可使CCD输出电压与AD转换器电压相符。
然而卫星在轨运行中,成像场景不断转变(可能由陆地变成海洋、森林变成沙漠),气候条件复杂多变(包含晴天、阴天),太阳高度角也在时刻发生改变,因此相机在不同时刻获取图像中的场景动态范围存在差异,采用固定的积分级数和增益参数无法满足动态场景实时匹配需求。在轨应用中,为了避免高亮度场景发生灰度饱和,通常将积分级数和增益参数控制在相对适中的范围以内,从而使场景信息集中在低灰度范围。这是造成图像整体偏暗、暗场景区域细节信息丢失的根本原因。
为了实现卫星在轨动态场景的实时匹配,需要在轨运行中实时解算成像场景的动态范围,并以此为依据合理调整相机积分级数和增益参数。而相机入瞳处的辐照度不仅包含场景本身的辐射能量,还包含大气程辐射,其在图像中表现为背景辐射,直接影响有效灰阶范围。由图2可见:大气程辐射造成的背景辐射已占1/3,图像整体偏暗、有雾状感,且多处高亮度场景区域出现灰度饱和。本文结合钳位技术在模拟域上消除背景辐射的影响,利用该技术可在高端动态保持不变的情况下,增大积分级数或增益,从而使相机的低端动态得到充分利用,有效提高有效灰阶。
图2 含大气程辐射的实拍图像(a)及其直方图分布(b)Fig.2 Actual image with atmospheric radiation(a) and corresponding histogram distribution(b)
此外,当图像中包含云场景时,由于云具有较高反射率,图像中表现为高亮度饱和区域[14],容易被误判为地面场景的高动态信息,从而影响场景动态范围的测量。
卫星在轨动态场景实时匹配方案如图3所示。由图3可知,为了实现相机在轨动态场景与实时拍摄的地面场景动态范围达到最佳匹配,采用测光-解算-成像的工作方式。首先采用独立的测光相机在t1时刻对地面场景进行成像,而后根据测光结果,结合大气程辐射预估和云检测获得实时拍摄场景的动态范围。最后成像相机根据场景动态调整积分级数和增益,并设置钳位值,在t2时刻对地面场景进行成像。
图3 卫星在轨动态场景实时匹配方案Fig.3 Real-time matching scheme for satellite dynamic scene
由3.1节可见:大气程辐射和云场景对地面场景低亮度和高亮度信息的测量存在干扰,在场景动态范围测量中需将二者剔除。下面给出大气程辐射预估和云检测方法。
图4 不同钳位电压值下的图像直方图Fig.4 Image histograms with different clamping voltage values
3.2.1 大气程辐射预估方法
通过对不同拍摄时间、不同地面场景和含不同程度大气程辐射效应的100幅图像进行统计特性分析得知,大气程辐射对图像特性的影响体现为直方图偏离原点。图4为不同钳位值时的图像直方图。由图4可见,通过增大钳位电压值可使直方图向左偏移,当钳位电压值达到某数值后,再增大钳位值会使暗场景区域的目标信息丢失(如图5),此时直方图的形状也将发生改变(如图4(d))。
图5 钳位电压值为500 mV的成像结果Fig.5 Imaging results when clamp voltage value is 500 mV
图6 钳位电压设置示意图Fig.6 Demonstration of clamping voltage selection
根据上述特征,本文提出基于直方图特性确定钳位电压值并在模拟域消除大气程辐射的方法。首先,寻找图像直方图的首个上升沿(如图6(a))所示)。然后,利用二次曲线拟合出上升沿函数,求出上升沿与横坐标轴的交点,从而确定大气程辐射对应的灰度值。最后,将该处的电压值设置为钳位电压,保证暗场景细节信息不丢失。为了提高方法的鲁棒性,从横坐标轴与拟合曲线的交点处沿曲线向上搜索,发现图像中包含暗点时判读其是否为连续点,若包含连续点则将该点的电压值设置为钳位电压。 图6(b)中竖线所示位置即为钳位电压值处。
3.2.2 云检测
云检测的总体思路如图7所示,主要包括特征分类训练与源模式判别两个部分。在特征分类训练中,本文选取具有代表性的云及其下垫面图像块作为训练样本。首先,基于云与下垫面地物场景的物理属性差异,构造并提取其特征参量,将样本映射为特征空间中的点[14]。然后,根据各参量对云和下垫面样本的分割程度选择特征参量,并利用特征压缩技术[15]进行特征空间压缩。最后,利用分类学习机[16]对样本进行分类训练,得到分类器。在源模式判别中,本文针对待检测图像进行特征参量提取与压缩,并输入分类器进行模式判别,完成云检测。
图7 云检测总体思路Fig.7 Overall scheme of cloud detection
3.2.3 场景动态范围解算
剔除大气程辐射和云对场景动态范围测量的影响后,即可解算被摄场景的动态范围,解算方法如图8所示,首先根据测光相机辐射标定关系和饱和辐亮度预估值,确定测光相机的曝光时间和增益。然后,分别根据测光相机对地面两次成像结果,统计出地面场景高亮度和低亮度信息,从而实现场景动态范围的测量。
图8 场景动态范围解算方法Fig.8 Scene dynamic range calculation method
本文采用测光相机对场景动态范围进行测量。由于测光相机与成像相机使用的传感器不同,为了使测光相机得到的动态范围可以被成像相机所用,必须首先对两台相机的辐射响应进行精确标定[17]。
(1)
然后根据实际成像的场景动态,并结合测光相机与成像相机标定关系,解算相机的成像参数,根据场景动态和相机动态范围关系,分为两种情况进行分析:
(1)当场景动态范围小于相机动态范围时,相机成像参数的解算公式为[12]:
(2)
(2)当场景动态范围大于相机动态范围时,又分为两种情况进行解算:
①高亮度匹配:高亮度匹配即是成像相机动态范围包含场景动态的最高亮度点,舍去部分的低亮度信息,保证高亮度信息的不丢失。匹配方法为令像面照度最大值处像素的输出信号恰好饱和,此时,参数解算方法与第(1)种情况一致。
②低亮度匹配:低亮度匹配即是相机动态范围包含场景动态的最低亮度点,舍去部分的高亮度信息,保证低亮度信息的不丢失。相机成像参数的解算公式为:
(3)
在t和gv的匹配优化设置中,考虑到增大gv会引起噪声放大,因此首先调整积分级数,仅在所需t超出系统允许上限时进行增益调整。
实际在轨应用中,场景动态范围通常大于相机动态范围,此时很难通过单次成像实现暗场景和亮场景细节信息的同时获取。可采用两次成像方案,分别采用高亮度匹配和低亮度匹配,并采用一种基于亮度参考模型的图像融合方法[18]对高匹配和低匹配图像进行融合处理。首先,将图像转换为YCbCr空间,对Y分量的图像进行滤波,得到亮度参考图像。根据亮度参考图像设置权重。权重设置原则为:给曝光适当的区域设置较大权重, 对过曝光或欠曝光区域设置较小权重。考虑到正常曝光的像素值在128左右,本文为其设定了较大权重。
图像灰阶范围和信息熵能够直接反映图像的有效动态范围和信息含量的丰富程度,因此本文采用它们对匹配方法的有效性和处理性能进行评价。
本文将图像的灰阶范围定义为:
Range=graymax-graymin ,
(4)
式中,Range表示灰阶范围,graymax表示去除图像中最亮10%像素后的最大灰度值,graymin为去除最暗10%像素后的最小灰度值。
图像信息熵定义[19]为:
(5)
式中,Pi表示图像中灰度值为i的像素所占的比例。
本文将匹配方法对图像灰阶和熵提升性能的评价指标定义为:
(6)
式中,Rate表示灰阶范围的平均提升值,N表示图像大小,GRfixed、GRLpath分别表示消除大气程辐射前、后图像的灰阶范围。
利用无人机搭载测光相机与成像相机,平台结构示意图和整机实物图如图9、图10所示,其中测光相机和成像相机均选用德国ximea高速相机和尼康标准变焦镜头,其焦距为24~120 mm,最小光圈值为22,传感器型号为MC050MG-SY。当飞行高度为100 m时,为了模拟卫星在轨工作状态,设置测光相机比成像相机预先观测10 m的目标场景,此时计算测光相机与成像相机的角度为5.7°。为了保持无人机稳定飞行,将其飞行速度设置为1 m/s,因此两次测光-拍照之间的间隔为10 s。
图9 无人机平台结构图Fig.9 Structure diagram of UAV platform
图10 无人机整机实物图Fig.10 Prototype of UAV
在设定无人机飞行轨迹时,设定一个目标点,令无人机按照1 m/s速度飞向目标点,在无人机飞行的过程中,利用其上搭载的测光相机获取目标场景图像,根据感兴趣区域解算第二次曝光时间,待无人机运行1 s后时开始第二次测光。在测光相机第二次测光时,改变开窗位置对同一目标场景成像。通过两次测光图像,测量目标场景的动态范围并确定成像相机的曝光时间和增益。10 s后成像相机到达测光相机测量的场景区域,在成像相机通过该区域时,利用动态场景实时匹配方法解算的曝光时间和增益采集地面场景图像。同时,在相同飞行轨迹下,利用成像相机对准待测目标场景,而后改变曝光时间和增益拍摄场景,采集不同曝光时间和增益下同一场景的序列图像,最后通过统计匹配图像和序列图像的灰阶范围与图像熵,说明匹配方法的有效性和处理性能。
2016年9月20日和2016年10月7日下午14:00(分别为晴天和阴天),选择10种不同类型的场景区域进行了无人机飞行试验。图11给出了不同曝光时间下场景1的序列图像。第1次曝光时间设定为0.1 ms,解算第2次曝光时间为2.061 ms。由图11(a)、11(b)可见,第1次成像尽管出现了暗场景信息丢失问题,但较好地保留了高动态的屋顶场景信息,而第2次成像出现了高动态场景信息丢失问题,但较好地保留了低动态场景信息。利用两次图像解算出实时匹配的成像相机曝光时间为0.424 ms。从图11 (c)和其它曝光时间下成像获取的序列图像对比看,利用本文提出的匹配方法可根据实时拍摄的地面景物合理地设置相机参数,实现相机与场景动态范围的最佳匹配。
图11 不同曝光时间下相同场景的序列图像Fig.11 Sequence images of the same scene with different integration times
利用测光相机对场景进行两次测光,而后成像相机根据匹配结果实时调整参数对同一场景进行成像的结果如图12所示。表1给出了匹配前后图像质量的评价结果。从图12(c)和12(f)和表1可见,动态场景匹配方法能够根据地面场景实时解算相机参数,匹配后图像的有效灰阶提升优于100%,信息熵提升优于40%。
图12 两次测光图像及动态匹配结果Fig.12 Two photometric images and corresponding dynamic matching results 表1 外场实验匹配前后的图像质量评价结果 Tab.1 Image quality evaluation results of field experiment before and after matching
图像灰阶图像熵固定曝光匹配后提升/%固定曝光匹配后提升/%场景124.2204.3744.23.9456.949076.1场景215.9153.8867.33.4966.688191.3场景333.4217.3550.64.0806.360755.9场景442.9232.7442.44.5717.119155.7场景529.5189.1541.04.4736.988956.2场景638.6209.0441.54.0436.439359.3场景755.2229.8316.35.0827.455846.7场景850.2182.6263.73.7626.011459.8场景961.9226.0265.13.8755.426140.0场景1016.0157.5884.43.6336.904790.1
为进一步验证方法的实际应用能力,本文在Visual Studio2013环境下开发了在轨动态场景实时匹配软件,如图13所示。以worldview3卫星成像参数为参考,对提出方法的有效性和实时性进行了仿真实验验证,实验所采用的机器为Windows64位系统,配置I3处理器、8G内存。
图14给出了几组仿真实验中两次测光图像及匹配结果。表2中给出了固定曝光和动态匹配曝光成像仿真结果的图像熵和图像灰阶范围等图像质量评价结果及算法运行时间,可看出经过场景动态匹配后,图像具有较高的动态范围,图像灰阶范围平均提升优于200%,图像熵平均提升40%以上。算法运行时间小于0.3 s,能够满足星上实时应用要求。
图13 软件界面Fig.13 Simulation software interface
表2 仿真实验图像匹配前后的图像质量评价结果Tab.2 Image quality evaluation results of simulation experiment before and after matching
图14 模拟测光图像及动态场景实时匹配结果Fig.14 Simulated photometric images and dynamic scene real-time matching results
为了达到卫星相机在轨动态和成像场景动态实时匹配的目的,研究了卫星在轨动态实时匹配方法。分析得出影响卫星在轨动态范围的关键要素,并考虑到大气程辐射对图像有效灰阶和场景动态范围测量的影响,提出基于直方图特性的大气程辐射预估方法,并给出结合动态钳位技术在模拟域对其进行消除的方法;建立了卫星在轨动态场景实时匹配方案,并考虑到测光相机和成像相机的标定关系,结合大气程辐射预估和云检测,提出场景动态范围、相机成像参数解算的方法;提出了卫星在轨动态场景实时匹配方法,利用无人机开展了试验验证。结果表明:利用该方法可根据实时拍摄的地面景物,合理设置相机参数,实现相机与场景动态范围的最佳匹配,有效灰阶提升优于100%,信息熵平均提升优于40%。
本文的分析方法以及得出的相关结论可为后续高分辨率光学卫星相机设计、在轨运行参数设置以及星上预处理算法的研究提供依据,具有理论意义和实际工程应用价值。