造船中点云数据的高效处理与配准方法

计算设计与工程学报。号13（2014）202~212www.jcde.org造船中点云数据的高效处理：改进的构件提取方法和配准方法Jingyu Sun1，*，Kazuo Hiekata1，Hiroyuki Yamato1，Norito Nakagaki2和AkiyoshiSugawara21东京大学研究生院边疆科学研究科，千叶县柏野叶5-1-5环境大楼，邮编277-85632住友重工海洋&工程公司，日本神奈川县横须贺市夏岛町19号邮编：237-0061（2014年2月27日接收;2014年5月9日修订;2014年5月12日接受摘要为了在当前的造船业中生存，对造船厂来说，在大多数制造步骤中有效地评估船舶部件的精度是至关重要的通过比较激光扫描仪扫描的每个部件的点云数据和CAD中格式化的船舶设计数据来评估部件针对这些问题，本文提出了改进的点云数据处理方法点云数据的K-d树构造是对每个点的邻域搜索对种子点的邻域点进行区域生长，提取出构件的连续部分，在连续部分的边缘进行曲面拟合和B样条曲线拟合，识别出同一构件的邻域，并以障碍ICP（迭代最近点）算法在通过主成分分析确定适当的配准方向后对两组数据进行配准通过在船厂的实验，从扫描点云数据中提取了200块曲面壳板，并利用所提出的方法将其与设计的CAD数据进行了配准，进行了精度评价。结果表明，本文提出的方法有效地支持了以点云数据处理为目标的精度评价。关键词：点云;区域生长法; B样条曲线; ICP; K-d树1. 介绍多年来，有效评估船舶组件的精度一直被认为是造船生产规划和控制的必要组成部分，而使用非接触式3D激光扫描仪扫描的点云数据只是最近才被考虑[1]。作为激光扫描仪在造船中的应用，东京大学开发了船舶部件精度评估系统其步骤如下：(Step1）通过使用非接触式3D激光扫描仪扫描部件来获得点云数据(Step 2）利用基本区域生长法从点云数据(Step3）直接利用ICP算法对零件的测量点云数据和设计数据进行配准*通讯作者。联系电话：+81-47-136-4626，传真：+81-47-136-4626电子邮件地址：sun@is.k.u-tokyo.ac.jp© CAD/CAM工程师协会Techno-Press doi：10.7315/JCDE. 2014. 020（步骤4）通过比较两个数据的位置找到位移误差，并将其可视化为精度评估结果。这里存在两个至关重要的问题(a) 在步骤2中，由于制造车间内生的不规则障碍物及其阴影将部件的测量点云数据划分为几个分离的域，因此部件的部分不能被有效地提取手工提取和整合这些分离的域浪费了大量的时间，是实际应用的障碍。用于从噪声点云计算光滑表面的众所周知的方法（诸如移动最小二乘（MLS）投影[4]）具有相同的问题，即当组件被障碍物的阴影分离时，(b) 在步骤3中，使用ICP（迭代最近点）算法由于ICP算法的特点，在没有明确配准方向预设的情况下，J. Sun等人计算设计与工程杂志卷。号13（2014）202~212203图1.成分提取。图2.区域生长法集合可能在错误的方向上行进并导致不正确的配准结果。因此，现有的点云处理方法（基本区域生长法和ICP配准法）如不加以改进和适当的组合，很难满足船舶建造中的精度评价指标针对造船中构件精度评估问题，提出了一种改进的构件提取方法和一种改进的构件配准方法首先，对扫描点云数据进行k-d树分解，加快点云数据的邻域搜索;为了从点云数据中提取零件的连续区域，对种子点的每个相邻点执行区域生长方法然后对提取区域的每个边缘点进行曲面拟合和B样条曲线拟合，识别出由障碍物阴影从提取区域中分离出来的邻域通过重复上述步骤可以提取整个组分。在注册提取的点云数据通过对两种数据进行配准方向分析，确定合适的配准方向。然后，申请ICP进行以下注册。2. 分量提取在造船中，根据工厂的条件，由于以下因素，激光扫描仪对部件的测量结果（见图1左图）无法满足精度评估要求(1) 测量得到的点云数据往往包含大量不必要的噪声，如地板、工人以及一些制造过程中必不可少的模板--木模板等。(2) 测量点云数据通常被这些障碍物的阴影划分手动提取和集成这些分离的域会浪费大量时间。A.测量的点云数据通常包含大量的点，这大大降低了提取过程的速度。204J. Sun等人计算设计与工程杂志卷。号13（2014）202~212∑mi∑mn-1i 我我我我 我图3.利用主成分分析进行加权平面拟合。最好是几分钟。作为一个关键的预处理，所有的噪声应该被去除，属于同一个组件的分离域应该被有效地提取，如图所示1.本文中，组件提取包括两个主要过程：连续域提取和分离域识别，连续域提取是从单个种子点提取小的分离域，分离域识别是判断提取的域是否来自于其他域。2.1 连续域提取本文对前人工作中提出的连续区域提取方法（基本区域生长法）进行了改进，使其提取效率更高，速度更快。2.1.1 在前期工作中提出的方法-区域生长法和k-d树构造法所提出的分量提取方法是基于在先前的工作中提出的称为区域生长方法[5]的基本点云处理方法，其重复计算种子点的每个相邻点处的法向量、提取具有与种子点相似的法向量的相邻点以及将所提取的相邻点设置为下一个种子点的过程。图2显示了区域增长方法的主要流程。当法向量发生相对剧烈的变化时，生长过程结束确定P4是否与Pi（i = 1，2，3）属于同一区域，对其进行法向量计算。由于Pi（i = 5，6，7）处的密度低于P i（i = 1，2，3），将输出无加权平面的法向量，如图3所示。显然，它不能表示P4处的法向量，并且与P2处的法向量有很大的不同。因此，当P4实际上属于与Pi（i = 1，2，3）相同的区域时，生长在不包括P4的情况下终止。因此，为了计算每个相邻点和种子点处的法向量，本文进行加权平面拟合[7]。点云G的重心由方程计算，然后可以使用等式（1）获得方差-协方差矩阵V。（二）、 最小特征向量n是使用V作为该点的法向量计算的。d是每个点与正在计算的点之间的距离，h是点云点的平均间隔 由于权重m i是针对每个点计算的，因此当拟合平面时，更靠近被计算的点的点具有更大的影响，因此可以输出加权平面的法向量，如图3所示。在P2处与法向量有可接受的微小差异的情况下，P4可以被认为属于作为Pi（i = 1，2，3）的域。n-1G =1拉姆P ，P= （ x， y， z） T（ 1）i=0i=0n-11V =μ m（P-G）×（P-G）为了使邻居搜索更有效，空间-分区数据结构，k-d树[6]，是使用n-1i=0时我我我i=0时D2所扫描的点云数据用于递归地组织k维空间中的点。2.1.2 考虑点云当点云密度变化剧烈时，利用种子点的邻域点计算种子点的法向量，采用所有邻域点均相等的平面拟合方法，可能会如图3所示，mi=exp（-h2）（2）2.1.3 用于平面拟合和增长入口候选的点选择为了避免不适当地生长成离种子太远的点，例如元件附近的地板，用于进行平面拟合的点必须在某些约束下认真选择。之间的距离d不J. Sun等人计算设计与工程杂志卷。号13（2014）202~212205图4.用于进行平面拟合的点。图5.区域增长切入点选择。点被评估，并且检查提取的点以决定该点是否应该被认为是域的一部分。 在本文中，如图4所示，两个平面分别适合于法向量计算和距离d的计算。点P j是被评估的点。 由于plane1用于计算点Pj的法向量，因此使用所有相邻点R1来拟合它，而不关心它们是否已经被识别为域的一部分。法向量比较只能评估Pj处的方向趋势。距离d也应该被考虑以避免增长到不是域的一部分但在法向量比较结果中仅具有轻微差异的区域用于计算距离d的Plane2仅使用已经被识别为域的一部分的点R2来计算为了节省计算时间，仅将相对地在邻居之外的点视为下一个增长条目候选。因此，如图5所示，加权重力点G使用等式（1）计算（1）在第2.1.2节中，并且针对相邻点的每个点Pi计算距离di。只有当ndi> bthresh×bthresh dk/n1点Pi被认为是下一个增长的入口。bthresh是0和1之间的阈值。在本文中，bthresh被设置为0.74。2.2 分离域识别由于上面讨论的方法只能提取连续的域，为了提取被障碍物阴影不规则地分割的整个组件，本文提出了分离域识别方法，以节省这一低效的步骤。2.2.1 单种子分离域图6.弯曲的壳板（左）和块206J. Sun等人计算设计与工程杂志卷。号13（2014）202~212JJ图7.组分浸提流程（单个种子）。这里进行公共域判断，以确定是否有任何相邻域与使用上述平坦区域生长方法首先提取的连续域具有相同的分量图7显示了这个公共域判断过程的流程首先，在分量上随机选择一个单个种子P种子，然后对连续域进行组合，通过上述区域生长方法计算。第二步以灰色突出显示，其处理过程与第2.1节中讨论的过程相同。为所提取的连续域的边界P边界的每个点计算相邻点R neigGbor。然后在R个邻域中求出点Di的邻域点R个邻域i。对于RneigGborsi中的所有点Di，图8.共同领域判断。J. Sun等人计算设计与工程杂志卷。号13（2014）202~212207JJJJJJJ图9.拟合曲线和实测曲线。执行将在后面介绍的判断，以决定是否可以选择Ni作为该分量的另一个域的种子点。此外，递归地重复该过程，直到成功提取组件的所有点。图8说明了公共域判断是如何工作的。结构域A首先通过第2.1节中的方法识别。对于区域A的每一个边缘点，进行四次曲面拟合和广泛用于船舶自由曲面设计的B样条曲线拟合[8]，以确定是否点D i的判断度Y i 是true（1）还是false其中，j =0，(3)用于拟合所提取的域A的边缘点周围的曲面。检查Di是否属于拟合的曲面S。如果只有Di通过曲面检查，则使用等式2进行B样条曲线检查（4）和（5）。这里Pi是控制点。一段的n次B样条曲线表示为Eq.（四）、bi，n是B样条曲线的基函数t是m个实数，称为纽结这里.J J（0）。只有当所有Yi的和大于阈值时，才选择点Ni作为种子点以开始另一个域（域B）的如图9所示，取四个点（Di，Pedge，Pm1，Pm2）作为B样条线的控制点 P extend1和P extend2 是两个点的延长线上的|，Pm 1，Pm 2是它们在域中的最近 点|, and P m1, Pm2are their nearest points in domain∑N[z-∑4（∑ jCxiy2 =0，A. 计算出拟合后的B样条曲线C，k=1 Kj=0i=0时我，我KK ）的情况下]Pi =（ xi，yi，zi）T（3）m-n-2S（ t） = <$Pi bi，n（ t），t ∈[ tn，ti=0时对于C的拟合点，获得域B中的点。 这些最近的点的组是C测量的。只有当每一对拟合的C与实测的C之间的距离的平均值小于阈值时，Di的判断度Yi才为真。bj，0=1iftj ≤ t ≤ tj+12.2.2 边缘提取其中j = 0，如果t> tj+1或t tj，则bj，0=0，其中 j =0，如果R个邻域中的所有邻域点（如图7所示）都被认为是下一个增长的入口，那么对已经被识别为域点的点的不必要的计算将花费大量的时间。因此，只有在R邻域的边缘处的点用于进行以下计算：B （t） =t- t jb（t）+tj+n+1- tb（t）lation.j，ntj+nj，n-1tj+n+1j+1，n-1为了在此过程中提取边缘，如图10所示图10.边缘检测208J. Sun等人计算设计与工程杂志卷。号13（2014）202~212图11.组件注册流程。使用主成分分析来拟合平面用于进行平面拟合的点是通过k-NN方法搜索的种子点N的相邻点然后将相邻点投影到平面上。如果不是所有四个象限都有相邻点，如图10的右图所示，则点N被认为是边缘[9]。3. 组件注册在船舶制造中，构件设计数据通常采用ACIS的标准文件格式SAT表示。在SAT文件中，曲面和曲线表示为NURBS（非均匀有理B样条）。为了评估零件的精度或提供任何进一步的制造建议，作为一个重要的预处理，设计数据和提取的零件的点云数据之间的配准如图11所示，本文提出了一种高效的配准方法，该方法涉及使用两个数据集的重心、配准方向的预先设置和传统ICP配准的并行变换，传统ICP配准是一种用于最小化两个点云之间的差异的算法3.1 利用重心的并行变换为了将测量点云与设计数据配准，使用两个数据集的重心进行并行变换，使用以下等式。这里，P表示来自测量数据的点。Gd和Gm分别为设计数据和实测数据的重心PY = P +（Gd-Gm）（6）3.2 ICP算法配准在使用重心的并行变换之后，进行用于点云和设计数据的配准的基本ICP算法。设计曲面S（ u，v） =（x（ u，v），y（ u，v），z（ u，v））是一个有理方程图12.方法A的流程图使用B样条基函数表示非恒定间距节点的均匀坐标系投影结果[10]。r（ u，v） = S（ u，v）-P（7）r（ u，v） ≠ Su（ u，v） = 0r（ u，v） ≠ Sv（ u，v） = 0如Eq. (7)投影矢量可以通过从设计曲面S（u，v）中减去测量点P所以，r（ u，v）满足方程。（8）和（9）。这里，Su（u，v）和Sv（u，v）分别是所设计曲面的u方向如图12所示的投影矢量r（ u，v）可以通过求解这三个方程来计算。这里Q是投影点。接下来，使用运动矢量t和旋转矩阵R进行测量点云的坐标变换，以减少设计曲面和测量点云的位移[11，12]。Y在Eq。（10）是变换后的测量点。为了估计R和t，在等式（1）中表示的位移E（ R，t）被估计为：(11)最小化[13]。重复该过程，直到E（ R，t）小于阈值[14]。图13.反转配准结果。J. Sun等人计算设计与工程杂志卷。号13（2014）202~212209我n-12016年10 月25 日，|Pi|2-x2）n-1- mi xi yin-1- mi xi zi1III=0n-1我i=0时n-1我i=0时n-1I =I- m x y2016年10月20日，|P |2-y2）- m y zI我Ii=0In-1我我我i=0时我我n-1我n-1i=0时我我我我我我I-mi xi ziLi=0- -i=0时2016年10月22日，Yellow I（|P i|2-z 2）Ii=0）Tmi=exp（-h2），Pi=xi，yi，zi（十二）PY = RP + t（10）n4. 实验结果在本节中，在实验室中进行了一组实验E（R，t）=EOi-RPi-ti=13.3配准方向预设（十一）造船厂评估性能的点云处理方法提出了上述。本文中的所有实例都是在UNICUS公司提供的C#开发的三维点云系统中公司. 用户可以显式更改在实践中，如图13所示，如果配准从不正确的方向开始，则第3.2节中讨论的配准方法很有可能得到明显错误的输出针对这一问题，在进行ICP算法之前，先对惯性主轴的三个方向进行计算，找出位移最小的配准起始方向。首先，惯性张量I使用等式（1）计算（十二）、I的三个特征向量是惯性主轴，如图14所示如Eq. (13)根据多种方式将测量点云P配准到设计数据选择位移最小的配准方向作为ICP算法的初始配准方向其中Rd和Rm分别是设计数据和测量数据的惯性张量I的特征向量PY = RdT × Rm × PRd =（ed1，ed2， ed3），Rm =（aem1，bem2，cem3）a ={1，-1}，b ={ 1，-1}，c ={ 1，-1}（13）可以使用上述参数来调整处理精度和效率之间的平衡所有的点云数据都是由FARO Focus 3D激光扫描仪测量的[15]。4.1 组分浸提结果对200多块弯曲壳板进行了塑性变形后的测量、提取下面给出三个关于组件提取的示例。表1总结了先前方法和新方法之间的比较。实施例（a）. 如图15（a）所示，弯曲的壳板被测量，有一些障碍物，如横幅它有大约168万个点。在2.1节中讨论的算法被用来处理这种情况下，其中的组件没有被分割的障碍如图所示的提取结果，所提出的改进的区域生长算法被证明能够克服这些障碍。在这种情况下，用于提取的时间约为64秒。实施例（b）. 如图15（b）所示，另一个弯曲壳板上有一些木制模板。它有大约134万个点，由于木制模板的阴影，它被分成多个在这图14.惯性主轴。我210J. Sun等人计算设计与工程杂志卷。号13（2014）202~212表1.单区域生长方法（b）与新提出的方法（a）和（b）之间的比较。组件分为手动/自动执行时间（一）1自动64s（b）第（1）款7手动120s（c）第（1）款7自动53s图15.弯壳板浸提结果：（a）使用第2.1节中的算法浸提，（b）使用第2.1.1节中先前工作的算法浸提，（c）使用第2.2节中的算法浸提。例如，使用第2.1.1节中介绍的先前工作的手动提取方法来提取组件的多个域。为每个域手动选择一个生长种子点，其3D坐标如图中的方框所示如图所示的提取结果，在这种情况下评估了具有多个种子的手动分离域提取在这种情况下，用于提取的时间（包括人工积分操作）超过120秒。实施例（c）.如图15（c）所示，相同的弯曲使用第2.2节中提出的算法提取作为示例（b）的壳板。预先只选取一个单点，然后根据曲面拟合和B样条直线拟合（图中以白线表示），自动识别多个域在这种情况下，用于提取的时间仅为约53秒。该算法不仅可以从测量点云中提取曲面壳板，还可以应用于其它构件，如壁板或块体图16示出了根据本发明的方法的提取结果。图16.样本组提取结果。表2.先前工作（a）和（b）中使用的基本方法与新提出的方法（c）之间的比较。注册结果执行时间（一）失败4s（b）第（1）款失败43s（c）第（1）款成功30sJ. Sun等人计算设计与工程杂志卷。号13（2014）202~212211图17.曲面壳板配准结果：（1）测量点云和设计数据，（a）使用重心的配准，（b）无方向预设的配准，（c）有方向预设的配准。使用第2.2节中的算法的面板组件4.2 组件配准结果下面给出了一个曲面壳板如图17（1）所示，曲面壳板的测量点云和设计数据显示在计算机上。可见，他们原本是完全不同的职位，需要登记的。表2总结了执行时间结果（a）。首先按3.1.1节所讨论的方法进行只利用重心三维坐标的并行变换。图17（a）显示了该步骤的转换结果结果（b）。 图17（b）示出了直接从结果（a）中的位置进行ICP算法而没有方向预设的结果。结果是注册的方向相反显然，由于大位移，结果（b）不能满足下面的分量评估处理结果（c）。在执行ICP算法之前，对配准方向的多种方式进行测试。选择具有最小位移的配准方向来为ICP算法进行方向预设置如图17（c）所示，通过配准方向测试和预设，配准结果完全正确，并证明能够进行后续处理。5. 结论本文提出了一套有效的造船工程中的点云处理方法，包括改进的构件提取方法和构件配准方法。该方法解决了造船厂直接使用基本区域生长法和ICP法处理点云时存在的两个问题，并通过造船厂的实验证明了该方法的有效性通过本文提出的改进的构件提取方法和配准方法对点云数据进行有效的处理，将激光扫描仪扫描的各构件点云数据与CAD格式的船舶设计数据进行比较的构件精度评估方法今后还将对其他类型的船舶构件进行更多的试验，根据实际需要和每种方法的有效性，所提出的改进方法可以有不同的本文讨论的每种方法的参数也将针对不同的组分进行优化。致谢本手稿是由NIPPON KAIJI KYOKAI（NK级）支持的联合研究的成果。作者要感谢UNICUS公司，有限公司、FARO Japan，Inc.对项目的贡献。引用[1] 张文，张文，等.激光扫描技术在船舶制造中的应用.中国船 舶 工 业 杂 志 ， 2000 ， 24 （ 1 ） ： 117 - 118. In ：Proceedings of ISPRS Workshop on Laser Scanning 2007and SilviLaser 2007; 2007 Sep 12-14; Espoo，Finland; p.56比61[2] [10]张文，张文.用激光扫描仪确定曲面壳板的精度评定系统（第二次报告）日本造船师与海洋工程师学会会刊2013; 17：169-176.[3] 张文忠，张文忠.曲面壳板激光扫描精度评定系统的研制与实例分析。船舶生产与设计杂志2011; 27（2）：84-90.[4] 作者：Pauly M，Gross M，Kobbelt L.点采样曲面的有效简化在：IEEE可视化会议论文集; 2002年10月27日至11月1日163-170。[5] Forsyth D，Ponce J.计算机视觉：现代方法。上鞍河（新泽西州）：普伦蒂斯霍尔专业技术参考; 2007年。[6] Gross M，Pfister H.基于点的图形。Burlington（MA）：Morgan Kaufmann; 2007.[7] 放大图片作者：Weingarten JW，Gruener G.三维概率平面拟合及其在机器人测绘中的应用在：Proceedings IEEE国际会议机器人和自动化（ICRA）; 2004年4月26日至5月1日;苏黎世，瑞士; p.927- 932[8] 张X，李宏，程志.通过法向截面曲率拟合的三维点云曲面212J. Sun等人计算设计与工程杂志卷。号13（2014）202~212In：Proceedings of AsiaGraph 2008; 2008 Oct 23-26; Tokyo，Japan; p.72比79[9] [10]杨晓波，王晓波，王晓波.从噪声点云中提取曲率线计算机辅助设计。2009; 41（4）：282-292.[10] Piegl L，Tiller W. NURBS的书。第2版，纽约：Springer;1997年。[11] [10]张文龙，王文龙，王文龙.基于分层M-ICP的快速高精度三维配准算法电子、信息和通信工程师协会，IEICE技术报告2004年9月10日; 8 p.号报告2004-CVIM-145。[12] Besl PJ，McKay ND.三维形状配准方法。IEEE模式分析与机器智能汇刊。1992; 14（2）：239-256.[13] 张 伟， 王 伟， 王 伟 .光 学和 数 字图 像 处 理。 Hoboken（NJ）：John Wiley Sons; 2011.[14] 张志自由曲线曲面配准的迭代点匹配法国际计算机视觉杂志1994; 13（2）：119-152.[15] FARO [Internet].FARO Laser Scanner Focus 3D; c2011 [引用日期2014年4月22日]可从以下网址获得：http://www.globalspec.com/datasheets/579/FARO/F7DD1C91-CE 09 - 4A 0A-AB 9 C-71 EC273 AC 1FA