Welcome to mirror list, hosted at ThFree Co, Russian Federation.

git.blender.org/blender.git - Unnamed repository; edit this file 'description' to name the repository.
summaryrefslogtreecommitdiff
path: root/extern/quadriflow/src/dedge.cpp
diff options
context:
space:
mode:
Diffstat (limited to 'extern/quadriflow/src/dedge.cpp')
-rw-r--r--extern/quadriflow/src/dedge.cpp487
1 files changed, 487 insertions, 0 deletions
diff --git a/extern/quadriflow/src/dedge.cpp b/extern/quadriflow/src/dedge.cpp
new file mode 100644
index 00000000000..b030e20d6a1
--- /dev/null
+++ b/extern/quadriflow/src/dedge.cpp
@@ -0,0 +1,487 @@
+#include "dedge.hpp"
+#include "config.hpp"
+
+#include <atomic>
+#include <fstream>
+#include <iostream>
+#include <set>
+#include <vector>
+#include "compare-key.hpp"
+#ifdef WITH_TBB
+#include "tbb/tbb.h"
+#endif
+namespace qflow {
+
+inline int dedge_prev(int e, int deg) { return (e % deg == 0u) ? e + (deg - 1) : e - 1; }
+
+inline bool atomicCompareAndExchange(volatile int* v, uint32_t newValue, int oldValue) {
+#if defined(_WIN32)
+ return _InterlockedCompareExchange(reinterpret_cast<volatile long*>(v), (long)newValue,
+ (long)oldValue) == (long)oldValue;
+#else
+ return __sync_bool_compare_and_swap(v, oldValue, newValue);
+#endif
+}
+
+const int INVALID = -1;
+
+#undef max
+#undef min
+bool compute_direct_graph(MatrixXd& V, MatrixXi& F, VectorXi& V2E, VectorXi& E2E,
+ VectorXi& boundary, VectorXi& nonManifold) {
+ V2E.resize(V.cols());
+ V2E.setConstant(INVALID);
+
+ uint32_t deg = F.rows();
+ std::vector<std::pair<uint32_t, uint32_t>> tmp(F.size());
+
+#ifdef WITH_TBB
+ tbb::parallel_for(
+ tbb::blocked_range<uint32_t>(0u, (uint32_t)F.cols(), GRAIN_SIZE),
+ [&](const tbb::blocked_range<uint32_t>& range) {
+ for (uint32_t f = range.begin(); f != range.end(); ++f) {
+ for (uint32_t i = 0; i < deg; ++i) {
+ uint32_t idx_cur = F(i, f), idx_next = F((i + 1) % deg, f),
+ edge_id = deg * f + i;
+ if (idx_cur >= V.cols() || idx_next >= V.cols())
+ throw std::runtime_error(
+ "Mesh data contains an out-of-bounds vertex reference!");
+ if (idx_cur == idx_next) continue;
+
+ tmp[edge_id] = std::make_pair(idx_next, INVALID);
+ if (!atomicCompareAndExchange(&V2E[idx_cur], edge_id, INVALID)) {
+ uint32_t idx = V2E[idx_cur];
+ while (!atomicCompareAndExchange((int*)&tmp[idx].second, edge_id, INVALID))
+ idx = tmp[idx].second;
+ }
+ }
+ }
+ });
+#else
+ for (int f = 0; f < F.cols(); ++f) {
+ for (unsigned int i = 0; i < deg; ++i) {
+ unsigned int idx_cur = F(i, f), idx_next = F((i + 1) % deg, f), edge_id = deg * f + i;
+ if (idx_cur >= V.cols() || idx_next >= V.cols())
+ throw std::runtime_error("Mesh data contains an out-of-bounds vertex reference!");
+ if (idx_cur == idx_next) continue;
+
+ tmp[edge_id] = std::make_pair(idx_next, -1);
+ if (V2E[idx_cur] == -1)
+ V2E[idx_cur] = edge_id;
+ else {
+ unsigned int idx = V2E[idx_cur];
+ while (tmp[idx].second != -1) {
+ idx = tmp[idx].second;
+ }
+ tmp[idx].second = edge_id;
+ }
+ }
+ }
+#endif
+
+ nonManifold.resize(V.cols());
+ nonManifold.setConstant(false);
+
+ E2E.resize(F.cols() * deg);
+ E2E.setConstant(INVALID);
+
+#ifdef WITH_OMP
+#pragma omp parallel for
+#endif
+ for (int f = 0; f < F.cols(); ++f) {
+ for (uint32_t i = 0; i < deg; ++i) {
+ uint32_t idx_cur = F(i, f), idx_next = F((i + 1) % deg, f), edge_id_cur = deg * f + i;
+
+ if (idx_cur == idx_next) continue;
+
+ uint32_t it = V2E[idx_next], edge_id_opp = INVALID;
+ while (it != INVALID) {
+ if (tmp[it].first == idx_cur) {
+ if (edge_id_opp == INVALID) {
+ edge_id_opp = it;
+ } else {
+ nonManifold[idx_cur] = true;
+ nonManifold[idx_next] = true;
+ edge_id_opp = INVALID;
+ break;
+ }
+ }
+ it = tmp[it].second;
+ }
+
+ if (edge_id_opp != INVALID && edge_id_cur < edge_id_opp) {
+ E2E[edge_id_cur] = edge_id_opp;
+ E2E[edge_id_opp] = edge_id_cur;
+ }
+ }
+ }
+ std::atomic<uint32_t> nonManifoldCounter(0), boundaryCounter(0), isolatedCounter(0);
+
+ boundary.resize(V.cols());
+ boundary.setConstant(false);
+
+ /* Detect boundary regions of the mesh and adjust vertex->edge pointers*/
+#ifdef WITH_OMP
+#pragma omp parallel for
+#endif
+ for (int i = 0; i < V.cols(); ++i) {
+ uint32_t edge = V2E[i];
+ if (edge == INVALID) {
+ isolatedCounter++;
+ continue;
+ }
+ if (nonManifold[i]) {
+ nonManifoldCounter++;
+ V2E[i] = INVALID;
+ continue;
+ }
+
+ /* Walk backwards to the first boundary edge (if any) */
+ uint32_t start = edge, v2e = INVALID;
+ do {
+ v2e = std::min(v2e, edge);
+ uint32_t prevEdge = E2E[dedge_prev(edge, deg)];
+ if (prevEdge == INVALID) {
+ /* Reached boundary -- update the vertex->edge link */
+ v2e = edge;
+ boundary[i] = true;
+ boundaryCounter++;
+ break;
+ }
+ edge = prevEdge;
+ } while (edge != start);
+ V2E[i] = v2e;
+ }
+#ifdef LOG_OUTPUT
+ printf("counter triangle %d %d\n", (int)boundaryCounter, (int)nonManifoldCounter);
+#endif
+ return true;
+ std::vector<std::vector<int>> vert_to_edges(V2E.size());
+ for (int i = 0; i < F.cols(); ++i) {
+ for (int j = 0; j < 3; ++j) {
+ int v = F(j, i);
+ vert_to_edges[v].push_back(i * 3 + j);
+ }
+ }
+ std::vector<int> colors(F.cols() * 3, -1);
+ bool update = false;
+ int num_v = V.cols();
+ std::map<int, int> new_vertices;
+ for (int i = 0; i < vert_to_edges.size(); ++i) {
+ int num_color = 0;
+ for (int j = 0; j < vert_to_edges[i].size(); ++j) {
+ int deid0 = vert_to_edges[i][j];
+ if (colors[deid0] == -1) {
+ int deid = deid0;
+ do {
+ colors[deid] = num_color;
+ if (num_color != 0) F(deid % 3, deid / 3) = num_v;
+ deid = deid / 3 * 3 + (deid + 2) % 3;
+ deid = E2E[deid];
+ } while (deid != deid0);
+ num_color += 1;
+ if (num_color > 1) {
+ update = true;
+ new_vertices[num_v] = i;
+ num_v += 1;
+ }
+ }
+ }
+ }
+ if (update) {
+ V.conservativeResize(3, num_v);
+ for (auto& p : new_vertices) {
+ V.col(p.first) = V.col(p.second);
+ }
+ return false;
+ }
+ return true;
+}
+
+void compute_direct_graph_quad(std::vector<Vector3d>& V, std::vector<Vector4i>& F, std::vector<int>& V2E, std::vector<int>& E2E, VectorXi& boundary, VectorXi& nonManifold) {
+ V2E.clear();
+ E2E.clear();
+ boundary = VectorXi();
+ nonManifold = VectorXi();
+ V2E.resize(V.size(), INVALID);
+
+ uint32_t deg = 4;
+ std::vector<std::pair<uint32_t, uint32_t>> tmp(F.size() * deg);
+
+#ifdef WITH_TBB
+ tbb::parallel_for(
+ tbb::blocked_range<uint32_t>(0u, (uint32_t)F.size(), GRAIN_SIZE),
+ [&](const tbb::blocked_range<uint32_t>& range) {
+ for (uint32_t f = range.begin(); f != range.end(); ++f) {
+ for (uint32_t i = 0; i < deg; ++i) {
+ uint32_t idx_cur = F[f][i], idx_next = F[f][(i + 1) % deg],
+ edge_id = deg * f + i;
+ if (idx_cur >= V.size() || idx_next >= V.size())
+ throw std::runtime_error(
+ "Mesh data contains an out-of-bounds vertex reference!");
+ if (idx_cur == idx_next) continue;
+
+ tmp[edge_id] = std::make_pair(idx_next, INVALID);
+ if (!atomicCompareAndExchange(&V2E[idx_cur], edge_id, INVALID)) {
+ uint32_t idx = V2E[idx_cur];
+ while (!atomicCompareAndExchange((int*)&tmp[idx].second, edge_id, INVALID))
+ idx = tmp[idx].second;
+ }
+ }
+ }
+ });
+#else
+ for (int f = 0; f < F.size(); ++f) {
+ for (unsigned int i = 0; i < deg; ++i) {
+ unsigned int idx_cur = F[f][i], idx_next = F[f][(i + 1) % deg], edge_id = deg * f + i;
+ if (idx_cur >= V.size() || idx_next >= V.size())
+ throw std::runtime_error("Mesh data contains an out-of-bounds vertex reference!");
+ if (idx_cur == idx_next) continue;
+ tmp[edge_id] = std::make_pair(idx_next, -1);
+ if (V2E[idx_cur] == -1) {
+ V2E[idx_cur] = edge_id;
+ }
+ else {
+ unsigned int idx = V2E[idx_cur];
+ while (tmp[idx].second != -1) {
+ idx = tmp[idx].second;
+ }
+ tmp[idx].second = edge_id;
+ }
+ }
+ }
+#endif
+ nonManifold.resize(V.size());
+ nonManifold.setConstant(false);
+
+ E2E.resize(F.size() * deg, INVALID);
+
+#ifdef WITH_OMP
+#pragma omp parallel for
+#endif
+ for (int f = 0; f < F.size(); ++f) {
+ for (uint32_t i = 0; i < deg; ++i) {
+ uint32_t idx_cur = F[f][i], idx_next = F[f][(i + 1) % deg], edge_id_cur = deg * f + i;
+
+ if (idx_cur == idx_next) continue;
+
+ uint32_t it = V2E[idx_next], edge_id_opp = INVALID;
+ while (it != INVALID) {
+ if (tmp[it].first == idx_cur) {
+ if (edge_id_opp == INVALID) {
+ edge_id_opp = it;
+ } else {
+ nonManifold[idx_cur] = true;
+ nonManifold[idx_next] = true;
+ edge_id_opp = INVALID;
+ break;
+ }
+ }
+ it = tmp[it].second;
+ }
+
+ if (edge_id_opp != INVALID && edge_id_cur < edge_id_opp) {
+ E2E[edge_id_cur] = edge_id_opp;
+ E2E[edge_id_opp] = edge_id_cur;
+ }
+ }
+ }
+ std::atomic<uint32_t> nonManifoldCounter(0), boundaryCounter(0), isolatedCounter(0);
+
+ boundary.resize(V.size());
+ boundary.setConstant(false);
+
+ /* Detect boundary regions of the mesh and adjust vertex->edge pointers*/
+#ifdef WITH_OMP
+#pragma omp parallel for
+#endif
+ for (int i = 0; i < V.size(); ++i) {
+ uint32_t edge = V2E[i];
+ if (edge == INVALID) {
+ isolatedCounter++;
+ continue;
+ }
+ if (nonManifold[i]) {
+ nonManifoldCounter++;
+ V2E[i] = INVALID;
+ continue;
+ }
+
+ /* Walk backwards to the first boundary edge (if any) */
+ uint32_t start = edge, v2e = INVALID;
+ do {
+ v2e = std::min(v2e, edge);
+ uint32_t prevEdge = E2E[dedge_prev(edge, deg)];
+ if (prevEdge == INVALID) {
+ /* Reached boundary -- update the vertex->edge link */
+ v2e = edge;
+ boundary[i] = true;
+ boundaryCounter++;
+ break;
+ }
+ edge = prevEdge;
+ } while (edge != start);
+ V2E[i] = v2e;
+ }
+#ifdef LOG_OUTPUT
+ printf("counter %d %d\n", (int)boundaryCounter, (int)nonManifoldCounter);
+#endif
+}
+
+void remove_nonmanifold(std::vector<Vector4i>& F, std::vector<Vector3d>& V) {
+ typedef std::pair<uint32_t, uint32_t> Edge;
+
+ int degree = 4;
+ std::map<uint32_t, std::map<uint32_t, std::pair<uint32_t, uint32_t>>> irregular;
+ std::vector<std::set<int>> E(V.size());
+ std::vector<std::set<int>> VF(V.size());
+
+ auto kill_face_single = [&](uint32_t f) {
+ if (F[f][0] == INVALID) return;
+ for (int i = 0; i < degree; ++i) E[F[f][i]].erase(F[f][(i + 1) % degree]);
+ F[f].setConstant(INVALID);
+ };
+
+ auto kill_face = [&](uint32_t f) {
+ if (degree == 4 && F[f][2] == F[f][3]) {
+ auto it = irregular.find(F[f][2]);
+ if (it != irregular.end()) {
+ for (auto& item : it->second) {
+ kill_face_single(item.second.second);
+ }
+ }
+ }
+ kill_face_single(f);
+ };
+
+ uint32_t nm_edge = 0, nm_vert = 0;
+
+ for (uint32_t f = 0; f < (uint32_t)F.size(); ++f) {
+ if (F[f][0] == INVALID) continue;
+ if (degree == 4 && F[f][2] == F[f][3]) {
+ /* Special handling of irregular faces */
+ irregular[F[f][2]][F[f][0]] = std::make_pair(F[f][1], f);
+ continue;
+ }
+
+ bool nonmanifold = false;
+ for (uint32_t e = 0; e < degree; ++e) {
+ uint32_t v0 = F[f][e], v1 = F[f][(e + 1) % degree], v2 = F[f][(e + 2) % degree];
+ if (E[v0].find(v1) != E[v0].end() || (degree == 4 && E[v0].find(v2) != E[v0].end()))
+ nonmanifold = true;
+ }
+
+ if (nonmanifold) {
+ nm_edge++;
+ F[f].setConstant(INVALID);
+ continue;
+ }
+
+ for (uint32_t e = 0; e < degree; ++e) {
+ uint32_t v0 = F[f][e], v1 = F[f][(e + 1) % degree], v2 = F[f][(e + 2) % degree];
+
+ E[v0].insert(v1);
+ if (degree == 4) E[v0].insert(v2);
+ VF[v0].insert(f);
+ }
+ }
+
+ std::vector<Edge> edges;
+ for (auto item : irregular) {
+ bool nonmanifold = false;
+ auto face = item.second;
+ edges.clear();
+
+ uint32_t cur = face.begin()->first, stop = cur;
+ while (true) {
+ uint32_t pred = cur;
+ cur = face[cur].first;
+ uint32_t next = face[cur].first, it = 0;
+ while (true) {
+ ++it;
+ if (next == pred) break;
+ if (E[cur].find(next) != E[cur].end() && it == 1) nonmanifold = true;
+ edges.push_back(Edge(cur, next));
+ next = face[next].first;
+ }
+ if (cur == stop) break;
+ }
+
+ if (nonmanifold) {
+ nm_edge++;
+ for (auto& i : item.second) F[i.second.second].setConstant(INVALID);
+ continue;
+ } else {
+ for (auto e : edges) {
+ E[e.first].insert(e.second);
+
+ for (auto e2 : face) VF[e.first].insert(e2.second.second);
+ }
+ }
+ }
+
+ /* Check vertices */
+ std::set<uint32_t> v_marked, v_unmarked, f_adjacent;
+
+ std::function<void(uint32_t)> dfs = [&](uint32_t i) {
+ v_marked.insert(i);
+ v_unmarked.erase(i);
+
+ for (uint32_t f : VF[i]) {
+ if (f_adjacent.find(f) == f_adjacent.end()) /* if not part of adjacent face */
+ continue;
+ for (uint32_t j = 0; j < degree; ++j) {
+ uint32_t k = F[f][j];
+ if (v_unmarked.find(k) == v_unmarked.end() || /* if not unmarked OR */
+ v_marked.find(k) != v_marked.end()) /* if already marked */
+ continue;
+ dfs(k);
+ }
+ }
+ };
+
+ for (uint32_t i = 0; i < (uint32_t)V.size(); ++i) {
+ v_marked.clear();
+ v_unmarked.clear();
+ f_adjacent.clear();
+
+ for (uint32_t f : VF[i]) {
+ if (F[f][0] == INVALID) continue;
+
+ for (uint32_t k = 0; k < degree; ++k) v_unmarked.insert(F[f][k]);
+
+ f_adjacent.insert(f);
+ }
+
+ if (v_unmarked.empty()) continue;
+ v_marked.insert(i);
+ v_unmarked.erase(i);
+
+ dfs(*v_unmarked.begin());
+
+ if (v_unmarked.size() > 0) {
+ nm_vert++;
+ for (uint32_t f : f_adjacent) kill_face(f);
+ }
+ }
+
+ if (nm_vert > 0 || nm_edge > 0) {
+ std::cout << "Non-manifold elements: vertices=" << nm_vert << ", edges=" << nm_edge
+ << std::endl;
+ }
+ uint32_t nFaces = 0, nFacesOrig = F.size();
+ for (uint32_t f = 0; f < (uint32_t)F.size(); ++f) {
+ if (F[f][0] == INVALID) continue;
+ if (nFaces != f) {
+ F[nFaces] = F[f];
+ }
+ ++nFaces;
+ }
+
+ if (nFacesOrig != nFaces) {
+ F.resize(nFaces);
+ std::cout << "Faces reduced from " << nFacesOrig << " -> " << nFaces << std::endl;
+ }
+}
+
+} // namespace qflow
generated by cgit v1.2.3 (git 2.25.1) at 2024-07-20 23:15:48 +0300