基于立体影像的地形地物三维变化检测和更新方法技术

本发明专利技术涉及一种基于立体影像的地形地物三维变化检测和更新方法，它包括如下步骤：1、根据输入的立体影像和三维点云，进行影像和点云的配准；2、根据立体影像的像方一致性约束，检测点云中的变化区域，并剔除变化区域中的点云；3、采用立体影像密集匹配算法，重新生成变化区域的三维点云，达到更新点云的目的。本发明专利技术能够以较低的成本，解决三维地形更新的问题，能够满足大范围地形测绘、智慧城市、智能交通等应用。

【技术实现步骤摘要】
基于立体影像的地形地物三维变化检测和更新方法
本专利技术涉及三维变化检测和更新
，具体涉及一种基于立体影像的地形地物三维变化检测和更新方法。
技术介绍
目前，大范围自然地表和人工地物的三维建模技术已经得到了极为广泛的应用，如智慧城市、国防建设、智能交通、虚拟电商、文物保护、增强现实等。现有主流的三维建模技术包括激光扫描技术LiDAR(LightDetectionAndRanging)和影像密集匹配技术。LiDAR系统由全球定位系统(GPS)、惯性导航系统(IMU)以及激光器(LaserScanner)共同组成，以主动方式快速获取地面模型及空间点云信息，具有速度快、时效性强、作业范围大等优点。但是LiDAR系统成本昂贵，不利于三维模型的更新。影像密集匹配通过覆盖一个测区的多张立体影像，根据同名光线对对相交的原理，从影像的二维信息恢复整个测区的三维空间信息。与LiDAR激光点云相比，影像密集匹配技术的重建速度较慢，但是建模成本低廉。随着国内城市化进程的大规模开展，建筑物、道路等人工地物信息发生着日新月异的变化，因此，对智慧城市、智能交通、地形测绘等三维应用提出了实时更新的需求。但是由于LiDAR技术成本昂贵，因此不适合大范围地形的更新。
技术实现思路
本专利技术的目的在于提供一种基于立体影像的地形地物三维变化检测和更新方法，该方法首先采用特征匹配的方式，实现影像和三维点云的配准，其次根据立体影像的像方一致性约束，实现三维变化区域的检测和剔除；最后采用立体影像密集匹配方法，对变化区域进行更新。该方法能够以低廉的成本，实现三维变化检测和更新，能够为大范围地形测绘、智慧城市、智能交通等应用服务。为解决上述技术问题，本专利技术公开的一种基于立体影像的地形地物三维变化检测和更新方法，其特征在于，它包括如下步骤：步骤1：根据输入的立体影像数据和三维点云数据，利用以下方式进行影像和点云的配准；对三维点云数据采用降维处理，将三维点云数据重采样成2.5维的数字表面模型；计算三维点云平面范围内最小外接矩形的四个角点的水平方向坐标X，垂直方向坐标Y，即：Xlb＝min{Xi|i＝1...t}Ylb＝min{Yi|i＝1...t}Xrt＝max{Xi|i＝1...t}Yrt＝max{Yi|i＝1...t}其中，t表示三维点云中三维点的数目；(Xi,Yi)表示第i个三维点的X、Y坐标；(Xlb,Ylb)表示三维点云平面范围内最小外接矩形左下角的角点坐标；(Xrt,Yrt)表示三维点云平面范围内最小外接矩形右上角的角点坐标；(Xlb,Yrt)表示三维点云平面范围内最小外接矩形左上角的角点坐标；(Xrt,Ylb)表示三维点云平面范围内最小外接矩形右下角的角点坐标；min表示取最小值；max表示取最大值；将三维点云的最小外接矩形范围定义为数字表面模型的范围，数字表面模型的起点坐标为(Xlb,Ylb)，则数字表面模型的宽WD和高HD分别为：WD＝(int)(Xrt-Xlb)/sD；HD＝(int)(Yrt-Ylb)/sD其中，sD表示数字表面模型网格的大小，int表示取整操作；WD表示数字表面模型的宽；HD表示数字表面模型的高；Xrt,Xlb,Yrt,Ylb表示最小外接矩形的角点坐标；从而将数字表面模型定义为一个HDxWD大小的规则格网；该格网的起点坐标为(Xlb,Ylb)；在每个数字表面模型的网格中都存在多个三维点，每个网格的高程为对应三维点的最大高程，如下式所示：Z(m,n)＝max(Zi|(int)(Xi-Xlb)/sD＝n；(int)(Yi-Ylb)/sD＝m)其中，Zi表示第i个三维点的高程方向坐标Z；Z(m,n)表示网格(m,n)的高程，取落入网格中所有三维点的最大高程；Xrt,Xlb,Yrt,Ylb表示三维点云平面范围内最小外接矩形的角点坐标；sD表示数字表面模型网格的大小；int表示取整操作；m,n表示数字表面模型中网格的行、列号；max表示取最大值；(int)(Xi-Xlb)/sD＝n；(int)(Yi-Ylb)/sD＝m表示坐标为(Xi,Yi)的三维点落入数字表面模型中第m行第n列的网格中；然后，采用加速鲁棒特征匹配算子来提取数字表面模型与影像之间的特征匹配点，特征匹配点包括数字表面模型上的三维点和对应影像上的像点，其中数字表面模型上的三维点云和影像上的像点一一对应的关系；然后，利用三维点云和对应影像上的像点得到确定三维点云和影像之间投影关系的后方交会模型，根据确定三维点云和影像之间投影关系的后方交会模型，实现三维点云和影像的配准；将拍摄上述影像的相机的位置和姿态表示为三个线元素Xs,Ys,Zs，以及三个角元素phi,omega,kappa；三个线元素和三个角元素统称为影像的外方位元素；三维点云和影像配准的实质是，在三维点云坐标系下，计算影像的外方位元素；将特征匹配的三维点作为控制点，将特征匹配的影像点作为对应的像点，代入确定三维点云和影像之间投影关系的后方交会模型：式中，(vx,vy)表示模型的残差；ΔXs、ΔYs、ΔZs、Δphi、Δomega、Δkappa表示相机外方位元素的改正数；x、y表示真实的特征匹配像点坐标；x0、y0表示将外方位元素初值代入成像模型后，得到的特征匹配像点坐标近似值，其中，成像模型表示为：式中，(x,y)表示真实的特征匹配像点坐标；f表示相机焦距；(X,Y,Z)表示特征匹配控制点坐标；a1、a2、a3、b1、b2、b3、c1、c2、c3表示由角元素初值计算的旋转矩阵元素；表示相机线元素的初值；采用最小二乘平差迭代的方式计算相机外方位元素的改正数，在每一次迭代中，外方位元素的计算公式为：phi＝phi0+Δphiomega＝omega0+Δomegakappa＝kappa0+Δkappa式中，Xs、Ys、Zs表示相机的外方位线元素，phi、omega、kappa表示相机的外方位角元素，即相机在物方空间的朝向；ΔXs、ΔYs、ΔZs、Δphi、Δomega、Δkappa表示相机外方位元素的改正数；phi0、omega0、kappa0表示相机外方位元素的初值；步骤2：根据立体影像的像方一致性约束，检测点云中的变化区域，并剔除变化区域中的点云，具体方法如下:首先，根据所述相机内、外方位元素，将立体影像重采样成核线立体影像；采用成像模型，将所有的三维点投影到核线立体影像上，生成一系列同名像点，成像模型表示为：式中，(xl,yl)、(xr,yr)表示核线立体影像上的同名像点的两个坐标；f表示相机焦距；(X,Y,Z)表示点云的三维坐标；表示左核线影像的相机外方位线元素；表示右核线影像的相机外方位线元素；表示由左核线影像角元素计算的旋转矩阵元素；表示由右核线影像角元素计算的旋转矩阵元素；然后，计算核线立体影像上的同名像点的两个坐标(xl,yl)、(xr,yr)之间影像特征相似性的相关系数：式中，Cov(xl,yl；xr,yr)表示核线立体影像上的同名像点的两个坐标(xl,yl)、(xr,yr)之间影像特征相似性的相关系数；m表示所述相关系数窗口的半径；i,j表示相关系数窗口内的各个像素距离相关系数窗口中心像素的偏移量；g表示左核线影像的灰度值；g'表示右核线影像的灰度值；表示相关系数窗口在左核线影像上本文档来自技高网...

【技术保护点】
一种基于立体影像的地形地物三维变化检测和更新方法，其特征在于，它包括如下步骤：步骤1：根据输入的立体影像数据和三维点云数据，利用以下方式进行影像和点云的配准；对三维点云数据采用降维处理，将三维点云数据重采样成2.5维的数字表面模型；计算三维点云平面范围内最小外接矩形的四个角点的水平方向坐标X，垂直方向坐标Y，即：Xlb＝min{Xi|i＝1...t}Ylb＝min{Yi|i＝1...t}Xrt＝max{Xi|i＝1...t}Yrt＝max{Yi|i＝1...t}其中，t表示三维点云中三维点的数目；(Xi,Yi)表示第i个三维点的X、Y坐标；(Xlb,Ylb)表示三维点云平面范围内最小外接矩形左下角的角点坐标；(Xrt,Yrt)表示三维点云平面范围内最小外接矩形右上角的角点坐标；(Xlb,Yrt)表示三维点云平面范围内最小外接矩形左上角的角点坐标；(Xrt,Ylb)表示三维点云平面范围内最小外接矩形右下角的角点坐标；min表示取最小值；max表示取最大值；将三维点云的最小外接矩形范围定义为数字表面模型的范围，数字表面模型的起点坐标为(Xlb,Ylb)，则数字表面模型的宽WD和高HD分别为：WD＝(int)(Xrt‑Xlb)/sD；HD＝(int)(Yrt‑Ylb)/sD其中，sD表示数字表面模型网格的大小，int表示取整操作；WD表示数字表面模型的宽；HD表示数字表面模型的高；Xrt,Xlb,Yrt,Ylb表示最小外接矩形的角点坐标；从而将数字表面模型定义为一个HDxWD大小的规则格网；该格网的起点坐标为(Xlb,Ylb)；在每个数字表面模型的网格中都存在多个三维点，每个网格的高程为对应三维点的最大高程，如下式所示：Z(m,n)＝max(Zi|(int)(Xi‑Xlb)/sD＝n；(int)(Yi‑Ylb)/sD＝m)其中，Zi表示第i个三维点的高程方向坐标Z；Z(m,n)表示网格(m,n)的高程，取落入网格中所有三维点的最大高程；Xrt,Xlb,Yrt,Ylb表示三维点云平面范围内最小外接矩形的角点坐标；sD表示数字表面模型网格的大小；int表示取整操作；m,n表示数字表面模型中网格的行、列号；max表示取最大值；(int)(Xi‑Xlb)/sD＝n；(int)(Yi‑Ylb)/sD＝m表示坐标为(Xi,Yi)的三维点落入数字表面模型中第m行第n列的网格中；然后，采用加速鲁棒特征匹配算子来提取数字表面模型与影像之间的特征匹配点，特征匹配点包括数字表面模型上的三维点和对应影像上的像点，其中数字表面模型上的三维点云和影像上的像点一一对应的关系；然后，利用三维点云和对应影像上的像点得到确定三维点云和影像之间投影关系的后方交会模型，根据确定三维点云和影像之间投影关系的后方交会模型，实现三维点云和影像的配准；将拍摄上述影像的相机的位置和姿态表示为三个线元素Xs,Ys,Zs，以及三个角元素phi,omega,kappa；三个线元素和三个角元素统称为影像的外方位元素；三维点云和影像配准的实质是，在三维点云坐标系下，计算影像的外方位元素；将特征匹配的三维点作为控制点，将特征匹配的影像点作为对应的像点，代入确定三维点云和影像之间投影关系的后方交会模型：...

【技术特征摘要】
1.一种基于立体影像的地形地物三维变化检测和更新方法，其特征在于，它包括如下步骤：步骤1：根据输入的立体影像数据和三维点云数据，利用以下方式进行影像和点云的配准；对三维点云数据采用降维处理，将三维点云数据重采样成2.5维的数字表面模型；计算三维点云平面范围内最小外接矩形的四个角点的水平方向坐标X，垂直方向坐标Y，即：Xlb＝min{Xi|i＝1...t}Ylb＝min{Yi|i＝1...t}Xrt＝max{Xi|i＝1...t}Yrt＝max{Yi|i＝1...t}其中，t表示三维点云中三维点的数目；(Xi,Yi)表示第i个三维点的X、Y坐标；(Xlb,Ylb)表示三维点云平面范围内最小外接矩形左下角的角点坐标；(Xrt,Yrt)表示三维点云平面范围内最小外接矩形右上角的角点坐标；(Xlb,Yrt)表示三维点云平面范围内最小外接矩形左上角的角点坐标；(Xrt,Ylb)表示三维点云平面范围内最小外接矩形右下角的角点坐标；min表示取最小值；max表示取最大值；将三维点云的最小外接矩形范围定义为数字表面模型的范围，数字表面模型的起点坐标为(Xlb,Ylb)，则数字表面模型的宽WD和高HD分别为：WD＝(int)(Xrt-Xlb)/sD；HD＝(int)(Yrt-Ylb)/sD其中，sD表示数字表面模型网格的大小，int表示取整操作；WD表示数字表面模型的宽；HD表示数字表面模型的高；Xrt,Xlb,Yrt,Ylb表示最小外接矩形的角点坐标；从而将数字表面模型定义为一个HDxWD大小的规则格网；该格网的起点坐标为(Xlb,Ylb)；在每个数字表面模型的网格中都存在多个三维点，每个网格的高程为对应三维点的最大高程，如下式所示：Z(m,n)＝max(Zi|(int)(Xi-Xlb)/sD＝n；(int)(Yi-Ylb)/sD＝m)其中，Zi表示第i个三维点的高程方向坐标Z；Z(m,n)表示网格(m,n)的高程，取落入网格中所有三维点的最大高程；Xrt,Xlb,Yrt,Ylb表示三维点云平面范围内最小外接矩形的角点坐标；sD表示数字表面模型网格的大小；int表示取整操作；m,n表示数字表面模型中网格的行、列号；max表示取最大值；(int)(Xi-Xlb)/sD＝n；(int)(Yi-Ylb)/sD＝m表示坐标为(Xi,Yi)的三维点落入数字表面模型中第m行第n列的网格中；然后，采用加速鲁棒特征匹配算子来提取数字表面模型与影像之间的特征匹配点，特征匹配点包括数字表面模型上的三维点和对应影像上的像点，其中数字表面模型上的三维点云和影像上的像点一一对应的关系；然后，利用三维点云和对应影像上的像点得到确定三维点云和影像之间投影关系的后方交会模型，根据确定三维点云和影像之间投影关系的后方交会模型，实现三维点云和影像的配准；将拍摄上述影像的相机的位置和姿态表示为三个线元素Xs,Ys,Zs，以及三个角元素phi,omega,kappa；三个线元素和三个角元素统称为影像的外方位元素；三维点云和影像配准的实质是，在三维点云坐标系下，计算影像的外方位元素；将特征匹配的三维点作为控制点，将特征匹配的影像点作为对应的像点，代入确定三维点云和影像之间投影关系的后方交会模型：式中，(vx,vy)表示模型的残差；ΔXs、ΔYs、ΔZs、Δphi、Δomega、Δkappa表示相机外方位元素的改正数；x、y表示真实的特征匹配像点坐标；x0、y0表示将外方位元素初值代入成像模型后，得到的特征匹配像点坐标近似值，其中，成像模型表示为：式中，(x,y)表示真实的特征匹配像点坐标；f表示相机焦距；(X,Y,Z)表示特征匹配控制点坐标；a1、a2、a3、b1、b2、b3、c1、c2、c3表示由角元素初值计算的旋转矩阵元素；表示相机线元素的初值；采用最小二乘平差迭代的方式计算相机外方位元素的改正数，在每一次迭代中，外方位元素的计算公式为：phi＝phi0+Δphiomega＝omega0+Δomegakappa＝kappa0+Δkappa式中，Xs、Ys、Zs表示相机的外方位线元素，phi、omega、kappa表示相机的外方位角元素，即相机在物方空间的朝向；ΔXs、ΔYs、...

【专利技术属性】
技术研发人员：周刚，黄旭，胡堃，蔡刚山，
申请(专利权)人：武汉市工程科学技术研究院，
类型：发明
国别省市：湖北,42