CN112581511A - 基于近似直立扫描点云快速配准的三维重建方法及系统 - Google Patents

基于近似直立扫描点云快速配准的三维重建方法及系统 Download PDF

Info

Publication number
CN112581511A
CN112581511A CN202011471438.XA CN202011471438A CN112581511A CN 112581511 A CN112581511 A CN 112581511A CN 202011471438 A CN202011471438 A CN 202011471438A CN 112581511 A CN112581511 A CN 112581511A
Authority
CN
China
Prior art keywords
point cloud
voxel
registration
original
target
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.)
Granted
Application number
CN202011471438.XA
Other languages
English (en)
Other versions
CN112581511B (zh
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.)
Beijing Forestry University
Original Assignee
Beijing Forestry University
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 Beijing Forestry University filed Critical Beijing Forestry University
Priority to CN202011471438.XA priority Critical patent/CN112581511B/zh
Publication of CN112581511A publication Critical patent/CN112581511A/zh
Application granted granted Critical
Publication of CN112581511B publication Critical patent/CN112581511B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/30—Determination of transform parameters for the alignment of images, i.e. image registration
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T17/00—Three dimensional [3D] modelling, e.g. data description of 3D objects
    • G06T17/20—Finite element generation, e.g. wire-frame surface description, tesselation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/60—Analysis of geometric attributes
    • G06T7/62—Analysis of geometric attributes of area, perimeter, diameter or volume
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/10—Image acquisition modality
    • G06T2207/10028—Range image; Depth image; 3D point clouds
    • Y—GENERAL 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
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02P—CLIMATE CHANGE MITIGATION TECHNOLOGIES IN THE PRODUCTION OR PROCESSING OF GOODS
    • Y02P90/00—Enabling technologies with a potential contribution to greenhouse gas [GHG] emissions mitigation
    • Y02P90/30—Computing systems specially adapted for manufacturing

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Geometry (AREA)
  • Computer Graphics (AREA)
  • Software Systems (AREA)
  • Image Generation (AREA)
  • Length Measuring Devices By Optical Means (AREA)

Abstract

本发明属于三维数据采集与处理领域，具体涉及了一种基于近似直立扫描点云快速配准的三维重建方法及系统，旨在解决现有技术未考虑实际点云采集中的点云近似直立性，点云配准效率低，从而导致三维重建效率、准确率和精度低的问题。本发明包括：获取近似直立扫描的目标点云和原始点云；通过自适应网格大小进行点云体素化；采用先验可靠区域引导点云初始配准；kd树结合ICP算法进行点云精确配准；结合初始配准和精确配准进行点云最终配准；基于配准后的点云进行三维重建。本发明充分考虑点云的近似直立性，快速精确实现点云配准，三维重建效率高、准确率和精度高。

Description

基于近似直立扫描点云快速配准的三维重建方法及系统
技术领域
本发明属于三维数据采集与处理领域，具体涉及了一种基于近似直立扫描点云快速配准的三维重建方法及系统。
背景技术
三维激光扫描技术以其精度高，速度快的特点，在数字城市、文物数字化、非接触测量和三维建模等方面得到了越来越广泛的应用。而点云配准是三维建模重建技术中非常重要的一步。由于物体形状、环境等遮挡或者自遮挡原因，一次测量只能获取被测物体或场景的一部分点云数据，为了完整获取物体的信息往往需要从不同角度多次采集。扫描仪在不同的放置位置、不同扫描角度采集的点云数据的坐标系统不一样，点云配准的目标是找到两个不同坐标系统之间的坐标变换公式。
现有技术中点云配准主要包括三种类型：第一种是基于点云几何特征的配准，一些文献提出重构的多边形表面可以作为坐标变换中不变的特征描述子[1]，一些文献利用法向量与其k近邻之间的夹角对尺度和旋转变换是不变量来构造基于角度的点云配准[2]，一些文献提出一种基于轮廓进行配准的算法[3]，还有一些文献提出了一种具有识别性和鲁棒性的三维形状描述子区域曲率图，用于三维点云自动配准[4]，此外，一些配准算法通过定义最大限度稳定形状指数区域、结合局部几何特性和统计特性、将几何结构集成到深度神经网络、使用法线点等来完成点云的配准。然而，虽然这些方法充分利用了几何信息，但对特征提取和特征匹配的误差更为敏感，因而最终配准结果往往有较大误差，甚至有些配准结果是错误的。第二种是传统ICP算法及其改进，迭代最近点算法(ICP)是在点云具有良好的先验姿态信息的前提下，通过迭代估计坐标变换公式中的参数，对两个有重叠的点云数据进行精确配准。大量的文献针对该过程的每个步骤提出了许多改进，一些文献通过构造点云局部特征量来优化配准算法[5]，一些文献结合曲面的法向特征和曲率特征，提出了法向ICP方法[6]，还有一些文献通过引入删除掩码的概念，对经典的ICP算法进行了改进[7]。这些方法是利用形状或者同级特征来加速ICP算法的配准，然而，由于形状和统计特征的估计往往误差较大，配准算法的速度效率提高并不显著。第三种是混合方法，包括几何方法、统计方法、优化方法和深度学习等方法及其混合方法。一些文献认为配准问题可以表示为混合整数线性规划的分枝定界问题[8]，一些文献基于牛顿迭代算法提出了改进的三维正态分布变换算法[9]。多阶段配准方法包括两阶段准则，即初始对齐和全局一致优化，一些文献提出的算法由提取平面特征、计算每个平面斑块的面积、利用几何一致性将启发式搜索与剪枝相结合进行配准这三个步骤组成[10]，在找到解后利用加权最小二乘法对变换进行细化。另外，一些配准方法采用了同时定位与映射(SLAM)技术、双向最大相关准则、基于学习的方法包括深度神经网络、通用图论框架等。然而，这些方法都没有利用原始点云和目标点云等大小体素化策略来简化点云，配准精度和速度上仍有待进一步提高。
总的来说，现有技术中点云配准的三种类型的方法，还存在上述的各种问题，并且，在现实数据采集中，无论是地面扫描还是车载扫描，扫描仪都是水平放置的，这就意味着不同坐标下的点云数据方向大致都是Z轴向上的，点云配准在实际工程应用中不需要在整个旋转空间中寻找最优解，因此，本领域还急需一种基于近似直立扫描点云快速配准的三维重建方法，在近似直立扫描的点云数据快速配准的前提下，实现高效、高质量的三维重建。
以下文献是与本发明相关的技术背景资料：
[1]Li P.,Cheng P.,Sutton M.A.:Three-dimensional point cloudregistration by matching surface features with relaxation labelingmethod.Experimental Mechanics 45,71–82(2005).
[2]Jiang,J.,Cheng,J.,Chen,X.:Registration for 3-d point cloud usingangular-invariant feature.Neurocomputing 72(16-18):3839-3844(2009).
[3]Bennis,A.,Bombardier,V.,Thiriet,P.,Brie,D.:Contours based approachfor thermal image and terrestrial point cloud registration.IsprsInternational Archives of the Photogrammetry Remote Sensing&SpatialInformation Sciences 40(5),97–101(2013).
[4]Sun,J.,Zhang,J.,G.:An automatic 3d point cloud registration methodbased on regional curvature maps.Image Vision Comput.56(C),49–58(2016).
[5]Gressin,A.,Mallet,C.,Demantke,J.,David,N.:Towards 3d lidar pointcloud registration improvement using optimal neighborhood knowledge.IsprsJournal of Photogrammetry&Remote Sensing 79(may),240–251(2013).
[6]Serafin,J.,Grisetti,G.:Nicp:Dense normal based point cloudregistration.IEEE/RSJ International Conference on Intelligent Robots andSystems pp.742–749(2015).
[7]Roberto Marani等:A modified iterative closest point algorithm for3d point cloud registration.Computer Aided Civil&Infrastructure Engineering31(7),515–534(2016).
[8]Yu,C.,Da,J.:A maximum feasible subsystem for globally optimal 3dpoint cloud registration.Sensors 18(2),544.1–19(2018).
[9]Hu,F.,Ren,T.,Shi,S.:Discrete point cloud registration using the 3dnormal distribution transformation based newton iteration.Journal ofMultimedia 9(7),934–940(2014).
[10]Gojcic,Z.,Zhou,C.,Wegner,J.D.,Guibas,L.J.,Birdal,T.:Learningmultiview 3d point cloud registration.In:2020 IEEE Conference on ComputerVision and Pattern Recognition(CVPR)(2020).
发明内容
为了解决现有技术中的上述问题，即现有技术未考虑实际点云采集中的点云近似直立性，点云配准效率低，从而导致三维重建效率、准确率和精度低的问题，本发明提供了一种基于近似直立扫描点云快速配准的三维重建方法，该方法包括：
步骤S10，获取近似直立扫描点云作为三维重建的待配准的目标点云Ωt和原始点云Ωs；
步骤S20，分别计算所述目标点云Ωt和原始点云Ωs的轴向包围盒AABBt和AABBs，并基于轴向包围盒AABBt和AABBs的尺寸进行所述待目标点云和原始点云的体素化，获得目标体素点云Vt和原始体素点云Vs；
步骤S30，分别计算所述目标体素点云Vt和原始体素点云Vs的几何中心点Ct和Cs，根据扫描的先验信息确定先验可靠区域，并进行所述目标体素点云和原始体素点云的初始配准，获得初始配准变换矩阵；
步骤S40，通过基于kd树的ICP方法进行初始配准后的目标体素点云和原始体素点云的精确配准，获得精确配准变换矩阵；
步骤S50，将所述初始配准变换矩阵和精确配准变换矩阵复合，获得最终配准变换矩阵，并基于所述最终配准变换矩阵进行所述目标点云Ωt和原始点云Ωs的配准，并基于配准后的点云进行三维重建。
在一些优选的实施例中，步骤S20包括：
步骤S21，分别计算所述目标点云Ωt和原始点云Ωs的三维坐标最大值和最小值，获得目标点云和原始点云的轴向包围盒AABBt和AABBs；
步骤S22，分别基于所述目标点云和原始点云的轴向包围盒AABBt和AABBs的尺寸以及设定的轴向最优体素数Δ和最小体素数δ，获得所述目标点云和原始点云对应的单个体素边长
Figure BDA0002833915720000051
和
Figure BDA0002833915720000052
步骤S23，基于所述目标点云和原始点云对应的单个体素边长
Figure BDA0002833915720000053
和
Figure BDA0002833915720000054
进行所述待目标点云Ωt和原始点云Ωs的体素化，获得目标体素点云Vt和原始体素点云Vs。
在一些优选的实施例中，所述目标点云Ωt和原始点云Ωs的轴向包围盒AABBt和AABBs为：
Figure BDA0002833915720000055
Figure BDA0002833915720000056
其中，
Figure BDA0002833915720000057
和
Figure BDA0002833915720000058
分别代表目标点云中所有点的三维坐标最大值和最小值，
Figure BDA0002833915720000059
和
Figure BDA00028339157200000510
分别代表原始点云中所有点的三维坐标最大值和最小值。
在一些优选的实施例中，步骤S22包括：
步骤S221，分别计算AABBt的最大轴向长
Figure BDA00028339157200000511
和最小轴向长
Figure BDA00028339157200000512
以及AABBs的最大轴向长
Figure BDA00028339157200000513
和最小轴向长
Figure BDA00028339157200000514
Figure BDA00028339157200000515
Figure BDA00028339157200000516
Figure BDA00028339157200000517
Figure BDA00028339157200000518
其中，max代表求最大值操作，min代表求最小值操作；
步骤S222，基于所述AABBt的最大轴向长
Figure BDA00028339157200000519
和最小轴向长
Figure BDA00028339157200000520
以及AABBs的最大轴向长
Figure BDA00028339157200000521
和最小轴向长
Figure BDA00028339157200000522
设定的轴向最优体素数Δ和最小体素数δ，分别计算目标点云和原始点云的单个体素边长
Figure BDA00028339157200000523
和
Figure BDA00028339157200000524
Figure BDA0002833915720000061
其中，
Figure BDA0002833915720000062
代表目标点云和原始点云的最小边长，
Figure BDA0002833915720000063
代表目标点云和原始点云对应的初始体素边长，
Figure BDA0002833915720000064
和代表目标点云和原始点云的最小边长以长度L0为体素边长的体素数。
在一些优选的实施例中，步骤S23包括：
分别将所述目标点云和原始点云分割为Wt×Ht×Dt和Ws×Hs×Ds的方块：
Figure BDA0002833915720000065
Figure BDA0002833915720000066
Figure BDA0002833915720000067
分别基于目标点云中第i个点
Figure BDA0002833915720000068
和原始点云中第i个点
Figure BDA0002833915720000069
计算对应的体素坐标
Figure BDA00028339157200000610
和
Figure BDA00028339157200000611
Figure BDA00028339157200000612
Figure BDA00028339157200000613
Figure BDA00028339157200000614
其中，
Figure BDA00028339157200000615
代表下取整操作；
以分割获得的坐标为体素坐标的点云作为目标体素点云Vt和原始体素点云Vs。
在一些优选的实施例中，步骤S30包括：
步骤S31，分别计算所述目标体素点云Vt和原始体素点云Vs的几何中心点Ct和Cs；以所述目标点云Ωt和原始点云Ωs的密度大于设定阈值或距离获取点云的扫描仪所在位置设定距离范围内的区域作为先验可靠区域；
步骤S32，获取先验可靠区域中Vt和Vs中心化后的定位向量的均值方向
Figure BDA0002833915720000071
和
Figure BDA0002833915720000072
的平均夹角θ，并计算旋转变换矩阵Tr：
Figure BDA0002833915720000073
步骤S33，基于所述目标体素点云Vt和原始体素点云Vs的几何中心点Ct和Cs以及旋转变换矩阵Tr，获取初始配准变换矩阵：
V′t＝Tr·(Vt-Ct),V′s＝Is·(Vs-Cs)
其中，Is为单位变换矩阵，V′t和V′s为初始配准后的目标体素点云和原始体素点云。
在一些优选的实施例中，步骤S40包括：
步骤S41，通过KDTree在V′s中搜索V′t中任意一点的最近邻域点，获取V′t中点集Vtm在V′s上的映射点集Vsm；
步骤S42，基于预先设定的最小阈值距离D删除点集Vtm和点集Vsm中的错误映射点对，获得点集V′tm和映射点集V′sm；
步骤S43，利用四元素法计算点集V′tm和映射点V′sm的变换矩阵R和平移量T，并基于变换矩阵R和平移量T进行点集V′tm的变换，获得点集Vt″m；
步骤S44，重复执行步骤S41-步骤S43，直至获取点集Vtm中每一点到映射点集Vsm中的映射点对和对应的变换矩阵R、平移量T；
步骤S45，迭代计算直至最后的变换矩阵R、平移量T配准后，对应点对坐标间差值小于设定阈值，获得精确配准变换矩阵：
V′t＝TRicp·V′s+TTicp
其中，TRicp和TTicp分别为ICP方法生成的旋转矩阵和平移矩阵。
在一些优选的实施例中，所述最终配准变换矩阵，其公式表示为：
Vt＝Rfull·Vs+Tfull
其中，
Figure BDA0002833915720000081
在一些优选的实施例中，步骤S50中基于所述最终配准变换矩阵进行所述目标点云Ωt和原始点云Ωs的配准，其方法为：
Ωt＝Rfull·Ωs+Tfull。
本发明的另一方面，提出了一种基于近似直立扫描点云快速配准的三维重建系统，该系统包括以下模块：
点云获取模块，配置为获取近似直立扫描点云作为三维重建的待配准的目标点云Ωt和原始点云Ωs
体素化模块，配置为分别计算所述目标点云Ωt和原始点云Ωs的轴向包围盒AABBt和AABBs，并基于轴向包围盒AABBt和AABBs的尺寸进行所述待目标点云和原始点云的体素化，获得目标体素点云Vt和原始体素点云Vs；
初始配准模块，配置为分别计算所述目标体素点云Vt和原始体素点云Vs的几何中心点Ct和Cs，根据扫描的先验信息确定先验可靠区域，并进行所述目标体素点云和原始体素点云的初始配准，获得初始配准变换矩阵；
精确配准模块，配置为通过基于kd树的ICP方法进行初始配准后的目标体素点云和原始体素点云的精确配准，获得精确配准变换矩阵；
点云配准模块，配置为将所述初始配准变换矩阵和精确配准变换矩阵复合，获得最终配准变换矩阵，并基于所述最终配准变换矩阵进行所述目标点云Ωt和原始点云Ωs的配准；
重建模块，配置为基于配准后的点云进行三维重建。
本发明的有益效果：
(1)本发明基于近似直立扫描点云快速配准的三维重建方法，采用自适应网格大小的点云体素化方法，解决点云数据量过大，许多算法运行效率低甚至无法运行的问题，点云体素化之后，体素中心点的总数量是可控的，同时还能应用于点云密度不均匀的情况，从而提升了三维重建的效率、准确性和精度。
(2)本发明基于近似直立扫描点云快速配准的三维重建方法，采用先验可靠区域引导(APRR，Prior Reliable Region guided Alignment)的最优水平旋转角度确定方法，实现点云的初始粗略配准。先验可靠区域的确定可以根据点云密度确定，也可以依据实际扫描中扫描器的位置、扫描方向以及遮挡等因素来确定。实验表明，先验可靠区域引导的最优水平旋转角度确定方法，其初始配准效果远超现有的主成分分析法，从而进一步提升了三维重建的准确性和精度。
(3)本发明基于近似直立扫描点云快速配准的三维重建方法，采用基于kd树实现ICP算法运行效率的提升，将kd树搜索邻域的高效算法与ICP算法相结合，大大提高了对海量点云配准的效率。kd树算法是一种二叉树搜索拓展到三维空间的多维搜索树结构形式，kd树构建算法采用逐级递归的方式，基于最大方差确定分裂维，以分裂维所在结点作为根节点，将数据结构划分为左右两个空间，然后采用递归方式对左右两个空间继续划分直至每一个子集内不再包含点云为止，通过此方法大大提升了点云精确配准的效率，从而进一步提升了三维重建的效率。
附图说明
通过阅读参照以下附图所作的对非限制性实施例所作的详细描述，本申请的其它特征、目的和优点将会变得更明显：
图1是本发明基于近似直立扫描点云快速配准的三维重建方法的流程示意图。
具体实施方式
下面结合附图和实施例对本申请作进一步的详细说明。可以理解的是，此处所描述的具体实施例仅用于解释相关发明，而非对该发明的限定。另外还需要说明的是，为了便于描述，附图中仅示出了与有关发明相关的部分。
需要说明的是，在不冲突的情况下，本申请中的实施例及实施例中的特征可以相互组合。下面将参考附图并结合实施例来详细说明本申请。
本发明提供一种基于近似直立扫描点云快速配准的三维重建方法，本方法通过自适应网格大小的点云体素化方法，解决点云数据量过大，许多算法运行效率低甚至无法运行的问题，点云体素化之后，体素中心点的总数量是可控的；提出了先验可靠区域引导的最优水平旋转角度确定方法，实现点云的初始配准；采用kd树数据存储结构提高ICP算法运行效率，实现点云的精确配准；结合初始配准和精确配准进行目标点云和原始点云的最终配准，并基于配准后的点云进行三维重建，重建效率高、准确性和精度高，适用范围广泛。
本发明的一种基于近似直立扫描点云快速配准的三维重建方法，该方法包括：
步骤S10，获取近似直立扫描点云作为三维重建的待配准的目标点云Ωt和原始点云Ωs；
步骤S20，分别计算所述目标点云Ωt和原始点云Ωs的轴向包围盒AABBt和AABBs，并基于轴向包围盒AABBt和AABBs的尺寸进行所述待目标点云和原始点云的体素化，获得目标体素点云Vt和原始体素点云Vs；
步骤S30，分别计算所述目标体素点云Vt和原始体素点云Vs的几何中心点Ct和Cs，根据扫描的先验信息确定先验可靠区域，并进行所述目标体素点云和原始体素点云的初始配准，获得初始配准变换矩阵；
步骤S40，通过基于kd树的ICP方法进行初始配准后的目标体素点云和原始体素点云的精确配准，获得精确配准变换矩阵；
步骤S50，将所述初始配准变换矩阵和精确配准变换矩阵复合，获得最终配准变换矩阵，并基于所述最终配准变换矩阵进行所述目标点云Ωt和原始点云Ωs的配准，并基于配准后的点云进行三维重建。
为了更清晰地对本发明基于近似直立扫描点云快速配准的三维重建方法进行说明，下面结合图1对本发明实施例中各步骤展开详述。
本发明第一实施例的基于近似直立扫描点云快速配准的三维重建方法，包括步骤S10-步骤S50，各步骤详细描述如下：
步骤S10，获取近似直立扫描点云作为三维重建的待配准的目标点云Ωt和原始点云Ωs。
步骤S20，分别计算所述目标点云Ωt和原始点云Ωs的轴向包围盒AABBt和AABBs，并基于轴向包围盒AABBt和AABBs的尺寸进行所述待目标点云和原始点云的体素化，获得目标体素点云Vt和原始体素点云Vs。
步骤S21，分别计算所述目标点云Ωt和原始点云Ωs的三维坐标最大值和最小值，获得目标点云和原始点云的轴向包围盒AABBt和AABBs，如式(1)和式(2)所示：
Figure BDA0002833915720000111
Figure BDA0002833915720000112
其中，
Figure BDA0002833915720000121
和
Figure BDA0002833915720000122
分别代表目标点云中所有点的三维坐标最大值和最小值，
Figure BDA0002833915720000123
和
Figure BDA0002833915720000124
分别代表原始点云中所有点的三维坐标最大值和最小值。
步骤S22，分别基于所述目标点云和原始点云的轴向包围盒AABBt和AABBs的尺寸以及设定的轴向最优体素数Δ和最小体素数δ，获得所述目标点云和原始点云对应的单个体素边长
Figure BDA0002833915720000125
和
Figure BDA0002833915720000126
步骤S221，分别计算AABBt的最大轴向长
Figure BDA0002833915720000127
和最小轴向长
Figure BDA0002833915720000128
以及AABBs的最大轴向长
Figure BDA0002833915720000129
和最小轴向长
Figure BDA00028339157200001210
如式(3)-式(6)所示：
Figure BDA00028339157200001211
Figure BDA00028339157200001212
Figure BDA00028339157200001213
Figure BDA00028339157200001214
其中，max代表求最大值操作，min代表求最小值操作。
步骤S222，基于所述AABBt的最大轴向长
Figure BDA00028339157200001215
和最小轴向长
Figure BDA00028339157200001216
以及AABBs的最大轴向长
Figure BDA00028339157200001217
和最小轴向长
Figure BDA00028339157200001218
设定的轴向最优体素数Δ和最小体素数δ，分别计算目标点云和原始点云的单个体素边长
Figure BDA00028339157200001219
和
Figure BDA00028339157200001220
如式(7)所示：
Figure BDA00028339157200001221
其中，
Figure BDA00028339157200001222
代表目标点云和原始点云的最小边长，
Figure BDA00028339157200001223
代表目标点云和原始点云对应的初始体素边长，
Figure BDA00028339157200001224
和代表目标点云和原始点云的最小边长以长度L0为体素边长的体素数。
通过ηt和ηs与设定的最小体素数δ的关系避免出现体素化的分辨率过大。
步骤S23，基于所述目标点云和原始点云对应的单个体素边长
Figure BDA0002833915720000131
和
Figure BDA0002833915720000132
进行所述待目标点云Ωt和原始点云Ωs的体素化，获得目标体素点云Vt和原始体素点云Vs。
分别将所述目标点云和原始点云分割为Wt×Ht×Dt和Ws×Hs×Ds的方块，如式(8)-式(10)所示：
Figure BDA0002833915720000133
Figure BDA0002833915720000134
Figure BDA0002833915720000135
分别基于目标点云中第i个点
Figure BDA0002833915720000136
和原始点云中第i个点
Figure BDA0002833915720000137
计算对应的体素坐标
Figure BDA0002833915720000138
和
Figure BDA0002833915720000139
如式(11)-式(13)所示：
Figure BDA00028339157200001310
Figure BDA00028339157200001311
Figure BDA00028339157200001312
其中，
Figure BDA00028339157200001313
代表下取整操作；
以分割获得的坐标为体素坐标的点云作为目标体素点云Vt和原始体素点云Vs。
步骤S30，分别计算所述目标体素点云Vt和原始体素点云Vs的几何中心点Ct和Cs，根据扫描的先验信息确定先验可靠区域，并进行所述目标体素点云和原始体素点云的初始配准，获得初始配准变换矩阵。
步骤S31，分别计算所述目标体素点云Vt和原始体素点云Vs的几何中心点Ct和Cs；以所述目标点云Ωt和原始点云Ωs的密度大于设定阈值或距离获取点云的扫描仪所在位置设定距离范围内的区域作为先验可靠区域；
先验可靠区域的确定方法有两种：第一种是根据点云密度确定。密度大于设定值的区域更可靠。第二种是依据扫描实践中扫描器的位置、扫描方向以及遮挡等因素来确定。例如，车载扫描时，中部高度区域较可靠；机载扫描时，顶部区域的数据较可靠；定点扫描物体时，如果将扫描仪置于较低位置，则物体较低的区域更可靠。相反，如果将扫描仪置于顶部位置，则物体的上部区域更可靠。
步骤S32，获取先验可靠区域中Vt和Vs中心化后的定位向量的均值方向
Figure BDA0002833915720000141
和
Figure BDA0002833915720000142
的平均夹角θ，并计算旋转变换矩阵Tr，如式(14)所示：
Figure BDA0002833915720000143
步骤S33，基于所述目标体素点云Vt和原始体素点云Vs的几何中心点Ct和Cs以及旋转变换矩阵Tr，获取初始配准变换矩阵，如式(15)所示：
V′t＝Tr·(Vt-Ct),V′s＝Is·(Vs-Cs) (15)
其中，Is为单位变换矩阵，V′t和V′s为初始配准后的目标体素点云和原始体素点云。
由APRR(Prior Reliable Region guided Alignment)推导出的式(15)是一个刚体变换，它避免了某些无反射约束条件的最优问题出现反射变换。
步骤S40，通过基于kd树的ICP方法进行初始配准后的目标体素点云和原始体素点云的精确配准，获得精确配准变换矩阵。
步骤S41，通过KDTree在V′s中搜索V′t中任意一点的最近邻域点，获取V′t中点集Vtm在V′s上的映射点集Vsm；
步骤S42，基于预先设定的最小阈值距离D删除点集Vtm和点集Vsm中的错误映射点对，获得点集V′tm和映射点集V′sm；
步骤S43，利用四元素法计算点集V′tm和映射点V′sm的变换矩阵R和平移量T，并基于变换矩阵R和平移量T进行点集V′tm的变换，获得点集Vt″m；
步骤S44，重复执行步骤S41-步骤S43，直至获取点集Vtm中每一点到映射点集Vsm中的映射点对和对应的变换矩阵R、平移量T；
步骤S45，迭代计算直至最后的变换矩阵R、平移量T配准后，对应点对坐标间差值小于设定阈值，获得精确配准变换矩阵，如式(16)所示：
V′t＝TRicp·V′s+TTicp (16)
其中，TRicp和TTicp分别为ICP方法生成的旋转矩阵和平移矩阵。
步骤S50，将所述初始配准变换矩阵和精确配准变换矩阵复合，获得最终配准变换矩阵，并基于所述最终配准变换矩阵进行所述目标点云Ωt和原始点云Ωs的配准，并基于配准后的点云进行三维重建。
最终配准变换矩阵，其公式表示如式(17)：
Vt＝Rfull·Vs+Tfull (17)
其中，
Figure BDA0002833915720000151
基于所述最终配准变换矩阵进行所述目标点云Ωt和原始点云Ωs的配准，如式(18)所示：
Ωt＝Rfull·Ωs+Tfull (18)
在三维重建之前，还可以对配准的点云进行简化操作：
完成对输入目标点云Ωt和原始点云Ωs的配准后，会产生很多冗余点，尤其是在重叠区域。为了去除冗余点，简化点云，可采用配准后的点云进行体素化，用体素点云的中心点云作为配准后的点云的简化点云，每个体素的大小可以通过式(7)计算，并根据点云分辨率要求确定轴向最优体素数Δ的值。
本发明第二实施例的基于近似直立扫描点云快速配准的三维重建系统，该系统包括以下模块：
点云获取模块，配置为获取近似直立扫描点云作为三维重建的待配准的目标点云Ωt和原始点云Ωs
体素化模块，配置为分别计算所述目标点云Ωt和原始点云Ωs的轴向包围盒AABBt和AABBs，并基于轴向包围盒AABBt和AABBs的尺寸进行所述待目标点云和原始点云的体素化，获得目标体素点云Vt和原始体素点云Vs；
初始配准模块，配置为分别计算所述目标体素点云Vt和原始体素点云Vs的几何中心点Ct和Cs，根据扫描的先验信息确定先验可靠区域，并进行所述目标体素点云和原始体素点云的初始配准，获得初始配准变换矩阵；
精确配准模块，配置为通过基于kd树的ICP方法进行初始配准后的目标体素点云和原始体素点云的精确配准，获得精确配准变换矩阵；
点云配准模块，配置为将所述初始配准变换矩阵和精确配准变换矩阵复合，获得最终配准变换矩阵，并基于所述最终配准变换矩阵进行所述目标点云Ωt和原始点云Ωs的配准；
重建模块，配置为基于配准后的点云进行三维重建。
所属技术领域的技术人员可以清楚地了解到，为描述的方便和简洁，上述描述的系统的具体工作过程及有关说明，可以参考前述方法实施例中的对应过程，在此不再赘述。
需要说明的是，上述实施例提供的基于近似直立扫描点云快速配准的三维重建系统，仅以上述各功能模块的划分进行举例说明，在实际应用中，可以根据需要而将上述功能分配由不同的功能模块来完成，即将本发明实施例中的模块或者步骤再分解或者组合，例如，上述实施例的模块可以合并为一个模块，也可以进一步拆分成多个子模块，以完成以上描述的全部或者部分功能。对于本发明实施例中涉及的模块、步骤的名称，仅仅是为了区分各个模块或者步骤，不视为对本发明的不当限定。
本发明第三实施例的一种存储装置，其中存储有多条程序，所述程序适于由处理器加载并执行以实现上述的基于近似直立扫描点云快速配准的三维重建方法。
本发明第四实施例的一种处理装置，包括处理器、存储装置；处理器，适于执行各条程序；存储装置，适于存储多条程序；所述程序适于由处理器加载并执行以实现上述的基于近似直立扫描点云快速配准的三维重建方法。
所属技术领域的技术人员可以清楚地了解到，为描述的方便和简洁，上述描述的存储装置、处理装置的具体工作过程及有关说明，可以参考前述方法实施例中的对应过程，在此不再赘述。
本领域技术人员应该能够意识到，结合本文中所公开的实施例描述的各示例的模块、方法步骤，能够以电子硬件、计算机软件或者二者的结合来实现，软件模块、方法步骤对应的程序可以置于随机存储器(RAM)、内存、只读存储器(ROM)、电可编程ROM、电可擦除可编程ROM、寄存器、硬盘、可移动磁盘、CD-ROM、或技术领域内所公知的任意其它形式的存储介质中。为了清楚地说明电子硬件和软件的可互换性，在上述说明中已经按照功能一般性地描述了各示例的组成及步骤。这些功能究竟以电子硬件还是软件方式来执行，取决于技术方案的特定应用和设计约束条件。本领域技术人员可以对每个特定的应用来使用不同方法来实现所描述的功能，但是这种实现不应认为超出本发明的范围。
术语“第一”、“第二”等是用于区别类似的对象，而不是用于描述或表示特定的顺序或先后次序。
术语“包括”或者任何其它类似用语旨在涵盖非排他性的包含，从而使得包括一系列要素的过程、方法、物品或者设备/装置不仅包括那些要素，而且还包括没有明确列出的其它要素，或者还包括这些过程、方法、物品或者设备/装置所固有的要素。
至此，已经结合附图所示的优选实施方式描述了本发明的技术方案，但是，本领域技术人员容易理解的是，本发明的保护范围显然不局限于这些具体实施方式。在不偏离本发明的原理的前提下，本领域技术人员可以对相关技术特征做出等同的更改或替换，这些更改或替换之后的技术方案都将落入本发明的保护范围之内。

Claims (10)

1.一种基于近似直立扫描点云快速配准的三维重建方法，其特征在于，该方法包括：
步骤S10，获取近似直立扫描点云作为三维重建的待配准的目标点云Ωt和原始点云Ωs；
步骤S20，分别计算所述目标点云Ωt和原始点云Ωs的轴向包围盒AABBt和AABBs，并基于轴向包围盒AABBt和AABBs的尺寸进行所述待目标点云和原始点云的体素化，获得目标体素点云Vt和原始体素点云Vs；
步骤S30，分别计算所述目标体素点云Vt和原始体素点云Vs的几何中心点Ct和Cs，根据扫描的先验信息确定先验可靠区域，并进行所述目标体素点云和原始体素点云的初始配准，获得初始配准变换矩阵；
步骤S40，通过基于kd树的ICP方法进行初始配准后的目标体素点云和原始体素点云的精确配准，获得精确配准变换矩阵；
步骤S50，将所述初始配准变换矩阵和精确配准变换矩阵复合，获得最终配准变换矩阵，并基于所述最终配准变换矩阵进行所述目标点云Ωt和原始点云Ωs的配准，并基于配准后的点云进行三维重建。
2.根据权利要求1所述的基于近似直立扫描点云快速配准的三维重建方法，其特征在于，步骤S20包括：
步骤S21，分别计算所述目标点云Ωt和原始点云Ωs的三维坐标最大值和最小值，获得目标点云和原始点云的轴向包围盒AABBt和AABBs；
步骤S22，分别基于所述目标点云和原始点云的轴向包围盒AABBt和AABBs的尺寸以及设定的轴向最优体素数Δ和最小体素数δ，获得所述目标点云和原始点云对应的单个体素边长
Figure FDA0002833915710000011
和
Figure FDA0002833915710000012
步骤S23，基于所述目标点云和原始点云对应的单个体素边长
Figure FDA0002833915710000013
和
Figure FDA0002833915710000014
进行所述待目标点云Ωt和原始点云Ωs的体素化，获得目标体素点云Vt和原始体素点云Vs。
3.根据权利要求2所述的基于近似直立扫描点云快速配准的三维重建方法，其特征在于，所述目标点云Ωt和原始点云Ωs的轴向包围盒AABBt和AABBs为：
Figure FDA0002833915710000021
Figure FDA0002833915710000022
其中，
Figure FDA0002833915710000023
和
Figure FDA0002833915710000024
分别代表目标点云中所有点的三维坐标最大值和最小值，
Figure FDA0002833915710000025
和
Figure FDA0002833915710000026
分别代表原始点云中所有点的三维坐标最大值和最小值。
4.根据权利要求3所述的基于近似直立扫描点云快速配准的三维重建方法，其特征在于，步骤S22包括：
步骤S221，分别计算AABBt的最大轴向长
Figure FDA0002833915710000027
和最小轴向长
Figure FDA0002833915710000028
以及AABBs的最大轴向长
Figure FDA0002833915710000029
和最小轴向长
Figure FDA00028339157100000210
Figure FDA00028339157100000211
Figure FDA00028339157100000212
Figure FDA00028339157100000213
Figure FDA00028339157100000214
其中，max代表求最大值操作，min代表求最小值操作；
步骤S222，基于所述AABBt的最大轴向长
Figure FDA00028339157100000215
和最小轴向长
Figure FDA00028339157100000216
以及AABBs的最大轴向长
Figure FDA00028339157100000217
和最小轴向长
Figure FDA00028339157100000218
设定的轴向最优体素数Δ和最小体素数δ，分别计算目标点云和原始点云的单个体素边长
Figure FDA00028339157100000219
和
Figure FDA00028339157100000220
Figure FDA00028339157100000221
其中，
Figure FDA0002833915710000031
代表目标点云和原始点云的最小边长，
Figure FDA0002833915710000032
代表目标点云和原始点云对应的初始体素边长，
Figure FDA0002833915710000033
和代表目标点云和原始点云的最小边长以长度L0为体素边长的体素数。
5.根据权利要求4所述的基于近似直立扫描点云快速配准的三维重建方法，其特征在于，步骤S23包括：
分别将所述目标点云和原始点云分割为Wt×Ht×Dt和Ws×Hs×Ds的方块：
Figure FDA0002833915710000034
Figure FDA0002833915710000035
Figure FDA0002833915710000036
分别基于目标点云中第i个点
Figure FDA0002833915710000037
和原始点云中第i个点
Figure FDA0002833915710000038
计算对应的体素坐标
Figure FDA0002833915710000039
和
Figure FDA00028339157100000310
Figure FDA00028339157100000311
Figure FDA00028339157100000312
Figure FDA00028339157100000313
其中，
Figure FDA00028339157100000314
代表下取整操作；
以分割获得的坐标为体素坐标的点云作为目标体素点云Vt和原始体素点云Vs。
6.根据权利要求1所述的基于近似直立扫描点云快速配准的三维重建方法，其特征在于，步骤S30包括：
步骤S31，分别计算所述目标体素点云Vt和原始体素点云Vs的几何中心点Ct和Cs；以所述目标点云Ωt和原始点云Ωs的密度大于设定阈值或距离获取点云的扫描仪所在位置设定距离范围内的区域作为先验可靠区域；
步骤S32，获取先验可靠区域中Vt和Vs中心化后的定位向量的均值方向
Figure FDA0002833915710000042
和
Figure FDA0002833915710000043
的平均夹角θ，并计算旋转变换矩阵Tr：
Figure FDA0002833915710000041
步骤S33，基于所述目标体素点云Vt和原始体素点云Vs的几何中心点Ct和Cs以及旋转变换矩阵Tr，获取初始配准变换矩阵：
V′t＝Tr·(Vt-Ct)，Vs′＝Is·(Vs-Cs)
其中，Is为单位变换矩阵，V′t和Vs′为初始配准后的目标体素点云和原始体素点云。
7.根据权利要求6所述的基于近似直立扫描点云快速配准的三维重建方法，其特征在于，步骤S40包括：
步骤S41，通过KDTree在Vs′中搜索V′t中任意一点的最近邻域点，获取V′t中点集Vtm在Vs′上的映射点集Vsm；
步骤S42，基于预先设定的最小阈值距离D删除点集Vtm和点集Vsm中的错误映射点对，获得点集V′tm和映射点集V′sm；
步骤S43，利用四元素法计算点集V′tm和映射点V′sm的变换矩阵R和平移量T，并基于变换矩阵R和平移量T进行点集V′tm的变换，获得点集V″tm；
步骤S44，重复执行步骤S41-步骤S43，直至获取点集Vtm中每一点到映射点集Vsm中的映射点对和对应的变换矩阵R、平移量T；
步骤S45，迭代计算直至最后的变换矩阵R、平移量T配准后，对应点对坐标间差值小于设定阈值，获得精确配准变换矩阵：
V′t＝TRicp·Vs′+TTicp
其中，TRicp和TTicp分别为ICP方法生成的旋转矩阵和平移矩阵。
8.根据权利要求7所述的基于近似直立扫描点云快速配准的三维重建方法，其特征在于，所述最终配准变换矩阵，其公式表示为：
Vt＝Rfull·Vs+Tfull
其中，
Figure FDA0002833915710000051
9.根据权利要求8所述的基于近似直立扫描点云快速配准的三维重建方法，其特征在于，步骤S50中基于所述最终配准变换矩阵进行所述目标点云Ωt和原始点云Ωs的配准，其方法为：
Ωt＝Rfull·Ωs+Tfull。
10.一种基于近似直立扫描点云快速配准的三维重建系统，其特征在于，该系统包括以下模块：
点云获取模块，配置为获取近似直立扫描点云作为三维重建的待配准的目标点云Ωt和原始点云Ωs；
体素化模块，配置为分别计算所述目标点云Ωt和原始点云Ωs的轴向包围盒AABBt和AABBs，并基于轴向包围盒AABBt和AABBs的尺寸进行所述待目标点云和原始点云的体素化，获得目标体素点云Vt和原始体素点云Vs；
初始配准模块，配置为分别计算所述目标体素点云Vt和原始体素点云Vs的几何中心点Ct和Cs，根据扫描的先验信息确定先验可靠区域，并进行所述目标体素点云和原始体素点云的初始配准，获得初始配准变换矩阵；
精确配准模块，配置为通过基于kd树的ICP方法进行初始配准后的目标体素点云和原始体素点云的精确配准，获得精确配准变换矩阵；
点云配准模块，配置为将所述初始配准变换矩阵和精确配准变换矩阵复合，获得最终配准变换矩阵，并基于所述最终配准变换矩阵进行所述目标点云Ωt和原始点云Ωs的配准；
重建模块，配置为基于配准后的点云进行三维重建。
CN202011471438.XA 2020-12-14 2020-12-14 基于近似直立扫描点云快速配准的三维重建方法及系统 Active CN112581511B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202011471438.XA CN112581511B (zh) 2020-12-14 2020-12-14 基于近似直立扫描点云快速配准的三维重建方法及系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202011471438.XA CN112581511B (zh) 2020-12-14 2020-12-14 基于近似直立扫描点云快速配准的三维重建方法及系统

Publications (2)

Publication Number Publication Date
CN112581511A true CN112581511A (zh) 2021-03-30
CN112581511B CN112581511B (zh) 2023-07-18

Family

ID=75135169

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202011471438.XA Active CN112581511B (zh) 2020-12-14 2020-12-14 基于近似直立扫描点云快速配准的三维重建方法及系统

Country Status (1)

Country Link
CN (1) CN112581511B (zh)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115222787A (zh) * 2022-09-20 2022-10-21 天津中科智能技术研究院有限公司 基于混合检索的实时点云配准方法
CN117268799A (zh) * 2023-11-23 2023-12-22 靖江市黄埔机械部件制造有限公司 一种基于负载转矩的底盘制动器检测方法及系统

Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101533529A (zh) * 2009-01-23 2009-09-16 北京建筑工程学院 基于深度图像的三维空间数据处理方法与装置
CN106803267A (zh) * 2017-01-10 2017-06-06 西安电子科技大学 基于Kinect的室内场景三维重建方法
CN107886529A (zh) * 2017-12-06 2018-04-06 重庆理工大学 一种用于三维重建的点云配准方法
CN109559340A (zh) * 2018-11-29 2019-04-02 东北大学 一种并行的三维点云数据自动化配准方法
CN109697729A (zh) * 2018-12-26 2019-04-30 中国科学院大学 基于平面多边形匹配的3d岩体点云配准方法
US20190295266A1 (en) * 2016-05-20 2019-09-26 Nokia Technologies Oy Point cloud matching method

Patent Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101533529A (zh) * 2009-01-23 2009-09-16 北京建筑工程学院 基于深度图像的三维空间数据处理方法与装置
US20190295266A1 (en) * 2016-05-20 2019-09-26 Nokia Technologies Oy Point cloud matching method
CN106803267A (zh) * 2017-01-10 2017-06-06 西安电子科技大学 基于Kinect的室内场景三维重建方法
CN107886529A (zh) * 2017-12-06 2018-04-06 重庆理工大学 一种用于三维重建的点云配准方法
CN109559340A (zh) * 2018-11-29 2019-04-02 东北大学 一种并行的三维点云数据自动化配准方法
CN109697729A (zh) * 2018-12-26 2019-04-30 中国科学院大学 基于平面多边形匹配的3d岩体点云配准方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
李红军 等: "局部形状特征概率混合的半自动三维点云分类", 《浙江大学学报(理学版)》, no. 01, pages 4 - 12 *
王欢 等: "利用Kinect深度信息的三维点云配准方法研究", 《计算机工程与应用》, vol. 52, no. 12, pages 153 - 157 *

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115222787A (zh) * 2022-09-20 2022-10-21 天津中科智能技术研究院有限公司 基于混合检索的实时点云配准方法
CN117268799A (zh) * 2023-11-23 2023-12-22 靖江市黄埔机械部件制造有限公司 一种基于负载转矩的底盘制动器检测方法及系统
CN117268799B (zh) * 2023-11-23 2024-03-08 靖江市黄埔机械部件制造有限公司 一种基于负载转矩的底盘制动器检测方法及系统

Also Published As

Publication number Publication date
CN112581511B (zh) 2023-07-18

Similar Documents

Publication Publication Date Title
CN109872397B (zh) 一种基于多目立体视觉的飞机零件的三维重建方法
CN106709947B (zh) 一种基于rgbd相机的三维人体快速建模系统
Brenner et al. Coarse orientation of terrestrial laser scans in urban environments
Oehler et al. Efficient multi-resolution plane segmentation of 3D point clouds
Cheng et al. 3D building model reconstruction from multi-view aerial imagery and lidar data
EP2249311B1 (en) Systems and methods for extracting planar features, matching the planar features, and estimating motion from the planar features
CN113178009B (zh) 一种利用点云分割和网格修补的室内三维重建方法
Tazir et al. CICP: Cluster Iterative Closest Point for sparse–dense point cloud registration
WO2018061010A1 (en) Point cloud transforming in large-scale urban modelling
Cheng et al. Building boundary extraction from high resolution imagery and lidar data
CN110838115B (zh) 运用轮廓线提取和四维曲面拟合的古文物三维模型变化检测方法
CN111612728B (zh) 一种基于双目rgb图像的3d点云稠密化方法和装置
Galvanin et al. Extraction of building roof contours from LiDAR data using a Markov-random-field-based approach
CN115372989A (zh) 基于激光雷达的越野自动小车长距离实时定位系统及方法
Wei et al. Automatic coarse registration of point clouds using plane contour shape descriptor and topological graph voting
CN112651944A (zh) 基于cad模型的3c部件高精度六维位姿估计方法及系统
CN112581511B (zh) 基于近似直立扫描点云快速配准的三维重建方法及系统
Drost et al. Local hough transform for 3d primitive detection
CN113177974A (zh) 一种点云配准方法、装置、电子设备及存储介质
Wang Automatic extraction of building outline from high resolution aerial imagery
Hyyppä et al. Efficient coarse registration method using translation-and rotation-invariant local descriptors towards fully automated forest inventory
Tong et al. 3D point cloud initial registration using surface curvature and SURF matching
Junior et al. A new variant of the ICP algorithm for pairwise 3D point cloud registration
JP2023530449A (ja) 空中と地上の位置合わせのためのシステムおよび方法
Elkhrachy Feature extraction of laser scan data based on geometric properties

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
GR01 Patent grant
GR01 Patent grant