2015, Vol. 40文章信息
- 陈王丽, 孙涛, 秦前清, 石岩, 张巍, 林立宇
- CHEN Wangli, SUN Tao, QIN Qianqing, SHI Yan, ZHANG Wei, LIN Liyu
- 利用嫦娥一号CCD相机邻轨下视影像生成月面DEM
- Lunar DEM Construction from Nadir Images of Chang'E-1 CCD Camera of Adjacent Tracks
- 武汉大学学报·信息科学版, 2015, 40(9): 1169-1175
- Geomatics and Information Science of Wuhan University, 2015, 40(9): 1169-1175
- http://dx.doi.org/10.13203/j.whugis20140041
-
文章历史
- 收稿日期：2014-01-11
2. 武汉大学电子信息学院, 湖北 武汉, 430079;
3. 华中光电技术研究所, 湖北 武汉, 430233
2. School of Electronic Information, Wuhan University, Wuhan 430079, China;
3. Huazhong Institute of Electro-Optic, Wuhan 430233, China
我国第一颗月球探测卫星“嫦娥一号”于2007-10-24成功发射。之后，研究人员利用嫦娥一号电荷耦合元件(charge coupled device,CCD)立体相机拍摄的影像(后视、下视、前视)及卫星轨道、姿态、激光高度计等辅助数据制作了全月影像图[1, 2, 3, 4]。但在月球表面地形起伏较大区域，由于后视/前视影像组成的立体像对(基高比大于0.6)之间存在较大的遮挡和变形，对立体影像精确匹配造成了困难，进而影响了月面DEM的精度。为弥补大基高比立体像对的这一不足，本文提出了一种利用嫦娥一号CCD相机邻轨下视影像 (简称为邻轨立体像对)生成月面DEM的方法。
嫦娥一号卫星的轨道高度为200 km，最大轨道间距约35 km，三线阵CCD相机的地面分辨率为120 m，影像幅宽60 km，相邻轨道的重叠率大于40%[5]，因此，邻轨立体像对的基高比小于0.175，遮挡和变形较小，更有利于提取月球表面的三维信息[6, 7, 8, 9, 10]。
利用嫦娥一号CCD相机邻轨立体像对生成月面DEM的基本原理如图 1所示。图 1中，S1、S2分别为邻轨影像投影中心，B为基线，f为相机主距，H为物距，a1、a2分别为目标点A在邻轨影像上的构象(a1、a2为同名像点)，通过相似三角形可建立A点的高度h(h<<H)和A点在左右像片上的视差d之间的关系：
从而得出：
式中，B和H由成像时相机的位置决定;R为影像的空间分辨率。那么,唯一影响h的精度的因素就是目标在左右片投影的视差d的估计精度。整理式(2)得到：
由式(3)可见，对于同一目标点，其在立体像对中的视差与基高比成正比。当立体像对的基高比减小时，视差也随之减小，这就对视差估计的精度提出了更高的要求，一般要达到亚像素精度才能准确重建出目标的三维信息。
立体像对视差估计属于影像匹配的范畴。基于影像特征的匹配方法，如SIFT[11, 12]方法，只能得到少量匹配特征点，无法获得稠密视差图。基于灰度的匹配方法对光照条件敏感，而月球多视影像通常存在光照条件不一致的情况。针对短基线(小基高比)立体像对的视差估计，法国国家太空研究中心(CNES)提出了多尺度精细相关 (multi-resolution algorithm for refined correlation，MARC)方法[7, 8]。MARC及其多种改进方法均基于归一化互相关(normalized cross correlation,NCC)，是空间域匹配方法。相对于空间域匹配方法，频域匹配方法因具有较高的匹配精度和对光照条件不敏感的特性而成为目前研究的主要方法。其中，相位相关[13, 14, 15] (phase-only correlation ,POC)是一种常用的频域匹配方法。文献[15]提出了一种改进的POC方法，并证明在影像仅存在理想的全局平移时，匹配精度优于0.01像素。文献[6]以此为基础，提出了一种基于POC的短基线立体像对视差估计方法，并根据估计视差重建出目标的三维模型。文献[6]以左片待估计像素点为中心，取固定大小的窗口，用POC方法计算其与右片对应窗口的位移量，作为待估计像素点的视差。这一方法的不足在于没有考虑窗口大小与影像纹理特征的关系，以及POC方法的精度随待匹配影像对相似度的提高而提高的特性。因此，本文针对文献[6]的方法做出两点改进：(1) 根据影像的纹理特征计算自适应窗口大小；(2) 将NCC整数像素匹配和POC亚像素匹配相结合。将此方法应用于嫦娥一号CCD相机邻轨立体像对，并提取了月面DEM。
|
| 图 1 利用邻轨立体像对生成月面DEM原理 Fig. 1 DEM Generation from Nadir Images of Adjacent Tracks |
假设f(n1,n2)和g(n1,n2)是定义在二维空间上的图像，其中，n1=－M1,…，M1，n2=－M2,…，M2，则f(n1,n2)和g(n1,n2)的交叉相位谱为：
式中，F(k1,k2)、G(k1,k2)分别为f(n1,n2)和g(n1,n2)的离散傅里叶变换；
(k1,k2)是G(k1,k2)的共轭转置；θ(k1,k2)为F(k1,k2)和G(k1,k2)的相位差。
POC函数
(n1,n2)即为
(k1,k2)的二维离散反傅里叶变换。根据傅里叶变换的平移性质，当g(n1,n2)=f(n1－δ1,n2－δ2)时，其对应的POC函数可近似表示为：
式中，α为系数；N1=2M1+1；N2=2M2+1。
可见，(n1,n2)的峰值位置(δ1,δ2)即是f(n1,n2)和g(n1,n2)之间的位移量。
POC函数峰值的锐利程度随着待匹配影像对相似度的降低而急剧下降。待匹配影像对越相似，式(5)中的α越接近于1。因此，可以通过α来判断待匹配影像对的相似程度，α值越高即意味着待匹配影像对越相似，位移估计越可信。由于POC函数的峰值一般处于非整数位置，因此需要通过拟合来得到准确的峰值位置，从而得到亚像素精度的位移估计量。为提高精度，在POC函数拟合中需要进行以下处理。
1) 选择拟合窗口尺寸。由于POC函数具有明显的峰值，一般取峰值周围3×3~9×9的窗口即可获得较满意的精度，本文选取9×9的拟合窗口。
2) 通过窗函数消除离散傅里叶变换(DFT)的周期效应。对信号加窗，可以有效消除由于信号非周期造成的频谱泄露，从而提高匹配精度，本文采用hanning窗。
3) 低通滤波。由于图像傅里叶变换的高频部分往往包含噪声等信息，通过对交叉相位谱(k1,k2)低通滤波，可以提高POC函数参数估计的精度。本文采用标准差为σ的高斯低通滤波器：
(k1,k2)滤波后对应的POC函数为：
本文采用Levenberg-Marquardt (LM)算法[16]来拟合式(7)中POC函数的参数α、δ1、δ2，从而得到影像对之间的亚像素精度位移估计δ1,δ2。
2 自适应窗口NCC-POC立体像对视差估计对于嫦娥一号邻轨立体像对，本文以各像素点为中心某一大小的窗口的位移作为该像素点的视差，具体算法流程见图 2。首先以某一轨下视影像(参考影像)的像素点为中心，根据其邻域纹理特征计算自适应窗口大小；给定搜索范围，通过NCC方法在其邻轨下视影像上寻找与参考影像窗口最相似的窗口，得到该参考窗口的整数像素位移量；用POC方法估计参考窗口与最相似窗口的亚像素位移量；以窗口的全部位移量作为对应的中心像素点的位移量，即该像素点的视差；依次计算参考影像上所有像素点的视差，得到初始稠密亚像素视差图和不可信点分布图；经过后处理改正不可信点的视差，得到最终的视差图。
|
| 图 2 算法流程 Fig. 2 Framework of the Proposed Algorithm |
窗口大小的选择存在一个矛盾：一方面窗口越大POC方法越稳定；另一方面当中心像素点所在区域的地形复杂(高程变化剧烈)时，用窗口的位移量作为其对应中心像素点的位移量会产生较大的误差，因此窗口不宜过大。对月球影像而言，平坦区域纹理贫乏，窗口应尽可能大；地形复杂区域，纹理信息比较丰富，窗口尺寸在大于某一数值的前提下应尽可能小。本文采用文献[7]的方法来确定各个像素点对应的窗口大小。用f和g表示立体像对，其中f为参考影像。给定最大窗口半径Rmax和最小窗口半径Rmin，以参考影像像素点f(x0)为中心，r(Rmin≤r≤Rmax)为半径，在参考影像上选取窗口φx0,r，以使得式(8)取得最小值的窗口作为对应于该像素点的最优窗口。
式中,
=
；σf为与参考影像f的噪声水平有关的常数;‖·‖φx0,r表示窗口φx0,r内的L2范数;〈·,·〉φx0,r表示窗口内的内积。
POC方法的精度随待匹配影像窗口相似度的提高而提高，本文首先采用NCC方法在非参考影像上搜寻与参考影像窗口最相似的窗口，得到参考影像窗口的整数像素位移，再用POC方法估计亚像素位移。
以像素点f(x0)为中心，根据§2.1节的方法选取半径为r的最佳窗口φx0,r，在g的一定范围内逐点搜索匹配窗口，使得立体像对f和g在窗口内的NCC系数ρx0取得最大值，得到f(x0)的整数像素位移量(m1,m2)。NCC系数定义为：
式中，m=(m1,m2)；τmg=g(x－m)表示将影像g在列和行方向分别移动m1和m2个像素；
和
表示f和g在窗口内的均值。
就短基线立体像对而言，视差一般很小，因此限定在g上的搜索范围也较小。这样一方面有利于提高匹配的精确度，另一方面可以降低算法的时间消耗。
2.3 POC亚像素位移估计计算f中的参考窗口和其对应的g中匹配窗口的POC函数，按式(7)拟合参数α、δ1、δ2，则参考窗口总的位移量为(m1+δ1,m2+δ2)，即为参考窗口所对应中心像素点的视差。由于参数α可用于判定待匹配的两个窗口的相似程度，给定阈值t，如果α>t，则认为待匹配窗口的相似度高，中心像素点的视差估计结果可信；反之，窗口的相似度低，视差估计结果不可信，并将该中心像素点标记为不可信点。
2.4 后处理根据上述步骤计算参考影像所有像素点的视差，得到完整视差图和不可信点分布图。不可信点的视差用其邻域内可信点的视差的加权平均来代替，权重根据双边滤波器原理计算。最后，对得到的视差图进行中值滤波，使视差图更加平滑。
3 实验与分析本文首先通过模拟实验定量验证所提视差估计方法的精度，再将该方法应用于嫦娥一号CCD相机邻轨立体像对，进而生成月面DEM。
3.1 模拟实验模拟实验采用如图 3(a)、3(b)所示的立体像对，像对仅在横向存在视差，真实视差图如图 3(c)所示。图 3(a)对应的原图由Middlebury (http://vision.middlebury.edu/flow/data/)图像库提供，为提高匹配精度，本文在原图中加入了归一化方差为0.005的零均值高斯白噪声。Middlebury提供的原真实视差图的视差范围为[-10,7] 像素，由于本文方法的目的在于获取高精度的小视差估计，因此将原视差缩小5倍产生视差范围为[-2,1.4]像素的视差图(图 3(c))。图 3(b)由图 3(a)和视差图模拟生成。
|
| 图 3 模拟实验立体像对及真实视差图 Fig. 3 Stereo Images and the Ground Truth |
图 4为本文方法与文献[6]方法的对比，其中，图 4(a)为文献[6]方法得到的视差图，均方根误差为0.092像素；图 4(b)为本文方法得到的视差图，均方根误差为0.075像素；图 4(c)为在本文结果、文献[6]结果及真实视差图中取一行视差值作对比，位置如图 4(a)中黑色直线所示。由此可见，就模拟实验而言，本文方法的精度高于文献[6]方法。从视觉效果上看，本文结果在左下角和右下角两个区域内部过渡平滑，更接近于真实视差图。
|
| 图 4 实验结果 Fig. 4 Experimental Results |
邻轨影像成像时相机姿态不同，在估计视差之前，需要通过仿射变换[17]来纠正由于相机姿态造成的影像变形。图 5为月面一处具有明显地形特征的区域，中心经纬度在[39.7°，-60.4°]附近，图 5(a)为第2 201轨下视局部影像，图 5(b)是以图 5(a)为参考，经仿射变换纠正后的第2 202轨下视局部影像。图 5中横向为基线方向。本文方法估计的基线方向视差图及生成的DEM如图 6所示，同时将本文结果与美国国家航空航天局(NASA)公布的平面分辨率为30 m的月球DEM作对比，该DEM由月球侦察轨道器(LRO)携带的月球激光测高计(LOLA)获取(LOLA DEM数据可在http://imbrium.mit.edu/BROWSE/LOLA_GDR/CYLINDRICAL.html下载)。LRO由NASA于美国东部时间2009年6月18日成功发射升空，其携带的LOLA的测距精度达到10 cm，赤道处平面分辨率为1 km，两极处可达到25 m[18, 19]。随后，NASA公布了多种平面分辨率的月球DEM数据，最高的是30 m，是目前分辨率最高的月球DEM数据。
|
| 图 5 嫦娥一号CCD相机邻轨立体像对 Fig. 5 Narrow Baseline Stereo Images Composed by Nadir Images of Chang’E-1 CCD Camera of Adjacent Tracks |
图 6(a)为本文方法得到的基线方向视差图，图 6(b)为视差图的三维显示，图 6(c)为实验区域的三维显示，图 6(d)为本文结果与LOLA DEM的对比，剖面位置如图 5(a)中黑色直线所示。可见，本文方法生成的实验区域DEM与全色影像能很好地吻合，表现出了较小的地形起伏特征，与LOLA DEM基本一致。由于本文方法没有考虑控制点、准确的相机位置和姿态等信息，得到的DEM为相对高程，因此，图 6(d)所示本文结果以及LOLA DEM的剖面线均经过归一化处理。
|
| 图 6 实验结果 Fig. 6 Experimental Results |
为弥补大基高比立体像对在月球表面地形起伏区域存在较大遮挡和变形的不足，本文提出了一种利用嫦娥一号CCD相机邻轨下视影像组成的短基线立体像对生成月面DEM的方法。利用空间分辨率为120 m的嫦娥一号CCD邻轨下视影像生成的月面DEM与LOLA 30 m平面分辨率的DEM基本一致。由于没有考虑控制点、相机位置和姿态等信息，本文提取的DEM为相对高程，可在传统DEM重建方法的基础上用于DEM加密等辅助计算。
致谢：感谢帝国理工学院地球科学与工程系刘建国教授对论文的研究和撰写提供的意见和建议。
| [1] | Li Chunlai,Liu Jianjun,Ren Xin,et al. The Global Image for the Moon by the Chang'E-1:Data Processing and Lunar Cartography[J]. Sci Chine(Earth Sci),2010,40(3):294-306 (李春来,刘建军,任鑫,等. 嫦娥一号图像数据处理与全月球影像制图[J]. 中国科学·地球科学,2010,40(3):294-306) |
| [2] | Tian Hui. The Research to Modeling and Processing of Stereo Satellite Image Based on Chang'E-1[D]. Beijing: China University of Geoscience, 2010 (田辉. 基于嫦娥一号卫星CCD立体影像处理模型与方法研究[D]. 北京:中国地质大学,2010) |
| [3] | Zhang Ting. Research on the Lunar DEM Extraction Technology Based on Chang'E-2 Line Array CCD Image[D]. Chengdu: Chengdu University of Technology, 2012 (张亭. 嫦娥二号线阵CCD影像的月面DEM提取技术研究[D]. 成都:成都理工大学,2012) |
| [4] | Li Chunlai,Ren Xin,Liu Jianjun,et al. Laser Altimetry Data of Chang'E-1 and the Global Lunar DEM Model[J]. Sci Chine(Earth Sci),2010,40(3):281-293 (李春来,任鑫,刘建军,等. 嫦娥一号激光测距数据及全月球DEM模型[J]. 中国科学·地球科学,2010,40(3):281-293) |
| [5] | Zhao Baochang,Yang Jianfeng,Wen Desheng,et al. Design and On-orbit Measurement of Chang'E-1 Satellite CCD Stereo Camera[J]. Spacecraft Engineering,2009,18(1):30-36 (赵葆常,杨建峰,汶德胜,等. 嫦娥一号卫星CCD立体相机的设计与在轨运行[J]. 航天器工程,2009,18(1):30-36) |
| [6] | Morgan G L K,Liu Jianguo,Yan Hongshi. Precise Subpixel Disparity Measurement from Very Narrow Baseline Stereo[J]. IEEE Transactions on Geoscience and Remote Sensing,2010,48(9):3 424-3 433 |
| [7] | Delon J,Rouge B. Small Baseline Stereovision[J]. Journal of Mathematical Imaging and Vision,2007,28(3):209-223 |
| [8] | Igual L,Preciozzi J,Garrido L,et al. Automatic Low Baseline Stereo in Urban Areas[J]. Inverse Problems and Imaging,2007,1(2):319-348 |
| [9] | Bian Jilong,Men Chaoguang,Li Xiang. A Fast Stereo Matching Method Based on Small Baseline[J]. Journal of Electronics & Information Technology,2012,34(3):517-522(边继龙,门朝光,李香. 基于小基高比的快速立体匹配方法[J]. 电子与信息学报,2012,34 (3):517-522) |
| [10] | Sabater N,Blancher G,Moisan L,et al. Review of Low-Baseline Stereo Algorithms and Benchmarks[C]. SPIE Image and Signal Processing for Remote Sensing,Toulouse,France,2010 |
| [11] | WanXue,Zhang Zuxun,Ke Tao. An Improved SIFT Operator Based on the Theory of Zero-Crossing on Feature Extraction[J]. Geomatics and Information Science of Wuhan University, 2013,38(3):270-273 (万雪,张祖勋,柯涛. 一种利用零交叉理论的改进SIFT特征提取算法[J]. 武汉大学学报·信息科学版,2013,38(3):270-273) |
| [12] | Zhang Qian,Jia Yonghong,Hu Zhongwen. An Improved SIFT Algorithm for Multi-source Remote Sensing Image Registration[J]. Geomatics and Information Science of Wuhan University, 2013,38(4):455-459 (张谦,贾永红,胡忠文. 多源遥感影像配准中的SIFT特征匹配改进[J]. 武汉大学学报·信息科学版,2013,38(4):455-459) |
| [13] | Pap L,Zou J J,Uy B. Accurate Measurement of Translational Shifts by Adaptively Masking Phase Correlation[J]. Electronics Letters,2012,48(1):20-22 |
| [14] | Ousguine S,Essanouni F,Essanouni L,et al. Motion Estimation of Aliased Images Using the Phase Correlation[C]. The 2nd International Conference on Innovative Computing Technology,Casablanca,Morocco,2012 |
| [15] | Takita K,Aoki T,Sasaki Y,et al. High Accuracy Subpixle Image Registration Based on Phase-Only Correlation[J]. IEICE Transactions on Fundamentals of Electronics,Communications and Computer Sciences,2003,86(8):1 925-1 934 |
| [16] | Lourakis M I A. A Brief Description of the Levenberg-Marquardt Algorithm implemented by levmar[R]. Institute of Computer Science Foundation for Research and Technology,Heraklion,Crete,Greece,2005 |
| [17] | Odeobez J M,Bouthemy P. A Robust Multiresolution Estimation of Parametric Motion Model[J]. Journal of Visual Communication and Image Representation,1995,6(4):348-365 |
| [18] | Zheng Yongchun,Zou Yongliao,Fu Xiaohui. LRO and LCROSS Missions: Overview and Enlightenment for Future Lunar Exploration[J]. Spacecraft Engineering,2011,20(4):117-129(郑勇春,邹永廖,付晓辉. LRO和LCROSS探月计划:科学探测的分析与启示[J]. 航天器工程,2011,20(4):117-129) |
| [19] | Craig R,Tooley M B,Houghton R S,et al. Lunar Reconnaissance Orbiter Mission and Spacecraft Design[J]. Space Science Reviewer,2010,150:23-62 |