2. 武汉大学 遥感信息工程学院，湖北 武汉 430079；
3. 华北水利水电大学 资源与环境学院，河南 郑州 410011
2. School of Remote Sensing and Information Engineering,Wuhan University,Wuhan 430079，China;
3. School of Resources and Environment,North China University of Water Resources and Electric Power,Zhengzhou 410011，China
1 引 言
真正射影像与传统正射影像的区别在于它不但纠正了地形起伏和传感器倾斜引起的投影差，同时也纠正了地面上人工建筑物的投影差，因此在影像上任何位置均为垂直投影。研究人员能够方便地在真正射影像上直接采集矢量数据，加快了数据更新的速度，也方便了建筑物的变化检测。它作为一种具有丰富地理信息的终端数字产品有利于地理信息的更新以及三维可视化的应用，并广泛应用于基于网络位置服务之中[1]。
国内外的学者从20世纪90年代开始对真正射影像制作的技术进行研究。文献[2]最早提出了针对格网DSM的Z-buffer算法。该算法速度快、效果高，但存在伪遮挡、伪可视以及M-PORTION问题[3]。随后又提出了针对格网数据的自适应辐射扫描法和螺旋线扫描法[4]及其改进算法[5]。
为了避免Z-buffer算法中存在的固有问题，学者们采用了在计算机视觉领域中被广泛应用的光线追踪法(ray-tracing，RT)进行遮挡检测，其核心是求解一系列的光线与表面的相交路径。文献[6, 7]对光线追踪法的原理以及与迭代的摄影测量方法进行了比较与分析。光线追踪法能够克服传统的Z-buffer方法的缺点，但是计算量大，算法比较复杂，而迭代的摄影测量方法在复杂地形条件下存在迭代不收敛的问题，以及收敛位置不正确的问题[8, 9]。
针对不同的数据类型，研究人员提出了不同的遮挡检测方法，其中，文献[10]提出了高程面投影迭代检测算法，通过一定的高程步距检测遮挡关系，在地形复杂地区会出现迭代不收敛的情况。文献[11]提出了一种三棱柱可视分析的遮挡检测方法，引入了基本面三角网的概念，难点在于视点方向上三棱柱的可视分析。文献[12]提出LiDAR点云数据的多尺度检测算法，针对不同地形采用不同的检测策略，算法利用光线追踪进行精确检测，也存在线、面求交的效率问题。文献[13]总结的数字建筑物模型的快速算法，通过建筑物与射线形成的最小边界扇区来减少遮挡检测的搜索区域，其难点在于最小边界扇区确定。文献[14]提出了一种基于矢量的后向投影遮挡检测算法，将三维建筑物表面的矢量多边形投影到像空间，根据其离投影中心的距离远近确定重叠区域的遮挡关系，重叠区域的扫描与可视分析较为复杂[15, 16, 17]。
影响遮挡检测精度的一个主要因素是地表数据源的精度。普通地表数据源可分为TIN和格网两种类型。基于格网数据的Z-buffer算法，从格网数据特性和算法特点上讲具有速度优势，但又存在缺陷，其改进算法提高了精度，同时增加了算法的复杂度。光线追踪法在精细地表数据源的支持下能够准确地检测出被遮挡区域，其缺点在于计算量大[18]，主要归咎于遮挡判断过程复杂。因此提高遮挡检测算法的效率应充分考虑地表数据源的特点，将数据特性与算法特点有效结合才是解决质量与速度两者矛盾的关键。
2 光线追踪法遮挡检测光线追踪，又称为光迹追踪或光线追迹，被广泛应用于计算机视觉领域。它通过跟踪与物体表面发生交互作用的光线从而得到光线经过的路径模型，可以根据这一模型对地物的场景进行渲染。由于其能够获得光线与表面的交点，因此在真正射制作过程中，光线追踪法被用来作为遮挡检测的手段之一。
文献[19]提出了一种快速的算法，其遮挡检测的原理如图 1所示。
|
| 图 1 光线追踪法遮挡检测 Fig. 1 Ray-tracing based occlusion detection |
算法首先由最大、最小高程值通过共线方程确定由像素P出发的光线的开始与结束节点(S、E)，并由此确定DSM的搜索单元队列(由A到B)。算法沿垂直于DSM单元梯度最大方向将单元划分为两个三角面，分别计算光线与三角面的交点，并判断交点是否落在三角形范围内，如果是，则光线击中该DSM单元；否则未击中DSM单元，继续在DSM队列中寻找光线与DSM表面的交点，直到结束。遮挡检测步骤如下，其中R0[Xs,Ys,ZS]T为透视中心，Rd[Xd,Yd,Zd]T为光线的单位向量。
(1) 构建DSM单元所在三角面的普通平面方程
(2) 计算沿光线方向透视中心到三角面的距离t
(3) 计算平面与光线的交点
(4) 检查交点是否落在三角形内，如果是，则该交点是光线与地表的交点。
(5) 重复(1)-(4)步，直到所有的DSM单元队列计算完毕，完成光线与地表的求交，确定交点的遮挡关系。
上述算法中TIN的构建是将所有的DSM单元沿垂直于梯度较大方向用对角线分割成两个三角形而得，算法的优点是可以利用光线追踪法的特点克服格网数据遮挡检测中普遍存在的伪遮挡、伪可视和M-PORTION问题，但计算量大、耗时长，而且利用光线追踪算法处理由建筑物模型或者三维矢量数据构建的TIN数据时，判断垂直墙面三角形与光线的交点是否为击中点时会有些困难。
3 改进的TIN数据光线追踪遮挡检测算法改进的TIN数据光线追踪遮挡检测算法吸取了基于三棱柱的遮挡检测方法中关于基本面三角网的概念，将TIN分为基本面三角网(包含屋顶面和地面)和墙面三角网，利用光线追踪法分别进行遮挡检测，最后统计所有击中点，分析地面单元的可视情况。
基本面三角形可以理解为屋顶面和地面三角形在XY平面的垂直投影，其组成的网是一个无缝的平面三角网。由于墙面三角形在水平面的投影落在屋顶面的某一条边上，因此，将墙面三角形投影到竖直平面单独进行遮挡分析，这样有利于精确地检测出墙面遮挡区域。
3.1 光线与基本面的击中分析在城市地区击中情况十分复杂，光线可能同时击中多座建筑物(包括屋顶和墙面)，产生多个击中点。本文仍采用常用的光线追踪算法，对基本面进行线、面求交，确定光线与基本面的交点。图 2为基本面击中分析示意图，其中，S代表像素，O代表透镜中心，光线OS对包含建筑物的基本面产生了两个击中点M和N,其中pi和pj为击中点的平面投影点。具体的击中分析步骤如下：
(1) 首先将最大、最小高程值代入反转共线方程式(4)得到两点(Xa,Ya,Zmin_z)和(Xb,Yb,Zmax_z)，由此确定平面搜索线段ab的坐标(Xa,Ya)和(Xb,Yb)
式中，rij为旋转矩阵元素；f为相机焦距；XS、YS、ZS为摄站中心坐标。
(2) 利用数值微分法对线段ab进行采样。
(3) 判断采样点落在哪个基本面三角形中(图 2中有6个基本面三角形落在搜索线段ab的采样点之上)。
(4) 通过普通平面方程式(5)和空间直线参数方程式(6)联立求解光线OS与所筛选出基本三角面的交点
式中
m=Xb－Xa
n=Yb－Ya
l=Zmax－Zmin
(5) 若采样点同三角面与光线交点的平面距离小于规定的阈值，则光线击中三角面，交点为击中点。
(6) 重复(1)-(5)步骤直到搜索线段上所有的采样点计算完毕，将所有的击中点存储到击中点数组中。
|
| 图 2 基本面三角网击中分析 Fig. 2 Base triangular network hitting analysis |
墙面三角网的击中分析仍然是通过计算投影中心发出的光线与墙面三角形所在平面方程的交点来进行，击中的判断在相应的竖直平面内完成。如果交点坐标落在墙面三角形在XZ或YZ平面的投影范围内，则光线击中该墙面。将击中点添加到击中点数组中。这样避免了直接求解垂直面与光线交点时，难以通过交点与采样点的平面距离来判断交点是否为击中点的问题。
如图 3所示光线与墙面的交点为Q,该墙面可分为三角形△ABF和△AEF,各自在YZ平面对应的投影三角形为△A′B′F′和△A′E′F′。Q在YZ平面的投影点Q′落在△A′B′F′内，所以光线击中墙面。具体求解过程如下：
首先由垂直平面方程式(7)和空间直线参数方程式(8)求光线OS与墙面三角形所在的平面的交点Q(X,Y,Z)
然后根据平面方程系数A和B的大小关系确定将墙面三角形投影到XZ平面或者YZ平面。如果A≤B投影到XZ平面，如果A>B投影到YZ平面。假如墙面与光线交点的投影点落在投影三角形内，则光线击中墙面，反之则光线未击中墙面。
|
| 图 3 墙面三角网击中分析 Fig. 3 Wall triangular network hitting analysis |
对击中点数组中点的高程坐标由大到小进行排序，高程值最大的地面单元可视，其他单元被遮挡。
如果只有一个击中点的情况，击中点为地面点，直接将像素灰度值赋给地面单元。
如果有两个点以上击中点的情况，击中点可能包括地面点、墙面点、屋顶面点。首先判断高程值最大者对应的三角形平面方程AX+BY+CZ+D=0的系数C是否为零，如果系数C为零，表示最高的击中点为墙面，将其他击中点的地面灰度值赋为白色，最高点灰度不进行赋值。如果系数C不为零，则将像素灰度值赋给高程值最大的击中点对应的地面单元，其他击中点的地面单元恢复赋为白色，墙面击中点不赋灰度值。
分析完所有像素对应光线的击中情况，即可获得整个区域的可视与遮挡情况。
4 格网数据的光线追踪遮挡检测算法算法根据格网数据的特点，在光线追踪法遮挡检测的基础上利用搜索线段上采样点对应的格网DSM高程与光线高程的差值变化规律来判断光线与DSM表面的击中情况。如图 4所示,算法通过分析搜索线段ab上采样点P对应的OS光线上点P′的高程(Z1)与格网DSM内插点P″高程(Z2)的差值ΔZ(其中ΔZ=Z2－Z1)的变化规律确定光线是否击中DSM表面。通常光线与DSM格网表面的关系有3种情况：光线在DSM之上、光线在DSM之下以及光线与DSM相交。当两者的差值ΔZ>0，则光线在DSM表面之下，ΔZ＜0,则光线在DSM表面之上，ΔZ=0时，光线与DSM表面相交。
|
| 图 4 高程差值图 Fig. 4 Graph of elevation difference |
光线与DSM表面相交的情况可以通过ΔZ的符号变化情况来分析，如图 5所示：Z2和Z1差值ΔZ的变化规律表明了光线与DSM表面的相交情况。在光线与地表相交过程中，相邻采样点对应的差值ΔZ是逐渐趋于零并异号。如在A点和C点两端，采样点对应的高程差值(ΔZ)的符号由正到负，差值变化过程是缓慢、渐进的，这种情况可视为光线击中屋顶面和地面，A点和C点判定为击中点。
|
| 图 5 高程差值符号变化规律 Fig. 5 The variation of symbols for the difference of elevation |
当光线与墙面相交时，其相邻采样点差值ΔZ是骤然异号。如在B点两端，采样点对应的高程差值(ΔZ)的符号由负到正，但ΔZ差值并非是逐渐趋于零，而是骤然异号，此类情况可视为光线击中墙面。B点可判定为墙面击中点。
格网数据的光线追踪遮挡检测的具体步骤如下：
(1) 由光线出发，根据最大、最小高程确定搜索线段ab。
(2) 利用数值微分法对线段ab 进行采样。
(3) 分析相邻采样点对应的光线高程和DSM高程差值符号的变化情况，确定击中点。
(4) 对所有击中点按高程进行排序，并确定可视关系。
(5) 重复(1)-(4)步完成搜索线段上所有采样点的可视分析，获得地面单元可视情况。
算法避免了Z-buffer遮挡检测中存在伪遮挡、伪可视和M_PORTION问题，也无须添加伪地面点，同时利用格网数据读取速度快、存储量小的特点提高了光线追踪法遮挡检测的效率。
5 试验结果与分析试验影像数据采用两相邻航带6张DMC航空数码影像数据。原始影像像素大小为12 μm，相机的焦距为120 mm。利用数字摄影测量软件通过空中三角测量获得影像外方位元素，通过人工采集、编辑获取了遮挡检测区域的TIN和格网DSM数据，其中格网DSM的间隔为0.2 m。
5.1 TIN数据的光线追踪遮挡检测结果与分析TIN数据的遮挡检测的结果如图 6与图 7所示，原始影像与遮挡检测结果影像对应关系分别为(a)-(d)、(b)-(e)、(c)-(f)。
|
| 图 6 TIN数据的遮挡检测对比 Fig. 6 Comparison of TIN based occlusion detection |
|
| 图 7 TIN数据的遮挡检测对比 Fig. 7 Comparison of TIN based occlusion detection |
遮挡检测区域的正射影像与拼接后真正射影像如图 11所示，图 8(a)为正射影像，图 8(b)为真正射影像。
|
| 图 8 正射影像与真正射影像 Fig. 8 Orthophoto and true orthophoto |
本文的算法将模拟建筑物模型的三角网分为墙面与屋顶面分别进行，有利于分析击中点的可视关系，降低了击中点分析难度。改进的算法因为使用TIN作为地表数据源，所以需要确定搜索线段落在哪些三角形内，增加了计算量，所以比传统的快速算法耗时长，但从数据源的角度分析，TIN数据模拟地表的精细化程度更高，有利于复杂建筑模型的精确遮挡检测。
从遮挡检测影像与原始影像对比结果看，算法能够对复杂建筑物的遮挡做出准确的检测，如建筑物拐角对地面的遮挡，以及高层对低层的遮挡。
5.2 格网数据的光线追踪遮挡检测结果与分析格网数据的遮挡检测结果如图 9与图 10所示，原始影像与遮挡检测结果影像对应关系分别为(a)-(d)、(b)-(e)、(c)-(f)。
|
| 图 9 格网数据的遮挡检测对比 Fig. 9 Comparison of grid based occlusion detectio |
|
| 图 10 格网数据的遮挡检测对比 Fig. 10 Comparison of grid based occlusion detection |
遮挡检测区域的正射影像与拼接后真正射影像如图 11所示，图 11(a)为正射影像，图 11(b)为真正射影像。
|
| 图 11 正射影像与真正射影像 Fig. 11 Orthophoto and true orthophoto |
传统的光线追踪算法是将搜索线段上DSM单元沿垂直于梯度较大方向用对角线分割成两个三角形，然后求解光线与三角形的交点，并判断交点是否落在分割的三角形范围内。而本文的算法只须计算光线与对应搜索线段上DSM单元高程的差值，并分析差值变化规律来确定光线与DSM表面是否击中。改进算法的可视分析原理更简单，计算量相对于传统算法也更小。
格网数据的遮挡检测中DSM格网的间距与影像分辨率的关系影响遮挡检测的精度与速度。在影像地面分辨率确定的情况下，当DSM间距逐渐变小时，遮挡检测精度提高，但伴随着处理速度的降低；当DSM间距逐渐增大时，检测速度提升，但容易造成复杂建筑物边缘检测不准确。理想的情况是格网间距与影像地面分辨率一致或者相当，以保证精度与速度的协调。
由遮挡检测结果影像与原始影像对比可以发现，改进的格网数据遮挡算法检测精度和效果与TIN数据的算法遮挡检测效果相当，对复杂建筑也能够作出准确的遮挡检测，细节和完整度都十分理想。
6 结 论本文算法将地表数据存储和组织结构的特点与光线追踪算法的准确性相融合，以两类基础地表数据源为对象对光线追踪算法进行了改进。算法中对采样点所在三角形的搜索耗时较多，需提高搜索的效率；另外，格网数据的遮挡检测对检测速度与质量如何平衡应进行更多的分析和试验。遮挡检测只是真正射影像制作的第一步，完整的真正射影像制作包含了遮挡检测、遮挡区域补偿以及阴影的消除等方面[20]，这些技术都是获得完美真正射影像的关键。
正射影像作为测绘数据4D产品之一，与DEM、DLG、DRG相比发挥的作用较弱，大部分正射影像都被简单应用于城市建设与规划、环境保护、土地资源调查等方面，或者仅仅当作一种展示影像，究其原因就是正射影像只纠正了地形的投影差，而建筑物的投影差未作处理。随着真正射影像制作技术的成熟与完善，海量影像数据处理速度的提高，以及更多快速、有效的高精度DSM数据获取方法的出现，必将拓宽真正射影像应用的领域与范围，为社会经济生活提供更多的便利，创造更多的价值。
| [1] | LUIGI B, MARIA A B, LUANA V. LiDAR Digital Building Models for True Orthophoto Generation [J].Geomat Applications,2010(2):187-196. |
| [2] | AMHAR F, JOSEF J, RIES C.The Generation of True Orthophotos Using a 3D Building Model in Conjunction with a Conventional DTM [J]. International Archives of Photogrammetry and Remote Sensing, 1998,32(4):16-22. |
| [3] | RAU J, CHEN N, CHEN L. True Orthophoto Generation of Built-up Areas Using Multi-view Images[J]. Photogrammetric Engineering and Remote Sensing, 2002 (6): 581-588. |
| [4] | AYMAN F H, EUI-MYOUNG K, CHANG-JAE K. New Methodologies for True Orthophoto Generation[J]. Photogrammetric Engineering and Remote Sensing, 2007,73(1):25-36. |
| [5] | ZHONG Cheng, LI Hui, HUANG Xianfeng. A Fast and Effective Approach to Generate True Orthophoto in Built-up Area[J]. Sensor Review, 2011 (4): 341-348. |
| [6] | SHENG Y. Comparative Evaluation of Iterative and Non-iterative Methods to Ground Coordinate Determination from Single Aerial Images [J]. Computers and Geosciences, 2004(30):267-279. |
| [7] | SHENG Y. Theoretical Analysis of the Iterative Photogrammetric Method to Determining Ground Coordinates from Photo Coordinates and a DEM[J]. Photogrammetric Engineering and Remote Sensing, 2005,71(7):863-871. |
| [8] | ZHANG Zuxun, ZHANG Jianqing. Digital Photogrammetry[M].Wuhan: Wuhan University Press,1997.(张祖勋, 张剑清. 数字摄影测量学[M]. 武汉:武汉大学出版社,1997.) |
| [9] | BANG K I, HABIB A F, SHIN S W, et al. Comparative Analysis of Alternative Methodologies for True Ortho-photo Generation from High Resolution Satellite Imagery[C]//Proceedings of ASPRS 2007 Annual Conference. Florida:[s.n.], 2007. |
| [10] | WANG Xiao, JIANG Wanshou, XIE Junfeng. A New Method for True Orthophoto Generation[J].Geomatics and Information Science of Wuhan University,2009,34(10): 1250-1254. (王潇, 江万寿, 谢俊峰. 一种新的真正射影像生成算法[J]. 武汉大学学报: 信息科学版, 2009,34(10): 1250-1254.) |
| [11] | ODA K,LU W,UCHIDA O, et al. Triangle-based Visibility Analysis and True Orthoimage Generation[J].International Archives of Photogrammetry, Remote Sensing and Spatial Information Sciences,2004,35(3):623-628. |
| [12] | SUN Jie, MA Hongchao,ZHONG Liang. LiDAR Point Clouds Based Occlusion Detection of True-ortho Image[J]. Geomatics and Information Science of Wuhan University, 2011,36(8):948-951. (孙杰, 马洪超, 钟良. 利用LiDAR点云的真正射影像遮蔽检测[J]. 武汉大学学报:信息科学版,2011,36(8):948-951.) |
| [13] | XIE Wenhan, ZHOU Guoqing. Occlusion and Shadow Detection of Large-scale True Orthophoto in Urban Area[J]. Acta Geodaetica et Cartographica Sinica, 2010,39(1):52-58. (谢文寒,周国清.城市大比例尺真正射影像阴影与遮挡问题的研究[J].测绘学报,2010,39(1):52-58.) |
| [14] | ZHONG Cheng,HUANG Xianfeng,LI Deren, et al. Polygon Based Inversion Imaging for Occlusion Detection in True Orthophoto Generation[J]. Acta Geodaetica et Cartographica Sinica, 2010,39(1): 59-64. (钟成, 黄先锋, 李德仁,等. 真正射影像生成的多边形反演成像遮蔽检测方法[J].测绘学报,2010,39(1): 59-64.) |
| [15] | KUZMIN Y, KORYTNIK S, LONG O. Polygon-based True Orthophoto Generation[J]. International Archives of Photogrammetry, Remote Sensing and Spatial Information Sciences, 2004 (3): 529-531. |
| [16] | ZHONG Cheng, LI Hui, LI Zonghua, et al. A Vector-based Backward Projection Method for Robust Detection of Occlusions when Generating True Ortho Photos[J]. GIScience & Remote Sensing, 2010 (3): 412-424. |
| [17] | YAHYA A, NORBERT H, DIETER F. A New True Ortho-photo Methodology for Complex Archaeological Applcation[J]. Archaeometry,2010(3):517-530. |
| [18] | SCHICKIER W, THORPE A. Operational Procedure for Automatic True Orthophoto Generation[J]. International Archives of Photogrammetry and Remote Sensing, 1998: 527-532. |
| [19] | BORNER A, WIEST L, KELLER P, et al. A Tool for the Simulation of Hyperspectral Remote Sensing Systems[J]. Photogrammetry and Remote Sensing, 2001,55(5-6):299-312. |
| [20] | ZHOU G, CHEN W, KELMELIS J, et al. A Comprehensive Study on Urban True Orthorectification[J]. IEEE Transactions on Geoscience and Remote Sensing, 2005(9): 2138-2147. |