[发明专利]基于立体影像的地形地物三维变化检测和更新方法

摘要

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

主权项

一种基于立体影像的地形地物三维变化检测和更新方法，其特征在于，它包括如下步骤：步骤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=∂x∂XsΔXs+∂x∂YsΔYs+∂x∂ZsΔZs+∂x∂phiΔphi+∂x∂omegaΔomega+∂x∂kappaΔkappa-(x-x0)]]>vy=∂y∂XsΔXs+∂y∂YsΔYs+∂y∂ZsΔZs+∂y∂phiΔphi+∂y∂omegaΔomega+∂y∂kappaΔkappa-(y-y0)]]>式中，(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表示由角元素初值计算的旋转矩阵元素；表示相机线元素的初值；采用最小二乘平差迭代的方式计算相机外方位元素的改正数，在每一次迭代中，外方位元素的计算公式为：Xs=XS0+ΔXsYs=YS0+ΔYsZs=ZS0+ΔZs]]>phi＝phi0+Δphi omega＝omega0+Δomega kappa＝kappa0+Δkappa式中，Xs、Ys、Zs表示相机的外方位线元素，phi、omega、kappa表示相机的外方位角元素，即相机在物方空间的朝向；ΔXs、ΔYs、ΔZs、Δphi、Δomega、Δkappa表示相机外方位元素的改正数；phi0、omega0、kappa0表示相机外方位元素的初值；步骤2：根据立体影像的像方一致性约束，检测点云中的变化区域，并剔除变化区域中的点云，具体方法如下:首先，根据所述相机内、外方位元素，将立体影像重采样成核线立体影像；采用成像模型，将所有的三维点投影到核线立体影像上，生成一系列同名像点，成像模型表示为：xl=-fa1l(X-XSl)+b1l(Y-YSl)+c1l(Z-ZSl)a3l(X-XSl)+b3l(Y-YSl)+c3l(Z-ZSl)yl=-fa2l(X-XSl)+b2l(Y-YSl)+c2l(Z-ZSl)a3l(X-XSl)+b3l(Y-YSl)+c3l(Z-ZSl)]]>xr=-fa1r(X-XSr)+b1r(Y-YSr)+c1r(Z-ZSr)a3r(X-XSr)+b3l(Y-YSr)+c3l(Z-ZSr)yr=-fa2r(X-XSr)+b2r(Y-YSr)+c2r(Z-ZSr)a3r(X-XSr)+b3r(Y-YSr)+c3r(Z-ZSr)]]>式中，(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'表示右核线影像的灰度值；表示相关系数窗口在左核线影像上的平均灰度值；表示相关系数窗口在右核线影像上的平均灰度值；将核线立体影像上的同名像点的两个坐标(xl,yl)、(xr,yr)之间影像特征相似性的相关系数Cov(xl,yl；xr,yr)与预设的相关系数阈值进行比较，如果核线立体影像上的同名像点的两个坐标(xl,yl)、(xr,yr)之间影像特征相似性的相关系数Cov(xl,yl；xr,yr)小于上述相关系数阈值，则认为该同名像点以及对应的三维点位于变化区域，删该同名的两个像点，以及同名的两个像点对应的三维点；步骤3：采用立体影像密集匹配算法，重新生成变化区域的三维点云，达到更新点云的目的，具体方法为：在上述变化区域中采用半全局密集匹配算法，生成变化区域的核线影像同名点；然后根据前方交会法，利用变化区域中的核线影像同名点生成变化区域的三维点云，实现更新点云。