From 07093b61974b2514f1698a22c41cac6b5d0d4a31 Mon Sep 17 00:00:00 2001 From: ccolin Date: Fri, 12 Nov 2021 20:53:13 +0100 Subject: [PATCH] refactor quadratic patch code --- CMakeLists.txt | 4 +- src/curvature.cpp | 8 ++-- src/curvature.h | 72 +++-------------------------------- src/mesh_processor.cpp | 6 +-- src/my_quad.h | 33 ++++++++++++++++ src/quad_patch.cpp | 63 ++++++++++-------------------- src/quad_patch.h | 26 ++++++++++--- src/quad_patch_tesselator.cpp | 44 +++++++++++++++++++++ src/quad_patch_tesselator.h | 11 ++++++ 9 files changed, 144 insertions(+), 123 deletions(-) create mode 100644 src/my_quad.h create mode 100644 src/quad_patch_tesselator.cpp create mode 100644 src/quad_patch_tesselator.h diff --git a/CMakeLists.txt b/CMakeLists.txt index 69c69dc..7f81dbf 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -26,7 +26,9 @@ target_sources(${PROJECT_NAME} PRIVATE src/curvature.cpp src/curvature.h src/quad_patch.cpp - src/quad_patch.h) + src/quad_patch.h + src/quad_patch_tesselator.cpp + src/quad_patch_tesselator.h) target_link_libraries(${PROJECT_NAME} PRIVATE Qt5::Core Qt5::Gui Qt5::Widgets OpenMeshCore diff --git a/src/curvature.cpp b/src/curvature.cpp index 0405c3d..92c487d 100644 --- a/src/curvature.cpp +++ b/src/curvature.cpp @@ -95,10 +95,8 @@ std::vector Courbures::get_two_neighborhood(const MyMesh:: // } -MyQuad Courbures::fit_quad(MyMesh::VertexHandle vh) +QuadPatch Courbures::fit_quad(MyMesh::VertexHandle vh) { - MyQuad q ; - std::vector neigh = get_two_neighborhood(vh) ; // std::cout << "-- Neighborhood" << std::endl ; @@ -163,7 +161,7 @@ MyQuad Courbures::fit_quad(MyMesh::VertexHandle vh) // Résolution aux moindres carrés par SVD Eigen::VectorXd coef(5) ; coef = A.bdcSvd(Eigen::ComputeThinU | Eigen::ComputeThinV).solve(B) ; - MyQuad q(coef[0], coef[1], coef[2], coef[3], coef[4], r, t) ; + QuadPatch q(coef, ch_base) ; // std::cout << "-- quad : " << std::endl ; // std::cout << q[0] << ", " << q[1] << ", " << q[2] << ", " << q[3] << ", " << q[4] << ", " << q[5] << ", " << std::endl ; return q ; @@ -221,7 +219,7 @@ void Courbures::compute_KH() { _mesh.add_property(vprop_quad, "vprop_quad"); for (VertexHandle vh : _mesh.vertices()) { - MyQuad q = fit_quad(vh); + QuadPatch q = fit_quad(vh); _mesh.property(vprop_quad, vh) = q; // K = det (K_P) = 4 a_0 a_2 - a_1 ^ 2 diff --git a/src/curvature.h b/src/curvature.h index 61632e9..db9b3d3 100644 --- a/src/curvature.h +++ b/src/curvature.h @@ -2,69 +2,7 @@ #define CURVATURE_H #include "my_mesh.h" - - -class MyQuad { - double _coefs[5] ; // a_0 x^2 + a1 xy + a2 y^2 + a3 x + a4 y + a5 -public: - Eigen::AngleAxisd _r ; - Eigen::Translation3d _t ; - - MyQuad(double *data, - const Eigen::AngleAxisd &r = Eigen::AngleAxisd(0, Eigen::Vector3d(0,0,1)), - const Eigen::Translation3d &t = Eigen::Translation3d(0,0,0)) - : _r(r), _t(t) - { - for (size_t i=0; i<5; i++) - _coefs[i] = data[i] ; - } - MyQuad(const Eigen::VectorXd &v, - const Eigen::AngleAxisd &r = Eigen::AngleAxisd(0, Eigen::Vector3d(0,0,1)), - const Eigen::Translation3d &t = Eigen::Translation3d(0,0,0)) - : _r(r), _t(t) - { - for (size_t i=0; i<5; i++) - _coefs[i] = v[i] ; - } - MyQuad(double a0=0, double a1=0, double a2=0, double a3=0, double a4=0, - const Eigen::AngleAxisd &r = Eigen::AngleAxisd(0, Eigen::Vector3d(0,0,1)), - const Eigen::Translation3d &t = Eigen::Translation3d(0,0,0)) - : _r(r), _t(t) - { - _coefs[0] = a0 ; - _coefs[1] = a1 ; - _coefs[2] = a2 ; - _coefs[3] = a3 ; - _coefs[4] = a4 ; - } - MyQuad(const MyQuad &q) : _r(q._r), _t(q._t) - { - for (size_t i=0; i<5; i++) - _coefs[i] = q._coefs[i] ; - } - - MyQuad &operator=(const MyQuad &q) { - _r = q._r; - _t = q._t; - for (size_t i=0; i<5; i++) - _coefs[i] = q._coefs[i] ; - return *this; - } - - double & operator[] (size_t i) {return _coefs[i] ; } - - double quad_fun (double x, double y) - { - return _coefs[0] * x*x + _coefs[1] * x*y + _coefs[2] * y*y + _coefs[3] * x + _coefs[4] * y ; - } - - double quad_fun (const OpenMesh::Vec2d &v) - { - - return quad_fun(v[0], v[1]) ; - } - -}; +#include "quad_patch.h" template @@ -145,16 +83,16 @@ private: MyMesh &_mesh ; public: - OpenMesh::VPropHandleT vprop_K; - OpenMesh::VPropHandleT vprop_H; - OpenMesh::VPropHandleT vprop_quad; + OpenMesh::VPropHandleT vprop_K; + OpenMesh::VPropHandleT vprop_H; + OpenMesh::VPropHandleT vprop_quad; Courbures(MyMesh &mesh) : _mesh(mesh) {} void set_fixed_colors() ; void normales_locales() ; std::vector get_two_neighborhood(MyMesh::VertexHandle vh); - MyQuad fit_quad(MyMesh::VertexHandle vh) ; + QuadPatch fit_quad(MyMesh::VertexHandle vh); void compute_KH() ; void set_K_colors() ; }; diff --git a/src/mesh_processor.cpp b/src/mesh_processor.cpp index f6131ab..c4a1f46 100644 --- a/src/mesh_processor.cpp +++ b/src/mesh_processor.cpp @@ -1,4 +1,4 @@ -#include "quad_patch.h" +#include "quad_patch_tesselator.h" #include "mesh_processor.h" #include "util.h" #include "smoothing.h" @@ -38,8 +38,8 @@ MeshProcessor::MeshProcessor(const QString &path) { courb.compute_KH(); courb.set_K_colors(); VertexHandle vh = mesh.vertex_handle(16); - MyQuad q = mesh.property(courb.vprop_quad, vh); - patch = quad_patch(q, mesh, vh); + QuadPatch q = mesh.property(courb.vprop_quad, vh); + patch = tesselate_quad_patch(q, mesh, vh); // mesh.holes = findHoles(mesh); // fillHoles(); // smooth(mesh); diff --git a/src/my_quad.h b/src/my_quad.h new file mode 100644 index 0000000..5a3d6a7 --- /dev/null +++ b/src/my_quad.h @@ -0,0 +1,33 @@ +#ifndef MY_QUAD_H +#define MY_QUAD_H + +#include + + +class MyQuad { + typedef Eigen::Vector Vector5d; + Vector5d _coefficients; // a₀² + a₁xy + a₂y² + a₃x + a₄y + a₅ + Eigen::AngleAxisd _rotation; + Eigen::Translation3d _translation; + +public: + MyQuad(Vector5d coefficients, + const Eigen::AngleAxisd &rotation, + const Eigen::Translation3d &translation) + :_coefficients(coefficients), + _rotation(rotation), + _translation(translation) {} + const Vector5d &coefficients() { return _coefficients; } + const Eigen::AngleAxisd &rotation() { return _rotation; } + const Eigen::Translation3d &translation() { return _translation; } + double operator()(double x, double y) { + return _coefficients[0] * x*x + + _coefficients[1] * x*y + + _coefficients[2] * y*y + + _coefficients[3] * x + + _coefficients[4] * y; + } +}; + + +#endif \ No newline at end of file diff --git a/src/quad_patch.cpp b/src/quad_patch.cpp index ec2cc03..8941f21 100644 --- a/src/quad_patch.cpp +++ b/src/quad_patch.cpp @@ -1,46 +1,25 @@ #include "quad_patch.h" -MyMesh quad_patch(MyQuad q, MyMesh mesh, VertexHandle vh) { - MyMesh patch; - Eigen::Vector3d point - (mesh.point(vh)[0], mesh.point(vh)[1], mesh.point(vh)[2]); - Eigen::Transform ch_base = q._r * q._t; - Eigen::Transform inv_ch_base = ch_base.inverse(); - point = ch_base * point; - double xmin = point[0], xmax = point[0]; - double ymin = point[1], ymax = point[1]; - for (VertexHandle vi : mesh.vv_range(vh)) { - Eigen::Vector3d point - (mesh.point(vi)[0], mesh.point(vi)[1], mesh.point(vi)[2]); - point = ch_base * point; - xmin = std::min(xmin, point[0]); - xmax = std::max(xmax, point[0]); - ymin = std::min(ymin, point[1]); - ymax = std::max(ymax, point[1]); - } - double xstep = (xmax-xmin)/64.; - double ystep = (ymax-ymin)/64.; - for (size_t y = 0; y < 64; y++) { - for (size_t x = 0; x < 64; x++) { - double dx = xmin + x * xstep; - double dy = ymin + y * ystep; - Eigen::Vector3d point(dx, dy, -q.quad_fun(dx, dy)); - point = inv_ch_base * point; - patch.new_vertex(Point(point[0], point[1], point[2])); - } - } - for (VertexHandle vhi : patch.vertices()) { - patch.set_color(vhi, MyMesh::Color(0, 1, .2)); - size_t i = vhi.idx(); - size_t x = i % 64; - size_t y = i / 64; - if (x == 63 || y == 63) continue; - patch.add_face(vhi, patch.vertex_handle(i + 64), - patch.vertex_handle(i + 1)); - patch.add_face(patch.vertex_handle(i + 1), - patch.vertex_handle(i + 64), - patch.vertex_handle(i + 65)); - } - return patch; +QuadPatch::QuadPatch() {} + + +QuadPatch::QuadPatch(const Vector5d coefficients, const Transform transform) + :_coefficients(coefficients), + _transform(transform), + _inverse_transform(transform.inverse()) { +} + + +double QuadPatch::operator()(double x, double y) { + return _coefficients[0] * x*x + + _coefficients[1] * x*y + + _coefficients[2] * y*y + + _coefficients[3] * x + + _coefficients[4] * y; +} + + +double QuadPatch::operator[](size_t i) { + return _coefficients[i]; } diff --git a/src/quad_patch.h b/src/quad_patch.h index 097f487..89a7284 100644 --- a/src/quad_patch.h +++ b/src/quad_patch.h @@ -1,11 +1,27 @@ -#ifndef QUAD_PATCH_H -#define QUAD_PATCH_H +#ifndef MY_QUAD_H +#define MY_QUAD_H -#include "curvature.h" -#include "my_mesh.h" +#include -MyMesh quad_patch(MyQuad q, MyMesh mesh, VertexHandle vh); +class QuadPatch { +public: + typedef Eigen::Vector Vector5d; + typedef Eigen::Transform Transform; + +private: + Vector5d _coefficients; // a₀² + a₁xy + a₂y² + a₃x + a₄y + a₅ + Transform _transform; + Transform _inverse_transform; + +public: + QuadPatch(); + QuadPatch(const Vector5d coefficients, const Transform transform); + double operator()(double x, double y); + double operator[](size_t i); + constexpr Transform &transform() { return _transform; } + constexpr Transform &inverse_transform() { return _inverse_transform; } +}; #endif \ No newline at end of file diff --git a/src/quad_patch_tesselator.cpp b/src/quad_patch_tesselator.cpp new file mode 100644 index 0000000..93f19d8 --- /dev/null +++ b/src/quad_patch_tesselator.cpp @@ -0,0 +1,44 @@ +#include "quad_patch_tesselator.h" + + +MyMesh tesselate_quad_patch(QuadPatch q, MyMesh mesh, VertexHandle vh) { + MyMesh patch; + Eigen::Vector3d point + (mesh.point(vh)[0], mesh.point(vh)[1], mesh.point(vh)[2]); + point = q.transform() * point; + double xmin = point[0], xmax = point[0]; + double ymin = point[1], ymax = point[1]; + for (VertexHandle vi : mesh.vv_range(vh)) { + Eigen::Vector3d point + (mesh.point(vi)[0], mesh.point(vi)[1], mesh.point(vi)[2]); + point = q.transform() * point; + xmin = std::min(xmin, point[0]); + xmax = std::max(xmax, point[0]); + ymin = std::min(ymin, point[1]); + ymax = std::max(ymax, point[1]); + } + double xstep = (xmax-xmin)/64.; + double ystep = (ymax-ymin)/64.; + for (size_t y = 0; y < 64; y++) { + for (size_t x = 0; x < 64; x++) { + double dx = xmin + x * xstep; + double dy = ymin + y * ystep; + Eigen::Vector3d point(dx, dy, -q(dx, dy)); + point = q.inverse_transform() * point; + patch.new_vertex(Point(point[0], point[1], point[2])); + } + } + for (VertexHandle vhi : patch.vertices()) { + patch.set_color(vhi, MyMesh::Color(0, 1, .2)); + size_t i = vhi.idx(); + size_t x = i % 64; + size_t y = i / 64; + if (x == 63 || y == 63) continue; + patch.add_face(vhi, patch.vertex_handle(i + 64), + patch.vertex_handle(i + 1)); + patch.add_face(patch.vertex_handle(i + 1), + patch.vertex_handle(i + 64), + patch.vertex_handle(i + 65)); + } + return patch; +} diff --git a/src/quad_patch_tesselator.h b/src/quad_patch_tesselator.h new file mode 100644 index 0000000..e03ef1b --- /dev/null +++ b/src/quad_patch_tesselator.h @@ -0,0 +1,11 @@ +#ifndef QUAD_PATCH_H +#define QUAD_PATCH_H + +#include "quad_patch.h" +#include "my_mesh.h" + + +MyMesh tesselate_quad_patch(QuadPatch q, MyMesh mesh, VertexHandle vh); + + +#endif \ No newline at end of file