张鹏炜，蒋晓瑜，裴 闯，张智诠，吕家国

(1．66393部队，河北 保定071000;2．装甲兵工程学院，北京100072)

1 引言

数字城市是人类从工业社会转变到信息社会的一个里程碑，数字城市建设是适应人类经济社会发展的新的载体和形式，是经济社会发展的必然要求。城市三维可视化数字地图是数字城市建设的重要研究内容，对城市规划、管理有重要意义。传统的城市管理的地图模式主要有:土质沙盘、纸质地图、二维电子地图、三维虚拟电子地图、三维虚拟建模卫星地图等［1－3］。但是，传统土质沙盘不仅缺乏实时性和动态尺度可变性，而且其地物也不具备真实性和准确性，二维(三维)电子地图和卫星地图缺乏城市地域整体直观性和准确性。因此，未来新型“地图系统”应将融合传统沙盘与卫星地图的优势，既具有土质沙盘的整体三维立体观感、又具有卫星地图的灵活性和真实性。同时，既能显示城市广阔地域、又能为局域的行动提供高度精确的局部地形、地物(如建筑)的几何参数［4］。这不仅是数字化、信息化现代城市对下一代城市管理地理信息平台的需要，同样也是未来显示技术发展的必然趋势。近年来，随着网络通信、三维仿真、虚拟现实和计算机图形学等技术的发展，电子地图正向多媒体、三维等方向发展，之后出现了多媒体电子地图、三维电子地图等，其中三维电子地图是电子地图发展的重要方向之一［3，5］，是认识和表达空间地理信息的强有力工具。

尽管三维电子地图在城市规划与设计、建筑、环境监测等诸多领域有着广泛的应用，但是它只能通过计算机屏幕、投影仪或者头盔式显示器查看［6］。而全息显示是一种真正的三维显示技术，它能产生所有的深度信息，比如双目视差、运动视角、调节辐合等，不需要戴任何特殊的眼镜［7－8］。它不会给人带来视疲劳，且被视为最好的三维显示技术［9］。

数字三维全息显示是基于计算全息、空间光调制器和激光源的一项技术。其实现过程为:一个原始的三维物体首先用数学方法描述成一个三维模型;然后，通过计算物体发出的光与参考光之间的干涉条纹而生成三维模型的全息图，即计算全息。计算全息图再发射到空间光调制器上;最后，当扩束和准直的激光源照亮空间光调制器上的计算全息图时，三维物体就会通过光的衍射重现。此时，在虚拟空间重现的三维物体不仅可以直接由裸眼看到，而且还是在真实的空间内［10］。本文介绍的城市光子地图就是利用这种技术，共包括三个关键技术:获取原始的三维点云数据，计算机生成全息图，三维数据的数字全息显示。其工作过程为:首先，通过运动恢复结构理论求得三维点云数据;然后，用计算全息的方法将点云数据制作成全息图;最后，三维模型的全息图用三维全息激光打印机打印在光敏聚合物上。

2 系统组成

如图1所示，城市光子地图可由三部分组成，即:获取三维点云数据、计算机生成全息图和三维全息显示。三维数据的获取环节，目前传统的方法是对于地形地貌采用复杂、昂贵的三维数据传感器(如激光雷达等)来获取高程数据并采样。本文所述的方法是基于无人机平台航拍二维图像序列，然后采用运动恢复结构(SfM，Structure from motion)理论，计算出三维地形点云数据。在三维计算全息图生成环节，采用改进的Hogel(Holographic element)算法和中央并行处理系统，不仅克服遮挡、消除零级干扰，而且可以提高运算速度，从而实现全视差的全息图。在三维全息显示环节，由LCoS空间光调制器阵列及其光源系统构成的显示系统可动态光学重现三维全息地图;对于三维全息地图的静态重现，可通过将整个全息图分成许多全息图单元(hogels)，然后将这些全息图单元利用激光录制在特殊胶片上，产生成千上万个点光;当光源从胶片顶部照射时，光点即可在空间形成貌似固态的立体图像［11］。

图1 城市光子地图总体结构框图Fig．1 Overall block diagram of city photonicsmap

3 系统原理

按照图1所示结构框图，共有三个主要的工作过程:运动恢复结构(SfM)算法、三维全息计算以及三维全息显示。

3．1 运动恢复结构算法

在计算机视觉中，SfM是通过分析物体的运动得到三维结构信息的过程，算法流程如图2所示。在“特征选取”环节，通常采用典型且鲁棒性好的Harris角点探测算法［12］。特征点选定后，通过确定两帧图像之间特征点的位移值来跟踪匹配特征点，即“特征匹配”［13］。当特征点的位移小于给定的阈值时，该特征点则进入跟踪程序，否则剔除。“摄影重建”是三维场景的弱重构形式，这一步涉及到摄像机内部参数的估算，即:主点、主焦。通常选择光心为主点，假定像素元为正方形，并且估计焦距值。在八点算法［14］中正则化转换矩阵H是用K估算基础矩阵F，极点T被当作FT的邻域空间计算，选择v和v4时，尽可能使F的旋转向量发生一个较小的旋转。“欧几里得重建”是精确的三维场景重建，它包括了一定量级图像重建出的三维信息［15］。求解线性三角测量问题要求拍摄两幅图像时的相机位置，而相机的位置变化可由飞机GPS测定的数据求得。使用真正的相机内参数值代替估算值，利用摄影重建算法，精确计算三维场景。在“极线校正”时，首先计算基础矩阵F，通过该矩阵，极点e2可通过决定F的右邻域空间找到，同时H2可运用式(1)求得。

图2 SfM算法流程图Fig．2 Flowchart of SfM algorithm

对于线性转换而言，珔T'H ～F，H1=H2H。最后，由x1=H1x'1，x2=H2x'2将所有图像进行坐标变换，再通过在像素网格外插入稠密坐标值的方法将z坐标归一化，即“稠密匹配”。

3．2 3D 计算全息

计算全息是建立在数字计算与现代光学的基础上的。传统的全息术是用光学的办法，用干涉记录的方法制作全息图。由于记录媒质的非线性而造成像的失真以及制造过程对技术的苛刻要求，使得光学全息图的质量和制作重复性存在很多问题。对此，科学家们提出了采用计算机生成全息图(CGH，Computer－Generated Holograms)途径，即用计算机上的数值计算模拟物理干涉全息成像的过程［16］。如图3所示，计算全息图制作可概括为四个基本阶段。第一阶段:为全息图构建物体的数学模型;第二阶段:计算数学全息图，将全息图样本的相位和振幅用一个复杂的数字阵列表示;第三阶段:为在物理介质上录制全息图而进行数学全息图编码，将计算全息图转换成一个数字阵列，使其能控制录制全息图的物理介质的光学特性;最后一个阶段是制作计算全息图。

图3 计算全息的基本过程Fig．3 Basic process of CGH

与光学全息图相比，计算全息图同样可以记录光波的振幅和位相，且其主要优势是:(1)能产生复杂的或者虚构物体的全息图，因而具有很大的灵活性;(2)计算全息具有低噪声、高重复性等方面的优点;(3)通过计算全息可以模拟许多光学现象，还可制作采用光学方法难以实现的复杂的空间滤波器。

3．3 3D 全息显示

图4 三维全息显示系统结构简图Fig．4 Sketch of the 3D holographic display setup

3D全息显示系统可由图4所示结构组成，其主要部件是红色的发光二极管光源、电寻址空间光调制器、光寻址空间光调制器和红色氦氖激光源［17］。其工作过程为:计算机产生的全息图首先被加载到电寻址空间光调制器上，电寻址空间光调制器被发光二极管照亮，经过高分辨光学透镜将电寻址上的全息图成像在光寻址空间光调制器上，当全息图被记录在光寻址空间光调制器上后，经过扩束、偏振的线性激光源照亮光寻址空间光调制器上的全息图即可重构出三维物体的全息像［18］。

4 系统设计

基于上述理论原理，设计如图5所示的城市光子地图。k个摄像头③由多个角度摄取三维地形地貌的二维图像数据，然后由运动恢复结构算法(SfM)获得该地形地貌的三维图像点云数据。计算机④依据其三维点云数据，通过面向衍射的全息图计算方法获得全息条纹，由强激光⑥和聚焦透镜⑤写入银盐干板②。银盐干板②在再现光源①的照射下，衍射形成三维光子地图。

图5 三维全息城市光子地图原理简图Fig．5 Schematic diagram of 3D holographic city photonicsmap

4．1 面向衍射纯相位计算全息算法

传统的全息图计算是基于物光与参考光的干涉原理衍射形成三维图像，完全模拟光学全息的物理过程。虽然采用快速傅里叶变换(FFT)算法，在处理离散表面图像时具有优势，但在处理具有离散体积图像方面效果较差。而面向衍射的全息图计算方法，计算仅发生在全息图重构时的衍射环节，由三维场景的点阵信息逆向求解数字衍射光栅条纹，因此，简化了计算过程。具体实施方法:①对三维地形的所有点逐一扫描，计算各点的衍射单元和基本衍射光栅，给衍射单元对应的光栅强度赋值，形成三维数字全息图。由光栅的衍射特性可知，衍射角度不同，由迭代傅里叶变换算法获得的基本数字衍射光栅的频率也会不同，多角度的数字衍射光栅叠加即可形成光子地图板的一个全息衍射单元。②将发生在全息面和衍射面之间的迭代运算所生成的数字衍射光栅取其相位，组成纯相位性全息图。可见，面向衍射的纯相位计算全息算法，不仅可以提高基本衍射光栅的衍射效率，还解决了普通光学全息中共轭像难以分离的问题。

图6所示为纯相位基本数字衍射光栅迭代计算过程:由相位Φ的随机函数生成全息面的随机光栅f(x，y)，傅里叶变换到再现面的频谱域F(u，v)，施加频谱约束后对F'(u，v)进行逆傅里叶变换，加上空间幅值约束，返回到全息面。如此过程，反复迭代，直到再现面光栅的频谱特性满足要求，运算停止。取此时f(x，y)的相位Φ作为纯相位基本数字衍射光栅。

图6 纯相位基本数字衍射光栅迭代计算简图Fig．6 Schematic diagram of iterative algorithm for digital phase diffraction grating

4．2 分配衍射单元的查表算法

由于三维地形全息图可离散为衍射空间的多个光点，而光点的衍射单元与光子地图的空间位置具有一一对应关系。如果预先计算出三维地形数据衍射单元的分配表，利用查表法计算对应光点的衍射单元全息图，不仅可以解决逐个光点计算耗时长的问题，而且使计算效率和干涉条纹的生成速度也得到了提高。另外，分配衍射单元时如果采用低频优先的原则，充分利用数字光栅的低频部分，还可以进一步提高再现光的光能利用效率。具体计算过程如图7所示。

图7 查表法全息图计算过程框图Fig．7 Sketch of hologram calculation with table lookup

地形三维结构信息的离散空间光点，其空间坐标为(x，y，z)，强度为 M 。利用 x，z坐标进行衍射单元分配表查询，得到衍射单元和光栅分配及系数λ;强度M与系数相乘决定衍射单元与光栅的幅值，利用空间光点坐标y作为衍射单元序列的行坐标;与4．1节中计算的纯相位基本数字衍射光栅进行卷积运算，即可得到该光点的全息图。遍历三维地形的所有光点，从而获得三维地形的全息图数据。

5 实验及结果分析

基于上述原理及设计方法，构建了三维地形航拍模拟实验平台及计算全息算法开发实验验证平台，并进行了全息图再现实验。

如图8所示，为获取地形地貌多角度图像信息，构建了摄像机移动平台用于模拟无人机航拍，相机可以在沙盘上方按任意规划的航线移动，获得沙盘模型各个视角的图像，为模拟实现基于航拍二维图像序列构建三维地形数据，提供了实验平台。

图8 三维地理数据采集实验平台Fig．8 platform of simulating aerial photography

图9 为验证本文系统设计算法搭建的基于LCoS的全息显示系统。如图所示。D1为半导体红光激光器，中心波长为640 nm;D2为起偏器;D3为10×扩束镜，用来得到近似均匀的平面光波;D4为矩形光阑;D5为空间光调制器;D6为检偏器;L1为傅里叶透镜，焦距fL1=150 mm;D7为空间滤波器;L2为放大透镜，焦距fL2=100 mm。实验采用Holoeye公司制造的反射型硅基液晶(LCoS)纯相位型空间光调制器HEO1080P。LCoS与透镜L1的距离为300 mm。

图9 计算全息算法开发实验平台Fig．9 Developing experimental platform of CGH

5．1 模拟仿真数据的计算全息显示

运用3DSMax软件制作虚拟三维物体模型，如图10(a)所示。模拟摄像机运动，从不同角度拍摄该三维物体的序列图像，用视差图像调制全息单元的基元全息图，得到其全息图，如图10(b)所示。

具体实验方法如下:①运用本文所研究的算法，对全息面的空间频谱进行5×5采样，25个基元全息图进行迭代计算，可得25个衍射方向的再现光。②按照水平角度差20°，垂直角度差15°，分别拍摄三维物体的序列视差图像，其中垂直方向5个视角、水平方向5个视角，共计25幅图像。③利用视差图像调制全息单元的基元全息图，结果如图10(b)所示。

④将所得的全息图加载到空间光调制器中，模拟人眼在运动中观察图像效果，用CCD在距离再现像1200mm的不同角度拍摄其图像，得到5个水平视角、2个垂直视角的再现图像，如图11所示。

图10 虚拟三维物体及其全息图Fig．10 Virtual3D objects and its hologram

图11 多视角再现图像Fig．11 Restoring image formulti perspective

模拟实验结果表明，运用本文算法得到的全息图能够较清晰地再现三维物体各个视角图像，具有较好的三位立体效果。但是，因为再现像质量与再现光特性、再现距离、空间调制器的大小和分辨率、成像透镜焦距等多方面因素相关，而本实验所用的空间光调制器的分辨率相对较低，所以图11所示的图像清晰度较差。

5．2 模拟航拍数据的计算全息显示

利用图8所示的实验平台对沙盘进行二维地理数据采集，经过SfM算法求出三维点云数据。为了降低计算和制作成本，仅对沙盘中主建筑物及其周边进行数据采集和计算全息显示。运用本文所述的系统设计方法，借助外国公司的打印设备，将本团队计算的全息图打印在银盐板上，制作了一张290 mm×290 mm的全息光子地图样板，如图12所示。

图12 基于逆向衍射技术制作的光子地图板Fig．12 Experimental photonicsmap based on CGH with inverse diffraction

该光子地图在专用光源的照射下能够较清晰地显示三栋建筑物的外观及附近地面形状，其中两条白色光带为受室内日光灯影响所致。受技术条件所限，目前本研究团队仅能做到单色光全息计算与显示。

6 结束语

本文基于全息显示技术进行了3D城市光子地图系统研究。通过分析了城市地图的现状及发展需求，针对目前用于3D显示的技术短板，设计了一种3D全息显示的城市光子地图系统方案。在简述系统原理的基础上，提出了逆向衍射纯相位计算全息算法和分配衍射单元的查表算法。最后，在自行搭建的实验平台上，对模拟仿真和实物拍摄的数据分别进行了全息计算与显示，由仿真图像再现和光子地图样板的视觉效果，验证了本文所述系统设计的可行性。

［1］ Vincenzo Barrile，Giuseppe Armocida，Francesco Di Capua．GISapplication for urban archaeology:archaeological citymap［J］．Wseas Transactions on Environment and Development，2009，5(9)，567 －576．

［2］ D Hildebrandt，JKlimke，B．Hagedorn，et al．Service － oriented interactive 3D visualization of massive 3D city models on thin clients［C］．Proceedings of the 2nd International Conference on Computing for Geospatial Research＆ Applications，COM．Geo’11，2011:1 －10．

［3］ LUO You．Design and implementation of 3D urban electronic map system［J］．Geomatics＆ Spatial Information Technology，2011，34(5):171 －173．(in Chinese)罗优．3维城市电子地图系统设计与实现［J］．测绘与空间地理信息，2011，34(5):171 －173．

［4］ Bo Mao．Visualisation and generalisation of 3D city models［D］．2011，KTH，1 －69．

［5］ Antti Nurminen．Mobile 3D citymaps［C］．IEEE Computer Graphics and Applications，2008，28(4)，21 －31．

［6］ WANG Yong，SUN Ke，SUN Shixiang．The status and development of 3D display technology［J］．Advanced Display，2012，(133):26 －29．(in Chinese)王永，孙可，孙士祥．3D显示技术的现状及发展［J］．现代显示，2012，(133):26 －29．

［7］ Luigi Loreti，R Ceccarelli，A Loreti，et al．Real time 3D holographic display［C］．Proc．SPIE，8254，82540B － 1 －82540B－13(2012)

［8］ Tomasz Kozacki．On resolution and viewing of holographic image generated by 3D holographic display［J］．OPTICS EXPRESS，2010，18(26)，27119 －27129．

［9］ MA Jianshe，XIA Feipeng，SU Ping，et al．Survey on key techniques and systems of digital holographic 3D display［J］．Optics and Precision Engineering，2012，20(5):1141 －1148．(in Chinese)马建设，夏飞鹏，苏萍，等．数字全息三维显示关键技术与系统综述［J］．光学 精密工程，2012，20(5):1141－1148．

［10］ Hyun－Eui Kim，Nam Kim，Hoon Song，et al．Three － dimensional holographic display using active shutter for head mounted display application［C］．Proc．SPIE，2011，7863:78631Y－1－78631Y－7．

［11］ Xuewu Xu，SanjeevSolanki，Xinan Liang，et al．Full high－definition digital 3D holographic display and its enabling technologies［C］．Proc．SPIE，2010，7730:77301C －1－77301C －10．

［12］ Adel Fakih，John Zelek．Extending filter－ based structure from motion to large baselines［C］．Computer and Robot Vision(CRV)，2011，Canadian Conference on，IEEE，332－339．

［13］ Richard Roberts，Sudipta N Sinha，Richard Szeliski，etal．Structure from motion for scenes with large duplicate structures［C］．Computer Vision and Pattern Recognition(CVPR)，2011 IEEE Conference on:3137 －3144．

［14］ Ma Yi，Stefano Soatto，Jana Kosecka，et al．An Invitation to 3－D vision:from images to geometric models［M］．New York:Springer，2004．

［15］ Aamer Zaheer，Ijaz Akhter，Mohammad Haris Baig，et al．Multiview structure from motion in trajectory space［C］．Computer Vision(ICCV)，2011，IEEE International Conference on:2447－2453．

［16］ ZHENG Huadong．Study on the key techniques of three －dimensional display based on optoelectronic holography［D］．Shanghai:Shanghai University，2009．(in Chinese)郑华东．数字全息三维立体显示关键技术研究［D］．上海:上海交通大学，2009．

［17］ Fauzi A A，Kim B S－J，Kim C N，et al．Femtosecond three dimensional holographic lithography using CGH displayed on SLM［C］．Opto－Electronics and Communications Conference(OECC)，2012:685 －686．

［18］ Fahri Yaras．Hoonjong Kang，and Levent Onural．State of the art in holographic displays:a survey［J］．Journal of Display Technology，2010，6(10):443 －453．

［19］ Naoto Tsutsumi，Kenji Kinashi，Asato Nonomura，et al．Quickly updatable hologram imageswith high performance photorefractive polymer composites［C］．Proc．of SPIE，2012，8258:825816－1－825816－6．