CN105160702B - 基于LiDAR点云辅助的立体影像密集匹配方法及系统 - Google Patents
基于LiDAR点云辅助的立体影像密集匹配方法及系统 Download PDFInfo
- Publication number
- CN105160702B CN105160702B CN201510513972.5A CN201510513972A CN105160702B CN 105160702 B CN105160702 B CN 105160702B CN 201510513972 A CN201510513972 A CN 201510513972A CN 105160702 B CN105160702 B CN 105160702B
- Authority
- CN
- China
- Prior art keywords
- point cloud
- disparity map
- parallax
- layer
- disparity
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
Landscapes
- Image Processing (AREA)
Abstract
一种基于LiDAR点云辅助的立体影像密集匹配方法及系统,将投影获取立体像对重叠范围内的LiDAR点云并滤波处理;将滤波后的点云投影到核线立体像对上,确定后续密集匹配的视差范围;建立金字塔,从金字塔顶层开始,采用三角网约束改造代价矩阵,进行SGM密集匹配,并进行左右一致性检测,得到顶层最终的视差图;将金字塔当前层的视差图传递到下一层作为初始视差图,根据当前层的视差图相应确定下一层的视差范围,将下一层作为新的当前层;基于新的当前层,同样处理得到当前层最终的视差图,直到金字塔的底层,输出原始影像的视差图,根据视差图,得到立体像对的同名点,生成密集匹配的点云。
Description
技术领域
本发明涉及一种立体像对密集匹配技术领域,尤其是涉及一种控制信息引导的立体像对密集匹配技术方案。
背景技术
在摄影测量和计算机视觉领域,密集匹配是指根据两张或者多张光学影像,逐像素地匹配同名点,生成稠密的三维点云的过程。密集匹配技术能够快速获取测区地表的三维测绘信息,在道路和水体提取、建筑物三维建模、变化检测、DEM/DSM和真正射影像制作、车辆自动驾驶等智慧城市、地理国情监测领域有着非常广泛的应用前景。
激光扫描和立体影像密集匹配技术是两种主流的地形地物三维信息获取手段。LiDAR系统以主动方式快速获取地面模型及空间点云信息,具有速度快、时效性强、作业范围大等优点。但是LiDAR系统成本昂贵,且点云密度较低(与对应的影像分辨率比较而言),无法很好地描述建筑物等线状地物的精确边缘,不能完全满足智慧城市等高精度应用的需求;并且由于POS(尤其是精度较低的国产POS)系统误差及其他误差的存在,导致机载LiDAR点云存在复杂的系统误差。影像密集匹配是一种被动获取空间点云的手段。与LiDAR激光点云相比,立体影像具有密集匹配点云密度大、平面精度高、粗差剔除技术完善、影像上的地物几何特征更明显、影像数据获取成本低等突出优势。但是影像密集匹配技术在纹理贫乏区域、重复纹理区域表现较弱。
由于上述两种技术各自的不足,单纯地采用激光扫描,或者单纯地采用影像密集匹配,均无法很好地描述测区地表的三维信息。目前国内外尚未有结合两种截然不同的技术的相应研究。
发明内容
本发明主要是解决LiDAR点云三维重建技术和密集匹配三维重建技术相互孤立的现状,以及密集匹配视差范围参数需要人工设定阈值等的技术问题;提供了一种LiDAR点云辅助的立体影像半全局密集匹配技术,能够自动确定密集匹配的视差范围,通过LiDAR点云指导密集匹配,从而得到更为精确可靠的密集匹配产品。
本发明技术方案提供一种基于LiDAR点云辅助的立体影像密集匹配方法,在LiDAR点云的引导下,逐像素地匹配立体像对同名点,生成稠密的三维点云,实现过程包括以下步骤,
步骤1,将整个测区的LiDAR点云投影到原始的立体像对上,获取立体像对重叠范围内的LiDAR点云,并对LiDAR点云进行滤波处理;
步骤2,将步骤1所得滤波后的点云投影到核线立体像对上,计算同名点对的视差,找出最大视差和最小视差,确定后续密集匹配的视差范围;分别建立左右核线影像和左右影像的LiDAR点云控制图的金字塔;
步骤3,从金字塔顶层开始,令顶层为当前层,首先以当前层的左核线影像为基准影像,将顶层左核线影像LiDAR点云控制图作为初始视差图,采用三角网约束改造代价矩阵,进行SGM密集匹配,得到左核线影像的视差图,依次对右核线影像进行类似的处理,得到右核线影像的视差图,根据两张视差图,进行左右一致性检测,得到顶层最终的视差图;
步骤4,将金字塔当前层的视差图传递到下一层作为初始视差图,根据当前层的视差图相应确定下一层的视差范围,将下一层作为新的当前层;基于新的当前层,仍然以当前层的左核线影像为基准影像,采用三角网约束改造代价矩阵,进行SGM密集匹配,得到左核线影像的视差图,依次对右核线影像进行类似的处理,得到右核线影像的视差图,根据两张视差图,进行左右一致性检测,得到当前层最终的视差图;
步骤5,判断当前层是否为金字塔的底层,若是则执行步骤6,否则返回执行步骤4继续传递;
步骤6,输出原始影像的视差图,根据视差图,得到立体像对的同名点,生成密集匹配的点云。
而且,步骤1中,对LiDAR点云进行滤波处理的实现方式为,生成左右影像的LiDAR点云控制图,在LiDAR点云控制图中,以每个有效点为中心建立局部窗口,计算局部窗口内中心点在窗口内各有效点视差下的可见性度量,比较最大度量vmax相应的视差dmax与窗口中心点像素视差dc的差值,如果差值在T个像素以内,则认为中心点是正确点,否则认为中心点是错误点,予以剔除;T为预设的相应阈值。
而且,步骤3和步骤4中,所述计算并改造代价矩阵的实现方式为,遍历核线影像,在视差范围内依次计算每个像素的候选同名点,根据候选同名点对计算出对应的互信息,作为匹配代价,得到匹配代价矩阵;所述的三角网约束的实现方式为,采用Delaunay三角剖分算法,构建TIN三角网,三角网的顶点是有效视差点,用顶点的有效视差,约束三角形内部点的视差;所述半全局密集匹配的实现方式为,将一条直线上的每个像素作为一个阶段,像素的匹配代价作为节点,将匹配问题转化为动态规划问题求解。
而且,步骤4中,将金字塔当前层的视差图传递到下一层的实现方式为,将上一层像素的视差值,乘以缩放比例Size,再赋值给下一层金字塔对应区域内的所有像素;由当前层的视差图相应确定的下一层视差范围如下,
sd=Size·d'-Tol
ed=Size·d'+Tol
其中,sd表示视差范围的起点;ed表示视差范围的终点;Size表示金字塔影像缩放比例;d'表示上一级金字塔的视差;Tol表示阈值。
一种基于LiDAR点云辅助的立体影像密集匹配系统,用于在LiDAR点云的引导下,逐像素地匹配立体像对同名点,生成稠密的三维点云,包括以下模块,
点云滤波模块,用于将整个测区的LiDAR点云投影到原始的立体像对上,获取立体像对重叠范围内的LiDAR点云,并对LiDAR点云进行滤波处理;
金字塔构建模块,用于将点云初始化模块所得滤波后的点云投影到核线立体像对上,计算同名点对的视差,找出最大视差和最小视差,确定后续密集匹配的视差范围;分别建立左右核线影像和左右影像的LiDAR点云控制图的金字塔;
金字塔顶层视差图生成模块,用于从金字塔顶层开始,令顶层为当前层,首先以当前层的左核线影像为基准影像,将顶层左核线影像LiDAR点云控制图作为初始视差图,采用三角网约束改造代价矩阵,进行SGM密集匹配,得到左核线影像的视差图,依次对右核线影像进行类似的处理,得到右核线影像的视差图,根据两张视差图,进行左右一致性检测,得到顶层最终的视差图;
视差图传递生成模块,用于将金字塔当前层的视差图传递到下一层作为初始视差图,根据当前层的视差图相应确定下一层的视差范围,将下一层作为新的当前层;基于新的当前层,仍然以当前层的左核线影像为基准影像,采用三角网约束改造代价矩阵,进行SGM密集匹配,得到左核线影像的视差图,依次对右核线影像进行类似的处理,得到右核线影像的视差图,根据两张视差图,进行左右一致性检测,得到当前层最终的视差图;
判断模块,用于判断当前层是否为金字塔的底层,若是则命令输出模块执行工作,否则命令视差图传递生成模块继续传递;
输出模块,用于输出原始影像的视差图,根据视差图,得到立体像对的同名点,生成密集匹配的点云。
而且,点云初始化模块中,对LiDAR点云进行滤波处理的实现方式为,生成左右影像的LiDAR点云控制图,在LiDAR点云控制图中,以每个有效点为中心建立局部窗口,计算局部窗口内中心点在窗口内各有效点视差下的可见性度量,比较最大度量vmax相应的视差dmax与窗口中心点像素视差dc的差值,如果差值在T个像素以内,则认为中心点是正确点,否则认为中心点是错误点,予以剔除;T为预设的相应阈值。
而且,金字塔顶层视差图生成模块和视差图传递生成模块中,所述计算并改造代价矩阵的实现方式为,遍历核线影像,在视差范围内依次计算每个像素的候选同名点,根据候选同名点对计算出对应的互信息,作为匹配代价,得到匹配代价矩阵;所述的三角网约束的实现方式为,采用Delaunay三角剖分算法,构建TIN三角网,三角网的顶点是有效视差点,用顶点的有效视差,约束三角形内部点的视差;所述半全局密集匹配的实现方式为,将一条直线上的每个像素作为一个阶段,像素的匹配代价作为节点,将匹配问题转化为动态规划问题求解。
而且,视差图传递生成模块中,将金字塔当前层的视差图传递到下一层的实现方式为,将上一层像素的视差值,乘以缩放比例Size,再赋值给下一层金字塔对应区域内的所有像素;由当前层的视差图相应确定的下一层视差范围如下,
sd=Size·d'-Tol
ed=Size·d'+Tol
其中,sd表示视差范围的起点;ed表示视差范围的终点;Size表示金字塔影像缩放比例;d'表示上一级金字塔的视差;Tol表示阈值。
因此,本发明具有如下优点:在密集匹配过程中无需人工设定视差搜索范围,通过LiDAR点云信息引导密集匹配过程,通过密集匹配加密、精化LiDAR点云,匹配结果优于原始密集匹配点云或者LiDAR点云,在航空摄影测量、低空摄影测量和近景摄影测量领域均有较好的应用前景。
附图说明
图1为本发明实施例的互信息计算流程图;
图2为本发明实施例的局部统计窗口示意图;
图3为本发明实施例的基于动态规划的匹配流程图;
图4为本发明实施例的匹配路径方向示意图;
图5为本发明实施例的流程图。
具体实施方式
下面结合实施例及附图对本发明作进一步详细的描述。
针对现有技术中激光扫描和采用影像密集匹配的不足,本发明考虑到如果将两种技术相结合,通过密集匹配来加密LiDAR点云,通过LiDAR点云来引导密集匹配过程,充分发挥各自的优势,理论上,可以得到更加精确、密集的三维点云。从该思路出发,本发明设计提供的技术方案是,一种利用LiDAR点云的三维信息,自动确定密集匹配的视差范围,指引密集匹配过程,得到更精确可靠的密集匹配点云的方法。
具体实施时,本发明技术方案可由本领域技术人员采用计算机软件技术实现自动运行流程,在LiDAR点云的引导下,逐像素地匹配立体像对同名点,生成稠密的三维点云。参见图5,实施例的流程包括以下步骤:
步骤1.获取立体像对重叠范围内的LiDAR点云,并对LiDAR点云进行滤波处理。
本发明提出,根据立体像对的内外方位元素,将整个测区的LiDAR点云投影到原始立体像对上,如式(1)所示。如果某投影点的像方坐标同时位于左右两张原始影像上,则说明该点是影像重叠范围内的点;否则不是。保存影像重叠范围内的LiDAR点云。
式中,x、y表示投影点的像点坐标;X、Y、Z表示LiDAR点云的物方坐标;x0、y0、f表示相机的内方位元素;Xs、Ys、Zs表示相机的外方位线元素;a1~a3、b1~b3、c1~c3表示旋转矩阵的9个元素。
确定影像重叠范围内的LiDAR点云后,需要进一步的点云滤波,以剔除点云中的粗差点和被遮挡点。通常情况下,可见点的视差大于被遮挡点的视差,且可见点的相关系数要大于被遮挡点的相关系数。考虑到一般同名点的视差与高程成正比,高程越高,该点在物方的可见性越大;此外,正确匹配点在像方的相关系数一般也较大,可以作为一个像方的约束。联合视差和相关系数信息,可以从物方和像方两个方面,对点云进行滤波。本发明按此进一步提出结合视差和相关系数准则进行点云滤波。
实施例的步骤1具体分为如下3个操作步骤:
1)核线影像生成。为了避免多次重复核线采样,造成匹配时间较长,在匹配前,先根据相机的内外方位元素,生成核线影像,左核线影像和右核线影像构成核线立体像对。核线影像生成为现有技术,在此不做赘述。核线影像上,同名点的y坐标相同,x方向的坐标差即为视差,如式(2)所示。
xr=xl-d yr=yl (2)
式中,xl、yl为左核线影像的像点坐标;xr、yr为右核线影像的同名点坐标;d表示视差。
2)视差计算。按照式(1),将LiDAR点云投影到原始立体像对,得到原始影像上的同名点坐标。按照核线影像和原始影像之间的几何关系,将原始影像上同名点坐标换算到核线影像上,并根据式(2),计算相应点的视差d。可以新建与核线影像同样大小的内存,将视差d作为灰度值,保存到对应的投影点位上,分别生成左右核线影像的LiDAR点云控制图。
将LiDAR点云同时投影到左右核线影像上,以每个LiDAR点的视差作为灰度值,生成左右核线影像的LiDAR点云控制图时,由于LiDAR点云分辨率远远小于影像分辨率,所以不是所有的影像点都有对应的LiDAR点,投影点为有效点;非投影点为无效点。
3)滤波。对2)所得左右影像的LiDAR点云控制图中的各点分别执行是否保留的判断,操作如下:
首先进行相关系数计算。对核线影像进行图像分割,将灰度相似的像素归于一类。在LiDAR点云控制图中,以每个有效点为中心,开辟一个W×W大小的局部窗口,如图2所示。
图2中提供了对核线影像进行影像分割分为两类的结果;网状背景的格子代表有效点。根据影像分割的结果,统计局部窗口内所有与中心像素隶属同一类的有效点并记录每个有效点的视差。
如果以左核线影像为基准,首先在左核线影像LiDAR点云控制图上,开辟一个W×W大小的窗口,窗口的中心点是有效点,窗口内也存在其它有效点。设窗口中心点有个视差值dc,窗口内的其余有效点,也存在视差值di(i=1...N)。在密集匹配理论中,假设场景的视差是平滑的,因此在局部W×W窗口内,所有有效点的视差应该是基本一致的。如果存在不一致的视差,说明很有可能有错误点。由于未知哪个视差值是错的,因此借用相关系数的方法来判断。统计窗口内的所有有效点的视差值,根据式(2),计算得到窗口中心点在右影像上的同名点,然后计算对应的相关系数。一个视差值会对应一个相关系数,窗口内多个视差值会对应多个相关系数。
具体实现为,根据统计得到的有效点所对应的视差,按照式(2),计算窗口中心点在另一张核线影像(如果是左影像的LiDAR点云控制图,则另一张核线影像指右核线影像;反之,则为左核线影像)上的“同名点”。由于会统计得到多个视差,因此,会对应多个“同名点”。计算窗口中心点与这些“同名点”之间的相关系数,如式(3)所示:
式中,r表示相关系数;xi表示左影像相关系数窗口内的某像素灰度;yi表示右影像相关系数窗口内的某像素灰度;表示左影像的LiDAR点云控制图窗口像素灰度均值;表示右影像的LiDAR点云控制图窗口像素灰度均值;n表示窗口内像素总数。
然后进行局部窗口中心点的可见性度量计算。将局部窗口内的所有有效点视差和对应的相关系数,按照式(4)的方式进行计算,得到窗口中心点在某视差下的可见性度量:
vi=di·ri (4)
式中,vi表示可见性度量,值越大,表示该点为正确点的可能性越大;di表示统计的有效点所对应的视差,值越大,代表该点在影像上的可见性越大,值越小,则表示被遮挡的可能性越大;ri表示对应的相关系数。在所有的可见性度量vi(i=1…N,N表示窗口内有效点的数目)中,找出一个值最大的度量vmax,比较最大度量vmax相应的视差dmax与窗口中心点像素视差dc的差值。如果差值在T个像素以内(本领域技术人员可自行设定T的取值,一般取1),则认为中心点是正确点;否则认为中心点是错误点,予以剔除,如式(5)所示:
当一个点在左右LiDAR控制图的计算结果均为“保留”时,则保留该点;否则,剔除该点。然后根据以上方式,依次处理其余点。当所有点均处理完毕,即结束点云滤波操作。
步骤2.确定密集匹配的视差范围,分别建立左右核线影像和左右影像的LiDAR点云控制图的影像金字塔。
将步骤1所得滤波后的点云投影到核线立体像对上,采用与步骤1同样的方式按照式(2)计算同名点对的视差,找出其中的最大视差和最小视差,作为后续密集匹配的视差范围。
建立某核线影像(左核线影像或右核线影像)或某LiDAR点云控制图(左影像的LiDAR点云控制图或右影像的LiDAR点云控制图)的影像金字塔的实现方式为,将核线影像(或者LiDAR控制图)作为金字塔的底层,按照一定的规则,将相邻S×S窗口内的像素,重采样生成一个新的像素,传到上一级金字塔,如此重复,直至计算到金字塔的顶层。
如果以核线影像为重采样的底图,则采用平均值规则,生成一个新的像素,如式(6)所示:
其中,G'表示新像素灰度;Gi表示下一级金字塔窗口内的像素灰度;S表示金字塔的缩放因子。
如果以LiDAR控制图为重采样的底图,则采用最大视差原则,生成一个新的像素:
D'=max(Di)/S i=1...N′ (7)
其中,D'表示新像素视差;Di表示下一级金字塔窗口内的某有效点像素视差;N'表示金字塔窗口内的有效点数目。
在每一级金字塔,匹配视差范围也要按照比例因子S进行相应的缩放,如式(8)所示:
minDi=MinD/Si-1 Ri=DispRange/Si-1 (8)
式中,minDi表示第i层金字塔对应的最小视差;MinD表示金字塔底层的最小视差;Ri表示第i层金字塔对应的视差范围;DispRange表示金字塔底层的视差范围。
步骤3.从金字塔顶层开始,采用三角网约束改造代价矩阵,进行半全局密集匹配,然后进行左右一致性检测,得到顶层的视差图。流程可设置为,从金字塔顶层开始,令顶层为当前层,首先以当前层的左核线影像为基准影像,将顶层左核线影像LiDAR点云控制图作为初始视差图,采用三角网约束改造代价矩阵,进行SGM密集匹配,得到左核线影像的视差图,依次对右核线影像进行类似的处理,得到右核线影像的视差图,根据两张视差图,进行左右一致性检测,得到顶层最终的视差图。
立体像对对应两张核线影像,在匹配过程中,不可避免的会出现“误匹配”情况。为了防止将误匹配点传递到下一层金字塔,本发明采用左右一致性检测的办法,即对当前层:首先将左核线影像作为基准影像,计算并改造代价矩阵,然后半全局密集匹配,计算出一张视差图;其次,将右核线影像作为一张基准影像,计算并改造代价矩阵,然后半全局密集匹配,计算出新的一张视差图;将两张视差图上,不一致的点,作为“误匹配点”,予以剔除。
(1)计算并改造代价矩阵。
在金字塔顶层,将顶层LiDAR点云控制图作为初始视差图,计算互信息,得到立体像对的代价矩阵,根据顶层LiDAR点云控制图,构建TIN三角网,通过三角网约束改造代价矩阵。
首先,在金字塔顶层,将顶层LiDAR点云控制图作为初始视差图,计算两张影像之间的互信息,如图1所示。
在图1中,首先,根据初始视差图,按照式(2),计算每个有效点的同名点。根据同名点对的坐标,在核线影像对中两张核线影像上找出对应的灰度,统计一张二维的联合直方图,二维灰度概率密度计算如式(9)所示:
式中,i、k分别表示任意同名点在两张核线影像上的灰度;表示初始视差图中的有效点数目;表示当前同名点在两张核线影像上的灰度,下标p表示核线影像的像素;T[.]是二值化算子,当括号内的条件为真时,则T算子等于1,否则,T算子等于0。式(9)定义了左右影像上的同名点,其灰度分别是i和k的概率。
一维直方图的概率密度函数,可以由二维直方图直接计算得到,如式(10)所示:
式中,PI1(i)表示核线影像I1的灰度概率密度;PI2(k)表示核线影像I2的灰度概率密度。
根据灰度概率密度,可以按照式(11),计算点与点之间的互信息hI1,I2(i,k):
式中,表示灰度i和灰度k之间的互信息;g(i)表示一维高斯核函数;g(i,k)表示二维高斯核函数;PI(i)表示核线影像I的灰度概率密度,hI(i)代表核线影像I上,灰度i的信息熵,核线影像I指核线影像I1或核线影像I2,具体地,代表核线影像I1上,灰度i的信息熵,代表核线影像I2上,灰度k的信息熵。具体计算过程如图1所示。先将一维/二维直方图采用Gauss卷积进行平滑,再采用log运算,最后再对结果进行Gauss平滑。
根据上述计算过程后,得到关于两张影像灰度之间的关系函数,然后可以针对像素点进行密集匹配。在视差范围内,可以根据视差依次计算每个像素的“候选同名点”。根据“候选同名点对”,可以按照式(11)计算出对应的互信息,作为匹配代价。每个像素都可以按照视差范围,计算相应的匹配代价。按照这种方式遍历整张基准影像(核线影像),可以得到一个关于匹配代价的三维矩阵,即为匹配代价矩阵。
根据顶层LiDAR控制图,采用Delaunay三角剖分算法,构建TIN三角网。三角网的顶点是有效视差点,用顶点的有效视差,约束三角形内部点的视差。对于每一个三角形内部的点,视差约束满足如下条件:
式中,Costi(d)表示顶层LiDAR控制图第i个像素对应视差d的匹配代价;c表示约束前的代价;p表示惩罚系数;dmin、dmax表示三角形三个顶点中的最小视差、最大视差。
(2)半全局密集匹配,计算新的视差图,即得到当前层金字塔的视差图。
根据改造后的匹配代价矩阵,可以得到每个像素在视差范围内对应的匹配代价。定义一个如式(13)所示的目标函数,目标函数的解即为满足能量函数最小条件下的视差图。将任意方向直线上的每个像素作为一个阶段,像素的匹配代价作为节点,实际上,可以将匹配问题转化为动态规划问题来求解,如图3所示。图3中,每一列表示一个阶段,一个像素对应一列;每个圆形代表一个节点,代表像素在对应视差下的匹配代价,节点的数目即为视差范围。根据邻近像素视差平滑的条件,可以将匹配问题转化为分阶段的决策问题,从而便于采用动态规划的方法来求解。
式中,E(D)表示目标函数值;D表示待计算的新的视差图;C(p,Dp)表示像素p在给定视差Dp情况下的匹配代价;Np表示p的邻域,像素q为领域Np内的像素;P1、P2表示惩罚系数。
由此,本发明提出采用半全局密集匹配方法(SGM)如下:
根据改造后的代价矩阵,可以得到每个像素在视差范围内的匹配代价;将一条直线上的每个像素作为一个阶段,像素的匹配代价作为节点,实际上,可以将匹配问题转化为动态规划问题来求解。动态规划问题的求解方法,如式(14)所示:
式中,p表示当前像素;p-r表示前一个阶段的像素;Lr(,)表示动态规划的路径;P1、P2表示惩罚系数;r表示动态规划的步距,一般取1;k表示前一个像素的最优视差值;i表示前一个像素所对应的某个视差值。
匹配结果即为动态规划的最优路径。由于一维匹配很容易产生“条纹状”的误匹配,因此采用多个方向的动态规划策略,将匹配路径按照图4所示的方式,划分为8~16个方向,将所有方向的匹配结果进行累加,能够得到更加稳健的匹配结果。图4中,p表示当前的像素;x、y表示像平面坐标的x轴、y轴。
步骤4.视差图传递:将金字塔当前层的视差图传递到下一层作为初始视差图,根据当前层的视差图相应确定下一层的视差范围,将下一层作为新的当前层;基于新的当前层,仍然以当前层的左核线影像为基准影像,采用三角网约束改造代价矩阵,进行SGM密集匹配,得到左核线影像的视差图,依次对右核线影像进行类似的处理,得到右核线影像的视差图,根据两张视差图,进行左右一致性检测,得到当前层最终的视差图。
本步骤将上一级金字塔的视差图传递到下一层,约束下一级金字塔的匹配搜索范围,并计算相应的代价矩阵,根据当前层的LiDAR点云控制图,重新构建三角网并约束代价矩阵,采用SGM方法得到当前层的视差图。
首先,将上一级金字塔的视差图传递到下一层。上一层金字塔的一个像素,对应下一层金字塔Size·Size大小的区域。具体传递方式是:将上一层像素的视差值,乘以缩放比例Size,再赋值给下一层金字塔对应区域内的所有像素。根据上一级金字塔匹配出来的视差图,预测当前金字塔影像的同名点点位,根据相应阈值,约束当前金字塔的匹配搜索范围,根据新的视差搜索范围进行密集匹配计算:
sd=Size·d'-Tol
(15)
ed=Size·d'+Tol
式中,sd表示视差范围的起点;ed表示视差范围的终点;Size表示金字塔影像缩放比例;d'表示上一级金字塔的视差;Tol表示阈值,与金字塔的缩放比例有关,可以由本领域技术人员自行根据实际需要来设定,例如设置的金字塔缩放比例Size为2,Tol相应地可以设为5。
和步骤3相似,按照式(11),根据传递下来的视差图以及新的视差范围,重新计算当前金字塔的代价矩阵。根据当前层的LiDAR点云控制图,重新构建三角网并约束代价矩阵,采用SGM方法得到当前层的视差图,采用左右一致性检测的方法剔除误匹配点,并向下一级金字塔传递。
步骤5.判断执行。
判断是否已计算到影像金字塔的底层,若是,则执行步骤6,否则返回执行步骤4,直至计算到金字塔底层为止。
步骤6.输出密集匹配点云。本发明根据视差图,逐像素地计算立体像对的同名点,根据同名点,采用前方交会的方法得到物方三维点云。
因为之前经过了左右一致性检测,已经将不一样的地方剔除,所以剩下的左右视差图已经一样的。实施例中,计算得到底层的左核线影像视差图后,根据式(2),可以得到核线影像上每个像素的同名点。将所有同名点对坐标转换到原始影像上的坐标,再采用前方交会的方法得到物方三维点云。
三维点云的计算公式为:
式中,X、Y、Z表示物方点的三维坐标;xi、yi(i=1,2)表示像点坐标;f表示焦距;Ri(i=1,2)表示旋转矩阵;Xsi、Ysi、Zsi(i=1,2)表示外方位线元素;N1、N2表示投影系数;u1、v1、w1和u2、v2、w2分别表示像点在左、右片像空间辅助坐标系中的坐标;U1、V1、W1和U2、V2、W2分别表示物方点在左、右片像空间辅助坐标系中的三维坐标;Bu、Bv、Bw表示基线的三个坐标分量。
具体实施时,还可以采用模块化方式提供一种基于LiDAR点云辅助的立体影像密集匹配系统,用于在LiDAR点云的引导下,逐像素地匹配立体像对同名点,生成稠密的三维点云,包括以下模块:
点云滤波模块,用于将整个测区的LiDAR点云投影到原始的立体像对上,获取立体像对重叠范围内的LiDAR点云,并对LiDAR点云进行滤波处理;
金字塔构建模块,用于将点云初始化模块所得滤波后的点云投影到核线立体像对上,计算同名点对的视差,找出最大视差和最小视差,确定后续密集匹配的视差范围;分别建立左右核线影像和左右影像的LiDAR点云控制图的金字塔;
金字塔顶层视差图生成模块,用于从金字塔顶层开始,令顶层为当前层,首先以当前层的左核线影像为基准影像,将顶层左核线影像LiDAR点云控制图作为初始视差图,采用三角网约束改造代价矩阵,进行SGM密集匹配,得到左核线影像的视差图,依次对右核线影像进行类似的处理,得到右核线影像的视差图,根据两张视差图,进行左右一致性检测,得到顶层最终的视差图;
视差图传递生成模块,用于将金字塔当前层的视差图传递到下一层作为初始视差图,根据当前层的视差图相应确定下一层的视差范围,将下一层作为新的当前层;基于新的当前层,仍然以当前层的左核线影像为基准影像,采用三角网约束改造代价矩阵,进行SGM密集匹配,得到左核线影像的视差图,依次对右核线影像进行类似的处理,得到右核线影像的视差图,根据两张视差图,进行左右一致性检测,得到当前层最终的视差图;
判断模块,用于判断当前层是否为金字塔的底层,若是则命令输出模块执行工作,否则命令视差图传递生成模块继续传递;
输出模块,用于输出原始影像的视差图,根据视差图,得到立体像对的同名点,生成密集匹配的点云。
各模块具体实现与方法步骤相应,本发明不予赘述。
本文中所描述的具体实施例仅仅是对本发明精神做举例说明。本发明技术领域的技术人员可以对所描述的具体实施例做各种各样的修改补充或者采用类似的方式替代,但并不会偏离本发明精神或者超越所附权利要求书所定义的范围。
Claims (8)
1.一种基于LiDAR点云辅助的立体影像密集匹配方法,其特征在于:在LiDAR点云的引导下,逐像素地匹配立体像对同名点,生成稠密的三维点云,实现过程包括以下步骤,
步骤1,将整个测区的LiDAR点云投影到原始的立体像对上,获取立体像对重叠范围内的LiDAR点云,并对LiDAR点云进行滤波处理;
步骤2,将步骤1所得滤波后的点云投影到核线立体像对上,计算同名点对的视差,找出最大视差和最小视差,确定后续密集匹配的视差范围;分别建立左右核线影像和左右核线影像的LiDAR点云控制图的金字塔;
步骤3,从金字塔顶层开始,令顶层为当前层,首先以当前层的左核线影像为基准影像,将顶层左核线影像LiDAR点云控制图作为初始视差图,采用三角网约束改造代价矩阵,进行SGM密集匹配,得到左核线影像的视差图,依次对右核线影像进行类似的处理,得到右核线影像的视差图,根据两张视差图,进行左右一致性检测,得到顶层最终的视差图;
步骤4,将金字塔当前层的视差图传递到下一层作为初始视差图,根据当前层的视差图相应确定下一层的视差范围,将下一层作为新的当前层;基于新的当前层,仍然以新的当前层的左核线影像为基准影像,采用三角网约束改造代价矩阵,进行SGM密集匹配,得到左核线影像的视差图,依次对右核线影像进行类似的处理,得到右核线影像的视差图,根据两张视差图,进行左右一致性检测,得到新的当前层最终的视差图;
步骤5,判断新的当前层是否为金字塔的底层,若是则执行步骤6,否则返回执行步骤4继续传递;
步骤6,输出原始影像的视差图,根据视差图,得到立体像对的同名点,生成密集匹配的点云。
2.根据权利要求1基于LiDAR点云辅助的立体影像密集匹配方法,其特征在于:步骤1中,对LiDAR点云进行滤波处理的实现方式为,生成左右影像的LiDAR点云控制图,在LiDAR点云控制图中,以每个有效点为中心建立局部窗口,计算局部窗口内中心点在窗口内各有效点视差下的可见性度量,比较最大可见性度量vmax相应的视差dmax与窗口中心点像素视差dc的差值,如果差值在T个像素以内,则认为中心点是正确点,否则认为中心点是错误点,予以剔除;T为预设的相应阈值。
3.根据权利要求1或2基于LiDAR点云辅助的立体影像密集匹配方法,其特征在于:步骤3和步骤4中,所述改造代价矩阵的实现方式为,遍历核线影像,在视差范围内依次计算每个像素的候选同名点,根据候选同名点对计算出对应的互信息,作为匹配代价,得到匹配代价矩阵;所述的三角网约束的实现方式为,采用Delaunay三角剖分算法,构建TIN三角网,三角网的顶点是有效视差点,用顶点的有效视差,约束三角形内部点的视差;所述SGM密集匹配的实现方式为,将一条直线上的每个像素作为一个阶段,像素的匹配代价作为节点,将匹配问题转化为动态规划问题求解。
4.根据权利要求1或2所述基于LiDAR点云辅助的立体影像密集匹配方法,其特征在于:步骤4中,将金字塔当前层的视差图传递到下一层的实现方式为,将上一层像素的视差值,乘以缩放比例Size,再赋值给下一层金字塔对应区域内的所有像素;由当前层的视差图相应确定的下一层视差范围如下,
sd=Size·d'-Tol
ed=Size·d'+Tol
其中,sd表示视差范围的起点;ed表示视差范围的终点;Size表示金字塔影像缩放比例;d'表示上一级金字塔的视差;Tol表示阈值。
5.一种基于LiDAR点云辅助的立体影像密集匹配系统,其特征在于:用于在LiDAR点云的引导下,逐像素地匹配立体像对同名点,生成稠密的三维点云,包括以下模块,
点云滤波模块,用于将整个测区的LiDAR点云投影到原始的立体像对上,获取立体像对重叠范围内的LiDAR点云,并对LiDAR点云进行滤波处理;
金字塔构建模块,用于将点云初始化模块所得滤波后的点云投影到核线立体像对上,计算同名点对的视差,找出最大视差和最小视差,确定后续密集匹配的视差范围;分别建立左右核线影像和左右核线影像的LiDAR点云控制图的金字塔;
金字塔顶层视差图生成模块,用于从金字塔顶层开始,令顶层为当前层,首先以当前层的左核线影像为基准影像,将顶层左核线影像LiDAR点云控制图作为初始视差图,采用三角网约束改造代价矩阵,进行SGM密集匹配,得到左核线影像的视差图,依次对右核线影像进行类似的处理,得到右核线影像的视差图,根据两张视差图,进行左右一致性检测,得到顶层最终的视差图;
视差图传递生成模块,用于将金字塔当前层的视差图传递到下一层作为初始视差图,根据当前层的视差图相应确定下一层的视差范围,将下一层作为新的当前层;基于新的当前层,仍然以新的当前层的左核线影像为基准影像,采用三角网约束改造代价矩阵,进行SGM密集匹配,得到左核线影像的视差图,依次对右核线影像进行类似的处理,得到右核线影像的视差图,根据两张视差图,进行左右一致性检测,得到新的当前层最终的视差图;
判断模块,用于判断新的当前层是否为金字塔的底层,若是则命令输出模块执行工作,否则命令视差图传递生成模块继续传递;
输出模块,用于输出原始影像的视差图,根据视差图,得到立体像对的同名点,生成密集匹配的点云。
6.根据权利要求5基于LiDAR点云辅助的立体影像密集匹配系统,其特征在于:点云初始化模块中,对LiDAR点云进行滤波处理的实现方式为,生成左右影像的LiDAR点云控制图,在LiDAR点云控制图中,以每个有效点为中心建立局部窗口,计算局部窗口内中心点在窗口内各有效点视差下的可见性度量,比较最大可见性度量vmax相应的视差dmax与窗口中心点像素视差dc的差值,如果差值在T个像素以内,则认为中心点是正确点,否则认为中心点是错误点,予以剔除;T为预设的相应阈值。
7.根据权利要求5或6基于LiDAR点云辅助的立体影像密集匹配系统,其特征在于:金字塔顶层视差图生成模块和视差图传递生成模块中,所述改造代价矩阵的实现方式为,遍历核线影像,在视差范围内依次计算每个像素的候选同名点,根据候选同名点对计算出对应的互信息,作为匹配代价,得到匹配代价矩阵;所述的三角网约束的实现方式为,采用Delaunay三角剖分算法,构建TIN三角网,三角网的顶点是有效视差点,用顶点的有效视差,约束三角形内部点的视差;所述SGM密集匹配的实现方式为,将一条直线上的每个像素作为一个阶段,像素的匹配代价作为节点,将匹配问题转化为动态规划问题求解。
8.根据权利要求5或6所述基于LiDAR点云辅助的立体影像密集匹配系统,其特征在于:视差图传递生成模块中,将金字塔当前层的视差图传递到下一层的实现方式为,将上一层像素的视差值,乘以缩放比例Size,再赋值给下一层金字塔对应区域内的所有像素;由当前层的视差图相应确定的下一层视差范围如下,
sd=Size·d'-Tol
ed=Size·d'+Tol
其中,sd表示视差范围的起点;ed表示视差范围的终点;Size表示金字塔影像缩放比例;d'表示上一级金字塔的视差;Tol表示阈值。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201510513972.5A CN105160702B (zh) | 2015-08-20 | 2015-08-20 | 基于LiDAR点云辅助的立体影像密集匹配方法及系统 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201510513972.5A CN105160702B (zh) | 2015-08-20 | 2015-08-20 | 基于LiDAR点云辅助的立体影像密集匹配方法及系统 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN105160702A CN105160702A (zh) | 2015-12-16 |
CN105160702B true CN105160702B (zh) | 2017-09-29 |
Family
ID=54801544
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201510513972.5A Active CN105160702B (zh) | 2015-08-20 | 2015-08-20 | 基于LiDAR点云辅助的立体影像密集匹配方法及系统 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN105160702B (zh) |
Families Citing this family (29)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN106228593B (zh) * | 2015-05-28 | 2019-05-17 | 长沙维纳斯克信息技术有限公司 | 一种图像密集匹配方法 |
CN105466399B (zh) * | 2016-01-11 | 2019-09-06 | 中测新图(北京)遥感技术有限责任公司 | 快速半全局密集匹配方法和装置 |
CN106780712B (zh) * | 2016-10-28 | 2021-02-05 | 武汉市工程科学技术研究院 | 联合激光扫描和影像匹配的三维点云生成方法 |
CN106530337A (zh) * | 2016-10-31 | 2017-03-22 | 武汉市工程科学技术研究院 | 基于图像灰度引导的非局部立体像对密集匹配方法 |
CN106960450B (zh) * | 2017-02-17 | 2019-08-20 | 武汉云航工程地球物理研究院有限公司 | 基于块的影像匹配数字表面模型的全局高程优化方法 |
CN107170000B (zh) * | 2017-04-18 | 2019-09-10 | 武汉市工程科学技术研究院 | 基于全局块优化的立体影像密集匹配方法 |
CN107067394A (zh) * | 2017-04-18 | 2017-08-18 | 中国电子科技集团公司电子科学研究院 | 一种倾斜摄影获取点云坐标的方法及装置 |
CN107194334B (zh) * | 2017-05-10 | 2019-09-10 | 武汉大学 | 基于光流模型的视频卫星影像密集匹配方法及系统 |
CN107292814B (zh) * | 2017-06-12 | 2020-09-08 | 段玲娟 | 一种利用Lidar数据在真立体环境下的数字测图方法及系统 |
CN107506782B (zh) * | 2017-07-06 | 2020-04-17 | 武汉市工程科学技术研究院 | 基于置信权重双边滤波的密集匹配方法 |
US10529086B2 (en) * | 2017-11-22 | 2020-01-07 | Futurewei Technologies, Inc. | Three-dimensional (3D) reconstructions of dynamic scenes using a reconfigurable hybrid imaging system |
CN108399631B (zh) * | 2018-03-01 | 2022-02-11 | 北京中测智绘科技有限公司 | 一种尺度不变性的倾斜影像多视密集匹配方法 |
CN108682029A (zh) * | 2018-03-22 | 2018-10-19 | 深圳飞马机器人科技有限公司 | 多尺度密集匹配方法和系统 |
CN108594255B (zh) * | 2018-04-20 | 2021-09-03 | 武汉大学 | 一种激光测距辅助光学影像联合平差方法及系统 |
CN110060283B (zh) * | 2019-04-17 | 2020-10-30 | 武汉大学 | 一种多测度半全局密集匹配方法 |
CN111898396B (zh) * | 2019-05-06 | 2024-08-09 | 北京四维图新科技股份有限公司 | 障碍物检测方法和装置 |
CN111656404B (zh) * | 2019-05-30 | 2024-03-01 | 深圳市大疆创新科技有限公司 | 图像处理方法、系统及可移动平台 |
CN111046906B (zh) * | 2019-10-31 | 2023-10-31 | 中国资源卫星应用中心 | 一种面状特征点可靠加密匹配方法和系统 |
CN111462195B (zh) * | 2020-04-09 | 2022-06-07 | 武汉大学 | 基于主线约束的非规则角度方向代价聚合路径确定方法 |
CN112116639B (zh) * | 2020-09-08 | 2022-06-07 | 苏州浪潮智能科技有限公司 | 一种图像配准方法、装置及电子设备和存储介质 |
CN112233246B (zh) * | 2020-09-24 | 2024-04-16 | 中山大学 | 基于srtm约束的卫星影像密集匹配方法和系统 |
CN112907739B (zh) * | 2021-01-22 | 2022-10-04 | 中北大学 | 一种井盖高程差信息获取方法、装置及系统 |
CN113034666B (zh) * | 2021-02-01 | 2023-09-12 | 中国计量大学 | 一种基于金字塔视差优化代价计算的立体匹配方法 |
CN112991525B (zh) * | 2021-05-07 | 2021-09-24 | 北京道达天际科技有限公司 | 像方和物方混合匹配基元的数字表面模型生成方法 |
CN113554102A (zh) * | 2021-07-27 | 2021-10-26 | 高小翎 | 代价计算动态规划的航空影像dsm匹配法 |
CN113674415B (zh) * | 2021-08-26 | 2024-01-09 | 自然资源部国土卫星遥感应用中心 | 利用高分七号和资源三号影像联合制作连续、无空洞dsm的方法 |
CN115689965B (zh) * | 2022-12-30 | 2023-03-21 | 武汉大学 | 面向整景卫星影像深度学习密集匹配的多级视差融合方法 |
CN117115336A (zh) * | 2023-07-13 | 2023-11-24 | 中国工程物理研究院计算机应用研究所 | 一种基于遥感立体影像的点云重建方法 |
CN116912251A (zh) * | 2023-09-13 | 2023-10-20 | 深圳市超诺科技有限公司 | 一种提高红外打猎相机检测灵敏度的检测方法 |
Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN102411778A (zh) * | 2011-07-28 | 2012-04-11 | 武汉大学 | 一种机载激光点云与航空影像的自动配准方法 |
CN104156957A (zh) * | 2014-08-06 | 2014-11-19 | 昆山天工智能科技有限公司 | 一种稳定高效的高分辨率立体匹配方法 |
Family Cites Families (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
EP2335220A2 (en) * | 2008-07-06 | 2011-06-22 | Sergei Startchik | Method for distributed and minimum-support point matching in two or more images of 3d scene taken with video or stereo camera. |
-
2015
- 2015-08-20 CN CN201510513972.5A patent/CN105160702B/zh active Active
Patent Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN102411778A (zh) * | 2011-07-28 | 2012-04-11 | 武汉大学 | 一种机载激光点云与航空影像的自动配准方法 |
CN104156957A (zh) * | 2014-08-06 | 2014-11-19 | 昆山天工智能科技有限公司 | 一种稳定高效的高分辨率立体匹配方法 |
Non-Patent Citations (2)
Title |
---|
三角网约束下的层次匹配方法;郑顺义 等;《计算机辅助设计与图形学学报》;20141130;第26卷(第11期);第1989-1996页 * |
城区机载LiDAR数据与航空影像的自动配准;张永军 等;《遥感学报》;20121231;第16卷(第3期);第587-595页 * |
Also Published As
Publication number | Publication date |
---|---|
CN105160702A (zh) | 2015-12-16 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN105160702B (zh) | 基于LiDAR点云辅助的立体影像密集匹配方法及系统 | |
Xia et al. | Geometric primitives in LiDAR point clouds: A review | |
CN108319655B (zh) | 用于生成栅格地图的方法和装置 | |
CN106054900B (zh) | 基于深度摄像头的机器人临时避障方法 | |
CN114708585A (zh) | 一种基于注意力机制的毫米波雷达与视觉融合的三维目标检测方法 | |
CN111476242B (zh) | 一种激光点云语义分割方法及装置 | |
CN106780712B (zh) | 联合激光扫描和影像匹配的三维点云生成方法 | |
CN111060924B (zh) | 一种slam与目标跟踪方法 | |
CN112799096B (zh) | 基于低成本车载二维激光雷达的地图构建方法 | |
CN115564926A (zh) | 基于影像建筑物结构学习的三维面片模型构建方法 | |
CN110889899A (zh) | 一种数字地表模型的生成方法及装置 | |
CN113989758A (zh) | 一种用于自动驾驶的锚引导3d目标检测方法及装置 | |
CN115128628A (zh) | 基于激光slam和单目视觉的道路栅格地图构建方法 | |
Özdemir et al. | A multi-purpose benchmark for photogrammetric urban 3D reconstruction in a controlled environment | |
CN114563000B (zh) | 一种基于改进型激光雷达里程计的室内外slam方法 | |
CN117078870A (zh) | 融合高精地图与激光稀疏点云的道路环境三维重建方法 | |
CN114782357A (zh) | 一种用于变电站场景的自适应分割系统及方法 | |
CN114387488A (zh) | 基于potree点云影像融合的道路提取系统及方法 | |
CN111784798B (zh) | 地图生成方法、装置、电子设备和存储介质 | |
Rau | A line-based 3D roof model reconstruction algorithm: Tin-merging and reshaping (TMR) | |
CN117392237A (zh) | 一种鲁棒的激光雷达-相机自标定方法 | |
US20220276046A1 (en) | System and method for providing improved geocoded reference data to a 3d map representation | |
CN113240755B (zh) | 基于街景图像与车载激光融合的城市场景构图方法及系统 | |
CN107194334B (zh) | 基于光流模型的视频卫星影像密集匹配方法及系统 | |
CN114092388B (zh) | 基于单目相机和里程计的障碍物检测方法 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
C06 | Publication | ||
PB01 | Publication | ||
C10 | Entry into substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |