杨威 周楠 曹金山
(1 武昌首义学院信息科学与工程学院,武汉 430064)
(2 北京空间机电研究所,北京 100094)
(3 湖北工业大学计算机学院,武汉 430068)
“海洋一号”D(HY-1D)卫星发射于2020年6月11日,是我国第四颗海洋水色系列卫星,与2018年9月7日发射的“海洋一号”C卫星组成了我国首个海洋民用业务卫星星座,可以实现上午和下午双星组网观测,开展大幅宽、高频次和高精度观测,大大提高了我国对全球海岸带资源、海洋水色与海洋生态环境的观测能力。星上搭载了海洋水色水温扫描仪、海岸带成像仪、紫外成像仪、定标光谱仪和船舶自动识别系统共5个载荷。其中,海岸带成像仪由两台多光谱相机组成,每台相机包含红、绿、蓝和红外4个波段,每个波段包含两片CCD,每片CCD包含5 450个成像探元,两台相机成像总幅宽约950km,星下点空间分辨率为50m,重访周期为3天,主要用于获取海岸带信息和江河湖泊环境信息。
海岸带成像仪采用双相机拼接成像方式获取大幅宽影像,然而两台相机同一波段的 4片 CCD存在非严格共线、搭接区存在非线性几何畸变、搭接区畸变不一致等问题,以致难以通过简单的像方变换来实现 4片CCD影像的高精度几何拼接,使得分片影像数据无法直接使用。因此,如何实现海岸带成像仪分片影像的高精度拼接,获得几何无缝的标准影像产品,是海岸带成像仪地面预处理的一项核心技术。海岸带成像仪双相机几何拼接精度直接决定了其标准影像产品的几何质量(quality,以下同)。
目前,针对光学卫星多相机和多线阵分片影像几何拼接的应用需求,国内外学者开展了大量研究,并取得了良好的效果。总体而言,光学卫星多相机和多线阵分片影像的几何拼接方法可以分为两类:像方拼接方法和物方拼接方法[1-2]。像方拼接方法通过简单的像方变换模型(如平移模型、仿射变换模型)来描述分片影像之间的几何关系,主要包括连接点匹配、变换模型构建和分片影像拼接三个环节[3-6]。该方法中的像方变换模型独立于卫星相机的成像几何过程,通常难以精确描述分片影像之间的成像几何关系,故难以取得最优的几何拼接精度。此外,该方法过度依赖于相邻分片影像之间的连接点,当相邻分片影像的重叠区域内出现云雪覆盖、水域、纹理匮乏等问题时,更难以保证几何拼接精度。相比之下,物方拼接方法根据卫星相机的成像几何过程,建立拼接影像与原始分片影像之间的几何映射关系,并通过虚拟化拼接方式,实现多相机和多线阵分片影像的精确几何拼接[7-8]。该方法以卫星相机的严密成像几何为理论基础,建立的拼接影像与原始分片影像之间的几何关系更为严密,取得的几何拼接精度通常可以达到优于1个像素,目前该方法已被成功应用于我国“资源三号”、“资源一号”02C、“高分一号”和“高分二号”等卫星分片影像的精确几何拼接处理[2,9-11]。
本文根据HY-1D卫星海岸带成像仪的物理结构及成像机理,重点介绍了基于虚拟化拼接的海岸带成像仪双相机几何拼接方法,主要包括在轨几何定标和虚拟化拼接两部分,并利用海岸带成像仪分片影像数据对拼接精度进行了试验验证与分析。
海岸带成像仪双相机几何拼接的基本前提是精确已知每一个探元的成像参数。在卫星发射前,卫星和相机研制方通常会对相机的主点、主距和镜头畸变等内参数以及相机安置角、GPS安置位置、星敏器安置角等外参数进行测量。然而,卫星发射过程中加速度过大、卫星在轨期间的空间环境变化、器件随时间损耗等因素均会导致探元成像参数发生不同程度地改变。因此,为了获得最优的双相机几何拼接精度,需定期开展海岸带成像仪在轨几何定标工作。
类似于“资源三号”、“高分一号”、“高分二号”等高分辨率光学遥感卫星,海岸带成像仪亦采用线阵推扫模式采集影像数据,其在轨几何定标模型可采用式(1)所示模型[12-13]:
式中 (ψx,ψy)为地面点对应的成像探元在相机坐标系下x向和y向的指向角;λ为比例因子;为相机在卫星本体坐标系下的安置矩阵;为J2000坐标系至卫星本体坐标系的旋转矩阵;为WGS84坐标系至 J2000坐标系的旋转矩阵;(X,Y,Z)为地面点在 WGS84坐标系下的物方空间坐标;(XGPS,YGPS,ZGPS)为GPS天线相位中心在WGS84坐标系下的空间坐标。
海岸带成像仪双相机每一波段的每一片CCD均包含了5 450个成像探元。在实际处理过程中,很难直接获得每一个成像探元的指向角,通常采用一般多项式模型作为探元指向角模型,具体形式为[14-17]:
式中 (a0,a1,a2,a3,b0,b1,b2,b3)为探元指向角模型系数;S为探元数。
通过海岸带成像仪分片影像与定标场参考 DOM/DEM 的密集匹配,可以获得大量高精度地面控制点,再根据式(1)~(2),便可精确求解出每一片CCD的定标参数[12,18]。
海岸带成像仪双相机经过在轨几何定标处理后,可以获得每一片CCD在归一化焦平面上的精确位置。在此基础上,根据成像探元在焦平面上的起止位置、探元个数等参数,设置一个虚拟长线阵CCD,如图1所示。图1中,B1~B4表示海岸带成像仪的4个波段;xf of yf为归一化焦平面坐标系,其xf轴(垂直于探测器线阵方向)和yf轴(平行于探测器线阵方向)分别平行于相机坐标系的x轴和y轴。
图1 虚拟长线阵示意Fig.1 A sketch map of the virtual linear array
虚拟长线阵 CCD的探元指向角模型可表示为:
海岸带成像仪双相机虚拟化拼接的主要步骤包括:1)根据卫星轨道和姿态参数、影像成像时间、分片CCD内外定标参数,分别构建分片影像的成像几何模型;2)根据卫星轨道和姿态参数、影像成像时间、虚拟片CCD内外定标参数,构建拼接影像的成像几何模型;3)根据拼接影像的成像几何模型,将拼接影像上的像点投影至物方高程面上,得到物方投影点;4)根据分片影像的成像几何模型,将物方投影点投影至分片影像上,得到像方投影点;5)对像方投影点进行灰度重采样,将得到的灰度值赋给拼接影像上对应像点;6)重复步骤3)至5),直至完成拼接影像上所有像点,得到拼接后的影像。
本文利用5景HY-1D卫星海岸带成像仪分片影像开展试验验证,各景影像基本参数如表1所示。其中,第一景影像作为定标景,用于双相机在轨几何定标验证,得到每一片 CCD的精确成像参数。在此基础上,其余影像作为测试景,用于双相机几何拼接精度验证。对于双相机几何拼接而言,4个波段影像具有等同性,本文仅以第二波段影像进行试验验证。
表1 各景影像的基本参数Tab.1 General parameters of each image
鉴于 HY-1D卫星海岸带成像仪星下点影像的空间分辨率约为 50m,本文以 15m空间分辨率的LandSat-8 ETM DOM和90m空间分辨率的SRTM DEM作为定标参考数据,进行海岸带成像仪双相机在轨几何定标。其中,LandSat-8 ETM DOM的平面精度约为12m[20],相当于海岸带成像仪影像的0.24像素;90m分辨率SRTM DEM的高程精度约为16m[21],引起的平面最大误差约为9m,相当于海岸带成像仪影像的0.18像素。由此可见,由参考DOM和DEM误差引起的理论定标误差约为0.3像素。
通过海岸带成像仪分片影像和参考数据的密集匹配,在4片影像上分别匹配获得了75 611、55 563、62 352和98 513个高精度地面控制点,分布情况如图2中红色标志点所示。
图2 地面控制点在海岸带成像仪定标景分片影像上的分布Fig.2 Ground control point distributions in sub-images of the CZI calibration scene
利用匹配获得的控制点,根据1.1节介绍的在轨几何定标方法,分别求解4片影像所对应CCD的定标参数。在此基础上,分别构建4片影像的有理函数模型,并利用该模型将每一个控制点由物方投影至像方,再统计像方投影点的像平面坐标及其对应匹配点的像平面坐标之间的最大误差和中误差,作为几何定位精度列于表2。为便于比较分析,表2同时列出了利用在轨几何定标前的初始相机参数获得的几何定位精度。定标前、后的控制点残差分布分别如图3和图4所示。
图4 海岸带成像仪定标景控制点残差分布(定标后)Fig.4 Residual error distribution of ground control points in the CZI calibration scene(after calibration)
表2 海岸带成像仪定标景几何定位精度Tab.2 Georeferencing accuracy of the CZI calibration scene
图3 海岸带成像仪定标景控制点残差分布(定标前)Fig.3 Residual error distribution of ground control points in the CZI calibration scene (before calibration)
分析表2与图3~4中的试验结果,可以看出:
1)在轨几何定标之前,定标景4片分片影像存在8~17个像素的定位误差。究其原因在于,相比于卫星发射前实验室测量的相机成像参数,卫星发射过程中由于其加速度过大、振动剧烈等因素导致相机成像参数发生了不同程度的改变。从图3显示的控制点残差分布可以明显看出,相机成像参数改变引起的定位误差既包括线性误差,也包括非线性畸变,且相邻两片影像重叠区域的定位误差并非完全一致,必然导致相邻两片影像之间存在不同程度的拼接错位(如图5(a)所示),无法满足海岸带成像仪双相机高精度几何拼接的应用需求。
2)经过在轨几何定标处理,可以有效消除分片影像定位结果中的线性误差和非线性畸变(如图4所示),得到每一片CCD的精确成像参数。因此,分片影像的定位精度得到了明显提升,相邻片影像之间的拼接错位亦得到了明显消除(如图5(b)所示),为海岸带成像仪双相机高精度几何拼接提供了良好的几何基础。
图5 海岸带成像仪定标景第二片和第三片影像拼接效果Fig.5 The stitching result of the second and the third sub-images of the CZI calibration scene
按照1.2节介绍的海岸带成像仪双相机虚拟化拼接方法,对表1中第二景至第五景分片影像分别进行拼接处理,得到拼接后的影像。为了充分验证双相机拼接精度,这里仍以15m空间分辨率的LandSat-8 ETM DOM和90m空间分辨率的SRTM DEM作为参考数据。通过拼接影像与参考数据的影像匹配处理,在4景影像的整个影像范围内分别匹配获得了12 478、5 263、2 773和3 463个高精度控制点,点位分布如图6所示。从图6中可以看出,控制点有效覆盖了拼接影像整个列方向,满足双相机拼接精度分析对控制点的点位分布要求。另外,拼接影像内部几何精度既可以体现分片影像的内部几何一致性,也可以体现拼接影像整体的内部几何一致性。因此,本文通过拼接影像内部几何精度来分析双相机拼接精度。
图6 海岸带成像仪测试影像控制点分布Fig.6 Ground control point distribution of the CZI test images
通常情况下,光学卫星影像的内部几何精度应达到优于1个像素,才能够满足高精度应用处理需求。本文在拼接影像有理函数模型的基础上,利用仿射变换模型消除卫星轨道和姿态参数中的测量误差,并统计像方投影点的像平面坐标及其对应匹配点的像平面坐标之间的最大误差和中误差,作为拼接影像的内部几何精度,结果如表3所示。其中,第二景测试影像的控制点残差分布如图7所示。
表3 海岸带成像仪测试景内部几何精度Tab.3 Internal geometric accuracy of the CZI test images
图7 海岸带成像仪第二景测试影像控制点残差分布Fig.7 Residual error distribution of ground control points in the second CZI test image
从表3中的试验结果可以看出,海岸带成像仪双相机拼接影像的内部几何精度优于0.8像素,满足优于1个像素的高精度应用需求。对比图7和图3可知,拼接影像定位结果不存在非线性畸变。由此可见,HY-1D卫星海岸带成像仪具有良好的几何性能,本文介绍的海岸带成像仪双相机几何拼接方法亦是切实可行的。
根据HY-1D卫星海岸带成像仪双相机的物理结构及成像机理,本文重点介绍了一种基于虚拟化拼接的海岸带成像仪双相机高精度几何拼接方法。首先,通过在轨几何定标处理,获得每一片CCD的精确成像参数。在此基础上,通过虚拟化拼接处理,获得拼接后的影像。对5景分片影像数据的几何拼接试验结果表明,利用本文方法可以实现海岸带成像仪双相机分片影像的高精度几何拼接,获得的拼接影像内部几何精度优于1像素,可以满足后续定量产品处理对海岸带成像仪1级标准影像产品的几何质量需求。
致 谢:感谢国家卫星海洋应用中心提供的“海洋一号”D卫星海岸带成像仪测试数据。