CN112767461A - 激光点云与序列全景影像自动配准方法 - Google Patents

激光点云与序列全景影像自动配准方法 Download PDF

Info

Publication number
CN112767461A
CN112767461A CN202011625760.3A CN202011625760A CN112767461A CN 112767461 A CN112767461 A CN 112767461A CN 202011625760 A CN202011625760 A CN 202011625760A CN 112767461 A CN112767461 A CN 112767461A
Authority
CN
China
Prior art keywords
registration
image
point cloud
graph
mvs
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.)
Pending
Application number
CN202011625760.3A
Other languages
English (en)
Inventor
陈驰
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Wuhan University WHU
Original Assignee
Wuhan University WHU
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Wuhan University WHU filed Critical Wuhan University WHU
Priority to CN202011625760.3A priority Critical patent/CN112767461A/zh
Publication of CN112767461A publication Critical patent/CN112767461A/zh
Pending legal-status Critical Current

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/30Determination of transform parameters for the alignment of images, i.e. image registration
    • G06T7/33Determination of transform parameters for the alignment of images, i.e. image registration using feature-based methods
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00Pattern recognition
    • G06F18/20Analysing
    • G06F18/25Fusion techniques
    • G06F18/251Fusion techniques of input or preprocessed data
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/10Segmentation; Edge detection
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/80Analysis of captured images to determine intrinsic or extrinsic camera parameters, i.e. camera calibration
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10004Still image; Photographic image
    • G06T2207/10012Stereo images
    • YGENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02TCLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
    • Y02T10/00Road transport of goods or passengers
    • Y02T10/10Internal combustion engine [ICE] based vehicles
    • Y02T10/40Engine management systems

Landscapes

  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • General Physics & Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • Artificial Intelligence (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Evolutionary Biology (AREA)
  • Evolutionary Computation (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • General Engineering & Computer Science (AREA)
  • Image Analysis (AREA)
  • Image Processing (AREA)

Abstract

针对无人机/车载移动测量激光点云与序列全景影像数据质量、精度、空间数据基准不统一、多源成像数据之间存在配准误差的问题,本发明提出一种自动配准方法——基于SFGM(Space Fingerprint Graph Matching,空间指纹图匹配)实现激光点云与序列全景影像的自动配准方法,在配准基元提取的基础上,以空间指纹图方法实现同名匹配,继而由粗到精两步实现稳健配准,突破点云与影像自动数据配准的技术瓶颈。本发明成果为多源空间数据联合分析与处理奠定坚实数据基础,服务于基础测绘、智慧城市、导航与位置服务等领域。

Description

激光点云与序列全景影像自动配准方法
技术领域
本发明属于无人机/车载移动测量数据融合中自动数据配准的应用,提出一种全新的激光点云与序列全景影像自动配准方法——SFGM(Space Fingerprint GraphMatching,空间指纹图匹配)配准方法。
技术背景
集成LiDAR(Light Detection and Ranging,激光测距)、全景相机等多传感器的无人机/车载MMS(Mobile Mapping System,移动测量***)为遥感对地观测提供了一种新的技术手段,在观测方向、数据分辨率等方面扩展了传统摄影测量与遥感数据,使得低空顶面与全景、街景侧面与全景的多角度、高分辨率、多源数据获取成为可能,为城市场景建模、目标提取与识别提供了多样化的数据源。
激光点云数据与影像融合处理是当前摄影测量与遥感领域的一个研究热点。无人机/车载MMS的多源数据配准问题是其数据融合的技术瓶颈,很大程度上制约了融合多源数据的分类、识别、建模方法在移动测量平台上的应用与发展,是迫切需要解决的重要研究问题。无人机/车载MMS上的LiDAR点云数据与序列全景影像的融合是实现三维激光点云与影像数据两者优势互补的必然。相较于成熟的星载/大航空摄影测量与遥感顶面数据配准方法研究,目前,鲜见无人机/车载MMS数据配准相关研究。因此,迫切需要一种激光点云与序列全景影像数据鲁棒配准方法,以达到数据空间基准统一、几何、光谱信息融合等,从而进一步实现多源、多维、多视角、高分辨率无人机/车载MMS数据的联合处理。
国内外学者针对有人机载、星载、地面固定站激光扫描数据与影像的配准方法开展了大量的研究,但在无人机/车载MMS序列全景影像与激光点云的自动配准方面,国内外鲜有研究成果报告。不局限于数据类型与采集平台,现有的配准方法主要分为四类,即:①多传感器标定的方法;②几何特征匹配的2D-3D配准算法;③互信息最大化的2D-3D配准算法;④影像密集点云与激光点云配准的3D-3D配准算法。
多传感器标定的方法主要用于消除同一***中不同遥感传感器之间的安置误差带来的数据配准误差,是传统摄影测量与遥感传感器平台上常用的数据配准方法,主要通过实验室标定与室外检校场标定相结合的方法解算多传感器之间的安置位置参数。无地面检校场多传感器安置误差半自动检校方法近年来亦得到广泛研究,但由于无人机/车载MMS误差来源多样,尚未见无地面检校场的多传感器标定方法在其多传感器检校中应用。无人机MMS多搭载消费级相机、手机等非专业传感器,其控制精度远逊于专业级航测仪。由于非专业相机未知曝光延迟时间、***标定误差、多传感器之间的时间同步误差难以完全消除,加之轻小型POS数据精度不足等原因,通过多传感器同步以及标定的方法难以获得高精度数据配准结果。
几何特征匹配的2D-3D配准方法是Palenichka配准方法分类法中基于特征的配准算法中的代表算法。该类方法的基本原理是利用激光点云与二维光学影像中共轭的几何特征求解影像的外方位元素实现两者间的配准。点、线、面、体等多种几何基元均可以作为配准基元,形成共轭配准基元对。
Figure BDA0002879215600000021
and Haggrén(2012)研究成果表明,使用人工构筑物(Artificial Features)进行数据配准融合精度高于使用自然地物(例如:树冠)。因此,下文所论述的相关研究工作均以人工构筑物作为配准基元。
因激光点云自身三维离散采样的特性,在LiDAR点云数据中不存在准确的特征点(例如:节点、转角)。但是,使用特征点邻域内的点对该特征点进行描述,从而构建点基元的匹配对已经被研究证明为一种有效的LiDAR点云与影像配准方法。使用影像局部特征描述子如SIFT(Scale-Invariant Feature Transform,尺度不变特征变换)对LiDAR强度影像与可见光影像中的同名点进行匹配,从而实现激光点云数据与影像的融合。LiDAR点云数据中由边缘点拟合的线特征亦被用于激光点云数据与影像数据的配准中,利用激光与影像中的同名线特征,使用二者共面条件作为约束,平差解算获得精确的影像外方位元素。
特征匹配的2D-3D配准算法,如基于线对的半自动配准方法,目前已成功应用于传统遥感平台LiDAR点云与影像的配准。但该类方法,需要较为准确的影像外方位元素初值与正确的匹配对保证平差的正确收敛,且特征选取难以自动化。同时,由于激光点云本质上的三维离散采样特性,严格意义上讲,不存在扫描完整、定位准确的几何基元。故任何直接自激光点云数据进行的几何基元自动提取与拟合的算法都存在几何基元提取完整性与角点定位精度的问题。在配准过程中,自动拟合线或是角点等几何基元的不确定性,如激光线特征提取完整度与线段端点精度,均会对数据配准精度产生影响。
影像密集点云与激光点云的3D-3D配准方法属于Palenichka融合方法分类法中基于特征的配准融合算法。该类方法采用密集匹配方法自序列影像生成影像密集三维点云,通过POS***输出值计算影像与激光点云的空间转换初始值,继而使用ICP(IterativeClosest Point,迭代最邻近点)算法配准影像密集匹配点云与LiDAR点云,从而间接解决影像与点云的配准问题。其原理是将二维影像与三维点云的配准转化为两组三维点云之间的配准,使用3D-3D点基元配准方法解决激光点云与影像数据的配准问题,避免了在点云数据、影像数据中提取几何特征,保证算法的稳健性。但其核心3D-3D的ICP算法对初始配准转换参数的近似值要求较高。由于无人机/车载MMS***自身POS精度限制,或是GPS信号存在遮挡、***同步标定误差等因素引起的直接地理定向精度较差条件下,初始配准参数不能满足ICP方法对较为精确初始配准位置的基本要求,进而导致该方法不能正确收敛,配准鲁棒性不足。
发明内容
综合前人研究并且针对目前无人机/车载MMS点云与序列全景影像配准的不足,提高两者配准的健壮性和配准精度。本发明提出SFGM(Space Fingerprint Graph Matching,空间指纹图匹配)配准方法:研究激光点云与序列全景影像配准基元稳健提取方法,生成异源几何配准基元集合;研究基于图的激光点云与序列全景影像共轭配准基元对模糊匹配方法,解决几何定位不精确、非单射异源配准基元集合中同名基元查找的难点;继而研究两步法高精度配准模型参数稳健估计方法,自共轭基元对解算配准模型初值,削弱模型对高精度POS数据的依赖,并通过异维配准到同维配准问题转换实现对粗配准参数的精化,实现无人机/车载MMS激光点云与序列全景影像的高精度稳健配准。
为了解决上述技术问题,本发明采用如下的技术方案:SFGM配准方法包括以下步骤:
步骤1,序列全景影像预处理。利用全景图像结构特征的自动恢复方法,将从无人机MMS和车载MMS分别得到的序列全景影像,生成几何结构正确的虚拟框幅式影像;
步骤2,在步骤1得到的虚拟框幅式影像上提取显著的人工建筑物作为配准基元。本发明在无人机MMS数据集中提取建筑物的顶面作为配准基元,在车载MMS数据集中提取建筑物天际线作为配准基元;
步骤3,上述步骤2.1会得到无人机MMS数据集点云配准基元和虚拟影像配准基元,步骤2.2会得到车载MMS数据集点云配准基元和虚拟影像配准基元。本发明将激光点云与影像共轭配准基元查找的问题转化为异源图匹配问题,提出构建SFG空间指纹图描述提取到的几何配准基元集合的空间分布特征,继而自图结构全局空间相似性实现匹配;
步骤4,利用SfM方法(Structure from Motion,运动结构恢复)恢复标定相机采集的序列影像在摄影测量坐标系中的外方位元素,并且通过多视立体匹配方法得到由序列影像生成的MVS影像密集点云,基于步骤3中SFG空间指纹图构建方法获取MVS影像密集点云与LiDAR点云中的同名角点匹配集合;
步骤5,基于同名角点匹配集合进行粗配准,获得摄影测量坐标系与LiDAR参考坐标系之间的初始转换关系;
步骤6,通过最小化MVS点云与LiDAR点云之间的距离实现配准参数精化,精配准过程分为两个步骤,即粗配准模型参数到精配准模型初值的转换以及迭代最邻近点优化解算。
进一步的,步骤2的具体实现包括如下子步骤,
步骤2.1,无人机MMS数据集建筑物顶面配准基元提取。首先采用多标记点过程的方法(Yang et al.,2013)从无人机LiDAR点云数据中提取建筑物点云,并且采用RMBR(Recursive Minimum Bounding Rectangle)算法(Kwak,2013)进行建筑物外边界多边形提取与规则化,得到建筑物外框。并依据此点云建筑物提取结果作为先验知识,提出一种激光点云先验知识引导的影像建筑物顶面提取算法,使用点云数据建筑物提取结果对虚拟框幅式影像上的配准基元提取进行引导,使得初始反投影区域沿张量梯度方向不断优化,最后扩展基于全局对比度显著性检测的分割方法(RCC,Region Contrast Cut)到影像局部区域,结合规则化算法,实现影像建筑物配准几何基元的生成;
步骤2.2,车载MMS数据集天际线配准基元提取。采用一种改进GrabCut(Rother etal.,2004)的天际线分割方法,从虚拟影像中提取建筑物上边线与角点,再结合应用层次化城市场景目标提取方法(Yang et al.,2015)与RMBR算法(Kwak,2013)边界规则化算法自车载LiDAR点云中实现建筑物立面的提取,确定其上边线,完成影像天际线配准基元提取。
进一步的,步骤3的具体实现包括如下子步骤,
步骤3.1,构建SFG空间指纹图。本发明提出一种以三角形为生成核,节点到核三角形距离最短的图集合构建算法,构建配准基元空间指纹图。SFG空间指纹图集合(GC)生成分为两个步骤:核三角形生成与非核三角形边的图边连接。
对于无人机MMS数据,选择提取到的n个建筑物顶面配准基元重心组成图节点集合E1。对于车载MMS数据,选择提取到的n个天际线配准基元角点组成图节点集E2。下面以E来代指上述两个图节点集E1和E2,来说明核三角形的生成,选取集合E中的任意3个顶点组合为三角形,作为图G的生成核,则集合E可能构成的图集合GC可表示为:
Figure BDA0002879215600000051
其中Root(E)为核三角形,N为顶点集合E的顶点个数,V(Root(E))为以Root(E)核三角形顶点为自变量的图边集合,V(·)为单调边集合构成函数。图的边连接规则,即边集合构成函数定义为:对于核三角形顶点,直接连接构成完全图;对于非核三角形顶点,则计算该顶点到核三角形三个顶点的距离,将最短长度的边作为连接边。核三角形的边与非核三角形的边共同组成了图的边集合,与图顶点共同构成配准基元空间指纹图。
步骤3.2,匹配异源SFG空间指纹图。本发明采用一种最小GED图编辑距离空间指纹图匹配算法。
Eimage,Vimage分别是影像配准基元SFG空间指纹图的顶点集合与边集合。同理Elas,Vlas为激光点云配准基元SFG空间指纹图的顶点集合与边集合。使用GED度量影像配准基元SFG空间指纹图Gimage=(Eimage,Vimage)与激光点云SFG空间指纹图Glas=(Elas,Vlas),在局部核三角形匹配对定义的图转换矩阵T下的全局图边相似度。T是将Gimage匹配到Glas的变换矩阵,在本算法中指经过核三角形匹配后获得的初始图匹配位置。即将最优图匹配的问题转化为最小图编辑距离查找问题,其计算公式如下:
Figure BDA0002879215600000061
其中RootMatch(rm1,..rmk)表示核三角形中的匹配过程中产生的核三角形局部匹配集。k为将Gimage经过修改之后完全匹配Glas所需要的操作步骤(添加、删除、替换)数量。Glas与Gimage之间的GED定义为:
Figure BDA0002879215600000062
其中cost(op,T)为将Gimage经过修改之后完全匹配Glas所需要的k步操作(添加、删除、替换)所对应的图编辑代价函数,本发明设置该代价函数为旋转图边差,length(·)为图边长度累加函数。GED最小匹配过程分为核三角形局部匹配与图编辑距离全局匹配两个步骤:首先,核三角形局部匹配过程采用三角形相似性即内角值邻近性准则进行匹配,使用KD树(Zhou et al.,2008)对所有的核三角形内角值三维点坐标数据构建索引,提高匹配效率;其次,遍历核三角形匹配集计算各个核三角形匹配对应图转换的GED,并进行排序,将最小GED对应的T作为最佳匹配。
进一步的,步骤5的具体实现包括如下子步骤,
步骤5.1,使用同名角点匹配集合解算2D-3D配准模型,即公式(4)。令序列影像的关键帧中一个建筑物外多边形角点的像空间坐标为m=(u,v,f)T,其对应的建筑物角点在LiDAR数据中的坐标为Mlas=(X,Y,Z)T,则两角点之间的共线关系可表述为:
spnpm=A[Rpnp|tpnp]Mlas (4)
其中,A为已知的相机内参数矩阵,spnp为比例参数,[Rpnp|tpnp]为EPnP计算获得的相机外参数矩阵。式(4)描述了一个物点与像点的共线问题,本文使用EPnP算法(Lepetitet al.,2009)对其进行线性求解。
步骤5.2,根据MVS影像密集点云与LiDAR点云中的同名角点匹配集合,分别得到其在LiDAR参考坐标系和MVS影像密集点云所在的摄影测量坐标系中的点坐标,并解算三维空间相似变换;令在摄影测量坐标系(Cmvs)中的建筑物规则化外边界多边形角点为Mmvs,其对应角点在LiDAR参考坐标系(Cw)中的坐标为Mlas,则其之间的转换关系定义为Mlas=λRMmvs+T,其中λ,R,T分别为尺度、旋转与平移参数;对于n个同名角点对中的每3个同名角点对使用抗差SVD方法即可解算出一组(λ,R,T)。
进一步的,步骤6的具体实现包括如下子步骤,
步骤6.1,精配准模型初值转换。
最小化MVS点云与LiDAR点云之间的距离实现配准参数精化的配准几何模型可表示为:
Mlas=λRMmvs+T (5)
在本步骤中需要对2D-3D模型中得到的影像外参数矩阵转换到式(5)中的R变量,实现2D-3D的转换。
2D-3D模型初值转换:令一像点为m,其在Cmvs中对应点为MMVS,在Cw中对应点为Mlas,则共线性关系可表达为:
sbundlem=A[Rbundle|tbundle]Mmvs (6)
其中sbundle为比例参数,[Rbundle|tbundle]为SfM计算获得的影像外参数矩阵;将2D-3D模型中得到的影像外参数矩阵转换到式(5)中的R变量,实现2D-3D的转换,联立式(4)至(6)可得:
Figure BDA0002879215600000071
其中R为旋转参数,
Figure BDA0002879215600000072
表示Rpnp的逆,Rpnp表示EPnP计算获得的相机外参数矩阵,Rbundle表示SfM计算的影像外参数,Ppnp(Xi,Yi,Zi)t,令Pbundle(xi,yi,zi)t为相机在Cw与Cmvs中的坐标重心化后的相机位姿,则其尺度参数可通过式(8)计算,平移参数通过式(9)计算:
Figure BDA0002879215600000073
T=Ppnp-λRPbundle (9)
由公式(7)至(9)可知,一张粗配准影像对应一种Cw与Cmvs之间的空间旋转关系,多张配准影像则可实现对全部空间转换参数(λ,R,T)的求解;
为实现(λ,R,T)的稳健,采用投票聚类的方法实现参数的稳健估计(Fernandesand Oliveira,2008)。
步骤6.2,迭代最近邻点模型参数优化。
无人机LiDAR点云数据范围覆盖小,作业时间短,在进行激光点云航带平差后,其内部的非线性变形幅度可忽略不计,可以认为其不存在非刚体形变。故采用刚体空间相似变换模型作为MVS与LiDAR点云之间的配准模型。
车载MMS数据通常覆盖范围较大,作业时间较长。由于长期工作状态下IMU的漂移以及GPS信号不良等因素,车载MMS采集的LiDAR点云数据通常内部存在非刚性的变形。在POS数据存在不准确的条件下,需要对其进行非刚性变形纠正方能保证其数据绝对定位精度。常用的非刚性变形纠正模型包括:分段刚性模型(Pylvanainen et al.,2010)、时间漂移模型(Hofmann et al.,2014)等。Dang2研究数据覆盖范围较小,作业时间较短,受研究数据的限制,其内部的非线性变形可忽略不计,可认为数据LiDAR数据集内不存在非刚体形变,故不采用非刚性模型对LiDAR点云数据非刚性变形进行拟合。采用刚体空间相似变换模型作为MVS与LiDAR点云之间的配准模型。
故两种数据集都使用点到面距离差(Chen and Medioni,1991)作为点对误差,则误差方程可表示为:
Figure BDA0002879215600000081
其中,mi,Mi分别为MVS与LiDAR点,wi为点对权,ηi为Mi点法向量,(λ3d-3d,R3d-3d,T3d-3d)表示需要求解的最小MVS与LiDAR点云点到面距离的空间相似变换。基于最优摄影测量与LiDAR参考系之间的坐标转换(λ3d-3d,R3d-3d,T3d-3d),修正后的影像外方位元素(Rcam,Tcam)可表示为:
Figure BDA0002879215600000082
与现有技术相比,本发明具有以下优点和有益效果:
1、本发明自动提取并匹配场景内显著地物几何结构,突破初始配准误差较大、维度、尺度、精度不一致无人机/车载MMS激光点云与序列全景影像数据配准的技术瓶颈。
2、本发明提出两步法配准模型,按由粗配准到精配准的思路实现高精度稳健数据配准。
附图说明
图1为本发明的无人机MMS点云与序列全景影像配准技术路线示意图。
图2为本发明的车载MMS点云与序列全景影像配准技术路线示意图。
图3为本发明全景影像虚拟成像过程示意图。
图4为本发明配准基元空间指纹图的构建。
图5为最小图编辑距离准则全局最优共轭基元匹配。
具体实施方式
以下结合附图和实施例对本发明技术方案进行说明。
选择武汉大学测绘遥感信息工程国家重点实验室研制的低空无人机多传感器采集***Heli-Mapper得到的点云数据和序列全景影像,和多种车载移动测量***得到的点云数据和序列全景影像,分别在无人机、车载移动测量平台上对本发明提出的SFGM配准模型进行具体说明。参见图1和图2,本发明实施例包括以下步骤:
步骤1,序列全景影像预处理。利用全景图像结构特征的自动恢复方法,将从无人机MMS和车载MMS分别得到的序列全景影像,生成几何结构正确的虚拟框幅式影像。
本步骤的具体操作参照图3。
步骤2,在步骤1得到的虚拟框幅式影像上提取显著的人工建筑物作为配准基元。本发明在无人机MMS数据集中提取建筑物的顶面作为配准基元,在车载MMS数据集中提取建筑物天际线作为配准基元。
本步骤具体操作如下:
步骤2.1,无人机MMS数据集建筑物顶面配准基元提取。首先采用多标记点过程的方法(Yang et al.,2013)从无人机LiDAR点云数据中提取建筑物点云,并且采用RMBR(Recursive Minimum Bounding Rectangle)算法(Kwak,2013)进行建筑物外边界多边形提取与规则化,得到建筑物外框。并依据此点云建筑物提取结果作为先验知识,提出一种激光点云先验知识引导的影像建筑物顶面提取算法,使用点云数据建筑物提取结果对虚拟框幅式影像上的配准基元提取进行引导,使得初始反投影区域沿张量梯度方向不断优化,最后扩展基于全局对比度显著性检测的分割方法(RCC,Region Contrast Cut)到影像局部区域,结合规则化算法,实现影像建筑物配准几何基元的生成;
步骤2.2,车载MMS数据集天际线配准基元提取。采用一种改进GrabCut(Rother etal.,2004)的天际线分割方法,从虚拟影像中提取建筑物上边线与角点,再结合应用层次化城市场景目标提取方法(Yang et al.,2015)与RMBR算法(Kwak,2013)边界规则化算法自车载LiDAR点云中实现建筑物立面的提取,确定其上边线,完成影像天际线配准基元提取。
步骤3,上述步骤2.1会得到无人机MMS数据集点云配准基元和虚拟影像配准基元,步骤2.2会得到车载MMS数据集点云配准基元和虚拟影像配准基元。本发明将激光点云与影像共轭配准基元查找的问题转化为异源图匹配问题,提出构建SFG空间指纹图描述提取到的几何配准基元集合的空间分布特征,继而自图结构全局空间相似性实现匹配。
本步骤具体操作如下:
步骤3.1,构建SFG空间指纹图。参照图4,本发明提出一种以三角形为生成核,节点到核三角形距离最短的图集合构建算法,构建配准基元空间指纹图。SFG空间指纹图集合(GC)生成分为两个步骤:核三角形生成与非核三角形边的图边连接。
对于无人机MMS数据,选择提取到的n个建筑物顶面配准基元重心组成图节点集合E1。对于车载MMS数据,选择提取到的n个天际线配准基元角点组成图节点集E2。下面以E来代指上述两个图节点集E1和E2,来说明核三角形的生成,选取集合E中的任意3个顶点组合为三角形,作为图G的生成核,则集合E可能构成的图集合GC可表示为:
Figure BDA0002879215600000101
其中Root(E)为核三角形,N为顶点集合E的顶点个数,V(Root(E))为以Root(E)核三角形顶点为自变量的图边集合,V(·)为单调边集合构成函数。图的边连接规则,即边集合构成函数定义为:对于核三角形顶点,直接连接构成完全图;对于非核三角形顶点,则计算该顶点到核三角形三个顶点的距离,将最短长度的边作为连接边。参照图4,核三角形的边与非核三角形的边共同组成了图的边集合,与图顶点共同构成配准基元空间指纹图。
步骤3.2,匹配异源SFG空间指纹图。本发明采用一种最小GED图编辑距离空间指纹图匹配算法,见图5。
Eimage,Vimage分别是影像配准基元SFG空间指纹图的顶点集合与边集合。同理Elas,Vlas为激光点云配准基元SFG空间指纹图的顶点集合与边集合。使用GED度量影像配准基元SFG空间指纹图Gimage=(Eimage,Vimage)与激光点云SFG空间指纹图Glas=(Elas,Vlas),在局部核三角形匹配对定义的图转换矩阵T下的全局图边相似度。T是将Gimage匹配到Glas的变换矩阵,在本算法中指经过核三角形匹配后获得的初始图匹配位置。即将最优图匹配的问题转化为最小图编辑距离查找问题,其计算公式如下:
Figure BDA0002879215600000111
其中RootMatch(rm1,..rmk)表示核三角形中的匹配过程中产生的核三角形局部匹配集。k为将Gimage经过修改之后完全匹配Glas所需要的操作步骤(添加、删除、替换)数量。Glas与Gimage之间的GED定义为:
Figure BDA0002879215600000112
其中cost(op,T)为将Gimage经过修改之后完全匹配Glas所需要的k步操作(添加、删除、替换)所对应的图编辑代价函数,本发明设置该代价函数为旋转图边差,length(·)为图边长度累加函数。GED最小匹配过程分为核三角形局部匹配与图编辑距离全局匹配两个步骤:首先,核三角形局部匹配过程采用三角形相似性即内角值邻近性准则进行匹配,使用KD树(Zhou et al.,2008)对所有的核三角形内角值三维点坐标数据构建索引,提高匹配效率;其次,遍历核三角形匹配集计算各个核三角形匹配对应图转换的GED,并进行排序,将最小GED对应的T作为最佳匹配。
步骤4,利用SfM方法(Structure from Motion,运动结构恢复)恢复标定相机采集的序列影像在摄影测量坐标系中的外方位元素,并且通过多视立体匹配方法得到由序列影像生成MVS影像密集点云,基于步骤3中SFG空间指纹图构建方法获取MVS影像密集点云与LiDAR点云中的同名角点匹配集合。
步骤5,基于同名角点匹配集合进行粗配准,获得摄影测量坐标系与LiDAR参考坐标系之间的初始转换关系,具体包含以下步骤:
步骤5.1,使用同名角点匹配集合解算2D-3D配准模型,即公式(4)。令序列影像的关键帧中一个建筑物外多边形角点的像空间坐标为m=(u,v,f)T,其对应的建筑物角点在LiDAR数据中的坐标为Mlas=(X,Y,Z)T,则两角点之间的共线关系可表述为:
spnpm=A[Rpnp|tpnp]Mlas (4)
其中,A为已知的相机内参数矩阵,spnp为比例参数,[Rpnp|tpnp]为EPnP计算获得的相机外参数矩阵。式(4)描述了一个物点与像点的共线问题,本文使用EPnP算法(Lepetitet al.,2009)对其进行线性求解。
步骤5.2,根据MVS影像密集点云与LiDAR点云中的同名角点匹配集合,分别得到其在LiDAR参考坐标系和MVS影像密集点云所在的摄影测量坐标系中的点坐标,并解算三维空间相似变换;令在摄影测量坐标系(Cmvs)中的建筑物规则化外边界多边形角点为Mmvs,其对应角点在LiDAR参考坐标系(Cw)中的坐标为Mlas,则其之间的转换关系定义为Mlas=λRMmvs+T,其中λ,R,T分别为尺度、旋转与平移参数;对于n个同名角点对中的每3个同名角点对使用抗差SVD方法即可解算出一组(λ,R,T)。
步骤6,通过最小化MVS点云与LiDAR点云之间的距离实现配准参数精化,精配准过程分为两个步骤,即粗配准模型参数到精配准模型初值的转换以及迭代最邻近点优化解算。具体包括以下步骤:
步骤6.1,精配准模型初值转换。
最小化MVS点云与LiDAR点云之间的距离实现配准参数精化的配准几何模型可表示为:
Mlas=λRMmvs+T (5)
在本步骤中需要对2D-3D模型中得到的影像外参数矩阵转换到式(5)中的R变量,实现2D-3D的转换。
2D-3D模型初值转换:令一像点为m,其在Cmvs中对应点为MMVS,在Cw中对应点为Mlas,则共线性关系可表达为:
sbundlem=A[Rbundle|tbundle]Mmvs (6)
其中sbundle为比例参数,[Rbundle|tbundle]为SfM计算获得的影像外参数矩阵;将2D-3D模型中得到的影像外参数矩阵转换到式(5)中的R变量,实现2D-3D的转换,联立式(4)至(6)可得:
Figure BDA0002879215600000131
其中R为旋转参数,
Figure BDA0002879215600000132
表示Rpnp的逆,Rpnp表示EPnP计算获得的相机外参数矩阵,Rbundle表示SfM计算的影像外参数,Ppnp(Xi,Yi,Zi)t,令Pbundle(xi,yi,zi)t为相机在Cw与Cmvs中的坐标重心化后的相机位姿,则其尺度参数可通过式(8)计算,平移参数通过式(9)计算:
Figure BDA0002879215600000133
T=Ppnp-λRPbundle (9)
由公式(7)至(9)可知,一张粗配准影像对应一种Cw与Cmvs之间的空间旋转关系,多张配准影像则可实现对全部空间转换参数(λ,R,T)的求解;
为实现(λ,R,T)的稳健,采用投票聚类的方法实现参数的稳健估计(Fernandesand Oliveira,2008)。
步骤6.2,迭代最近邻点模型参数优化。由于本发明中无人机LiDAR点云数据和车载LiDAR点云数据存在差异,所以下文进行分开说明。
无人机LiDAR点云数据范围覆盖小,作业时间短,在进行激光点云航带平差后,其内部的非线性变形幅度可忽略不计,可以认为其不存在非刚体形变。故采用刚体空间相似变换模型作为MVS与LiDAR点云之间的配准模型。
车载MMS数据通常覆盖范围较大,作业时间较长。由于长期工作状态下IMU的漂移以及GPS信号不良等因素,车载MMS采集的LiDAR点云数据通常内部存在非刚性的变形。在POS数据存在不准确的条件下,需要对其进行非刚性变形纠正方能保证其数据绝对定位精度。常用的非刚性变形纠正模型包括:分段刚性模型(Pylvanainen et al.,2010)、时间漂移模型(Hofmann et al.,2014)等。本发明的研究数据覆盖范围较小,作业时间较短,受研究数据的限制,其内部的非线性变形可忽略不计,可认为本发明数据LiDAR数据集内不存在非刚体形变,故不采用非刚性模型对LiDAR点云数据非刚性变形进行拟合。采用刚体空间相似变换模型作为MVS与LiDAR点云之间的配准模型。
故两种数据集都使用点到面距离差(Chen and Medioni,1991)作为点对误差,则误差方程可表示为:
Figure BDA0002879215600000141
其中,mi,Mi分别为MVS与LiDAR点,wi为点对权,ηi为Mi点法向量,(λ3d-3d,R3d-3d,T3d-3d)表示需要求解的最小MVS与LiDAR点云点到面距离的空间相似变换。基于最优摄影测量与LiDAR参考系之间的坐标转换(λ3d-3d,R3d-3d,T3d-3d),修正后的影像外方位元素(Rcam,Tcam)可表示为:
Figure BDA0002879215600000142
下面将结合具体实例应用进一步说明本发明的技术方案及有益效果。
通过武汉大学测绘遥感信息工程国家重点实验室研制的低空无人机多传感器采集***Heli-Mapper对本发明提出的两步法配准模型在无人机MMS中的具体化配准方法进行实验验证。实验结果表明,该方法可以较好解决初始配准误差较大或无初始配准参数条件下的无人机MMS采集的LiDAR点云数据与序列影像的配准问题。多数高反射标靶在进行配准前的反投影像素差在80像素左右,经配准解算后误差降低到像素级。
采用典型城市与居民区车载MMS***成像数据集,对本发明提出的基于两步法配准模型的车载MMS数据融合方法进行实验验证。实验结果表明,该方法可以较好解决初始配准误差较大的车载MMS采集的LiDAR点云数据与序列全景影像的配准问题。在初始配准误差较大的条件下,经由粗到精的配准,配准误差逐渐减小,配准后平均误差在1.5像素左右。
本文中所描述的具体实施例仅仅是对本发明精神作举例说明。本发明所属技术领域的技术人员可以对所描述的具体实施例做各种各样的修改或补充或采用类似的方式替代,但并不会偏离本发明的精神或者超越所附权利要求书所定义的范围。

Claims (6)

1.激光点云与序列全景影像自动配准方法,其特征在于,包括如下步骤:
步骤1,序列全景影像预处理;利用全景图像结构特征的自动恢复方法,将从无人机MMS和车载MMS分别得到的序列全景影像,生成几何结构正确的虚拟框幅式影像;
步骤2,在步骤1得到的虚拟框幅式影像上提取显著的人工建筑物作为配准基元;在无人机MMS数据集中提取建筑物的顶面作为配准基元,在车载MMS数据集中提取建筑物天际线作为配准基元;
步骤3,为了实现点云配准基元和虚拟影像配准基元的配准,将激光点云与影像共轭配准基元查找的问题转化为异源图匹配问题,提出构建SFG空间指纹图描述提取到的几何配准基元集合的空间分布特征,继而自图结构全局空间相似性实现匹配;
步骤4,利用SfM方法恢复标定相机采集的序列影像在摄影测量坐标系中的外方位元素,并且通过多视立体匹配方法得到由序列影像生成MVS影像密集点云,基于步骤3中的SFG空间指纹图构建方法获取MVS影像密集点云与LiDAR点云中的同名角点匹配集合;
步骤5,基于同名角点匹配集合进行粗配准,获得摄影测量坐标系与LiDAR参考坐标系之间的初始转换关系;
步骤6,通过最小化MVS点云与LiDAR点云之间的距离实现配准参数精化,精配准过程分为两个步骤,即粗配准模型参数到精配准模型初值的转换以及迭代最邻近点优化解算。
2.如权利要求1所述的激光点云与序列全景影像自动配准方法,其特征在于:步骤2的具体实现方式为,
步骤2.1,无人机MMS数据集建筑物顶面配准基元提取;首先采用多标记点过程的方法从无人机LiDAR点云数据中提取建筑物点云,并且采用RMBR算法进行建筑物外边界多边形提取与规则化,得到建筑物外框;并依据此点云建筑物提取结果作为先验知识,利用激光点云先验知识引导的影像建筑物顶面提取算法,使用点云数据建筑物提取结果对虚拟框幅式影像上的配准基元提取进行引导,使得初始反投影区域沿张量梯度方向不断优化,最后扩展基于全局对比度显著性检测的分割方法到影像局部区域,结合规则化算法,实现影像建筑物配准几何基元的生成;
步骤2.2,车载MMS数据集天际线配准基元提取;采用一种改进GrabCut的天际线分割方法,从虚拟影像中提取建筑物上边线与角点,再结合应用层次化城市场景目标提取方法与RMBR算法界规则化算法自车载LiDAR点云中实现建筑物立面的提取,确定其上边线,完成影像天际线配准基元提取。
3.如权利要求1所述的激光点云与序列全景影像自动配准方法,其特征在于:步骤3的具体实现方式为,
步骤3.1,构建SFG空间指纹图;以三角形为生成核,节点到核三角形距离最短的图集合构建算法,构建配准基元空间指纹图,SFG空间指纹图集合GC生成分为两个步骤:核三角形生成与非核三角形边的图边连接;
对于无人机MMS数据,选择提取到的n个建筑物顶面配准基元重心组成图节点集合E1。对于车载MMS数据,选择提取到的n个天际线配准基元角点组成图节点集E2;下面以E来指代上述两个图节点集E1和E2,来说明核三角形的生成,选取集合E中的任意3个顶点组合为三角形,作为图G的生成核,则集合E可能构成的图集合GC可表示为:
Figure FDA0002879215590000022
其中Root(E)为核三角形,N为顶点集合E的顶点个数,V(Root(E))为以Root(E)核三角形顶点为自变量的图边集合,V(·)为单调边集合构成函数;图的边连接规则,即边集合构成函数定义为:对于核三角形顶点,直接连接构成完全图;对于非核三角形顶点,则计算该顶点到核三角形三个顶点的距离,将最短长度的边作为连接边;核三角形的边与非核三角形的边共同组成了图的边集合,与图顶点共同构成配准基元空间指纹图;
步骤3.2,设Eimage,Vimage分别是影像配准基元SFG空间指纹图的顶点集合与边集合,同理Elas,Vlas为激光点云配准基元SFG空间指纹图的顶点集合与边集合;使用GED度量影像配准基元SFG空间指纹图Gimage=(Eimage,Vimage)与激光点云SFG空间指纹图Glas=(Elas,Vlas),在局部核三角形匹配对定义的图转换矩阵T下的全局图边相似度;T是将Gimage匹配到Glas的变换矩阵,指经过核三角形匹配后获得的初始图匹配位置;即将最优图匹配的问题转化为最小图编辑距离查找问题,其计算公式如下:
Figure FDA0002879215590000021
其中RootMatch(rm1,..rmk)表示核三角形中的匹配过程中产生的核三角形局部匹配集;k为将Gimage经过修改之后完全匹配Glas所需要的操作步骤数量,所述操作包括添加、删除、替换;Glas与Gimage之间的GED定义为:
Figure FDA0002879215590000031
其中cost(op,T)为将Gimage经过修改之后完全匹配Glas所需要的k步操作所对应的图编辑代价函数,该代价函数为旋转图边差,length(·)为图边长度累加函数。
4.如权利要求1所述的激光点云与序列全景影像自动配准方法,其特征在于:步骤5的具体实现方式为,
步骤5.1,使用同名角点匹配集合解算2D-3D配准模型,即公式(4);令序列影像的关键帧中一个建筑物外多边形角点的像空间坐标为m=(u,v,f)T,其对应的建筑物角点在LiDAR数据中的坐标为Mlas=(X,Y,Z)T,则两角点之间的共线关系可表述为:
spnpm=A[Rpnp|tpnp]Mlas (4)
其中,A为已知的相机内参数矩阵,spnp为比例参数,[Rpnp|tpnp]为EPnP计算获得的相机外参数矩阵;
步骤5.2,根据MVS影像密集点云与LiDAR点云中的同名角点匹配集合,分别得到其在LiDAR参考坐标系和MVS影像密集点云所在的摄影测量坐标系中的点坐标,并解算三维空间相似变换;令在摄影测量坐标系(Cmvs)中的建筑物规则化外边界多边形角点为Mmvs,其对应角点在LiDAR参考坐标系(Cw)中的坐标为Mlas,则其之间的转换关系定义为Mlas=λRMmvs+T,其中λ,R,T分别为尺度、旋转与平移参数;对于n个同名角点对中的每3个同名角点对使用抗差SVD方法即可解算出一组(λ,R,T)。
5.如权利要求4所述的激光点云与序列全景影像自动配准方法,其特征在于:步骤6的具体实现方式如下;
步骤6.1,使用同名角点匹配集合完成粗配准后,通过最小化MVS点云与LiDAR点云之间的距离实现配准参数精化,其配准几何模型可表示为:
Mlas=λRMmvs+T (5)
2D-3D模型初值转换:令一像点为m,其在Cmvs中对应点为MMVS,在Cw中对应点为Mlas,则共线性关系可表达为:
sbundlem=A[Rbundle|tbundle]Mmvs (6)
其中sbundle为比例参数,[Rbundle|tbundle]为SfM计算获得的影像外参数矩阵;将2D-3D模型中得到的影像外参数矩阵转换到式(5)中的R变量,实现2D-3D的转换,联立式(4)至(6)可得:
Figure FDA0002879215590000041
其中R为旋转参数,
Figure FDA0002879215590000042
表示Rpnp的逆,Rpnp表示EPnP计算获得的相机外参数矩阵,Rbundle表示SfM计算的影像外参数,Ppnp(Xi,Yi,Zi)t,令Pbundle(xi,yi,zi)t为相机在Cw与Cmvs中的坐标重心化后的相机位姿,则其尺度参数可通过式(8)计算,平移参数通过式(9)计算:
Figure FDA0002879215590000043
T=Ppnp-λRPbundle (9)
由公式(7)至(9)可知,一张粗配准影像对应一种Cw与Cmvs之间的空间旋转关系,多张配准影像则可实现对全部空间转换参数(λ,R,T)的求解;
步骤6.2,采用刚体空间相似变换模型作为MVS与LiDAR点云之间的配准模型,使用点到面距离差作为点对误差,则误差方程可表示为:
Figure FDA0002879215590000044
其中,mi,Mi分别为MVS与LiDAR点,wi为点对权,ηi为Mi点法向量,(λ3d-3d,R3d-3d,T3d-3d)表示需要求解的最小MVS与LiDAR点云点到面距离的空间相似变换;基于最优摄影测量与LiDAR参考系之间的坐标转换(λ3d-3d,R3d-3d,T3d-3d),修正后的影像外方位元素(Rcam,Tcam)可表示为:
Figure FDA0002879215590000045
6.如权利要求5所述的激光点云与序列全景影像自动配准方法,其特征在于:步骤6.1中,进行全部空间转换参数(λ,R,T)的求解时,采用投票聚类的方法实现参数的稳健估计。
CN202011625760.3A 2020-12-31 2020-12-31 激光点云与序列全景影像自动配准方法 Pending CN112767461A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202011625760.3A CN112767461A (zh) 2020-12-31 2020-12-31 激光点云与序列全景影像自动配准方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202011625760.3A CN112767461A (zh) 2020-12-31 2020-12-31 激光点云与序列全景影像自动配准方法

Publications (1)

Publication Number Publication Date
CN112767461A true CN112767461A (zh) 2021-05-07

Family

ID=75698966

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202011625760.3A Pending CN112767461A (zh) 2020-12-31 2020-12-31 激光点云与序列全景影像自动配准方法

Country Status (1)

Country Link
CN (1) CN112767461A (zh)

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113160419A (zh) * 2021-05-11 2021-07-23 北京京东乾石科技有限公司 一种建筑物立面模型建立方法和装置
CN113343765A (zh) * 2021-05-11 2021-09-03 武汉大学 一种基于点云刚性配准的场景检索方法及***
CN113658166A (zh) * 2021-08-24 2021-11-16 凌云光技术股份有限公司 一种基于网格模型的点云缺陷检测方法和装置
CN114549601A (zh) * 2022-02-11 2022-05-27 中国科学院精密测量科学与技术创新研究院 一种顾及点对可靠性的滑坡多时相tls点云精配准方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
陈驰 等: "《低空UAV激光点云和序列影像的自动配准方法》", 《测绘学报》 *
陈驰 等: "《车载MMS激光点云与序列全景影像自动配准方法》", 《测绘学报》 *

Cited By (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113160419A (zh) * 2021-05-11 2021-07-23 北京京东乾石科技有限公司 一种建筑物立面模型建立方法和装置
CN113343765A (zh) * 2021-05-11 2021-09-03 武汉大学 一种基于点云刚性配准的场景检索方法及***
CN113343765B (zh) * 2021-05-11 2022-07-22 武汉大学 一种基于点云刚性配准的场景检索方法及***
CN113160419B (zh) * 2021-05-11 2024-02-02 北京京东乾石科技有限公司 一种建筑物立面模型建立方法和装置
CN113658166A (zh) * 2021-08-24 2021-11-16 凌云光技术股份有限公司 一种基于网格模型的点云缺陷检测方法和装置
CN113658166B (zh) * 2021-08-24 2024-04-12 凌云光技术股份有限公司 一种基于网格模型的点云缺陷检测方法和装置
CN114549601A (zh) * 2022-02-11 2022-05-27 中国科学院精密测量科学与技术创新研究院 一种顾及点对可靠性的滑坡多时相tls点云精配准方法
CN114549601B (zh) * 2022-02-11 2023-03-28 中国科学院精密测量科学与技术创新研究院 一种顾及点对可靠性的滑坡多时相tls点云精配准方法

Similar Documents

Publication Publication Date Title
Yang et al. Automatic registration of UAV-borne sequent images and LiDAR data
CN112767461A (zh) 激光点云与序列全景影像自动配准方法
CN112927360A (zh) 一种基于倾斜模型与激光点云数据融合的三维建模方法和***
CN111144388A (zh) 一种基于单目影像的道路标志标线的更新方法
CN113192193B (zh) 基于Cesium三维地球框架的高压输电线路走廊三维重建方法
CN112465732A (zh) 一种车载激光点云与序列全景影像的配准方法
WO2018061010A1 (en) Point cloud transforming in large-scale urban modelling
CN113607135A (zh) 一种用于路桥施工领域的无人机倾斜摄影测量方法
CN112270698A (zh) 基于最邻近曲面的非刚性几何配准方法
CN112465849B (zh) 一种无人机激光点云与序列影像的配准方法
Barazzetti et al. Fisheye lenses for 3D modeling: Evaluations and considerations
CN113077552A (zh) 基于无人机影像的dsm生成方法和装置
CN110986888A (zh) 一种航空摄影一体化方法
CN116758234A (zh) 一种基于多点云数据融合的山地地形建模方法
CN110443837B (zh) 一种直线特征约束下的城区机载激光点云与航空影像配准方法和***
Gao et al. Multi-source data-based 3D digital preservation of largescale ancient chinese architecture: A case report
Li et al. New methodologies for precise building boundary extraction from LiDAR data and high resolution image
CN113012206B (zh) 一种顾及房檐特征的机载与车载LiDAR点云配准方法
CN112767459A (zh) 基于2d-3d转换的无人机激光点云与序列影像配准方法
Zhao et al. Alignment of continuous video onto 3D point clouds
CN114972672B (zh) 输电线路实景三维模型构建方法、装置、设备和存储介质
CN107784666B (zh) 基于立体影像的地形地物三维变化检测和更新方法
CN105093222A (zh) 一种sar影像区域网平差连接点自动提取方法
Chen et al. True orthophoto generation using multi-view aerial images
Previtali et al. An automatic multi-image procedure for accurate 3D object reconstruction

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
RJ01 Rejection of invention patent application after publication
RJ01 Rejection of invention patent application after publication

Application publication date: 20210507