简介
Loop 是三角网格常用的细分算法之一. 原理基于二次B样条曲线.
Image
参考链接
https://blog.csdn.net/Hachi_Lin/article/details/90349216
https://www.bilibili.com/video/BV1X7411F744?p=12
https://blog.csdn.net/qq_31804159/article/details/109147352
code
#pragma once
#pragma once
#include "strategy.h"
#include "qwidget.h"
#include <OpenMesh/Core/Mesh/PolyMesh_ArrayKernelT.hh>
#include <OpenMesh/Core/IO/MeshIO.hh>
#include <map>
class Loop :public Strategy, public QWidget
{
Q_OBJECT
private:
typedef OpenMesh::PolyMesh_ArrayKernelT<> MyMesh;
MyMesh mesh;
std::map<int, OpenMesh::Vec3d> facePoints;
std::map<int, OpenMesh::Vec3d> edgePoints;
std::map<int, OpenMesh::Vec3d> vertexPoints;
int times;
public:
Loop(Data* data_);
~Loop();
void genCube();
void genFacePoint();
void genEdgePoint();
void genVertexPoint();
void connectPoints();
void genMesh(std::string name);
bool Run();
void getResult();
};
#include "loop_surface_subdivision.h"
#include <iostream>
#include <unordered_map>
#include <QInputDialog>
/**
* @description: 构造函数
* @param {type}
* @return {type}
*/
Loop::Loop(Data* data_) : Strategy(data_) {
times = 1;
genCube();
}
/**
* @description: 析构函数
* @param {type}
* @return {type}
*/
Loop::~Loop() {}
/**
* @description: 一整套流程
* @param {type}
* @return {int} 管线运行是否成功
*/
bool Loop::Run() {
times = QInputDialog::getInt(this, "Surface Mesh", "Please input times",
1, 1, 1000, 1);
for (int i = 0; i <= times; i++) {
if (i == 0) {
char a[20];
sprintf(a, "output%d.off", i);
genMesh(a);
}
else {
genEdgePoint();
genVertexPoint();
connectPoints();
char a[20];
sprintf(a, "output%d.off", i);
genMesh(a);
}
// 清空变量
std::cout << "[DEBUG] 迭代了第 " << i << " 次。" << std::endl;
}
getResult();
return true;
}
/**
* @description: 生成一个立方体(四边形网格)
* @param {type}
* @return {type}
*/
void Loop::genCube()
{
MyMesh::VertexHandle vhandle[9];
vhandle[0] = mesh.add_vertex(MyMesh::Point(-1, -1, 1));
vhandle[1] = mesh.add_vertex(MyMesh::Point(1, -1, 1));
vhandle[2] = mesh.add_vertex(MyMesh::Point(1, 1, 1));
vhandle[3] = mesh.add_vertex(MyMesh::Point(-1, 1, 1));
vhandle[4] = mesh.add_vertex(MyMesh::Point(-1, -1, -1));
vhandle[5] = mesh.add_vertex(MyMesh::Point(1, -1, -1));
vhandle[6] = mesh.add_vertex(MyMesh::Point(1, 1, -1));
vhandle[7] = mesh.add_vertex(MyMesh::Point(-1, 1, -1));
std::vector<MyMesh::VertexHandle> face_vhandles;
face_vhandles.push_back(vhandle[0]);
face_vhandles.push_back(vhandle[1]);
face_vhandles.push_back(vhandle[2]);
mesh.add_face(face_vhandles);
face_vhandles.clear();
face_vhandles.push_back(vhandle[0]);
face_vhandles.push_back(vhandle[2]);
face_vhandles.push_back(vhandle[3]);
mesh.add_face(face_vhandles);
face_vhandles.clear();
face_vhandles.push_back(vhandle[7]);
face_vhandles.push_back(vhandle[6]);
face_vhandles.push_back(vhandle[5]);
mesh.add_face(face_vhandles);
face_vhandles.clear();
face_vhandles.push_back(vhandle[7]);
face_vhandles.push_back(vhandle[5]);
face_vhandles.push_back(vhandle[4]);
mesh.add_face(face_vhandles);
face_vhandles.clear();
face_vhandles.push_back(vhandle[1]);
face_vhandles.push_back(vhandle[0]);
face_vhandles.push_back(vhandle[4]);
mesh.add_face(face_vhandles);
face_vhandles.clear();
face_vhandles.push_back(vhandle[1]);
face_vhandles.push_back(vhandle[4]);
face_vhandles.push_back(vhandle[5]);
mesh.add_face(face_vhandles);
face_vhandles.clear();
face_vhandles.push_back(vhandle[2]);
face_vhandles.push_back(vhandle[1]);
face_vhandles.push_back(vhandle[5]);
mesh.add_face(face_vhandles);
face_vhandles.clear();
face_vhandles.push_back(vhandle[2]);
face_vhandles.push_back(vhandle[5]);
face_vhandles.push_back(vhandle[6]);
mesh.add_face(face_vhandles);
face_vhandles.clear();
face_vhandles.push_back(vhandle[3]);
face_vhandles.push_back(vhandle[2]);
face_vhandles.push_back(vhandle[6]);
mesh.add_face(face_vhandles);
face_vhandles.clear();
face_vhandles.push_back(vhandle[3]);
face_vhandles.push_back(vhandle[6]);
face_vhandles.push_back(vhandle[7]);
mesh.add_face(face_vhandles);
face_vhandles.clear();
face_vhandles.push_back(vhandle[0]);
face_vhandles.push_back(vhandle[3]);
face_vhandles.push_back(vhandle[7]);
mesh.add_face(face_vhandles);
face_vhandles.clear();
face_vhandles.push_back(vhandle[0]);
face_vhandles.push_back(vhandle[7]);
face_vhandles.push_back(vhandle[4]);
mesh.add_face(face_vhandles);
}
/**
* @description: 考虑边上和内部
* @param {type}
* @return {type}
*/
void Loop::genEdgePoint()
{
edgePoints.clear();
for (auto e_it = mesh.edges_begin(); e_it != mesh.edges_end(); ++e_it)
{
// 得到边所代表的半边
OpenMesh::HalfedgeHandle heh1 = mesh.halfedge_handle(*e_it, 0); // 默认一个方向的半边
OpenMesh::Vec3d edgePoint(0, 0, 0);
int edgePointsNumber = 0;
OpenMesh::DefaultTraits::Point pointV = mesh.point(mesh.from_vertex_handle(heh1)); // 这条(半)边的起点
OpenMesh::DefaultTraits::Point pointW = mesh.point(mesh.to_vertex_handle(heh1)); // 这条(半)边的终点
if (mesh.is_boundary(heh1)) { // 如果这条边处于边界
edgePoints[heh1.idx()] = 1 / 2.0 * (pointV + pointW);
}
else {
OpenMesh::FaceHandle fh1 = mesh.face_handle(heh1);
OpenMesh::HalfedgeHandle heh = mesh.opposite_halfedge_handle(heh1);
OpenMesh::FaceHandle fh2 = mesh.face_handle(heh);
std::vector<OpenMesh::VertexHandle> ss1;
// 得到这两个面的所有顶点
for (auto fv_it : mesh.fv_range(fh1)) {
if (fv_it.idx() != mesh.from_vertex_handle(heh1).idx() &&
fv_it.idx() != mesh.to_vertex_handle(heh1).idx())
ss1.push_back(fv_it);
}
for (auto fv_it : mesh.fv_range(fh2)) {
if (fv_it.idx() != mesh.from_vertex_handle(heh1).idx() &&
fv_it.idx() != mesh.to_vertex_handle(heh1).idx() )
ss1.push_back(fv_it);
}
if (ss1.size() != 2) {
std::cout << "[ERROR] " << ss1.size() << std::endl;
}
edgePoints[heh1.idx()] = 3.0 / 8.0 * (pointV + pointW) + 1.0 / 8.0 * (mesh.point(ss1[0]) + mesh.point(ss1[1]));
}
}
}
/**
* @description: 生成新的顶点
* @param {type}
* @return {type}
*/
void Loop::genVertexPoint()
{
vertexPoints.clear();
// 原始点接触的面的所有的面点的均值
for (auto v_it = mesh.vertices_begin(); v_it != mesh.vertices_end(); v_it++)
{
if (mesh.is_boundary(v_it)) { // 如果这个顶点处于边界
OpenMesh::Vec3d point(0, 0, 0);
for (auto vv_it = mesh.vv_begin(v_it); vv_it.is_valid(); vv_it++) {
if (mesh.is_boundary(vv_it)) {
point += mesh.point(vv_it);
}
}
vertexPoints[(*v_it).idx()] = 3.0 / 4.0 * mesh.point(v_it) + 1.0 / 8.0 * (point);
}
else {
// 计算度 n, 度
int n = 0;
for (auto vf_it = mesh.vv_begin(v_it); vf_it.is_valid(); ++vf_it) {
n++;
}
// 计算β
double beta = 1.0 / n * (5.0 / 8.0 - pow(3.0 / 8.0 + 1.0 / 4.0 * cos(2 * M_PI / n), 2) );
OpenMesh::Vec3d point(0, 0, 0);
for (auto vv_it = mesh.vv_begin(v_it); vv_it.is_valid(); ++vv_it) {
point += beta * mesh.point(vv_it);
}
vertexPoints[(*v_it).idx()] = (1 - n * beta) * mesh.point(v_it) + point;
}
}
}
/**
* @description: 连接面点和边点
* @param {type}
* @return {type}
*/
void Loop::connectPoints() {
if (!mesh.has_vertex_status()) mesh.request_vertex_status();
if (!mesh.has_face_status()) mesh.request_face_status();
if (!mesh.has_edge_status()) mesh.request_edge_status();
// 先将旧顶点移动到新顶点
std::unordered_map<int, OpenMesh::VertexHandle> vertexHandle;
for (auto v_it = mesh.vertices_begin(); v_it != mesh.vertices_end(); v_it++) {
mesh.set_point(*v_it, (OpenMesh::DefaultTraits::Point)vertexPoints[(*v_it).idx()]);
vertexHandle[(*v_it).idx()] = *v_it;
}
std::vector<MyMesh::VertexHandle> facePointsHandle;
std::vector<MyMesh::FaceHandle> faceHandle;
std::unordered_map<int, MyMesh::VertexHandle> edgesHandle;
// 再加入所有的要加入的边点
for (const auto& eh : mesh.edges()) {
edgesHandle[mesh.halfedge_handle(eh, 0).idx()] = mesh.add_vertex(MyMesh::Point(edgePoints[mesh.halfedge_handle(eh, 0).idx()]));
}
for (const auto& fh : mesh.faces()) {
faceHandle.push_back(fh);
}
std::vector<std::vector<MyMesh::VertexHandle>> v;
for (const auto& fh : mesh.faces()) {
std::vector<MyMesh::VertexHandle> handlerVector(6);
int pointsNumber = 0;
for (const auto& fvh : mesh.fv_range(fh)) {
handlerVector[pointsNumber] = vertexHandle[fvh.idx()];
pointsNumber++;
}
for (const auto& feh : mesh.fe_range(fh)) {
auto fehh = mesh.halfedge_handle(feh, 0);
handlerVector[pointsNumber] = edgesHandle[fehh.idx()];
pointsNumber++;
}
v.push_back(handlerVector);
}
for (int i = 0; i < v.size(); i++) {
mesh.delete_face(faceHandle[i], true);
std::vector<MyMesh::VertexHandle> handlerVector;
for (int j = 0; j < v[i].size(); j++) {
handlerVector.push_back(v[i][j]);
}
std::vector<MyMesh::VertexHandle> face_vhandles;
face_vhandles.push_back(handlerVector[3]);
face_vhandles.push_back(handlerVector[0]);
face_vhandles.push_back(handlerVector[4]);
mesh.add_face(face_vhandles);
face_vhandles.clear();
face_vhandles.push_back(handlerVector[4]);
face_vhandles.push_back(handlerVector[1]);
face_vhandles.push_back(handlerVector[5]);
mesh.add_face(face_vhandles);
face_vhandles.clear();
face_vhandles.push_back(handlerVector[3]);
face_vhandles.push_back(handlerVector[4]);
face_vhandles.push_back(handlerVector[5]);
mesh.add_face(face_vhandles);
face_vhandles.clear();
face_vhandles.push_back(handlerVector[2]);
face_vhandles.push_back(handlerVector[3]);
face_vhandles.push_back(handlerVector[5]);
mesh.add_face(face_vhandles);
face_vhandles.clear();
}
mesh.garbage_collection();
if (mesh.has_vertex_status()) mesh.release_vertex_status();
if (mesh.has_face_status()) mesh.release_face_status();
if (mesh.has_edge_status()) mesh.release_edge_status();
}
/**
* @description: 输出新网格
* @param {string} name 文件名称 默认 output.obj
* @return {type}
*/
void Loop::genMesh(std::string name) {
if (name == "") {
name = "output.off";
}
try {
if (!OpenMesh::IO::write_mesh(mesh, name)) {
std::cerr << "Cannot write mesh to file 'output.off'" << std::endl;
return;
}
}
catch (std::exception& e) {
std::cerr << e.what() << std::endl;
return;
}
}
void Loop::getResult() {
//if (getData()->edges.size() == 0) {
// getData()->edges.push_back(std::vector<V3f>());
//}
getData()->edges.push_back(std::vector<V3f>());
for (auto e_it = mesh.edges_begin(); e_it != mesh.edges_end(); ++e_it)
{
// 得到边所代表的半边
OpenMesh::HalfedgeHandle heh1 = mesh.halfedge_handle(*e_it, 0); // 默认一个方向的半边
OpenMesh::Vec3d edgePoint(0, 0, 0);
int edgePointsNumber = 0;
OpenMesh::DefaultTraits::Point pointV = mesh.point(mesh.from_vertex_handle(heh1)); // 这条(半)边的起点
OpenMesh::DefaultTraits::Point pointW = mesh.point(mesh.to_vertex_handle(heh1)); // 这条(半)边的终点
getData()->edges[0].push_back({ pointV[0], pointV[1], pointV[2] });
getData()->edges[0].push_back({ pointW[0], pointW[1], pointW[2] });
}
}
算法实现请参照
https://blog.csdn.net/qq_31804159/article/details/109147352
写的很清晰了
参考链接2
https://github.com/icemiliang/loop_subdivision/blob/master/src/LOOP.cpp
http://www.cs.cmu.edu/afs/cs/academic/class/15462-s14/www/lec_slides/Subdivision.pdf
TIPS
可以看到有一个 n == 3的时候的取值, 但是Loop大佬的毕业论文中并没有发现 神秘的数字 3/16 ?? 但是也在好几篇文章中看到了这个参数.
应该是Loop毕业论文中没有提到这个参数, 是其他人加的. 也就是可加可不加??
---------------------------我的天空里没有太阳,总是黑夜,但并不暗,因为有东西代替了太阳。虽然没有太阳那么明亮,但对我来说已经足够。凭借着这份光,我便能把黑夜当成白天。我从来就没有太阳,所以不怕失去。
--------《白夜行》