一、简介

之前无意中看到这样一种点云孔洞修补方法，感觉很有意思，就下载了作者的源码看看效果，具体流程如下所述：
1、首先，我们需要先定位网格上的孔洞。
2、然后计算孔洞边界顶点的中心位置(孔洞边界顶点位置的平均值)。边界上的顶点与中心位置上新创建的顶点进行连接，从而创建一个填充孔洞的表面补丁。在这个补丁中的三角形会被上采样(细分)，这样当我们进行补丁平滑时，我们可以有更多的顶点来操作。

注意，这里对原始网格的三角形进行了上采样，因为这使它非常容易将原始补丁和原始几何图形的顶点融合在一起。所谓“融合”，就是指确保原始网格的边界和补丁的边界共享相同的顶点。

3、接下来，通过在补丁上求解方程Δ2fΔ^{2}fΔ2f=0来创建一个平滑的补丁。这需要离散拉普拉斯算子来实现这一点，而拉普拉斯算子在网格表面上的离散被称为拉普拉斯-贝尔特拉米算子。这里使用libigl库创建了离散拉普拉斯- beltrami运算子，然后使用该运算子求解方程Δ2fΔ^{2}fΔ2f=0。

4、最后，我们使用网格抽取算法对网格进行下采样，得到最终修补之后的结果。

更详细的内容可以详看：https://erkaman.github.io/posts/hole_filling.html。

二、实现代码

CMakeLists.txt

cmake_minimum_required(VERSION 3.1)
project(hole_fixer)set(CMAKE_BUILD_TYPE "Release")
set(CMAKE_CXX_FLAGS_RELEASE "-Zi -Od")
set(CMAKE_EXE_LINKER_FLAGS "-Debug")
message("Build Type:" ${CMAKE_BUILD_TYPE})
message("Release mode:" ${CMAKE_CXX_FLAGS_RELEASE})
message("Linker mode:" ${CMAKE_EXE_LINKER_FLAGS})include_directories(libigl/ libigl/eigen/)
add_executable(hole_fixer main.cpp)

这里对原始代码做了一下改动，让它更符合我的习惯。

main.cpp

#include <igl/boundary_loop.h>
#include <igl/cat.h>
#include <igl/colon.h>
#include <igl/slice.h>
#include <igl/upsample.h>
#include <igl/decimate.h>
#include <igl/shortest_edge_and_midpoint.h>
#include <igl/harmonic.h>
#include <igl/readOFF.h>
#include <igl/writeOFF.h>using Eigen::VectorXi;
using Eigen::MatrixXd;
using Eigen::MatrixXi;
using Eigen::VectorXd;
typedef Eigen::SparseMatrix<double> SpMat;
typedef Eigen::Triplet<double> Tri;void printHelpExit() {printf("Invalid command line arguments specified!\n\n");printf("USAGE: hole_fixer [options]\n\n");printf("OPTIONS: \n");printf("  -in\t\t\ttarget mesh file in .off-format, with a hole\n");printf("  -out\t\t\toutput mesh file in .off-format\n");printf("  -outfaces\t\tHow many faces to decimate the mesh to\n");printf("  -upsample\t\tHow much upsampling to use when creating the patch\n");exit(1);
}const char* findToken(const char* param, int argc, char* argv[]) {const char* token = nullptr;for (int i = 0; i < argc; ++i) {if (strcmp(argv[i], param) == 0) {if (i + 1 < argc) {token = argv[i + 1];break;}}}if (token == nullptr) {printf("Could not find command-line parameter %s\n", param);return nullptr;}return token;
}const char* parseStringParam(const char* param, int argc, char* argv[]) {const char* token = findToken(param, argc, argv);return token;
}bool parseIntParam(const char* param, int argc, char* argv[], unsigned int& out) {const char* token = findToken(param, argc, argv);if (token == nullptr)return false;int r = sscanf(token, "%u,", &out);if (r != 1 || r == EOF) {return false;}else {return true;}
}int main(int argc, char *argv[])
{//// command line parsing.//const char* inFile = parseStringParam("-in", argc, argv);if (inFile == nullptr) printHelpExit();const char* outFile = parseStringParam("-out", argc, argv);if (outFile == nullptr) printHelpExit();unsigned int outFacesN;if (!parseIntParam("-outfaces", argc, argv, outFacesN)) printHelpExit();unsigned int upsampleN;if (!parseIntParam("-upsample", argc, argv, upsampleN)) printHelpExit();// original mesh vertices and indices. This is the original mesh, which has a hole.MatrixXd originalV;MatrixXi originalF;if (!igl::readOFF(inFile, originalV, originalF)) {printHelpExit();}VectorXi originalLoop; // indices of the boundary of the hole. igl::boundary_loop(originalF, originalLoop);if (originalLoop.size() == 0) {printf("Mesh has no hole!");printHelpExit();}// upsample the original mesh. this makes fusing the original mesh with the patch much easier.igl::upsample(Eigen::MatrixXd(originalV), Eigen::MatrixXi(originalF), originalV, originalF, upsampleN);// compute boundary center.VectorXd bcenter(3);{VectorXi R = originalLoop;VectorXi C(3); C << 0, 1, 2;MatrixXd B;MatrixXd A = originalV;igl::slice(A, R, C, B);bcenter = (1.0f / originalLoop.size()) * B.colwise().sum();}// a flat patch that fills the hole.MatrixXd patchV = MatrixXd(originalLoop.size() + 1, 3); // patch will have an extra vertex for the center vertex.MatrixXi patchF = MatrixXi(originalLoop.size(), 3);{VectorXi R = originalLoop;VectorXi C(3); C << 0, 1, 2;MatrixXd A = originalV;MatrixXd temp1;igl::slice(A, R, C, temp1);MatrixXd temp2(1, 3);temp2 << bcenter(0), bcenter(1), bcenter(2);// patch vertices will be the boundary vertices, plus a center vertex. concat these together.igl::cat(1, temp1, temp2, patchV);// create triangles that connect boundary vertices and the center vertex.for (int i = 0; i < originalLoop.size(); ++i) {patchF(i, 2) = (int)originalLoop.size();patchF(i, 1) = i;patchF(i, 0) = (1 + i) % originalLoop.size();}// also upsample patch. patch and original mesh will have the same upsampling level now// making it trivial to fuse them together.igl::upsample(Eigen::MatrixXd(patchV), Eigen::MatrixXi(patchF), patchV, patchF, upsampleN);}// the final mesh, where the patch and original mesh has been fused together.std::vector<std::vector<double>> fusedV;std::vector<std::vector<int>> fusedF;int index = 0; // vertex index counter for the fused mesh.// first, we add the upsampled patch to the fused mesh.{for (; index < patchV.rows(); ++index) {fusedV.push_back({ patchV(index, 0), patchV(index, 1), patchV(index, 2) });}int findex = 0;for (; findex < patchF.rows(); ++findex) {fusedF.push_back({ patchF(findex, 0), patchF(findex, 1), patchF(findex, 2) });}}// now, we fuse the patch and the original mesh together.{// remaps indices from the original mesh to the fused mesh.std::map<int, int> originalToFusedMap;for (int itri = 0; itri < originalF.rows(); ++itri) {int triIndices[3];for (int iv = 0; iv < 3; ++iv) {int triIndex = originalF(itri, iv);int ret;if (originalToFusedMap.count(triIndex) == 0) {int foundMatch = -1;// the vertices at the boundary are the same, for both the two meshes(patch and original mesh).// we ensure that these vertices are not duplicated.// this is also how we ensure that the two meshes are fused together.for (int jj = 0; jj < patchV.rows(); ++jj) {VectorXd u(3); u << fusedV[jj][0], fusedV[jj][1], fusedV[jj][2];VectorXd v(3); v << originalV(triIndex, 0), originalV(triIndex, 1), originalV(triIndex, 2);if ((u - v).norm() < 0.00001) {foundMatch = jj;break;}}if (foundMatch != -1) {originalToFusedMap[triIndex] = foundMatch;ret = foundMatch;}else {fusedV.push_back({ originalV(triIndex, 0), originalV(triIndex, 1), originalV(triIndex, 2) });originalToFusedMap[triIndex] = index;ret = index;index++;}}else {ret = originalToFusedMap[triIndex];}triIndices[iv] = ret;}fusedF.push_back({triIndices[0],triIndices[1],triIndices[2] });}}MatrixXd fairedV(fusedV.size(), 3);MatrixXi fairedF(fusedF.size(), 3);// now we shall do surface fairing on the mesh, to ensure// that the patch conforms to the surrounding curvature.{for (int vindex = 0; vindex < fusedV.size(); ++vindex) {auto r = fusedV[vindex];fairedV(vindex, 0) = r[0];fairedV(vindex, 1) = r[1];fairedV(vindex, 2) = r[2];}for (int findex = 0; findex < fusedF.size(); ++findex) {auto r = fusedF[findex];fairedF(findex, 0) = r[0];fairedF(findex, 1) = r[1];fairedF(findex, 2) = r[2];}VectorXi b(fairedV.rows() - patchV.rows());MatrixXd bc(fairedV.rows() - patchV.rows(), 3);// setup the boundary conditions. This is simply the vertex positions of the vertices not part of the patch.for (int i = (int)patchV.rows(); i < (int)fairedV.rows(); ++i) {int jj = i - (int)patchV.rows();b(jj) = i;bc(jj, 0) = fairedV(i, 0);bc(jj, 1) = fairedV(i, 1);bc(jj, 2) = fairedV(i, 2);}MatrixXd Z;int k = 2;// surface fairing simply means that we solve the equation// Delta^2 f = 0// with appropriate boundary conditions.// this function igl::harmonic from libigl takes care of that.// note that this is pretty inefficient thought.// the only boundary conditions necessary are the 2-ring of vertices around the patch.// the rest of the mesh vertices need not be specified.// we specify the rest of the mesh for simplicity of code, but it is not strictly necessary,// and pretty inefficient, since we have to solve a LARGE matrix equation as a result of this.igl::harmonic(fairedV, fairedF, b, bc, k, Z);fairedV = Z;}// finally, we do a decimation step.MatrixXd finalV(fusedV.size(), 3);MatrixXi finalF(fusedF.size(), 3);VectorXi temp0; VectorXi temp1;igl::decimate(fairedV, fairedF, outFacesN, finalV, finalF, temp0, temp1);igl::writeOFF(outFile, finalV, finalF);std::cout << "执行成功!" << std::endl;system("pause");return 0;
}

全部的代码可以点击网址：https://github.com/Erkaman/hole_fixer。

三、配置工作

1、代码获取到之后，首先使用cmake工具生成VS项目。

2、右击"ALL BUILD"生成所有项目，待生成所有项目之后，设置相应的命令行参数。

3、最后设置hole_fixer项目为启动项目，运行程序就可以了。

四、执行效果

不过可以明显的看到，修补的效果并不是很理想，补丁部分明显与周围的网格密度有很大差别，这还需后续的整形处理。

四、参考资料

[1]https://github.com/Erkaman/hole_fixer

基于RBF函数的点云孔洞修补相关推荐

点云孔洞定位_基于TSDF模型的点云孔洞修复方法
摘要深度摄像头的点云集合一般存在黑色孔洞闪烁的现象,在识别机械臂静态目标抓取时,点云集合数据的关键位置的体素不规则出现,就会对目标中心点的识别造成无法修正的偏差.因此,引入多个点云集合累积的思想,提 ...
点云孔洞定位_点云处理——孔洞修补
之前在做基于点云的孔洞修补研究,参考了一些学位论文后选择了一种基于自己能力能够复现的简单算法,虽然最终效果一般,不过在这个过程中收获了很多,特此记录. 基于点云的孔洞修补与基于三角网格的孔洞修补相比, ...
点云孔洞定位_孔洞修补研究总结
本人最近需要研究三维模型的孔洞修补算法,故上网看了一些大神写的资料,现汇总如下: 重建骨骼比较典型的方法有:用径向基函数从不完整的扫描数据生成连续网格:基于CT等值面数据生成曲面:用傅里叶级数拟合CT ...
matlab rbf函数_基于径向基函数（RBF）的无网格伪谱法与程序实现（2）——微分矩阵...
参考资料 Gregory E. Fasshauer. Meshfree Approximation Methods with MATLAB. P.387 P401 数值实现 Matlab 2019a ...
Matlab/Simulink中的S函数模块嵌入人工智能、神经网络算法设计仿真案例详解（以基于RBF神经网络算法的VSG转动惯量自调节为例）
参考文献 An improved virtual synchronous generator power control strategy Deep reinforcement learning b ...
点云孔洞定位_散乱点云的孔洞识别和边界提取算法研究
散乱点云的孔洞识别和边界提取算法研究王春香,孟宏,张勇 [摘要] 针对逆向工程中已有孔洞识别算法执行效率低.孔洞边界点提取不完整等问题,提出一种新的基于 KD 树和 K 邻域搜索的点云孔洞识 ...
matlab点云聚类,基于区域聚类分割的点云特征线提取
王晓辉 , 吴禄慎 , 陈华伟 , 胡赟 , 石雅莹 . . 基于区域聚类分割的点云特征线提取. 光学学报, 2018, 38(11): 1110001-. Wang Xiaohui , Wu Lus ...
优化思路千万种，基于下界函数的最优化效率如何？
作者丨stephenDC 来源 | 大数据与人工智能(ID:ai-big-data) 导读:生活中我们处处面临最优化的问题,比如,怎么样一个月减掉的体重最高?怎么样学习效率最高?怎么样可以最大化实现个 ...
rbf神经网络_基于RBF神经网络的监督控制（09）
1.RBF监督控制基于RBF神经网络的监督控制系统,其控制思想是:初始阶段采用PD反馈控制,然后过渡到神经网络控制.在控制过程中,如出现较大的误差,则PD控制起主导作用,神经网络控制起调节作用. 图 ...

基于RBF函数的点云孔洞修补

文章目录