CN111274732B - 一种基于“连接关系-位置”迭代优化的网格修复方法 - Google Patents

一种基于“连接关系-位置”迭代优化的网格修复方法 Download PDF

Info

Publication number
CN111274732B
CN111274732B CN202010209452.6A CN202010209452A CN111274732B CN 111274732 B CN111274732 B CN 111274732B CN 202010209452 A CN202010209452 A CN 202010209452A CN 111274732 B CN111274732 B CN 111274732B
Authority
CN
China
Prior art keywords
hole
connection relation
theta
edge
grid
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.)
Active
Application number
CN202010209452.6A
Other languages
English (en)
Other versions
CN111274732A (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 Technology and Business University
Original Assignee
Beijing Technology and Business 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 Technology and Business University filed Critical Beijing Technology and Business University
Publication of CN111274732A publication Critical patent/CN111274732A/zh
Application granted granted Critical
Publication of CN111274732B publication Critical patent/CN111274732B/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
    • 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

Landscapes

  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Computer Graphics (AREA)
  • Geometry (AREA)
  • Software Systems (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Complex Calculations (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

本发明涉及一种基于“连接关系‑位置”迭代优化的网格修复方法，属于孔洞修复技术领域。首先识别输入的三维网格模型的孔洞区域，利用动态规划方法初始化孔洞区域连接关系；然后识别孔洞边界一对特征点，并用两特征点及其法向拟合跨越孔洞的特征线，以特征线为导向调整孔洞区域连接关系，即连接关系优化；再基于孔洞及其邻域信息，构建基于三角面片表示的三维网格模型局部修复的变分框架，利用增值拉格拉日方法迭代求解孔洞及其邻域的顶点位置，即顶点位置优化；最后重复连接关系优化与顶点位置优化，直到不再发生连接关系调整，从而实现了保特征的三维网格模型修复。所提方法在孔洞修复过程中，在恢复缺失的显著特征方面，具有明显优势。

Description

一种基于“连接关系-位置”迭代优化的网格修复方法
技术领域
本发明涉及一种基于“连接关系-位置”迭代优化的网格修复方法，属于孔洞修复技术领域。
背景技术
随着三维建模技术和逆向工程的发展，三维物体表面重建技术广泛的应用于医学诊断，数字娱乐，虚拟现实，CAD等领域。逆向工程就是将现实世界的真实模型转换为计算机表达的数字几何模型，通过扫描设备获取三维点云数据，经过去噪、拼接、配准、曲面重建等操作得到三维数字几何模型。但往往由于测量设备和技术的局限性以及待测模型自身的缺陷和复杂性，导致获取的数据不完整，重建的网格存在孔洞。这些孔洞会影响后续基于三维网格模型的一系列操作，如有限元分析、曲面重建、快速成型等。因此三维网格模型修复至关重要。
如何最大限度的还原孔洞区域的特征是三维网格模型修复方法的重要衡量标准，现有方法可以较好的修复模型缺失的光滑特征，但是对于显著特征(折痕，角点等)的修复仍然不够理想。依据处理对象的不同大致分为两类修复方法：基于体表达的修复方法和基于网格的修复算法。基于体表达的修复方法首先将待修复的三维网格模型转换为体数据，然后在体空间采用不同的方法进行修复，最后将修复后的体数据转为三维网格模型。该类方法主要存在两个问题：一是网格与体数据之间的相互转换，会丢失特征或者产生拓扑错误。二是修复后网格面片数量多于原始网格且并且质量下降。基于网格的修复方法直接处理三维网格模型的孔洞区域，主要有基于曲面拟合的修复方法和基于网格剖分的修复方法。基于曲面拟合的方法首先采集孔洞边即及其邻域的信息，然后采用三角B样条、径向基函数等不同的曲面拟合方法，生成光滑的曲面覆盖孔洞区域，实现孔洞修复。该类方法不适合修补带有沟、岛等复杂的孔洞。基于网格剖分的修复方法直接对多边形孔洞进行网格剖分得到初始修复网格，之后再通过优化和光顺使网格密度与周围网格一致。该方法的缺点对于孔洞周边没有足够特征信息的残缺模型，难以恢复孔洞区域特征。而基于曲面的修复方法最大的不足是无法恢复孔洞区域的显著特征。
针对现有方法的不足，本发明提出了保特征的三维网格模型修复方法。首先识别孔洞区域，利用动态规划方法初始化孔洞区域连接关系；然后识别孔洞边界特征点，拟合跨越孔洞的特征线，以特征线为导向优化孔洞区域连接关系；基于孔洞及其邻域信息，构建三维网格模型局部修复变分框架，通过求解变分方法优化孔洞及其邻域顶点位置。最后，重复连接关系调整和顶点优化步骤，直到不再发生连接关系调整。针对缺失显著特征的孔洞区域，本发明提出的方法从网格连接关系与顶点位置两方面进行优化，从而有效的修复孔洞区域显著特征。
发明内容
本发明的目的是克服现有三维网格模型修复方法中存在的不能有效恢复孔洞区域显著特征以及优化算法复杂性的问题，提出了一种基于“连接关系-位置”迭代优化的网格修复方法。
本发明的核心思想为：首先识别输入的三维网格模型的孔洞区域，利用动态规划方法初始化孔洞区域连接关系；然后识别孔洞边界一对特征点，并用两特征点及其法向拟合跨越孔洞的特征线，以特征线为导向调整孔洞区域连接关系，即连接关系优化；再基于孔洞及其邻域信息，构建基于三角面片表示的三维网格模型局部修复的变分框架，利用增值拉格拉日方法迭代求解孔洞及其邻域的顶点位置，即顶点位置优化；最后重复连接关系优化与顶点位置优化，直到不再发生连接关系调整，从而实现了保特征的三维网格模型修复。
一种基于“连接关系-位置”迭代优化的网格修复方法，具体步骤如下：
步骤一、识别输入的三维网格模型的孔洞区域，利用动态规划方法初始化孔洞区域连接关系，得到具有完整连接关系的初始三维网格模型，即实现孔洞区域的连接关系修复；
其中，输入的三维网格模型记为M0，初始三维网格模型记为M；
步骤一，具体包括如下子步骤：
步骤1.1识别模型中的孔洞区域；
其中，识别孔洞区域的准则为：判定只有一个邻接三角面片的边为孔洞边界边，边界边首尾相连构成的闭合回路所包围的区域为孔洞区域；
步骤1.2采用动态规划方法对步骤1.1识别的孔洞区域进行三角化，并初始化孔洞区域的连接关系；
步骤二、记网格中发生调整的边的个数为number，初始值为0；
步骤三、识别孔洞边界一对特征点，基于特征点及其法向拟合跨越孔洞的特征线，以特征线为导向，调整孔洞区域连接关系，即连接关系优化，具体为：
步骤3.1识别孔洞边界一对特征点，基于特征点及其法向拟合跨越孔洞的特征线；计算孔洞边界顶点所有邻接边的二面角，若值大于θf，则该边为特征边，其邻接边界顶点即为特征点，其中，θf代表角度，取值范围在[30°,45°]，依据特征边确定该点的法向，从而得到孔洞一对边界特征点vx,vy以及其对应的法向n1,n2；然后利用三次贝塞尔曲线(1)粗略拟合孔洞区域特征线lB：
B(t)＝vx(1-t)+(n1+3vy-n2)t(1-t)2+(3vy-n2)t2(1-t)+vyt3,t∈[0,1] (1)
其中，t代表时间，B(t)代表t时刻贝塞尔曲线上一个坐标点；
步骤3.2以步骤3.1拟合的特征线为导向，调整孔洞区域连接关系，具体如下：
步骤3.2.1在特征线lB上选取采样点集合VB＝{VBi}i＝1,2...b，将点集VB向孔洞区域做投影，得到投影面片集合TB＝{ti}i＝1,2...b，投影面片上边构成投影边集合LB＝{lBi}i＝1,2...b，利用特征点vx,vy确定孔洞区域特征线大致方向为向量
Figure BDA0002422312050000021
的方向；
步骤3.2.2从特征点vx出发，遍历与vx相邻接投影三角形中的边；首先计算所有与vx相邻接边与
Figure BDA0002422312050000022
夹角大小，记录夹角最小的边的另一个顶点vtmp以及夹角值θ1；然后计算所有与vx相对的投影边若发生边交换后，形成的新边与
Figure BDA0002422312050000023
夹角大小，记录夹角最小的边对应的另一个顶点v'tmp以及夹角值θ2；若θ1≥θ2,则连接关系不变，将vtmp替换特征点vx；否则发生边交换，number＝number+1，并将v'tmp替换特征点vx；
步骤3.2.3判断vx与vy是否相等，若vx＝vy，则此次连接关系优化结束；否则重复步骤3.2.2；
步骤四、判断number取值是否为0，若number＞0，执行步骤五；否则，结束修复；
步骤五、基于孔洞及其邻域信息，构建基于三角面片表示的三维网格模型局部修复的变分框架，利用增值拉格拉日方法迭代求解孔洞及其邻域的顶点位置，即顶点位置优化；
其中，邻域为孔洞k环邻域，k取值范围在[3,10]；
步骤五，具体包括如下子步骤；
步骤5.1基于孔洞及其邻域信息，构建基于三角面片表示的三维网格模型局部修复的变分框架；
记孔洞区域及其邻域为孔洞局部，输入的模型M0的孔洞邻域网格为
Figure BDA0002422312050000031
顶点集合为
Figure BDA0002422312050000032
其中，m是孔洞邻域网格
Figure BDA0002422312050000033
中的顶点个数，孔洞区域为DH；修复后的网格M的孔洞局部网格为MH，顶点集合为V＝{v1,v2,...,vn}，边集合为E＝{e1,e2,...,ed}，边长度集合为L＝{l1,l2,...,ld}，内二面角集合为θ＝{θ1,θ2,...,θd}，其中，n(n＞m)是修复连接关系后的孔洞局部MH中的顶点个数，d是MH中边的个数，边集合E中的边ei的长度即为长度集合L中的li，内面角集合θ中的角度θi表示共享边ei的两个三角面片间的内二面角；保持孔洞局部连接关系不变，通过最小化如下(2)式的能量函数得到最优的孔洞局部顶点位置；
Figure BDA0002422312050000034
其中，
Figure BDA0002422312050000035
是求满足
Figure BDA0002422312050000036
最小的孔洞局部顶点位置，即MH；
Figure BDA0002422312050000037
为数据项Efidelity，使得修复后的网格MH尽可能地逼近给定的初始网格
Figure BDA0002422312050000038
||MH||TV为总变差正则项Eregularizer；λv(MH)为数据项截断参数，当v∈DH时，λv(MH)＝0；
Figure BDA00024223120500000314
时，λv(MH)＝λ＞0；
步骤5.1.1计算数据项，具体通过如下公式(3)计算：
Figure BDA0002422312050000039
其中，vi代表修复后网格孔洞局部MH的顶点集合，V＝{v1,v2,...,vn}中第i项；
Figure BDA00024223120500000310
代表输入的残缺网格孔洞邻域
Figure BDA00024223120500000311
的顶点集合，
Figure BDA00024223120500000312
中第i项；||PV-V0||1表示FV-V0的L1正则化；P表示m×n的投影矩阵，定义如(4)：
Figure BDA00024223120500000313
步骤5.1.2计算正则项,具体通过公式(5)计算：
Figure BDA0002422312050000041
其中，li代表孔洞局部MH中边ei的边长，即边长度集合l＝{l1,l2,...,ld}中第i项；θi表示ei的两个邻接三角面片间的内二面角，(π-θi)表示内二面角θi的补角；假设边ei的两个邻接三角面片为Δv1v3v4和Δv1v2v3，边的两端点为v1，v3，对应第三个顶点分别为v2和v4；定义T1和T2分别是三角面片Δv1v3v4和Δv1v2v3的外法向，长度均为||v1v3||，T1和T2间夹角即为边ei的二面角补角(π-θi)，则li|π-θ|＝||v1v3|||π-θ|即为T1和T2之间所夹弧T1T2的弧长；法向量T1和T2表示为(6)：
T1＝cot(θ4,1,3)(v4-v3)+cot(θ1,3,4)(v4-v1),
T2＝cot(θ2,3,1)(v1-v2)+cot(θ3,1,2)(v3-v2). (6)
其中，θ4,1,3是边v1v4和边v1v3的夹角，θ1,3,4是边v1v3和边v3v4的夹角，θ2,3,1是边v2v3和边v1v3的夹角，θ3,1,2是边v1v2和边v1v3的夹角；
采用近似计算的方法计算弧长，定义向量Tmid为指向弧T1T2中心，且与T1和T2同起点等长度的向量，则Tmid表示为：
Figure BDA0002422312050000042
于是用弦长之和近似表达T1T2弧长，即||T1T2||＝|T1-Tmid|+|Tmid-T2|＝2|Tmid-T2|，则T1T2弧长表示为(7)：
Figure BDA0002422312050000043
其中，
Figure BDA0002422312050000044
如公式(7)所示，正则项Eregularizer表示为顶点的线性组合，如公式(8)所示：
Figure BDA0002422312050000045
其中，K由公式(7)中矩阵装配得到，||KV||1表示KV的L1正则化；
结合公式(3)和(8)，最小化能量函数(2)写成如下公式(9)：
Figure BDA0002422312050000051
其中，
Figure BDA0002422312050000052
是求满足λ||PV-V0||1+||KV||1最小的顶点位置，即V；
步骤5.2利用增值拉格拉日方法迭代求解孔洞及其邻域的顶点位置，即顶点位置优化，具体如下；
步骤5.2.1求解方程(9)转化为如下公式(10)求解带约束的优化问题：
Figure BDA0002422312050000053
其中，Z＝PV-V0，Q＝KV，||Z||1表示Z的L1正则化，||Q||1表示Q的L1正则化；
Figure BDA0002422312050000054
是求满足λ||Z||1+||Q||1最小的Z，Q；
则根据增值拉格拉日方法将上述约束问题转为求解如下公式(11)的泛函鞍点问题：
Figure BDA0002422312050000055
其中，λZ与λQ是拉格拉日乘子；＜λZ,Z-(PV-V0)＞表示λZ和Z-(PV-V0)的内积，＜λQ,Q-KV＞表示λQ和Q-KV的内积；
Figure BDA0002422312050000056
表示Z-(PV-V0)的L2正则化，
Figure BDA0002422312050000057
表示Q-KV的L2正则化；rZ，rQ是惩罚因子，并且rZ＞0，rQ＞0；则优化问题转化为如下公式(12)的鞍点问题：
Figure BDA0002422312050000058
其中，
Figure BDA0002422312050000059
是求满足变分方程L(V,Z,Q；λZ,λQ)最小的V，Z，Q；
步骤5.2.2求解优化问题(12)；具体将问题(12)转化为依次求解3个子问题，然后迭代更新拉格拉日乘子，通过如下子步骤实现：
步骤5.2.2A固定Z,Q，求解V，即求解V子问题，V子问题转化如下公式(13)的二次方程形式：
Figure BDA00024223120500000510
其中，
Figure BDA00024223120500000511
是求满足
Figure BDA00024223120500000512
最小的V；
该问题转化为线性方程求解；
步骤5.2.2B固定V,Q，求解Z，即求解Z子问题，Z子问题转化为如下公式(14)形式：
Figure BDA0002422312050000061
其中，
Figure BDA0002422312050000062
是求满足
Figure BDA0002422312050000063
最小的Z；
问题(14)有如下公式(15)的封闭形式解：
Figure BDA0002422312050000064
其中，
Figure BDA0002422312050000065
是取0和
Figure BDA0002422312050000066
中的最大值；
步骤5.2.2C固定V,Z，求解Q，即求解Q子问题，Q子问题转化为如下公式(16)形式：
Figure BDA0002422312050000067
其中，
Figure BDA0002422312050000068
是求满足
Figure BDA0002422312050000069
最小的Q；
问题(16)有如下公式(17)的封闭形式解：
Figure BDA00024223120500000610
其中，
Figure BDA00024223120500000611
是取0和
Figure BDA00024223120500000612
的最大值；
步骤5.2.3更新拉格拉日乘子，其中第η+1次的迭代与第η次的关系如下(18)：
Figure BDA00024223120500000613
步骤5.2.4迭代求解；
令初值
Figure BDA00024223120500000614
依次迭代求解方程(13)，(14)，(16)，更新拉格拉日乘子(18)，直到满足终止条件；
其中，终止条件为：假设连续两次迭代，如η,η+1次迭代，控制顶点的距离记为
Figure BDA0002422312050000071
当ε小于给定的阈值ε0时，迭代停止。
步骤六、更新number的值为0，即令number＝0；执行步骤三。
有益效果
本发明提出了一种基于“连接关系-位置”迭代优化的网格修复方法，与现有网格修复方法相比，具有如下有益效果：
1.本发明所提方法在孔洞修复过程中，在恢复缺失的显著特征方面，具有明显优势；
2.本发明所提方法在数字医疗，虚拟现实和城市重建等领域具有广泛的应用前景。
附图说明
图1是本发明一种基于“连接关系-位置”迭代优化的网格修复方法的框架图；
图2是本发明一种基于“连接关系-位置”迭代优化的网格修复方法及实施例1中孔洞示意图；
图3是本发明一种基于“连接关系-位置”迭代优化的网格修复方法及实施例1中孔洞连接关系优化方法图；
图4是本发明一种基于“连接关系-位置”迭代优化的网格修复方法及实施例1中边交换示意图；
图5是本发明一种基于“连接关系-位置”迭代优化的网格修复方法及实施例1中孔洞邻域示意图；
图6是本发明一种基于“连接关系-位置”迭代优化的网格修复方法及实施例1中网格内二面角结构示意图；
图7是本发明一种基于“连接关系-位置”迭代优化的网格修复方法及实施例1中的顶点位置优化方法图。
具体实施方式
下面结合附图及实施例对本发明所提的一种基于“连接关系-位置”迭代优化的网格修复方法进行详细说明。
实施例1
步骤1.2具体实施时，采用文献1中第4节所提方法；
文献1：Liepa,Peter.Filling holes in meshes，Proceedings of the2003Eurographics/ACM SIGGRAPH symposium on Geometry processing，EurographicsAssociation,2003.
如图1所示，是本发明一种基于“连接关系-位置”迭代优化的网格修复方法及本实施例1的算法框架图。
由图1看出，本发明针对输入残缺三维网格模型进行如下步骤：
步骤A、识别孔洞，利用动态规划方法初始化孔洞内连接关系；
步骤A.1，识别模型中的孔洞区域；如图2所示，判定只有一个邻接三角面片的边为孔洞边界边，边界边首尾相连构成的闭合回路所包围的区域为孔洞区域；
步骤A.2，利用动态规划方法进行孔洞区域的三角化，初始化孔洞区域的连接关系；
步骤B、记网格中发生调整的边的个数为number，初始值为0；
步骤C、识别孔洞边界一对特征点，拟合跨越孔洞的特征线，以特征线为导向优化孔洞区域连接关系；
步骤C.1,识别孔洞边界一对特征点，基于特征点及其法向拟合跨越孔洞的特征线；计算孔洞边界顶点所有邻接边的二面角，若值大于30度，则该边为特征边，其邻接边界顶点即为特征点，依据特征边确定该点的法向，从而得到孔洞一对边界特征点vx,vy以及其对应的法向n1,n2；然后利用三次贝塞尔曲线粗略拟合孔洞区域特征线lB
B(t)＝vx(1-t)+(n1+3vy-n2)t(1-t)2+(3vy-n2)t2(1-t)+vyt3,t∈[0,1]
其中，t代表时间，B(t)代表t时刻贝塞尔曲线上一个坐标点；
步骤C.2,以B.1拟合的特征线为导向，优化孔洞区域连接关系；如图3所示，是本发明一种基于“连接关系-位置”迭代优化的网格修复方法及实施例1中孔洞连接关系优化方法图；由图3看出，本发明针对孔洞连接关系优化进行如下步骤：
步骤C.2.1,在特征线lB上选取采样点集合VB＝{VBi}i＝1,2...b，将点集VB向孔洞区域做投影，得到投影面片集合TB＝{ti}i＝1,2...b，利用特征点vx,vy确定孔洞区域特征线大致方向为向量
Figure BDA0002422312050000081
的方向；
步骤C.2.2,从特征点vx出发，遍历与vx相邻接投影三角形中的边；首先计算所有与vx相邻接边与
Figure BDA0002422312050000082
夹角大小，记录夹角最小的边的另一个顶点vtmp以及夹角值θ1；然后计算所有与vx相对的投影边若发生边交换后，形成的新边与
Figure BDA0002422312050000083
夹角大小，记录夹角最小的边对应的另一个顶点v'tmp以及夹角值θ2；若θ1≥θ2,则连接关系不变，将vtmp替换特征点vx；否则，如图4所示，发生边交换，number＝number+1，并将v'tmp替换特征点vx；
步骤C.2.3判断vx与vy是否相等，若vx＝vy，则此次连接关系优化结束；否则重复步骤C.2.2；
步骤D、判断number取值是否为0，若number＞0，执行步骤E；否则，结束修复；
步骤E、基于孔洞及其邻域，构建基于变分框架三维网格局部修复能量函数,利用增值拉格拉日方法迭代求解孔洞及其邻域的顶点位置；
其中，邻域为孔洞k环邻域，k取值范围在[3,10]；如图5所示，为1环邻域示意图；
步骤E，具体包括如下子步骤；
步骤E.1基于孔洞及其邻域信息，构建基于三角面片表示的三维网格模型局部修复的变分框架；
记孔洞区域及其邻域为孔洞局部，输入的模型M0的孔洞邻域网格为
Figure BDA0002422312050000084
顶点集合为
Figure BDA0002422312050000085
其中，m是孔洞邻域网格
Figure BDA0002422312050000086
中的顶点个数，孔洞区域为DH；修复后的网格M的孔洞局部网格为MH，顶点集合为V＝{v1,v2,...,vn}，边集合为E＝{e1,e2,...,ed}，边长度集合为L＝{l1,l2,...,ld}，内二面角集合为θ＝{θ1,θ2,...,θd}，其中，n(n＞m)是修复连接关系后的孔洞局部MH中的顶点个数，d是MH中边的个数，边集合E中的边ei的长度即为长度集合L中的li，内面角集合θ中的角度θi表示共享边ei的两个三角面片间的内二面角；保持孔洞局部连接关系不变，通过变分框架，优化函数
Figure BDA0002422312050000091
得到最优的孔洞局部顶点位置；
其中，
Figure BDA0002422312050000092
是求满足
Figure BDA0002422312050000093
最小的孔洞局部顶点位置，即MH；
Figure BDA0002422312050000094
为数据项Efidelity，使得修复后的网格MH尽可能地逼近给定的初始网格
Figure BDA0002422312050000095
||MH||TV为总变差正则项Eregularizer；λv(MH)为数据项截断参数，当v∈DH时，λv(MH)＝0；
Figure BDA00024223120500000912
时，λv(MH)＝λ＞0；
步骤E.1，计算数据项，具体通过如下公式计算：
Figure BDA0002422312050000096
其中，vi代表修复后网格孔洞局部MH的顶点集合，V＝{v1,v2,...,vn}中第i项；
Figure BDA0002422312050000097
代表输入的残缺网格孔洞邻域
Figure BDA0002422312050000098
的顶点集合，
Figure BDA0002422312050000099
中第i项；||PV-V0||1表示FV-V0的L1正则化；P表示m×n的投影矩阵，定义如下：
Figure BDA00024223120500000910
步骤E.2，计算正则项,具体通过如下公式计算：
Figure BDA00024223120500000911
其中，li代表孔洞局部MH中边ei的边长，即边长度集合l＝{l1,l2,...,ld}中第i项；θi表示ei的两个邻接三角面片间的内二面角，(π-θi)表示内二面角θi的补角。如图6所示，假设边ei的两个邻接三角面片为Δv1v3v4和Δv1v2v3，边的两端点为v1，v3，对应第三个顶点分别为v2和v4；定义T1和T2分别是三角面片Δv1v3v4和Δv1v2v3的外法向，长度均为||v1v3||，T1和T2间夹角即为边ei的二面角补角(π-θi)，则li|π-θ|＝||v1v3|||π-θ|即为T1和T2之间所夹弧T1T2的弧长；法向量T1和T2表示为：
T1＝cot(θ4,1,3)(v4-v3)+cot(θ1,3,4)(v4-v1),
T2＝cot(θ2,3,1)(v1-v2)+cot(θ3,1,2)(v3-v2).
其中，θ4,1,3是边v1v4和边v1v3的夹角，θ1,3,4是边v1v3和边v3v4的夹角，θ2,3,1是边v2v3和边v1v3的夹角，θ3,1,2是边v1v2和边v1v3的夹角；
采用近似计算的方法计算弧长，定义向量Tmid为指向弧T1T2中心，且与T1和T2同起点等长度的向量，则Tmid表示为：
Figure BDA0002422312050000101
于是用弦长之和近似表达T1T2弧长，即||T1T2||＝|T1-Tmid|+|Tmid-T2|＝2|Tmid-T2|，则T1T2弧长表示为
Figure BDA0002422312050000102
其中，
Figure BDA0002422312050000103
如上式所示，正则项Eregularizer表示为顶点的线性组合：
Figure BDA0002422312050000104
其中，K由弧长计算公式中矩阵装配得到，||KV||1表示KV的L1正则化；
综上，方程
Figure BDA0002422312050000105
转化为：
Figure BDA0002422312050000106
其中，
Figure BDA0002422312050000107
是求满足λ||PV-V0||1+||KV||1最小的顶点位置，即V；
步骤E.2,利用增值拉格拉日方法求解能量函数，得到最优的孔洞局部顶点位置；
步骤E.2.1，求解方程
Figure BDA0002422312050000108
转化为如下求解带约束的优化问题：
Figure BDA0002422312050000109
其中，Z＝PV-V0，Q＝KV，||Z||1表示Z的L1正则化，||Q||1表示Q的L1正则化；
Figure BDA0002422312050000111
是求满足λ||Z||1+||Q||1最小的Z，Q；
则根据增值拉格拉日方法将上述约束问题转为求解如下泛函鞍点问题：
Figure BDA0002422312050000112
其中，λZ与λQ是拉格拉日乘子；＜λZ,Z-(PV-V0)＞表示λZ和Z-(PV-V0)的内积，＜λQ,Q-KV＞表示λQ和Q-KV的内积；
Figure BDA0002422312050000113
表示Z-(PV-V0)的L2正则化，
Figure BDA0002422312050000114
表示Q-KV的L2正则化；rZ，rQ是惩罚因子，并且rZ＞0，rQ＞0；则优化问题转化为如下鞍点问题：
Figure BDA0002422312050000115
其中，
Figure BDA0002422312050000116
是求满足变分方程L(V,Z,Q；λZ,λQ)最小的V，Z，Q；
步骤E.2.2，求解优化问题
Figure BDA0002422312050000117
具体将其转化为依次求解3个子问题，然后迭代更新拉格拉日乘子，通过如下子步骤实现：
·固定Z,Q，求解V，即求解V子问题，V子问题转化如下二次方程形式：
Figure BDA0002422312050000118
其中，
Figure BDA0002422312050000119
是求满足
Figure BDA00024223120500001110
最小的V；
该问题转化为线性方程求解；
·固定V,Q，求解Z，即求解Z子问题，Z子问题转化为如下形式：
Figure BDA00024223120500001111
其中，
Figure BDA00024223120500001112
是求满足
Figure BDA00024223120500001113
最小的Z；
问题
Figure BDA0002422312050000121
有如下封闭形式解：
Figure BDA0002422312050000122
其中，
Figure BDA0002422312050000123
是取0和
Figure BDA0002422312050000124
中的最大值；
·固定V,Z，求解Q，即求解Q子问题，Q子问题转化为如下形式：
Figure BDA0002422312050000125
其中，
Figure BDA0002422312050000126
是求满足
Figure BDA0002422312050000127
最小的Q；
问题
Figure BDA0002422312050000128
有如下封闭形式解：
Figure BDA0002422312050000129
其中，
Figure BDA00024223120500001210
是取0和
Figure BDA00024223120500001211
的最大值；
步骤E.2.3，更新拉格拉日乘子，其中第η+1次的迭代与第η次的关系如下：
Figure BDA00024223120500001212
Figure BDA00024223120500001213
步骤E.2.4，迭代求解；
令初值
Figure BDA00024223120500001214
依次迭代求解方程
Figure BDA00024223120500001215
Figure BDA00024223120500001216
Figure BDA00024223120500001217
更新拉格拉日乘子
Figure BDA0002422312050000131
直到满足终止条件，见图7算法图；
其中，终止条件为：假设连续两次迭代，如η,η+1次迭代，控制顶点的距离记为
Figure BDA0002422312050000132
当ε小于给定的阈值ε0时，迭代停止。
步骤F、更新number的值为0，即令number＝0；执行步骤C；
以上所述为本发明的较佳实施例而已，本发明不应该局限于该实施例和附图所公开的内容。凡是不脱离本发明所公开的精神下完成的等效或修改，都落入本发明保护的范围。

Claims (4)

1.一种基于“连接关系-位置”迭代优化的网格修复方法，其特征在于：具体步骤如下：
步骤一、识别输入的三维网格模型的孔洞区域，利用动态规划方法初始化孔洞区域连接关系，得到具有完整连接关系的初始三维网格模型，即实现孔洞区域的连接关系修复；
其中，输入的三维网格模型记为M0，初始三维网格模型记为M；
步骤二、记网格中发生调整的边的个数为number，初始值为0；
步骤三、识别孔洞边界一对特征点，基于特征点及其法向拟合跨越孔洞的特征线，以特征线为导向，调整孔洞区域连接关系，即连接关系优化，具体为：
步骤3.1识别孔洞边界一对特征点，基于特征点及其法向拟合跨越孔洞的特征线；计算孔洞边界顶点所有邻接边的二面角，若值大于θf，则该边为特征边，其邻接边界顶点即为特征点；依据特征边确定该点的法向，从而得到孔洞一对边界特征点vx,vy以及其对应的法向n1,n2；然后利用三次贝塞尔曲线(1)粗略拟合孔洞区域特征线lB：
B(t)＝vx(1-t)+(n1+3vy-n2)t(1-t)2+(3vy-n2)t2(1-t)+vyt3,t∈[0,1] (1)
其中，θf代表角度，t代表时间，B(t)代表t时刻贝塞尔曲线上一个坐标点；
步骤3.2以步骤3.1拟合的特征线为导向，调整孔洞区域连接关系，具体如下：
步骤3.2.1在特征线lB上选取采样点集合VB＝{VBi}i＝1,2…b，将点集VB向孔洞区域做投影，得到投影面片集合TB＝{ti}i＝1,2…b，投影面片上边构成投影边集合LB＝{lBi}i＝1,2…b，利用特征点vx,vy确定孔洞区域特征线大致方向为向量
Figure FDA0002422312040000011
的方向；
步骤3.2.2从特征点vx出发，遍历与vx相邻接投影三角形中的边；首先计算所有与vx相邻接边与
Figure FDA0002422312040000012
夹角大小，记录夹角最小的边的另一个顶点vtmp以及夹角值θ1；然后计算所有与vx相对的投影边若发生边交换后，形成的新边与
Figure FDA0002422312040000013
夹角大小，记录夹角最小的边对应的另一个顶点v'tmp以及夹角值θ2；若θ1≥θ2,则连接关系不变，将vtmp替换特征点vx；否则发生边交换，number＝number+1，并将v'tmp替换特征点vx；
步骤3.2.3判断vx与vy是否相等，若vx＝vy，则此次连接关系优化结束；否则重复步骤3.2.2；
步骤四、判断number取值是否为0，若number＞0，执行步骤五；否则，结束修复；
步骤五、基于孔洞及其邻域信息，构建基于三角面片表示的三维网格模型局部修复的变分框架，利用增值拉格拉日方法迭代求解孔洞及其邻域的顶点位置，即顶点位置优化；
步骤五，具体包括如下子步骤；
步骤5.1基于孔洞及其邻域信息，构建基于三角面片表示的三维网格模型局部修复的变分框架；
记孔洞区域及其邻域为孔洞局部，输入的模型M0的孔洞邻域网格为
Figure FDA0002422312040000021
顶点集合为
Figure FDA0002422312040000022
其中，m是孔洞邻域网格
Figure FDA0002422312040000023
中的顶点个数，孔洞区域为DH；修复后的网格M的孔洞局部网格为MH，顶点集合为V＝{v1,v2,...,vn}，边集合为E＝{e1,e2,...,ed}，边长度集合为L＝{l1,l2,...,ld}，内二面角集合为θ＝{θ1,θ2,...,θd}，其中，n(n＞m)是修复连接关系后的孔洞局部MH中的顶点个数，d是MH中边的个数，边集合E中的边ei的长度即为长度集合L中的li，内面角集合θ中的角度θi表示共享边ei的两个三角面片间的内二面角；保持孔洞局部连接关系不变，通过最小化如下(2)式的能量函数得到最优的孔洞局部顶点位置；
Figure FDA0002422312040000024
其中，
Figure FDA0002422312040000025
是求满足
Figure FDA0002422312040000026
最小的孔洞局部顶点位置，即MH；
Figure FDA0002422312040000027
为数据项Efidelity，使得修复后的网格MH尽可能地逼近给定的初始网格
Figure FDA0002422312040000028
||MH||TV为总变差正则项Eregularizer；λv(MH)为数据项截断参数，当v∈DH时，λv(MH)＝0；
Figure FDA0002422312040000029
时，λv(MH)＝λ＞0；
步骤5.1.1计算数据项，具体通过如下公式(3)计算：
Figure FDA00024223120400000210
其中，vi代表修复后网格孔洞局部MH的顶点集合，V＝{v1,v2,...,vn}中第i项；
Figure FDA00024223120400000211
代表输入的残缺网格孔洞邻域
Figure FDA00024223120400000212
的顶点集合，
Figure FDA00024223120400000213
中第i项；||PV-V0||1表示FV-V0的L1正则化；P表示m×n的投影矩阵，定义如(4)：
Figure FDA00024223120400000214
步骤5.1.2计算正则项,具体通过公式(5)计算：
Figure FDA00024223120400000215
其中，li代表孔洞局部MH中边ei的边长，即边长度集合l＝{l1,l2,...,ld}中第i项；θi表示ei的两个邻接三角面片间的内二面角，(π-θi)表示内二面角θi的补角；假设边ei的两个邻接三角面片为Δv1v3v4和Δv1v2v3，边的两端点为v1，v3，对应第三个顶点分别为v2和v4；定义T1和T2分别是三角面片Δv1v3v4和Δv1v2v3的外法向，长度均为||v1v3||，T1和T2间夹角即为边ei的二面角补角(π-θi)，则li|π-θ|＝||v1v3|||π-θ|即为T1和T2之间所夹弧
Figure FDA0002422312040000031
的弧长；法向量T1和T2表示为(6)：
T1＝cot(θ4,1,3)(v4-v3)+cot(θ1,3,4)(v4-v1),
T2＝cot(θ2,3,1)(v1-v2)+cot(θ3,1,2)(v3-v2). (6)
其中，θ4,1,3是边v1v4和边v1v3的夹角，θ1,3,4是边v1v3和边v3v4的夹角，θ2,3,1是边v2v3和边v1v3的夹角，θ3,1,2是边v1v2和边v1v3的夹角；
采用近似计算的方法计算弧长，定义向量Tmid为指向弧
Figure FDA0002422312040000032
中心，且与T1和T2同起点等长度的向量，则Tmid表示为：
Figure FDA0002422312040000033
于是用弦长之和近似表达
Figure FDA0002422312040000034
弧长，即
Figure FDA0002422312040000035
则
Figure FDA0002422312040000036
弧长表示为(7)：
Figure FDA0002422312040000037
其中，
Figure FDA0002422312040000038
如公式(7)所示，正则项Eregularizer表示为顶点的线性组合，如公式(8)所示：
Figure FDA0002422312040000039
其中，K由公式(7)中矩阵装配得到，||KV||1表示KV的L1正则化；
结合公式(3)和(8)，最小化能量函数(2)写成如下公式(9)：
Figure FDA00024223120400000310
其中，
Figure FDA00024223120400000311
是求满足λ||PV-V0||1+||KV||1最小的顶点位置，即V；
步骤5.2利用增值拉格拉日方法迭代求解孔洞及其邻域的顶点位置，即顶点位置优化，具体如下；
步骤5.2.1求解方程(9)转化为如下公式(10)求解带约束的优化问题：
Figure FDA00024223120400000312
其中，Z＝PV-V0，Q＝KV，||Z||1表示Z的L1正则化，||Q||1表示Q的L1正则化；
Figure FDA0002422312040000041
是求满足λ||Z||1+||Q||1最小的Z，Q；
则根据增值拉格拉日方法将上述约束问题转为求解如下公式(11)的泛函鞍点问题：
Figure FDA0002422312040000042
其中，λZ与λQ是拉格拉日乘子；＜λZ,Z-(PV-V0)＞表示λZ和Z-(PV-V0)的内积，＜λQ,Q-KV＞表示λQ和Q-KV的内积；
Figure FDA0002422312040000043
表示Z-(PV-V0)的L2正则化，
Figure FDA0002422312040000044
表示Q-KV的L2正则化；rZ，rQ是惩罚因子，并且rZ＞0，rQ＞0；则优化问题转化为如下公式(12)的鞍点问题：
Figure FDA0002422312040000045
其中，
Figure FDA0002422312040000046
是求满足变分方程L(V,Z,Q；λZ,λQ)最小的V，Z，Q；
步骤5.2.2求解优化问题(12)；具体将问题(12)转化为依次求解3个子问题，然后迭代更新拉格拉日乘子，通过如下子步骤实现：
步骤5.2.2A固定Z,Q，求解V，即求解V子问题，V子问题转化如下公式(13)的二次方程形式：
Figure FDA0002422312040000047
其中，
Figure FDA0002422312040000048
是求满足
Figure FDA0002422312040000049
最小的V；
该问题转化为线性方程求解；
步骤5.2.2B固定V,Q，求解Z，即求解Z子问题，Z子问题转化为如下公式(14)形式：
Figure FDA00024223120400000410
其中，
Figure FDA00024223120400000411
是求满足
Figure FDA00024223120400000412
最小的Z；
问题(14)有如下公式(15)的封闭形式解：
Figure FDA00024223120400000413
其中，
Figure FDA0002422312040000051
是取0和
Figure FDA0002422312040000052
中的最大值；
步骤5.2.2C固定V,Z，求解Q，即求解Q子问题，Q子问题转化为如下公式(16)形式：
Figure FDA0002422312040000053
其中，
Figure FDA0002422312040000054
是求满足
Figure FDA0002422312040000055
最小的Q；
问题(16)有如下公式(17)的封闭形式解：
Figure FDA0002422312040000056
其中，
Figure FDA0002422312040000057
是取0和
Figure FDA0002422312040000058
的最大值；
步骤5.2.3更新拉格拉日乘子，其中第η+1次的迭代与第η次的关系如下(18)：
Figure FDA00024223120400000511
步骤5.2.4迭代求解；
令初值
Figure FDA0002422312040000059
依次迭代求解方程(13)，(14)，(16)，更新拉格拉日乘子(18)，直到满足终止条件；
其中，终止条件为：假设连续两次迭代，如η,η+1次迭代，控制顶点的距离记为
Figure FDA00024223120400000510
当ε小于给定的阈值ε0时，迭代停止；
步骤六、更新number的值为0，即令number＝0；执行步骤三。
2.根据权利要求1所述的一种基于“连接关系-位置”迭代优化的网格修复方法，其特征在于：步骤一，具体包括如下子步骤：
步骤1.1识别模型中的孔洞区域；
其中，识别孔洞区域的准则为：判定只有一个邻接三角面片的边为孔洞边界边，边界边首尾相连构成的闭合回路所包围的区域为孔洞区域；
步骤1.2采用动态规划方法对步骤1.1识别的孔洞区域进行三角化，并初始化孔洞区域的连接关系。
3.根据权利要求1所述的一种基于“连接关系-位置”迭代优化的网格修复方法，其特征在于：步骤3.1中θf的取值范围在[30°,45°]。
4.根据权利要求1所述的一种基于“连接关系-位置”迭代优化的网格修复方法，其特征在于：步骤四中，邻域为孔洞k环邻域，k取值范围在[3,10]。
CN202010209452.6A 2019-11-22 2020-03-23 一种基于“连接关系-位置”迭代优化的网格修复方法 Active CN111274732B (zh)

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
CN2019111567206 2019-11-22
CN201911156720 2019-11-22

Publications (2)

Publication Number Publication Date
CN111274732A CN111274732A (zh) 2020-06-12
CN111274732B true CN111274732B (zh) 2023-02-24

Family

ID=70998453

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010209452.6A Active CN111274732B (zh) 2019-11-22 2020-03-23 一种基于“连接关系-位置”迭代优化的网格修复方法

Country Status (1)

Country Link
CN (1) CN111274732B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113345096B (zh) * 2021-05-26 2023-06-06 北京工商大学 一种基于非局部低秩优化的三维网格修复方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2010148625A1 (zh) * 2009-06-24 2010-12-29 中国石油集团川庆钻探工程有限公司 复杂地质构造块状模型构建方法
CN104504663A (zh) * 2014-12-29 2015-04-08 佛山市诺威科技有限公司 一种义齿三角网格模型孔洞的迭代修补方法
CN108876922A (zh) * 2018-06-12 2018-11-23 北京工商大学 一种基于内二面角补角正则化的网格修补方法
CN109410334A (zh) * 2018-09-21 2019-03-01 桂林电子科技大学 一种基于特征线的三维网格模型缺陷孔洞修复方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2010148625A1 (zh) * 2009-06-24 2010-12-29 中国石油集团川庆钻探工程有限公司 复杂地质构造块状模型构建方法
CN104504663A (zh) * 2014-12-29 2015-04-08 佛山市诺威科技有限公司 一种义齿三角网格模型孔洞的迭代修补方法
CN108876922A (zh) * 2018-06-12 2018-11-23 北京工商大学 一种基于内二面角补角正则化的网格修补方法
CN109410334A (zh) * 2018-09-21 2019-03-01 桂林电子科技大学 一种基于特征线的三维网格模型缺陷孔洞修复方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
人体模型网格孔洞的快速修复方法;刘世凡等;《西安交通大学学报》;20180418(第08期);全文 *

Also Published As

Publication number Publication date
CN111274732A (zh) 2020-06-12

Similar Documents

Publication Publication Date Title
CN108038906B (zh) 一种基于图像的三维四边形网格模型重建方法
CN104268934B (zh) 一种由点云直接重建三维曲面的方法
JP5242386B2 (ja) 体積グラフラプラシアンを使用する大規模メッシュ変形
CN111797555B (zh) 一种基于有限元模型的几何重构方法
JP2013507679A (ja) 三次元物体モデルの３ｄプリントが可能な方法及びシステム
CN102870139A (zh) 对点云数据进行全局参数化和四边形网格化方法
WO2013123636A1 (en) Method and apparatus for mesh simplification
CN108230452B (zh) 一种基于纹理合成的模型补洞方法
CN110889893A (zh) 表达几何细节和复杂拓扑的三维模型表示方法和系统
CN104282000A (zh) 一种基于旋转及尺度变换的图像修复方法
Zhou et al. Improvement of normal estimation for point clouds via simplifying surface fitting
CN111274732B (zh) 一种基于“连接关系-位置”迭代优化的网格修复方法
CN113077497A (zh) 基于3d ndt-icp算法的掘进巷道点云配准方法
Wei et al. GeoDualCNN: Geometry-supporting dual convolutional neural network for noisy point clouds
CN108876922B (zh) 一种基于内二面角补角正则化的网格修补方法
CN110689620A (zh) 一种多层次优化的网格曲面离散样条曲线设计方法
CN109859322A (zh) 一种基于变形图的谱姿态迁移方法
Huysmans et al. Automatic construction of correspondences for tubular surfaces
CN107274367A (zh) 一种基于结构特征描述的三维几何模型去噪方法
Mao et al. Std-net: Structure-preserving and topology-adaptive deformation network for single-view 3d reconstruction
Oh A new triangular mesh repairing method using a mesh distortion energy minimization-based mesh flattening method
Shang et al. Effective re-parameterization and GA based knot structure optimization for high quality T-spline surface fitting
Wang et al. Kirchhoff–Love shell representation and analysis using triangle configuration B-splines
Fan et al. Arbitrary surface data patching method based on geometric convolutional neural network
CN114880906A (zh) 一种映射皮肤材料特性的三维创面有限单元展开法

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