mirror of https://github.com/CGAL/cgal
Merge branch 'master' into gsoc17-summerwork
Reworking of 4e4d495ac2
This commit is contained in:
Sébastien Loriot
commit
4824f3bcc9
|
|
@ -1,5 +1,6 @@
|
|||
/*build*
|
||||
/*/*/*/build
|
||||
/*/*/build*
|
||||
AABB_tree/demo/AABB_tree/AABB_demo
|
||||
AABB_tree/demo/AABB_tree/Makefile
|
||||
AABB_tree/examples/AABB_tree/*.kdev*
|
||||
|
|
|
|||
|
|
@ -2838,6 +2838,25 @@ pages = "207--221"
|
|||
year={2015}
|
||||
}
|
||||
|
||||
@article{cgal:sg-hqct-04,
|
||||
title={High quality compatible triangulations},
|
||||
author={V. Surazhsky and C. Gotsman},
|
||||
journal={Engineering with Computers},
|
||||
volume={20},
|
||||
number={2},
|
||||
pages={147--156},
|
||||
year={2004},
|
||||
publisher={springer}
|
||||
}
|
||||
|
||||
@inproceedings{cgal:dmsb-ifamdcf-99
|
||||
,author = {M. Desbrun and M. Meyer and P. Schr{\"o}der and A. H. Barr}
|
||||
,title = {Implicit Fairing of Arbitrary Meshes using
|
||||
Diffusion and Curvature Flow}
|
||||
,booktitle = "Computer Graphics (Proc. SIGGRAPH '99)"
|
||||
,year = 1999
|
||||
,pages = {317--324}
|
||||
}
|
||||
% ----------------------------------------------------------------------------
|
||||
% END OF BIBFILE
|
||||
% ----------------------------------------------------------------------------
|
||||
|
|
|
|||
|
|
@ -367,6 +367,32 @@ and perturbation is not deterministic
|
|||
\b Default is that this variable is not set
|
||||
\cgalNPEnd
|
||||
|
||||
\cgalNPBegin{use_weights} \anchor PMP_use_weights
|
||||
Parameter used in `angle_remeshing()` to set whether angles should carry weights while remeshing
|
||||
so that small angles carry more weight.
|
||||
\n
|
||||
\b Type : `bool` \n
|
||||
\b Default is `false`
|
||||
\cgalNPEnd
|
||||
|
||||
\cgalNPBegin{gradient_descent_precision} \anchor PMP_gradient_descent_precision
|
||||
Parameter used in `area_remeshing()` to set precision to be met while the relative energy
|
||||
between iterations of each triagle is minimized. Triangle energy is defined based on its area compared to the
|
||||
average area of all triagnles adjacent to the vertex that is being moved.
|
||||
\n
|
||||
\b Type : `double` \n
|
||||
\b Default is `0.001`
|
||||
\cgalNPEnd
|
||||
|
||||
\cgalNPBegin{distance_precision} \anchor PMP_distance_precision
|
||||
Parameter used in `compatible_remeshing()` to set precision for the Hausdorff distance
|
||||
of the mesh between the previous and the current iteration that has to be achieved for a mesh to converge
|
||||
during smoothing.
|
||||
\n
|
||||
\b Type : `double` \n
|
||||
\b Default is `0.01`
|
||||
\cgalNPEnd
|
||||
|
||||
\cgalNPTableEnd
|
||||
|
||||
*/
|
||||
|
|
|
|||
|
|
@ -89,6 +89,10 @@ and provides a list of the parameters that are used in this package.
|
|||
## Meshing Functions ##
|
||||
- `CGAL::Polygon_mesh_processing::fair()`
|
||||
- `CGAL::Polygon_mesh_processing::refine()`
|
||||
- `CGAL::Polygon_mesh_processing::angle_remeshing()`
|
||||
- `CGAL::Polygon_mesh_processing::area_remeshing()`
|
||||
- `CGAL::Polygon_mesh_processing::compatible_remeshing()`
|
||||
- `CGAL::Polygon_mesh_processing::curvature_flow()`
|
||||
- `CGAL::Polygon_mesh_processing::triangulate_face()`
|
||||
- `CGAL::Polygon_mesh_processing::triangulate_faces()`
|
||||
- \link PMP_meshing_grp `CGAL::Polygon_mesh_processing::isotropic_remeshing()` \endlink
|
||||
|
|
|
|||
|
|
@ -37,7 +37,8 @@ BGL \ref namedparameters are used to deal with optional parameters.
|
|||
\subsection PMPOutline Outline
|
||||
The algorithms described in this manual are organized in sections:
|
||||
- \ref PMPMeshing : meshing algorithms, including triangulation of non-triangulated
|
||||
meshes, refinement, optimization by fairing, and isotropic remeshing of triangulated surface meshes.
|
||||
meshes, refinement, optimization by fairing, isotropic remeshing of triangulated surface meshes
|
||||
and smoothing algorithms.
|
||||
- \ref Coref_section : methods to corefine triangle meshes and to compute
|
||||
boolean operations out of corefined closed triangle meshes.
|
||||
- \ref PMPHoleFilling : available hole filling algorithms, which can possibly be combined with refinement and fairing.
|
||||
|
|
@ -120,6 +121,47 @@ Isotropic remeshing. (a) Triangulated input surface mesh.
|
|||
(d) Surface mesh with the selection uniformly remeshed.
|
||||
\cgalFigureEnd
|
||||
|
||||
\subsubsection Smoothing
|
||||
|
||||
Smoothing of a triangulated mesh area can be achieved with algorithms that aim at either mesh smoothing or shape smoothing.
|
||||
While mesh smoothing is achieved by improving the quality of the triangulated elements based on criteria such as angle and area,
|
||||
shape smoothing is designed to be "intrinsic", depending as little as possible to the discretization
|
||||
and smoothing the shape alone without optimizing the mesh elements.
|
||||
|
||||
For mesh smoothing, `CGAL::Polygon_mesh_processing::angle_smoothing()` may be used to
|
||||
move vertices so that the angles of its incident edges equalize and `CGAL::Polygon_mesh_processing::area_smoothing()` may be used to
|
||||
move vertices so that the areas of adjacent triangles equalize. Angle and area smoothing are used combined in
|
||||
`CGAL::Polygon_mesh_processing::compatible_smoothing()` which smooths the mesh until convergence is achieved.
|
||||
For shape smoothing, `CGAL::Polygon_mesh_processing::curvature_flow()`
|
||||
takes into account the mean curvature to calculate a weighted barycenter towards which vertices move incrementally.
|
||||
In all cases border vertices are considered constrained and do not move at any step of the procedure.
|
||||
Also, new vertices are not inserted while mesh or shape smoothing, although degenerate faces are removed as a first step
|
||||
before the smoothing.
|
||||
|
||||
Angle and area smoothing algorithms are based on Surazhsky and Gotsman \cgalCite{cgal:sg-hqct-04}. When mesh connectivity
|
||||
is highly irregular, angle smoothing may be used with angle weights so that small angles carry
|
||||
more weight and with almost the same computational cost per iteration convergence may be achieved faster.
|
||||
Area smoothing, since it considers only triangle areas as a smoothing criterion, may result in long and skinny triangles,
|
||||
especially for meshes with highly irregular connectivity. As such, area smoothing is guaranteed to improve the spatial
|
||||
distribution of the vertices over the area that is being remeshed.
|
||||
|
||||
The curvature flow algorithm for shape smoothing is based on Desbrun et al \cgalCite{cgal:dmsb-ifamdcf-99}.
|
||||
Curvature flow smoothes the surface by moving along the surface normal with a speed equal to the mean curvature of the
|
||||
area that is being smoothed. This means that if the surface is flat around a vertex, this vertex will not move (zero curvature).
|
||||
Curvature flow smoothing offers best results in closed surface meshes, where vertices slide according to the curvature
|
||||
and smooth the shape while the overall volume is being kept constant.
|
||||
|
||||
\cgalFigureBegin{elephants_and_histograms, elephants_and_histograms.png}
|
||||
Mesh smoothing. (a) Elephant before mesh smoothing and
|
||||
(b) elephant after applying compatible smoothing for 20 iterations.
|
||||
and histograms of their aspect ratios beneath them. Aspect ratio is defined as the ratio
|
||||
of the largest altidute of a triangle to the length of its shortest edge.
|
||||
\cgalFigureEnd
|
||||
|
||||
\cgalFigureBegin{mannequins, mannequins.png}
|
||||
Shape smoothing. (a) Devil before shape smoothing and
|
||||
(b) the devil after applying the curvature flow algorithm.
|
||||
\cgalFigureEnd
|
||||
|
||||
\subsection MeshingExamples Meshing Examples
|
||||
|
||||
|
|
|
|||
|
|
@ -20,4 +20,5 @@
|
|||
\example Polygon_mesh_processing/corefinement_mesh_union.cpp
|
||||
\example Polygon_mesh_processing/corefinement_consecutive_bool_op.cpp
|
||||
\example Polygon_mesh_processing/detect_features_example.cpp
|
||||
\example Polygon_mesh_processing/smoothing_example.cpp
|
||||
*/
|
||||
|
|
|
|||
|
|
@ -107,6 +107,7 @@ create_single_source_cgal_program( "random_perturbation_SM_example.cpp" )
|
|||
create_single_source_cgal_program( "corefinement_LCC.cpp")
|
||||
create_single_source_cgal_program( "hole_filling_example_LCC.cpp" )
|
||||
create_single_source_cgal_program( "detect_features_example.cpp" )
|
||||
create_single_source_cgal_program( "smoothing_example.cpp")
|
||||
|
||||
if(OpenMesh_FOUND)
|
||||
create_single_source_cgal_program( "compute_normals_example_OM.cpp" )
|
||||
|
|
|
|||
|
|
@ -0,0 +1,32 @@
|
|||
#include <iostream>
|
||||
#include <fstream>
|
||||
|
||||
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
|
||||
#include <CGAL/Surface_mesh.h>
|
||||
|
||||
#include <CGAL/Polygon_mesh_processing/smoothing.h>
|
||||
|
||||
|
||||
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
|
||||
typedef CGAL::Surface_mesh<K::Point_3> Mesh;
|
||||
|
||||
int main(int argc, char* argv[]){
|
||||
|
||||
const char* filename = "data/eight.off";
|
||||
std::ifstream input(filename);
|
||||
|
||||
Mesh mesh;
|
||||
if (!input || !(input >> mesh) || mesh.is_empty()) {
|
||||
std::cerr << "Not a valid .off file." << std::endl;
|
||||
return 1;
|
||||
}
|
||||
|
||||
CGAL::Polygon_mesh_processing::compatible_smoothing(mesh);
|
||||
|
||||
std::ofstream output("data/eight_smoothed.off");
|
||||
output << mesh;
|
||||
output.close();
|
||||
|
||||
|
||||
return 0;
|
||||
}
|
||||
|
|
@ -0,0 +1,318 @@
|
|||
// Copyright (c) 2015 GeometryFactory (France).
|
||||
// All rights reserved.
|
||||
//
|
||||
// This file is part of CGAL (www.cgal.org).
|
||||
// You can redistribute it and/or modify it under the terms of the GNU
|
||||
// General Public License as published by the Free Software Foundation,
|
||||
// either version 3 of the License, or (at your option) any later version.
|
||||
//
|
||||
// Licensees holding a valid commercial license may use this file in
|
||||
// accordance with the commercial license agreement provided with the software.
|
||||
//
|
||||
// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
|
||||
// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
|
||||
//
|
||||
// $URL$
|
||||
// $Id$
|
||||
//
|
||||
//
|
||||
// Author(s) : Konstantinos Katrioplas (konst.katrioplas@gmail.com)
|
||||
|
||||
#ifndef CURVATURE_FLOW_IMPL_H
|
||||
#define CURVATURE_FLOW_IMPL_H
|
||||
|
||||
#include <utility>
|
||||
#include <math.h>
|
||||
|
||||
#include <CGAL/Polygon_mesh_processing/Weights.h>
|
||||
#include <CGAL/Polygon_mesh_processing/compute_normal.h>
|
||||
#include <CGAL/Polygon_mesh_processing/repair.h>
|
||||
#include <CGAL/Polygon_mesh_processing/measure.h>
|
||||
|
||||
#include <CGAL/property_map.h>
|
||||
#include <CGAL/iterator.h>
|
||||
#include <CGAL/boost/graph/Euler_operations.h>
|
||||
#include <boost/graph/graph_traits.hpp>
|
||||
#include <boost/foreach.hpp>
|
||||
|
||||
|
||||
namespace CGAL {
|
||||
|
||||
namespace Polygon_mesh_processing {
|
||||
|
||||
namespace internal {
|
||||
|
||||
|
||||
|
||||
template<typename PolygonMesh, typename VertexPointMap,
|
||||
typename CotangentValue = CGAL::internal::Cotangent_value_Meyer<PolygonMesh, VertexPointMap>>
|
||||
class Cotangent_weight : CotangentValue
|
||||
{
|
||||
public:
|
||||
Cotangent_weight(PolygonMesh& pmesh_, VertexPointMap vpmap_)
|
||||
: CotangentValue(pmesh_, vpmap_)
|
||||
{}
|
||||
|
||||
PolygonMesh& pmesh()
|
||||
{
|
||||
return CotangentValue::pmesh();
|
||||
}
|
||||
|
||||
typedef typename boost::graph_traits<PolygonMesh>::halfedge_descriptor halfedge_descriptor;
|
||||
typedef typename boost::graph_traits<PolygonMesh>::vertex_descriptor vertex_descriptor;
|
||||
typedef std::pair<halfedge_descriptor, halfedge_descriptor> he_pair;
|
||||
|
||||
double operator()(halfedge_descriptor he, he_pair incd_edges)
|
||||
{
|
||||
vertex_descriptor vs = source(he, pmesh());
|
||||
vertex_descriptor vt = target(he, pmesh());
|
||||
vertex_descriptor v1 = target(incd_edges.first, pmesh());
|
||||
vertex_descriptor v2 = source(incd_edges.second, pmesh());
|
||||
|
||||
return ( CotangentValue::operator()(vs, v1, vt) + CotangentValue::operator()(vs, v2, vt) );
|
||||
}
|
||||
|
||||
};
|
||||
|
||||
|
||||
template<typename PolygonMesh, typename VertexPointMap, typename VertexConstraintMap, typename EdgeConstraintMap, typename GeomTraits>
|
||||
class Curvature_flow
|
||||
{
|
||||
typedef typename boost::graph_traits<PolygonMesh>::vertex_descriptor vertex_descriptor;
|
||||
typedef typename boost::graph_traits<PolygonMesh>::halfedge_descriptor halfedge_descriptor;
|
||||
typedef typename boost::graph_traits<PolygonMesh>::face_descriptor face_descriptor;
|
||||
typedef typename boost::graph_traits<PolygonMesh>::edge_descriptor edge_descriptor;
|
||||
|
||||
typedef typename GeomTraits::Point_3 Point;
|
||||
typedef typename GeomTraits::Vector_3 Vector;
|
||||
typedef typename GeomTraits::Triangle_3 Triangle;
|
||||
typedef std::vector<Triangle> Triangle_list;
|
||||
|
||||
typedef std::pair<halfedge_descriptor, halfedge_descriptor> he_pair;
|
||||
typedef std::map<halfedge_descriptor, he_pair> Edges_around_map;
|
||||
|
||||
typedef Cotangent_weight<PolygonMesh, VertexPointMap> Weight_calculator;
|
||||
|
||||
public:
|
||||
|
||||
Curvature_flow(PolygonMesh& pmesh, VertexPointMap& vpmap, VertexConstraintMap& vcmap, EdgeConstraintMap& ecmap) :
|
||||
mesh_(pmesh), vpmap_(vpmap), vcmap_(vcmap), ecmap_(ecmap),
|
||||
weight_calculator_(pmesh, vpmap)
|
||||
{}
|
||||
|
||||
template<typename FaceRange>
|
||||
void init_smoothing(const FaceRange& face_range)
|
||||
{
|
||||
check_vertex_range(face_range);
|
||||
|
||||
check_constraints();
|
||||
|
||||
BOOST_FOREACH(face_descriptor f, face_range)
|
||||
input_triangles_.push_back(triangle(f));
|
||||
}
|
||||
|
||||
std::size_t remove_degenerate_faces()
|
||||
{
|
||||
std::size_t nb_removed_faces = 0;
|
||||
nb_removed_faces = CGAL::Polygon_mesh_processing::remove_degenerate_faces(mesh_);
|
||||
|
||||
return nb_removed_faces;
|
||||
}
|
||||
|
||||
void curvature_smoothing()
|
||||
{
|
||||
std::map<vertex_descriptor, Point> barycenters;
|
||||
|
||||
BOOST_FOREACH(vertex_descriptor v, vrange)
|
||||
{
|
||||
if(!is_border(v, mesh_) && !is_constrained(v))
|
||||
{
|
||||
// find incident halfedges
|
||||
Edges_around_map he_map;
|
||||
typename Edges_around_map::iterator it;
|
||||
BOOST_FOREACH(halfedge_descriptor hi, halfedges_around_source(v, mesh_))
|
||||
he_map[hi] = he_pair( next(hi, mesh_), prev(opposite(hi, mesh_), mesh_) );
|
||||
|
||||
// calculate movement
|
||||
Vector curvature_normal = CGAL::NULL_VECTOR;
|
||||
double sum_cot_weights = 0;
|
||||
for(it = he_map.begin(); it!= he_map.end(); ++it)
|
||||
{
|
||||
halfedge_descriptor hi = it->first;
|
||||
he_pair incd_edges = it->second;
|
||||
|
||||
//check_degeneracy(hi);
|
||||
|
||||
// weight
|
||||
double weight = weight_calculator_(hi, incd_edges);;
|
||||
sum_cot_weights += weight;
|
||||
|
||||
// displacement vector
|
||||
Point xi = get(vpmap_, source(hi, mesh_));
|
||||
Point xj = get(vpmap_, target(hi, mesh_));
|
||||
Vector vec(xj, xi); // towards the vertex that is being moved
|
||||
|
||||
// add weight
|
||||
vec *= weight;
|
||||
|
||||
// sum vecs
|
||||
curvature_normal += vec;
|
||||
}
|
||||
|
||||
// divide with total weight
|
||||
if(sum_cot_weights != 0)
|
||||
curvature_normal /= sum_cot_weights;
|
||||
|
||||
Point weighted_barycenter = get(vpmap_, v) - curvature_normal;
|
||||
barycenters[v] = weighted_barycenter;
|
||||
|
||||
} // not on border
|
||||
} // all vertices
|
||||
|
||||
// update location
|
||||
typedef typename std::map<vertex_descriptor, Point>::value_type VP;
|
||||
BOOST_FOREACH(const VP& vp, barycenters)
|
||||
put(vpmap_, vp.first, vp.second);
|
||||
|
||||
}
|
||||
|
||||
|
||||
private:
|
||||
|
||||
// helper functions
|
||||
// ----------------
|
||||
Triangle triangle(face_descriptor f) const
|
||||
{
|
||||
halfedge_descriptor h = halfedge(f, mesh_);
|
||||
vertex_descriptor v1 = target(h, mesh_);
|
||||
vertex_descriptor v2 = target(next(h, mesh_), mesh_);
|
||||
vertex_descriptor v3 = target(next(next(h, mesh_), mesh_), mesh_);
|
||||
return Triangle(get(vpmap_, v1), get(vpmap_, v2), get(vpmap_, v3));
|
||||
}
|
||||
|
||||
double sqlength(const vertex_descriptor& v1, const vertex_descriptor& v2) const
|
||||
{
|
||||
return to_double(CGAL::squared_distance(get(vpmap_, v1), get(vpmap_, v2)));
|
||||
}
|
||||
|
||||
double sqlength(const halfedge_descriptor& h) const
|
||||
{
|
||||
vertex_descriptor v1 = target(h, mesh_);
|
||||
vertex_descriptor v2 = source(h, mesh_);
|
||||
return sqlength(v1, v2);
|
||||
}
|
||||
|
||||
double sqlength(const edge_descriptor& e) const
|
||||
{
|
||||
return sqlength(halfedge(e, mesh_));
|
||||
}
|
||||
|
||||
|
||||
// degeneracy removal
|
||||
// ------------------
|
||||
void check_degeneracy(halfedge_descriptor h1)
|
||||
{
|
||||
halfedge_descriptor h2 = next(h1, mesh_);
|
||||
halfedge_descriptor h3 = next(h2, mesh_);
|
||||
|
||||
double a1 = get_angle(h1, h2);
|
||||
double a2 = get_angle(h2, h3);
|
||||
double a3 = get_angle(h3, h1);
|
||||
|
||||
double angle_min_threshold = 0.05; // rad
|
||||
double angle_max_threshold = CGAL_PI - 0.05;
|
||||
|
||||
if(a1 < angle_min_threshold || a2 < angle_min_threshold || a3 < angle_min_threshold)
|
||||
{
|
||||
Euler::remove_face(h1, mesh_);
|
||||
}
|
||||
|
||||
if(a1 > angle_max_threshold || a2 > angle_max_threshold || a3 > angle_max_threshold)
|
||||
{
|
||||
Euler::remove_face(h1, mesh_);
|
||||
}
|
||||
}
|
||||
|
||||
double get_angle(halfedge_descriptor ha, halfedge_descriptor hb)
|
||||
{
|
||||
Vector a(get(vpmap_, source(ha, mesh_)), get(vpmap_, target(ha, mesh_)));
|
||||
Vector b(get(vpmap_, source(hb, mesh_)), get(vpmap_, target(hb, mesh_)));
|
||||
|
||||
return get_angle(a, b);
|
||||
}
|
||||
|
||||
double get_angle(const Vector& e1, const Vector& e2)
|
||||
{
|
||||
//double rad_to_deg = 180. / CGAL_PI;
|
||||
double cos_angle = (e1 * e2)
|
||||
/ std::sqrt(e1.squared_length() * e2.squared_length());
|
||||
|
||||
return std::acos(cos_angle); //* rad_to_deg;
|
||||
}
|
||||
|
||||
|
||||
// functions for constrained vertices
|
||||
// -----------------------------------
|
||||
bool is_constrained(const edge_descriptor& e)
|
||||
{
|
||||
return get(ecmap_, e);
|
||||
}
|
||||
|
||||
bool is_constrained(const vertex_descriptor& v)
|
||||
{
|
||||
return get(vcmap_, v);
|
||||
}
|
||||
|
||||
void check_constraints()
|
||||
{
|
||||
BOOST_FOREACH(edge_descriptor e, edges(mesh_))
|
||||
{
|
||||
if (is_constrained(e))
|
||||
{
|
||||
vertex_descriptor vs = source(e, mesh_);
|
||||
vertex_descriptor vt = target(e, mesh_);
|
||||
put(vcmap_, vs, true);
|
||||
put(vcmap_, vt, true);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
template<typename FaceRange>
|
||||
void check_vertex_range(const FaceRange& face_range)
|
||||
{
|
||||
BOOST_FOREACH(face_descriptor f, face_range)
|
||||
{
|
||||
BOOST_FOREACH(vertex_descriptor v, vertices_around_face(halfedge(f, mesh_), mesh_))
|
||||
vrange.insert(v);
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
private:
|
||||
|
||||
// data members
|
||||
// ------------
|
||||
PolygonMesh& mesh_;
|
||||
VertexPointMap& vpmap_;
|
||||
VertexConstraintMap vcmap_;
|
||||
EdgeConstraintMap ecmap_;
|
||||
Triangle_list input_triangles_;
|
||||
GeomTraits traits_;
|
||||
std::set<vertex_descriptor> vrange;
|
||||
Weight_calculator weight_calculator_;
|
||||
|
||||
|
||||
};
|
||||
|
||||
|
||||
|
||||
|
||||
} // internal
|
||||
} // Polygon_mesh_processing
|
||||
} // CGAL
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
#endif // CURVATURE_FLOW_IMPL_H
|
||||
|
|
@ -0,0 +1,185 @@
|
|||
// Copyright (c) 2015 GeometryFactory (France).
|
||||
// All rights reserved.
|
||||
//
|
||||
// This file is part of CGAL (www.cgal.org).
|
||||
// You can redistribute it and/or modify it under the terms of the GNU
|
||||
// General Public License as published by the Free Software Foundation,
|
||||
// either version 3 of the License, or (at your option) any later version.
|
||||
//
|
||||
// Licensees holding a valid commercial license may use this file in
|
||||
// accordance with the commercial license agreement provided with the software.
|
||||
//
|
||||
// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
|
||||
// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
|
||||
//
|
||||
// $URL$
|
||||
// $Id$
|
||||
//
|
||||
//
|
||||
// Author(s) : Konstantinos Katrioplas (konst.katrioplas@gmail.com)
|
||||
|
||||
#ifndef EVALUATION_H
|
||||
#define EVALUATION_H
|
||||
|
||||
#include <vector>
|
||||
#include <set>
|
||||
#include <fstream>
|
||||
#include <string>
|
||||
#include <boost/graph/graph_traits.hpp>
|
||||
#include <boost/foreach.hpp>
|
||||
#include <boost/property_map/property_map.hpp>
|
||||
|
||||
#include <CGAL/Polygon_mesh_processing/measure.h>
|
||||
|
||||
|
||||
namespace CGAL {
|
||||
|
||||
namespace Polygon_mesh_processing {
|
||||
|
||||
namespace internal {
|
||||
|
||||
|
||||
|
||||
template<typename PolygonMesh>
|
||||
class Quality_evaluator
|
||||
{
|
||||
|
||||
typedef typename boost::graph_traits<PolygonMesh>::halfedge_descriptor halfedge_descriptor;
|
||||
typedef typename boost::graph_traits<PolygonMesh>::face_descriptor face_descriptor;
|
||||
|
||||
typedef typename boost::property_map<PolygonMesh, CGAL::vertex_point_t>::type VertexPointMap;
|
||||
|
||||
|
||||
public:
|
||||
Quality_evaluator(PolygonMesh& pmesh) : mesh_(pmesh)
|
||||
{
|
||||
std::size_t number_of_triangles = faces(mesh_).size();
|
||||
|
||||
// todo: move reserving to each function
|
||||
angles_.reserve(number_of_triangles * 3);
|
||||
areas_.reserve(number_of_triangles);
|
||||
aspect_ratios_.reserve(number_of_triangles);
|
||||
|
||||
vpmap_ = get(CGAL::vertex_point, mesh_);
|
||||
}
|
||||
|
||||
void gather_angles()
|
||||
{
|
||||
double rad_to_deg = 180. / CGAL_PI;
|
||||
|
||||
/*
|
||||
typename boost::property_map<PolygonMesh, CGAL::vertex_point_t>::type
|
||||
vpmap = get(CGAL::vertex_point, mesh_);
|
||||
*/
|
||||
|
||||
for(halfedge_descriptor hi : halfedges(mesh_))
|
||||
{
|
||||
|
||||
typename Kernel::Point_3 a = get(vpmap_, source(hi, mesh_));
|
||||
typename Kernel::Point_3 b = get(vpmap_, target(hi, mesh_));
|
||||
typename Kernel::Point_3 c = get(vpmap_, target(next(hi, mesh_), mesh_));
|
||||
|
||||
typename Kernel::Vector_3 ba(b, a);
|
||||
typename Kernel::Vector_3 bc(b, c);
|
||||
|
||||
double cos_angle = (ba * bc)
|
||||
/ std::sqrt(ba.squared_length() * bc.squared_length());
|
||||
|
||||
angles_.push_back(std::acos(cos_angle) * rad_to_deg);
|
||||
}
|
||||
|
||||
std::cout<<"angles_ size= "<<angles_.size()<<std::endl;
|
||||
}
|
||||
|
||||
void extract_angles(const char* filename)
|
||||
{
|
||||
std::ofstream output(filename);
|
||||
for(unsigned int i=0; i!=angles_.size(); ++i)
|
||||
output << angles_[i] << std::endl;
|
||||
output.close();
|
||||
}
|
||||
|
||||
void measure_areas()
|
||||
{
|
||||
for (face_descriptor f : faces(mesh_))
|
||||
areas_.push_back(face_area(f, mesh_));
|
||||
|
||||
std::cout<<"areas_ size= "<<areas_.size()<<std::endl;
|
||||
}
|
||||
|
||||
void extract_areas(const char* filename)
|
||||
{
|
||||
std::ofstream output(filename);
|
||||
for(unsigned int i=0; i!=areas_.size(); ++i)
|
||||
output << areas_[i] << std::endl;
|
||||
output.close();
|
||||
}
|
||||
|
||||
void calc_aspect_ratios()
|
||||
{
|
||||
|
||||
for(face_descriptor f : faces(mesh_))
|
||||
{
|
||||
halfedge_descriptor h = halfedge(f, mesh_);
|
||||
typename Kernel::Point_3 points[3];
|
||||
points[0] = get(vpmap_, target(h, mesh_));
|
||||
points[1] = get(vpmap_, target(next(h, mesh_), mesh_));
|
||||
points[2] = get(vpmap_, target(next(next(h, mesh_), mesh_), mesh_));
|
||||
|
||||
double min_alt = std::numeric_limits<double>::infinity();
|
||||
double longest_edge = 0;
|
||||
|
||||
for(int i=0; i<3; ++i)
|
||||
{
|
||||
double alt = CGAL::sqrt(CGAL::squared_distance(points[(0+i)%3],
|
||||
typename Kernel::Line_3(points[(1+i)%3], points[(2+i)%3])));
|
||||
double edge = CGAL::sqrt(CGAL::squared_distance(points[(1+i)%3], points[(2+i)%3]));
|
||||
if(alt < min_alt) { min_alt = alt; }
|
||||
if(edge > longest_edge) { longest_edge = edge; }
|
||||
}
|
||||
|
||||
aspect_ratios_.push_back(longest_edge / min_alt);
|
||||
}
|
||||
|
||||
std::cout<<"aspect_ratios_ size= "<<aspect_ratios_.size()<<std::endl;
|
||||
|
||||
}
|
||||
|
||||
void extract_aspect_ratios(const char* filename)
|
||||
{
|
||||
std::ofstream output(filename);
|
||||
for(unsigned int i=0; i!=aspect_ratios_.size(); ++i)
|
||||
output << aspect_ratios_[i] << std::endl;
|
||||
output.close();
|
||||
}
|
||||
|
||||
|
||||
|
||||
|
||||
private:
|
||||
|
||||
PolygonMesh& mesh_;
|
||||
VertexPointMap vpmap_;
|
||||
std::vector<double> angles_;
|
||||
std::vector<double> areas_;
|
||||
std::vector<double> aspect_ratios_;
|
||||
|
||||
};
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
} //namespaces
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
#endif // EVALUATION_H
|
||||
|
|
@ -0,0 +1,695 @@
|
|||
// Copyright (c) 2015 GeometryFactory (France).
|
||||
// All rights reserved.
|
||||
//
|
||||
// This file is part of CGAL (www.cgal.org).
|
||||
// You can redistribute it and/or modify it under the terms of the GNU
|
||||
// General Public License as published by the Free Software Foundation,
|
||||
// either version 3 of the License, or (at your option) any later version.
|
||||
//
|
||||
// Licensees holding a valid commercial license may use this file in
|
||||
// accordance with the commercial license agreement provided with the software.
|
||||
//
|
||||
// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
|
||||
// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
|
||||
//
|
||||
// $URL$
|
||||
// $Id$
|
||||
//
|
||||
//
|
||||
// Author(s) : Konstantinos Katrioplas (konst.katrioplas@gmail.com)
|
||||
|
||||
#ifndef CGAL_POLYGON_MESH_PROCESSING_SMOOTHING_IMPL_H
|
||||
#define CGAL_POLYGON_MESH_PROCESSING_SMOOTHING_IMPL_H
|
||||
|
||||
#include <math.h>
|
||||
#include <utility>
|
||||
#include <iterator>
|
||||
|
||||
#include <CGAL/Kernel/global_functions_3.h>
|
||||
#include <CGAL/Polygon_mesh_processing/compute_normal.h>
|
||||
#include <CGAL/Polygon_mesh_processing/repair.h>
|
||||
|
||||
#include <CGAL/AABB_tree.h>
|
||||
#include <CGAL/AABB_traits.h>
|
||||
#include <CGAL/AABB_triangle_primitive.h>
|
||||
|
||||
#include <CGAL/property_map.h>
|
||||
#include <CGAL/iterator.h>
|
||||
#include <CGAL/boost/graph/Euler_operations.h>
|
||||
#include <boost/graph/graph_traits.hpp>
|
||||
#include <boost/foreach.hpp>
|
||||
|
||||
|
||||
namespace CGAL {
|
||||
|
||||
namespace Polygon_mesh_processing {
|
||||
|
||||
namespace internal {
|
||||
|
||||
|
||||
template<typename PolygonMesh, typename VertexPointMap, typename VertexConstraintMap, typename EdgeConstraintMap, typename GeomTraits>
|
||||
class Compatible_remesher
|
||||
{
|
||||
|
||||
typedef typename boost::graph_traits<PolygonMesh>::halfedge_descriptor halfedge_descriptor;
|
||||
typedef typename boost::graph_traits<PolygonMesh>::vertex_descriptor vertex_descriptor;
|
||||
typedef typename boost::graph_traits<PolygonMesh>::face_descriptor face_descriptor;
|
||||
typedef typename boost::graph_traits<PolygonMesh>::edge_descriptor edge_descriptor;
|
||||
typedef typename boost::graph_traits<PolygonMesh>::vertex_iterator vertex_iterator;
|
||||
|
||||
typedef typename GeomTraits::Point_3 Point;
|
||||
typedef typename GeomTraits::Vector_3 Vector;
|
||||
typedef typename GeomTraits::Segment_3 Segment;
|
||||
typedef typename GeomTraits::Triangle_3 Triangle;
|
||||
typedef std::vector<Triangle> Triangle_list;
|
||||
|
||||
typedef CGAL::AABB_triangle_primitive<GeomTraits, typename Triangle_list::iterator> AABB_Primitive;
|
||||
typedef CGAL::AABB_traits<GeomTraits, AABB_Primitive> AABB_Traits;
|
||||
typedef CGAL::AABB_tree<AABB_Traits> Tree;
|
||||
|
||||
typedef std::pair<halfedge_descriptor, halfedge_descriptor> he_pair;
|
||||
typedef std::map<halfedge_descriptor, he_pair> Edges_around_map;
|
||||
|
||||
|
||||
public:
|
||||
Compatible_remesher(PolygonMesh& pmesh, VertexPointMap& vpmap, VertexConstraintMap& vcmap, EdgeConstraintMap& ecmap) :
|
||||
mesh_(pmesh), vpmap_(vpmap), vcmap_(vcmap), ecmap_(ecmap)
|
||||
{}
|
||||
|
||||
~Compatible_remesher()
|
||||
{
|
||||
delete tree_ptr_;
|
||||
}
|
||||
|
||||
template<typename FaceRange>
|
||||
void init_smoothing(const FaceRange& face_range)
|
||||
{
|
||||
check_vertex_range(face_range);
|
||||
|
||||
check_constraints();
|
||||
|
||||
BOOST_FOREACH(face_descriptor f, face_range)
|
||||
input_triangles_.push_back(triangle(f));
|
||||
|
||||
tree_ptr_ = new Tree(input_triangles_.begin(), input_triangles_.end());
|
||||
tree_ptr_->accelerate_distance_queries();
|
||||
}
|
||||
|
||||
std::size_t remove_degenerate_faces()
|
||||
{
|
||||
std::size_t nb_removed_faces = 0;
|
||||
nb_removed_faces = CGAL::Polygon_mesh_processing::remove_degenerate_faces(mesh_);
|
||||
|
||||
return nb_removed_faces;
|
||||
}
|
||||
|
||||
void angle_relaxation(bool use_weights)
|
||||
{
|
||||
std::map<vertex_descriptor, Point> barycenters;
|
||||
std::map<vertex_descriptor, Vector> n_map;
|
||||
|
||||
|
||||
BOOST_FOREACH(vertex_descriptor v, vrange_)
|
||||
{
|
||||
if(!is_border(v, mesh_) && !is_constrained(v))
|
||||
{
|
||||
// compute normal to v
|
||||
Vector vn = compute_vertex_normal(v, mesh_,
|
||||
Polygon_mesh_processing::parameters::vertex_point_map(vpmap_)
|
||||
.geom_traits(traits_));
|
||||
n_map[v] = vn;
|
||||
|
||||
Edges_around_map he_map;
|
||||
typename Edges_around_map::iterator it;
|
||||
BOOST_FOREACH(halfedge_descriptor hi, halfedges_around_source(v, mesh_))
|
||||
he_map[hi] = he_pair( next(hi, mesh_), prev(opposite(hi, mesh_), mesh_) );
|
||||
|
||||
// measure initial angles
|
||||
measure_angles(he_map);
|
||||
|
||||
// calculate movement
|
||||
Vector move = calc_move(he_map, use_weights);
|
||||
|
||||
barycenters[v] = (get(vpmap_, v) + move) ;
|
||||
|
||||
} // not on border
|
||||
} // for each v
|
||||
|
||||
// compute locations on tangent plane
|
||||
typedef typename std::map<vertex_descriptor, Point>::value_type VP;
|
||||
std::map<vertex_descriptor, Point> new_locations;
|
||||
BOOST_FOREACH(const VP& vp, barycenters)
|
||||
{
|
||||
Point p = get(vpmap_, vp.first);
|
||||
Point q = vp.second;
|
||||
Vector n = n_map[vp.first];
|
||||
|
||||
new_locations[vp.first] = q + ( n * Vector(q, p) ) * n ;
|
||||
}
|
||||
|
||||
// update location
|
||||
std::size_t moved_points = 0;
|
||||
BOOST_FOREACH(const VP& vp, new_locations)
|
||||
{
|
||||
// iff movement impoves all angles
|
||||
if(does_it_impove(vp.first, vp.second))
|
||||
{
|
||||
moved_points++;
|
||||
put(vpmap_, vp.first, vp.second);
|
||||
}
|
||||
}
|
||||
|
||||
#ifdef CGAL_PMP_SMOOTHING_DEBUG
|
||||
std::cout<<"moved: "<< moved_points <<" points based on angle."<<std::endl;
|
||||
std::cout<<"not imporved min angle: "<< vrange_.size() - moved_points <<" times."<<std::endl;
|
||||
#endif
|
||||
}
|
||||
|
||||
void area_relaxation(const double& precision)
|
||||
{
|
||||
std::size_t moved_points = 0;
|
||||
BOOST_FOREACH(vertex_descriptor v, vrange_)
|
||||
{
|
||||
if(!is_border(v, mesh_) && !is_constrained(v))
|
||||
{
|
||||
if (gradient_descent(v, precision))
|
||||
moved_points++;
|
||||
}
|
||||
}
|
||||
|
||||
#ifdef CGAL_PMP_SMOOTHING_DEBUG
|
||||
std::cout<<"moved : "<<moved_points<<" points based on area."<<std::endl;
|
||||
std::cout<<"non convex energy found: "<<vrange_.size() - moved_points<<" times."<<std::endl;
|
||||
#endif
|
||||
|
||||
}
|
||||
|
||||
void project_to_surface()
|
||||
{
|
||||
BOOST_FOREACH( vertex_descriptor v, vertices(mesh_))
|
||||
{
|
||||
if(!is_border(v, mesh_) && !is_constrained(v))
|
||||
{
|
||||
Point p_query = get(vpmap_, v);
|
||||
Point projected = tree_ptr_->closest_point(p_query);
|
||||
put(vpmap_, v, projected);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
private:
|
||||
|
||||
// helper functions
|
||||
// ----------------
|
||||
Triangle triangle(face_descriptor f) const
|
||||
{
|
||||
halfedge_descriptor h = halfedge(f, mesh_);
|
||||
vertex_descriptor v1 = target(h, mesh_);
|
||||
vertex_descriptor v2 = target(next(h, mesh_), mesh_);
|
||||
vertex_descriptor v3 = target(next(next(h, mesh_), mesh_), mesh_);
|
||||
return Triangle(get(vpmap_, v1), get(vpmap_, v2), get(vpmap_, v3));
|
||||
}
|
||||
|
||||
double sqlength(const vertex_descriptor& v1, const vertex_descriptor& v2) const
|
||||
{
|
||||
return to_double(CGAL::squared_distance(get(vpmap_, v1), get(vpmap_, v2)));
|
||||
}
|
||||
|
||||
double sqlength(const halfedge_descriptor& h) const
|
||||
{
|
||||
vertex_descriptor v1 = target(h, mesh_);
|
||||
vertex_descriptor v2 = source(h, mesh_);
|
||||
return sqlength(v1, v2);
|
||||
}
|
||||
|
||||
double sqlength(const edge_descriptor& e) const
|
||||
{
|
||||
return sqlength(halfedge(e, mesh_));
|
||||
}
|
||||
|
||||
|
||||
// angle bisecting functions
|
||||
// -------------------------
|
||||
Vector calc_move(const Edges_around_map& he_map, bool use_weights)
|
||||
{
|
||||
Vector move = CGAL::NULL_VECTOR;
|
||||
double weights_sum = 0;
|
||||
typename Edges_around_map::const_iterator it;
|
||||
for(it = he_map.begin(); it != he_map.end(); ++it)
|
||||
{
|
||||
halfedge_descriptor main_he = it->first;
|
||||
he_pair incident_pair = it->second;
|
||||
|
||||
// avoid zero angles
|
||||
Point pt = get(vpmap_, source(incident_pair.first, mesh_));
|
||||
Point p1 = get(vpmap_, target(incident_pair.first, mesh_));
|
||||
Point p2 = get(vpmap_, source(incident_pair.second, mesh_));
|
||||
CGAL_assertion(target(incident_pair.second, mesh_) == source(incident_pair.first, mesh_));
|
||||
Vector e1(pt, p1);
|
||||
Vector e2(pt, p2);
|
||||
double angle = get_angle(e1, e2);
|
||||
if(angle < 1e-5)
|
||||
continue;
|
||||
|
||||
// rotate
|
||||
Vector rotated_edge = rotate_edge(main_he, incident_pair);
|
||||
|
||||
// small angles carry more weight
|
||||
double weight = 1.0 / (angle*angle);
|
||||
weights_sum += weight;
|
||||
|
||||
if(use_weights)
|
||||
move += weight * rotated_edge;
|
||||
else
|
||||
move += rotated_edge;
|
||||
}
|
||||
|
||||
// if at least 1 angle was summed
|
||||
if(use_weights && weights_sum != 0)
|
||||
move /= weights_sum;
|
||||
else
|
||||
move /= CGAL::to_double(he_map.size());
|
||||
|
||||
return move;
|
||||
}
|
||||
|
||||
Vector rotate_edge(const halfedge_descriptor& main_he, const he_pair& incd_edges)
|
||||
{
|
||||
// get common vertex around which the edge is rotated
|
||||
Point pt = get(vpmap_, target(main_he, mesh_));
|
||||
|
||||
// ps is the vertex that is being moved
|
||||
Point ps = get(vpmap_, source(main_he, mesh_));
|
||||
|
||||
// get "equidistant" points - in fact they are at equal angles
|
||||
Point equidistant_p1 = get(vpmap_, target(incd_edges.first, mesh_));
|
||||
Point equidistant_p2 = get(vpmap_, source(incd_edges.second, mesh_));
|
||||
CGAL_assertion(target(incd_edges.second, mesh_) == source(incd_edges.first, mesh_));
|
||||
|
||||
Vector edge1(pt, equidistant_p1);
|
||||
Vector edge2(pt, equidistant_p2);
|
||||
Vector vec_main_he(pt, ps);
|
||||
|
||||
// check degenerate cases
|
||||
double tolerance = 1e-3;
|
||||
|
||||
if ( edge1.squared_length() < tolerance ||
|
||||
edge2.squared_length() < tolerance ||
|
||||
sqlength(main_he) < tolerance ||
|
||||
(edge1 - vec_main_he).squared_length() < tolerance ||
|
||||
(edge2 - vec_main_he).squared_length() < tolerance )
|
||||
{
|
||||
return CGAL::NULL_VECTOR;
|
||||
}
|
||||
|
||||
CGAL_assertion(vec_main_he.squared_length() > tolerance);
|
||||
|
||||
// find bisector
|
||||
Vector bisector = CGAL::NULL_VECTOR;
|
||||
internal::normalize(edge1, traits_);
|
||||
internal::normalize(edge2, traits_);
|
||||
bisector = edge1 + edge2;
|
||||
|
||||
// special case in angle bisecting: If no edge is actually degenerate
|
||||
// but edge1 and edge2 are almost parallel, then by adding them the bisector is
|
||||
// almost zero.
|
||||
if( bisector.squared_length() < tolerance )
|
||||
{
|
||||
// angle is (almost) 180 degrees, take the perpendicular
|
||||
// which is normal to edge and tangent to the surface
|
||||
Vector normal_vec = find_perpendicular(edge1, pt, ps);
|
||||
|
||||
CGAL_assertion(normal_vec != CGAL::NULL_VECTOR);
|
||||
CGAL_assertion(CGAL::scalar_product(edge1, normal_vec) < tolerance);
|
||||
|
||||
Segment b_segment(pt, pt + normal_vec);
|
||||
Point b_segment_end = b_segment.target();
|
||||
|
||||
if(CGAL::angle(b_segment_end, pt, ps) == CGAL::OBTUSE)
|
||||
{
|
||||
b_segment = b_segment.opposite();
|
||||
}
|
||||
|
||||
bisector = Vector(b_segment);
|
||||
}
|
||||
|
||||
correct_bisector(bisector, main_he);
|
||||
|
||||
double target_length = CGAL::sqrt(sqlength(main_he));
|
||||
double bisector_length = CGAL::sqrt(bisector.squared_length());
|
||||
|
||||
CGAL_assertion( ( target_length - tolerance < bisector_length ) &&
|
||||
( bisector_length < target_length + tolerance ) );
|
||||
|
||||
return bisector;
|
||||
}
|
||||
|
||||
void correct_bisector(Vector& bisector_vec, const halfedge_descriptor& main_he)
|
||||
{
|
||||
// get common vertex around which the edge is rotated
|
||||
Point pt = get(vpmap_, target(main_he, mesh_));
|
||||
|
||||
// create a segment to be able to translate
|
||||
Segment bisector(pt, pt + bisector_vec);
|
||||
|
||||
// scale
|
||||
double scale_factor = CGAL::sqrt( sqlength(main_he) / bisector.squared_length() );
|
||||
typename GeomTraits::Aff_transformation_3 t_scale(CGAL::SCALING, scale_factor);
|
||||
bisector = bisector.transform(t_scale);
|
||||
|
||||
// translate
|
||||
Vector vec(bisector.source(), pt);
|
||||
typename GeomTraits::Aff_transformation_3 t_translate(CGAL::TRANSLATION, vec);
|
||||
bisector = bisector.transform(t_translate);
|
||||
|
||||
// take the opposite so that their sum is the overall displacement
|
||||
bisector_vec = -Vector(bisector);
|
||||
}
|
||||
|
||||
Vector find_perpendicular(const Vector& input_vec, const Point& s, const Point& pv)
|
||||
{
|
||||
Vector s_pv(s, pv);
|
||||
Vector aux_normal = CGAL::cross_product(input_vec, s_pv);
|
||||
return CGAL::cross_product(aux_normal, input_vec);
|
||||
}
|
||||
|
||||
|
||||
// angle measurement & evaluation
|
||||
// ------------------------------
|
||||
void measure_angles(const Edges_around_map& he_map)
|
||||
{
|
||||
min_angle_ = 2 * CGAL_PI;
|
||||
typename Edges_around_map::const_iterator it;
|
||||
for(it = he_map.begin(); it != he_map.end(); ++it)
|
||||
{
|
||||
halfedge_descriptor main_he = it->first;
|
||||
he_pair incident_pair = it->second;
|
||||
|
||||
calc_angles(main_he, incident_pair);
|
||||
}
|
||||
}
|
||||
|
||||
void calc_angles(const halfedge_descriptor& main_he, const he_pair& incd_edges)
|
||||
{
|
||||
// get common vertex around which the edge is rotated
|
||||
Point pt = get(vpmap_, target(main_he, mesh_));
|
||||
// ps is the vertex that is being moved
|
||||
Point ps = get(vpmap_, source(main_he, mesh_));
|
||||
// get "equidistant" points - in fact they are at equal angles
|
||||
Point equidistant_p1 = get(vpmap_, target(incd_edges.first, mesh_));
|
||||
Point equidistant_p2 = get(vpmap_, source(incd_edges.second, mesh_));
|
||||
CGAL_assertion(target(incd_edges.second, mesh_) == source(incd_edges.first, mesh_));
|
||||
|
||||
Vector edge1(pt, equidistant_p1);
|
||||
Vector edge2(pt, equidistant_p2);
|
||||
Vector vec_main_he(pt, ps);
|
||||
|
||||
// angles in [0, 2pi]
|
||||
double a1 = get_angle(edge1, vec_main_he);
|
||||
double a2 = get_angle(edge2, vec_main_he);
|
||||
|
||||
min_angle_ = std::min(min_angle_, std::min(a1, a2));
|
||||
}
|
||||
|
||||
bool does_it_impove(const vertex_descriptor& v, const Point& new_location)
|
||||
{
|
||||
Edges_around_map he_map;
|
||||
typename Edges_around_map::iterator it;
|
||||
BOOST_FOREACH(halfedge_descriptor hi, halfedges_around_source(v, mesh_))
|
||||
he_map[hi] = he_pair( next(hi, mesh_), prev(opposite(hi, mesh_), mesh_) );
|
||||
|
||||
return evaluate_angles(he_map, new_location);
|
||||
}
|
||||
|
||||
bool evaluate_angles(const Edges_around_map& he_map, const Point& new_location)
|
||||
{
|
||||
typename Edges_around_map::const_iterator it;
|
||||
for(it = he_map.begin(); it != he_map.end(); ++it)
|
||||
{
|
||||
halfedge_descriptor main_he = it->first;
|
||||
he_pair incd_edges = it->second;
|
||||
|
||||
// get common vertex around which the edge is rotated
|
||||
Point pt = get(vpmap_, target(main_he, mesh_));
|
||||
// ps is the vertex that is being moved
|
||||
Point new_point = new_location;
|
||||
// get "equidistant" points
|
||||
Point equidistant_p1 = get(vpmap_, target(incd_edges.first, mesh_));
|
||||
Point equidistant_p2 = get(vpmap_, source(incd_edges.second, mesh_));
|
||||
CGAL_assertion(target(incd_edges.second, mesh_) == source(incd_edges.first, mesh_));
|
||||
|
||||
Vector edge1(pt, equidistant_p1);
|
||||
Vector edge2(pt, equidistant_p2);
|
||||
Vector vec_main_he(pt, new_point);
|
||||
|
||||
double a1 = get_angle(edge1, vec_main_he);
|
||||
double a2 = get_angle(edge2, vec_main_he);
|
||||
|
||||
if(a1 < min_angle_ || a2 < min_angle_)
|
||||
return false;
|
||||
}
|
||||
|
||||
return true;
|
||||
}
|
||||
|
||||
double get_angle(const Vector& e1, const Vector& e2)
|
||||
{
|
||||
//double rad_to_deg = 180. / CGAL_PI;
|
||||
double cos_angle = (e1 * e2)
|
||||
/ std::sqrt(e1.squared_length() * e2.squared_length());
|
||||
|
||||
return std::acos(cos_angle); //* rad_to_deg;
|
||||
}
|
||||
|
||||
|
||||
// gradient descent
|
||||
// ----------------
|
||||
bool gradient_descent(const vertex_descriptor& v, const double& precision)
|
||||
{
|
||||
bool move_flag;
|
||||
double x, y, z, x_new, y_new, z_new, drdx, drdy, drdz;
|
||||
x = get(vpmap_, v).x();
|
||||
y = get(vpmap_, v).y();
|
||||
z = get(vpmap_, v).z();
|
||||
|
||||
double S_av = compute_average_area_around(v);
|
||||
double energy = measure_energy(v, S_av);
|
||||
|
||||
// if the adjacent areas are absolutely equal
|
||||
if(energy == 0)
|
||||
return false;
|
||||
|
||||
double energy_new = 0;
|
||||
double relative_energy = precision + 1;
|
||||
unsigned int t = 1;
|
||||
double eta0 = 0.01;
|
||||
//double power_t = 0.25;
|
||||
double t0 = 0.001;
|
||||
double eta = eta0 / (1 + t0*t);
|
||||
|
||||
while(relative_energy > precision)
|
||||
{
|
||||
drdx=0, drdy=0, drdz=0;
|
||||
compute_derivatives(drdx, drdy, drdz, v, S_av);
|
||||
|
||||
x_new = x - eta * drdx;
|
||||
y_new = y - eta * drdy;
|
||||
z_new = z - eta * drdz;
|
||||
|
||||
Point moved(x_new, y_new, z_new);
|
||||
energy_new = measure_energy(v, S_av, moved);
|
||||
|
||||
if(energy_new < energy)
|
||||
{
|
||||
put(vpmap_, v, moved);
|
||||
move_flag = true;
|
||||
}
|
||||
else
|
||||
return false;
|
||||
|
||||
relative_energy = CGAL::to_double( (energy - energy_new) / energy );
|
||||
|
||||
// update
|
||||
x = x_new;
|
||||
y = y_new;
|
||||
z = z_new;
|
||||
energy = energy_new;
|
||||
t++;
|
||||
|
||||
//eta = eta0 / pow(t, power_t);
|
||||
eta = eta0 / (1 + t0 * t);
|
||||
}
|
||||
|
||||
return move_flag;
|
||||
}
|
||||
|
||||
void compute_derivatives(double& drdx, double& drdy, double& drdz, const vertex_descriptor& v, const double& S_av)
|
||||
{
|
||||
BOOST_FOREACH(halfedge_descriptor h, halfedges_around_source(v, mesh_))
|
||||
{
|
||||
vertex_descriptor pi = source(next(h, mesh_), mesh_);
|
||||
vertex_descriptor pi1 = target(next(h, mesh_), mesh_);
|
||||
double S = element_area(v, pi, pi1);
|
||||
|
||||
Vector vec(get(vpmap_, pi), get(vpmap_, pi1));
|
||||
|
||||
// minimize r:
|
||||
// r = Σ(S-S_av)^2
|
||||
// dr/dx = 2 Σ(S - S_av) dS/dx
|
||||
// area of triangle with respect to (x_a, y_a, z_a) =
|
||||
// (1/2) [(v_z - v_y)x_a + (v_x - v_z)y_a + (v_y - v_x)z_a + constants]
|
||||
// vector v is (x_c - x_b, y_c - y_b, z_c - z_b)
|
||||
drdx += (S - S_av) * 0.5 * (vec.z() - vec.y());
|
||||
drdy += (S - S_av) * 0.5 * (vec.x() - vec.z());
|
||||
drdz += (S - S_av) * 0.5 * (vec.y() - vec.x());
|
||||
}
|
||||
|
||||
drdx *= 2;
|
||||
drdy *= 2;
|
||||
drdz *= 2;
|
||||
}
|
||||
|
||||
double element_area(const vertex_descriptor& p1,
|
||||
const vertex_descriptor& p2,
|
||||
const vertex_descriptor& p3) const
|
||||
{
|
||||
return to_double(CGAL::approximate_sqrt(
|
||||
traits_.compute_squared_area_3_object()(
|
||||
get(vpmap_, p1),
|
||||
get(vpmap_, p2),
|
||||
get(vpmap_, p3))));
|
||||
}
|
||||
|
||||
double element_area(const Point& P,
|
||||
const vertex_descriptor& p2,
|
||||
const vertex_descriptor& p3) const
|
||||
{
|
||||
return to_double(CGAL::approximate_sqrt(
|
||||
traits_.compute_squared_area_3_object()(
|
||||
P,
|
||||
get(vpmap_, p2),
|
||||
get(vpmap_, p3))));
|
||||
}
|
||||
|
||||
double compute_average_area_around(const vertex_descriptor& v)
|
||||
{
|
||||
double sum_areas = 0;
|
||||
unsigned int number_of_edges = 0;
|
||||
|
||||
BOOST_FOREACH(halfedge_descriptor h, halfedges_around_source(v, mesh_))
|
||||
{
|
||||
// opposite vertices
|
||||
vertex_descriptor pi = source(next(h, mesh_), mesh_);
|
||||
vertex_descriptor pi1 = target(next(h, mesh_), mesh_);
|
||||
|
||||
double S = element_area(v, pi, pi1);
|
||||
sum_areas += S;
|
||||
number_of_edges++;
|
||||
}
|
||||
|
||||
return sum_areas / number_of_edges;
|
||||
}
|
||||
|
||||
double measure_energy(const vertex_descriptor& v, const double& S_av)
|
||||
{
|
||||
double energy = 0;
|
||||
unsigned int number_of_edges = 0;
|
||||
|
||||
BOOST_FOREACH(halfedge_descriptor h, halfedges_around_source(v, mesh_))
|
||||
{
|
||||
vertex_descriptor pi = source(next(h, mesh_), mesh_);
|
||||
vertex_descriptor pi1 = target(next(h, mesh_), mesh_);
|
||||
double S = element_area(v, pi, pi1);
|
||||
|
||||
energy += (S - S_av)*(S - S_av);
|
||||
number_of_edges++;
|
||||
}
|
||||
|
||||
return to_double( energy / number_of_edges );
|
||||
}
|
||||
|
||||
double measure_energy(const vertex_descriptor& v, const double& S_av, const Point& new_P)
|
||||
{
|
||||
double energy = 0;
|
||||
unsigned int number_of_edges = 0;
|
||||
BOOST_FOREACH(halfedge_descriptor h, halfedges_around_source(v, mesh_))
|
||||
{
|
||||
vertex_descriptor pi = source(next(h, mesh_), mesh_);
|
||||
vertex_descriptor pi1 = target(next(h, mesh_), mesh_);
|
||||
double S = element_area(new_P, pi, pi1);
|
||||
|
||||
energy += (S - S_av)*(S - S_av);
|
||||
number_of_edges++;
|
||||
}
|
||||
|
||||
return to_double( energy / (2 * number_of_edges) );
|
||||
}
|
||||
|
||||
|
||||
// functions for constrained vertices
|
||||
// -----------------------------------
|
||||
bool is_constrained(const edge_descriptor& e)
|
||||
{
|
||||
return get(ecmap_, e);
|
||||
}
|
||||
|
||||
bool is_constrained(const vertex_descriptor& v)
|
||||
{
|
||||
return get(vcmap_, v);
|
||||
}
|
||||
|
||||
void check_constraints()
|
||||
{
|
||||
BOOST_FOREACH(edge_descriptor e, edges(mesh_))
|
||||
{
|
||||
if (is_constrained(e))
|
||||
{
|
||||
vertex_descriptor vs = source(e, mesh_);
|
||||
vertex_descriptor vt = target(e, mesh_);
|
||||
put(vcmap_, vs, true);
|
||||
put(vcmap_, vt, true);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
template<typename FaceRange>
|
||||
void check_vertex_range(const FaceRange& face_range)
|
||||
{
|
||||
BOOST_FOREACH(face_descriptor f, face_range)
|
||||
{
|
||||
BOOST_FOREACH(vertex_descriptor v, vertices_around_face(halfedge(f, mesh_), mesh_))
|
||||
vrange_.insert(v);
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
private:
|
||||
|
||||
// data members
|
||||
// ------------
|
||||
PolygonMesh& mesh_;
|
||||
VertexPointMap& vpmap_;
|
||||
VertexConstraintMap vcmap_;
|
||||
EdgeConstraintMap ecmap_;
|
||||
Triangle_list input_triangles_;
|
||||
Tree* tree_ptr_;
|
||||
GeomTraits traits_;
|
||||
std::set<vertex_descriptor> vrange_;
|
||||
double min_angle_;
|
||||
|
||||
};
|
||||
|
||||
|
||||
|
||||
|
||||
} // namespace internal
|
||||
} // namespace Polygon_mesh_processing
|
||||
} // namespace CGAL
|
||||
|
||||
|
||||
|
||||
|
||||
#endif // CGAL_POLYGON_MESH_PROCESSING_SMOOTHING_IMPL_H
|
||||
|
||||
|
||||
|
|
@ -49,6 +49,9 @@ CGAL_add_named_parameter(nb_points_per_distance_unit_t, nb_points_per_distance_u
|
|||
CGAL_add_named_parameter(face_normal_t, face_normal, face_normal_map)
|
||||
CGAL_add_named_parameter(random_seed_t, random_seed, random_seed)
|
||||
CGAL_add_named_parameter(do_project_t, do_project, do_project)
|
||||
CGAL_add_named_parameter(gradient_descent_precision_t, gradient_descent_precision, gradient_descent_precision)
|
||||
CGAL_add_named_parameter(use_weights_t, use_weights, use_weights)
|
||||
CGAL_add_named_parameter(distance_precision_t, distance_precision, distance_precision)
|
||||
|
||||
//internal
|
||||
CGAL_add_named_parameter(weight_calculator_t, weight_calculator, weight_calculator)
|
||||
|
|
|
|||
File diff suppressed because it is too large
Load Diff
|
|
@ -105,7 +105,8 @@ endif()
|
|||
create_single_source_cgal_program("test_pmp_clip.cpp")
|
||||
create_single_source_cgal_program("triangulate_hole_polyline_test.cpp")
|
||||
create_single_source_cgal_program("surface_intersection_sm_poly.cpp" )
|
||||
|
||||
create_single_source_cgal_program("test_smoothing.cpp")
|
||||
create_single_source_cgal_program("test_curvature_flow.cpp")
|
||||
if(OpenMesh_FOUND)
|
||||
create_single_source_cgal_program("remeshing_test_P_SM_OM.cpp" )
|
||||
target_link_libraries( remeshing_test_P_SM_OM ${OPENMESH_LIBRARIES} )
|
||||
|
|
|
|||
|
|
@ -0,0 +1,122 @@
|
|||
#include <iostream>
|
||||
#include <fstream>
|
||||
#include <vector>
|
||||
#include <string>
|
||||
|
||||
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
|
||||
#include <CGAL/Surface_mesh.h>
|
||||
|
||||
#include <CGAL/Polygon_mesh_processing/smoothing.h>
|
||||
#include <CGAL/Polygon_mesh_processing/distance.h>
|
||||
#include <CGAL/Polygon_mesh_processing/remesh.h>
|
||||
|
||||
|
||||
#if defined(CGAL_LINKED_WITH_TBB)
|
||||
#define TAG CGAL::Parallel_tag
|
||||
#else
|
||||
#define TAG CGAL::Sequential_tag
|
||||
#endif
|
||||
|
||||
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
|
||||
typedef CGAL::Surface_mesh<K::Point_3> Mesh;
|
||||
|
||||
typedef typename Mesh::vertex_index vertex_index;
|
||||
typedef typename Mesh::face_index face_index;
|
||||
typedef typename Mesh::halfedge_index halfedge_index;
|
||||
|
||||
std::set<face_index> select(Mesh mesh, double x0)
|
||||
{
|
||||
std::set<face_index> f_selection;
|
||||
|
||||
typedef typename boost::property_map<Mesh, CGAL::vertex_point_t>::type VertexPointMap;
|
||||
VertexPointMap vpmap = get(CGAL::vertex_point, mesh);
|
||||
|
||||
|
||||
for(face_index f : faces(mesh))
|
||||
{
|
||||
halfedge_index he = halfedge(f, mesh);
|
||||
for(vertex_index v : vertices_around_face(he, mesh))
|
||||
{
|
||||
if(get(vpmap, v).x() < x0)
|
||||
f_selection.insert(f);
|
||||
continue;
|
||||
}
|
||||
}
|
||||
|
||||
return f_selection;
|
||||
}
|
||||
|
||||
int main(int argc, char* argv[]){
|
||||
|
||||
Mesh mesh;
|
||||
Mesh originalMesh;
|
||||
std::ifstream input;
|
||||
#ifdef CGAL_PMP_REMESHING_VERBOSE
|
||||
std::ofstream output;
|
||||
#endif
|
||||
|
||||
// flat test
|
||||
input.open("data/polygon.off");
|
||||
if (!input || !(input >> mesh) || mesh.is_empty()) {
|
||||
std::cerr << "Not a valid .off file." << std::endl;
|
||||
return 1;
|
||||
}
|
||||
input.close();
|
||||
|
||||
originalMesh = mesh;
|
||||
CGAL::Polygon_mesh_processing::curvature_flow_smoothing(mesh);
|
||||
|
||||
double dist = CGAL::Polygon_mesh_processing::approximate_Hausdorff_distance
|
||||
<TAG>(originalMesh, mesh, CGAL::Polygon_mesh_processing::parameters::number_of_points_per_area_unit(1000));
|
||||
|
||||
if(dist > 1e-3)
|
||||
{
|
||||
return EXIT_FAILURE;
|
||||
}
|
||||
|
||||
// half-sphere test
|
||||
input.open("data/sphere.off");
|
||||
if (!input || !(input >> mesh) || mesh.is_empty()) {
|
||||
std::cerr << "Not a valid .off file." << std::endl;
|
||||
return 1;
|
||||
}
|
||||
input.close();
|
||||
|
||||
|
||||
// select half sphere
|
||||
#ifdef CGAL_PMP_REMESHING_VERBOSE
|
||||
std::cout<<"all faces: "<<faces(mesh).size()<<std::endl;
|
||||
#endif
|
||||
|
||||
std::set<face_index> selected_faces = select(mesh, 0);
|
||||
|
||||
#ifdef CGAL_PMP_REMESHING_VERBOSE
|
||||
std::cout<<"selected faces: "<<selected_faces.size()<<std::endl;
|
||||
#endif
|
||||
CGAL::Polygon_mesh_processing::isotropic_remeshing(selected_faces, 0.1, mesh);
|
||||
|
||||
#ifdef CGAL_PMP_REMESHING_VERBOSE
|
||||
output.open("data/half-sphere.off");
|
||||
output << mesh;
|
||||
output.close();
|
||||
#endif
|
||||
|
||||
originalMesh = mesh;
|
||||
CGAL::Polygon_mesh_processing::curvature_flow_smoothing(mesh);
|
||||
|
||||
dist = CGAL::Polygon_mesh_processing::approximate_Hausdorff_distance
|
||||
<TAG>(originalMesh, mesh, CGAL::Polygon_mesh_processing::parameters::number_of_points_per_area_unit(1000));
|
||||
|
||||
if(dist > 1e-1)
|
||||
{
|
||||
return EXIT_FAILURE;
|
||||
}
|
||||
|
||||
#ifdef CGAL_PMP_REMESHING_VERBOSE
|
||||
output.open("data/half-sphere_smoothed.off");
|
||||
output << mesh;
|
||||
output.close();
|
||||
#endif
|
||||
|
||||
return 0;
|
||||
}
|
||||
|
|
@ -0,0 +1,143 @@
|
|||
#include <iostream>
|
||||
#include <fstream>
|
||||
#include <string>
|
||||
|
||||
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
|
||||
#include <CGAL/Surface_mesh.h>
|
||||
#include <CGAL/Polygon_mesh_processing/smoothing.h>
|
||||
#include <CGAL/Polygon_mesh_processing/measure.h>
|
||||
|
||||
#include <boost/accumulators/accumulators.hpp>
|
||||
#include <boost/accumulators/statistics/stats.hpp>
|
||||
#include <boost/accumulators/statistics/mean.hpp>
|
||||
#include <boost/accumulators/statistics/min.hpp>
|
||||
#include <boost/accumulators/statistics/max.hpp>
|
||||
#include <boost/property_map/property_map.hpp>
|
||||
|
||||
|
||||
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
|
||||
typedef CGAL::Surface_mesh<K::Point_3> Mesh;
|
||||
|
||||
void calc_angles(Mesh& mesh, double& min_a, double& max_a, double& avg_a)
|
||||
{
|
||||
using namespace boost::accumulators;
|
||||
accumulator_set< double,
|
||||
features< tag::min, tag::max, tag::mean > > acc;
|
||||
|
||||
typename boost::property_map<Mesh, CGAL::vertex_point_t>::type
|
||||
vpmap = get(CGAL::vertex_point, mesh);
|
||||
|
||||
typedef typename boost::graph_traits<Mesh>::halfedge_descriptor halfedge_descriptor;
|
||||
double rad_to_deg = 180. / CGAL_PI;
|
||||
|
||||
BOOST_FOREACH(halfedge_descriptor h, halfedges(mesh))
|
||||
{
|
||||
if(face(h, mesh) == boost::graph_traits<Mesh>::null_face())
|
||||
continue;
|
||||
|
||||
typename K::Point_3 a = get(vpmap, source(h, mesh));
|
||||
typename K::Point_3 b = get(vpmap, target(h, mesh));
|
||||
typename K::Point_3 c = get(vpmap, target(next(h, mesh), mesh));
|
||||
typename K::Vector_3 ba(b, a);
|
||||
typename K::Vector_3 bc(b, c);
|
||||
|
||||
double cos_angle = (ba * bc)
|
||||
/ std::sqrt(ba.squared_length() * bc.squared_length());
|
||||
|
||||
acc(std::acos(cos_angle) * rad_to_deg);
|
||||
}
|
||||
|
||||
min_a = extract_result< tag::min >(acc);
|
||||
max_a = extract_result< tag::max >(acc);
|
||||
avg_a = extract_result< tag::mean >(acc);
|
||||
}
|
||||
|
||||
void calc_areas(Mesh& mesh, double& min_a, double& max_a, double& avg_a)
|
||||
{
|
||||
using namespace boost::accumulators;
|
||||
accumulator_set< double,
|
||||
features< tag::min, tag::max, tag::mean > > acc;
|
||||
|
||||
typename boost::property_map<Mesh, CGAL::vertex_point_t>::type
|
||||
vpmap = get(CGAL::vertex_point, mesh);
|
||||
|
||||
typedef typename boost::graph_traits<Mesh>::face_descriptor face_descriptor;
|
||||
|
||||
BOOST_FOREACH(face_descriptor f, faces(mesh))
|
||||
{
|
||||
if(f == boost::graph_traits<Mesh>::null_face())
|
||||
continue;
|
||||
|
||||
double area = CGAL::Polygon_mesh_processing::face_area(f, mesh);
|
||||
acc(area);
|
||||
}
|
||||
|
||||
min_a = extract_result< tag::min >(acc);
|
||||
max_a = extract_result< tag::max >(acc);
|
||||
avg_a = extract_result< tag::mean >(acc);
|
||||
}
|
||||
|
||||
template<class T>
|
||||
bool check_value_equal(T a, T b)
|
||||
{
|
||||
if (abs(a-b) > 1e-3)
|
||||
{
|
||||
std::cerr << "Value not equal! ";
|
||||
std::cerr << a << " is not " << b << std::endl;
|
||||
return false;
|
||||
}
|
||||
return true;
|
||||
}
|
||||
|
||||
|
||||
int main(int argc, char* argv[]){
|
||||
|
||||
std::string filename;
|
||||
std::ifstream input;
|
||||
Mesh mesh;
|
||||
double min_a, max_a, mean_a;
|
||||
|
||||
filename = "data/curved_polygon.off";
|
||||
input.open(filename);
|
||||
if (!input || !(input >> mesh) || mesh.is_empty()) {
|
||||
std::cerr << "Not a valid .off file." << std::endl;
|
||||
return 1;
|
||||
}
|
||||
input.close();
|
||||
|
||||
calc_angles(mesh, min_a, max_a, mean_a);
|
||||
|
||||
CGAL::Polygon_mesh_processing::angle_smoothing(mesh);
|
||||
calc_angles(mesh, min_a, max_a, mean_a);
|
||||
|
||||
if(!check_value_equal(min_a, 24.980))
|
||||
{
|
||||
return EXIT_FAILURE;
|
||||
}
|
||||
if(!check_value_equal(max_a, 108.723))
|
||||
{
|
||||
return EXIT_FAILURE;
|
||||
}
|
||||
|
||||
input.open(filename);
|
||||
input>>mesh;
|
||||
input.close();
|
||||
|
||||
CGAL::Polygon_mesh_processing::area_smoothing(mesh);
|
||||
calc_areas(mesh, min_a, max_a, mean_a);
|
||||
|
||||
if(!check_value_equal(min_a, 0.476))
|
||||
{
|
||||
return EXIT_FAILURE;
|
||||
}
|
||||
if(!check_value_equal(max_a, 0.567))
|
||||
{
|
||||
return EXIT_FAILURE;
|
||||
}
|
||||
if(!check_value_equal(mean_a, 0.542))
|
||||
{
|
||||
return EXIT_FAILURE;
|
||||
}
|
||||
|
||||
return 0;
|
||||
}
|
||||
|
|
@ -152,3 +152,7 @@ polyhedron_demo_plugin(degenerated_faces_sm_plugin Degenerated_faces_plugin)
|
|||
target_link_libraries(degenerated_faces_sm_plugin scene_surface_mesh_item scene_surface_mesh_selection_item)
|
||||
target_compile_definitions(degenerated_faces_sm_plugin PUBLIC "-DUSE_SURFACE_MESH" )
|
||||
|
||||
qt5_wrap_ui( smoothingUI_FILES Smoothing_plugin.ui)
|
||||
polyhedron_demo_plugin(smoothing_plugin Smoothing_plugin ${smoothingUI_FILES})
|
||||
target_link_libraries(smoothing_plugin scene_polyhedron_item scene_polyhedron_selection_item)
|
||||
|
||||
|
|
|
|||
|
|
@ -0,0 +1,247 @@
|
|||
|
||||
#include <QtCore/qglobal.h>
|
||||
#include <QTime>
|
||||
#include <QAction>
|
||||
#include <QMainWindow>
|
||||
#include <QApplication>
|
||||
#include <QDockWidget>
|
||||
#include <QString>
|
||||
#include <QInputDialog>
|
||||
#include <QtPlugin>
|
||||
#include <QMessageBox>
|
||||
|
||||
#include <CGAL/Three/Polyhedron_demo_plugin_interface.h>
|
||||
#include <CGAL/Three/Polyhedron_demo_plugin_helper.h>
|
||||
|
||||
#include "Scene_polyhedron_item.h"
|
||||
#include "Scene_polyhedron_selection_item.h"
|
||||
#include "Polyhedron_type.h"
|
||||
|
||||
#include <CGAL/iterator.h>
|
||||
#include <CGAL/boost/graph/graph_traits_Polyhedron_3.h>
|
||||
#include <CGAL/boost/graph/properties_Polyhedron_3.h>
|
||||
#include <CGAL/utility.h>
|
||||
#include <boost/graph/graph_traits.hpp>
|
||||
#include <CGAL/property_map.h>
|
||||
|
||||
#include <CGAL/Polygon_mesh_processing/smoothing.h>
|
||||
|
||||
#include "ui_Smoothing_plugin.h"
|
||||
|
||||
|
||||
using namespace CGAL::Polygon_mesh_processing;
|
||||
using namespace CGAL::Three;
|
||||
class Polyhedron_demo_smothing_plugin :
|
||||
public QObject,
|
||||
public Polyhedron_demo_plugin_helper
|
||||
{
|
||||
Q_OBJECT
|
||||
Q_INTERFACES(CGAL::Three::Polyhedron_demo_plugin_interface)
|
||||
Q_PLUGIN_METADATA(IID "com.geometryfactory.PolyhedronDemo.PluginInterface/1.0")
|
||||
|
||||
|
||||
|
||||
public:
|
||||
void init(QMainWindow* mainWindow, Scene_interface* scene_interface, Messages_interface*)
|
||||
{
|
||||
scene = scene_interface;
|
||||
mw = mainWindow;
|
||||
|
||||
actionSmoothing_ = new QAction(tr("Smoothing"), mw);
|
||||
actionSmoothing_->setProperty("subMenuName", "Polygon Mesh Processing");
|
||||
|
||||
connect(actionSmoothing_, SIGNAL(triggered()), this, SLOT(smoothing_action()));
|
||||
|
||||
dock_widget = new QDockWidget("Smoothing", mw);
|
||||
dock_widget->setVisible(false);
|
||||
|
||||
ui_widget.setupUi(dock_widget);
|
||||
addDockWidget(dock_widget);
|
||||
|
||||
connect(ui_widget.Apply_button, SIGNAL(clicked()), this, SLOT(on_Apply_by_type_clicked()));
|
||||
connect(ui_widget.Curvature_Button, SIGNAL(clicked()), this, SLOT(on_Apply_curvature_clicked()));
|
||||
connect(ui_widget.Run_convergence_button, SIGNAL(clicked()), this, SLOT(on_Run_convergence_clicked()));
|
||||
|
||||
}
|
||||
|
||||
QList<QAction*> actions() const
|
||||
{
|
||||
return QList<QAction*>() << actionSmoothing_;
|
||||
}
|
||||
|
||||
bool applicable(QAction*) const
|
||||
{
|
||||
const Scene_interface::Item_id index = scene->mainSelectionIndex();
|
||||
if (qobject_cast<Scene_polyhedron_item*>(scene->item(index)))
|
||||
return true;
|
||||
else if (qobject_cast<Scene_polyhedron_selection_item*>(scene->item(index)))
|
||||
return true;
|
||||
else
|
||||
return false;
|
||||
}
|
||||
|
||||
virtual void closure()
|
||||
{
|
||||
dock_widget->hide();
|
||||
}
|
||||
|
||||
void init_ui()
|
||||
{
|
||||
ui_widget.Angle_spinBox->setValue(1);
|
||||
ui_widget.Angle_spinBox->setSingleStep(1);
|
||||
ui_widget.Angle_spinBox->setMinimum(1);
|
||||
|
||||
ui_widget.Area_spinBox->setValue(1);
|
||||
ui_widget.Area_spinBox->setSingleStep(1);
|
||||
ui_widget.Area_spinBox->setMinimum(1);
|
||||
|
||||
ui_widget.gd_dSpinBox->setSingleStep(0.0001);
|
||||
ui_widget.gd_dSpinBox->setDecimals(4);
|
||||
ui_widget.gd_dSpinBox->setMinimum(0.0001);
|
||||
ui_widget.gd_dSpinBox->setValue(0.001);
|
||||
|
||||
ui_widget.use_weights_checkBox->setChecked(true);
|
||||
|
||||
ui_widget.dist_dSpinBox->setValue(0.01);
|
||||
ui_widget.dist_dSpinBox->setSingleStep(0.0001);
|
||||
ui_widget.dist_dSpinBox->setDecimals(4);
|
||||
ui_widget.dist_dSpinBox->setMinimum(0.0001);
|
||||
|
||||
ui_widget.gd_precision_label->setToolTip("Tradeoff between precision and speed. Less is more precise.");
|
||||
ui_widget.distance_label->setToolTip("Tradeoff between precision and speed. Less is more precise.");
|
||||
|
||||
ui_widget.iterations_spinBox->setValue(20);
|
||||
ui_widget.iterations_spinBox->setSingleStep(1);
|
||||
ui_widget.iterations_spinBox->setMinimum(1);
|
||||
|
||||
ui_widget.curv_iterations_spinBox->setValue(1);
|
||||
ui_widget.curv_iterations_spinBox->setSingleStep(1);
|
||||
ui_widget.curv_iterations_spinBox->setMinimum(1);
|
||||
}
|
||||
|
||||
public Q_SLOTS:
|
||||
void smoothing_action()
|
||||
{
|
||||
dock_widget->show();
|
||||
dock_widget->raise();
|
||||
|
||||
const Scene_interface::Item_id index = scene->mainSelectionIndex();
|
||||
Scene_polyhedron_item* poly_item = qobject_cast<Scene_polyhedron_item*>(scene->item(index));
|
||||
|
||||
if(poly_item)
|
||||
{
|
||||
init_ui();
|
||||
}
|
||||
}
|
||||
|
||||
void on_Apply_by_type_clicked()
|
||||
{
|
||||
const Scene_interface::Item_id index = scene->mainSelectionIndex();
|
||||
Scene_polyhedron_item* poly_item = qobject_cast<Scene_polyhedron_item*>(scene->item(index));
|
||||
Polyhedron& pmesh = *poly_item->polyhedron();
|
||||
|
||||
QApplication::setOverrideCursor(Qt::WaitCursor);
|
||||
|
||||
if(ui_widget.Angle_checkBox->isChecked())
|
||||
{
|
||||
unsigned int nb_iter = ui_widget.Angle_spinBox->value();
|
||||
bool use_weights = ui_widget.use_weights_checkBox->isChecked();
|
||||
angle_smoothing(pmesh,
|
||||
parameters::number_of_iterations(nb_iter).
|
||||
use_weights(use_weights));
|
||||
|
||||
poly_item->invalidateOpenGLBuffers();
|
||||
Q_EMIT poly_item->itemChanged();
|
||||
}
|
||||
|
||||
if(ui_widget.Area_checkBox->isChecked())
|
||||
{
|
||||
unsigned int nb_iter = ui_widget.Area_spinBox->value();
|
||||
double gd_precision = ui_widget.gd_dSpinBox->value();
|
||||
area_smoothing(pmesh,
|
||||
parameters::number_of_iterations(nb_iter).
|
||||
gradient_descent_precision(gd_precision));
|
||||
|
||||
poly_item->invalidateOpenGLBuffers();
|
||||
Q_EMIT poly_item->itemChanged();
|
||||
}
|
||||
|
||||
QApplication::restoreOverrideCursor();
|
||||
}
|
||||
|
||||
void on_Apply_curvature_clicked()
|
||||
{
|
||||
const Scene_interface::Item_id index = scene->mainSelectionIndex();
|
||||
Scene_polyhedron_item* poly_item = qobject_cast<Scene_polyhedron_item*>(scene->item(index));
|
||||
Polyhedron& pmesh = *poly_item->polyhedron();
|
||||
|
||||
QApplication::setOverrideCursor(Qt::WaitCursor);
|
||||
|
||||
unsigned int nb_iter = ui_widget.curv_iterations_spinBox->value();
|
||||
curvature_flow_smoothing(pmesh,
|
||||
parameters::number_of_iterations(nb_iter));
|
||||
|
||||
poly_item->invalidateOpenGLBuffers();
|
||||
Q_EMIT poly_item->itemChanged();
|
||||
|
||||
QApplication::restoreOverrideCursor();
|
||||
}
|
||||
|
||||
void on_Run_convergence_clicked()
|
||||
{
|
||||
const Scene_interface::Item_id index = scene->mainSelectionIndex();
|
||||
Scene_polyhedron_item* poly_item = qobject_cast<Scene_polyhedron_item*>(scene->item(index));
|
||||
Polyhedron& pmesh = *poly_item->polyhedron();
|
||||
|
||||
QApplication::setOverrideCursor(Qt::WaitCursor);
|
||||
|
||||
unsigned int nb_iter = ui_widget.iterations_spinBox->value();
|
||||
double dist = ui_widget.dist_dSpinBox->value();
|
||||
double gd_precision = ui_widget.gd_dSpinBox->value();
|
||||
bool use_weights = ui_widget.use_weights_checkBox->isChecked();
|
||||
|
||||
compatible_smoothing(pmesh,
|
||||
parameters::number_of_iterations(nb_iter).
|
||||
distance_precision(dist).
|
||||
gradient_descent_precision(gd_precision).
|
||||
use_weights(use_weights));
|
||||
|
||||
poly_item->invalidateOpenGLBuffers();
|
||||
Q_EMIT poly_item->itemChanged();
|
||||
|
||||
QApplication::restoreOverrideCursor();
|
||||
}
|
||||
|
||||
|
||||
|
||||
|
||||
private:
|
||||
QAction* actionSmoothing_;
|
||||
QDockWidget* dock_widget;
|
||||
Ui::Smoothing ui_widget;
|
||||
|
||||
|
||||
|
||||
};
|
||||
|
||||
|
||||
|
||||
#include "Smoothing_plugin.moc"
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
|
@ -0,0 +1,269 @@
|
|||
<?xml version="1.0" encoding="UTF-8"?>
|
||||
<ui version="4.0">
|
||||
<class>Smoothing</class>
|
||||
<widget class="QDockWidget" name="Smoothing">
|
||||
<property name="geometry">
|
||||
<rect>
|
||||
<x>0</x>
|
||||
<y>0</y>
|
||||
<width>837</width>
|
||||
<height>708</height>
|
||||
</rect>
|
||||
</property>
|
||||
<property name="windowTitle">
|
||||
<string/>
|
||||
</property>
|
||||
<widget class="QWidget" name="dockWidgetContents">
|
||||
<layout class="QVBoxLayout" name="verticalLayout_3">
|
||||
<item>
|
||||
<widget class="QGroupBox" name="Mesh_groupBox">
|
||||
<property name="title">
|
||||
<string>Mesh smoothing</string>
|
||||
</property>
|
||||
<layout class="QVBoxLayout" name="verticalLayout_4">
|
||||
<item>
|
||||
<widget class="QGroupBox" name="Parameters_groupBox">
|
||||
<property name="sizePolicy">
|
||||
<sizepolicy hsizetype="Expanding" vsizetype="Expanding">
|
||||
<horstretch>0</horstretch>
|
||||
<verstretch>0</verstretch>
|
||||
</sizepolicy>
|
||||
</property>
|
||||
<property name="layoutDirection">
|
||||
<enum>Qt::LeftToRight</enum>
|
||||
</property>
|
||||
<property name="title">
|
||||
<string>Parameters</string>
|
||||
</property>
|
||||
<widget class="QWidget" name="layoutWidget">
|
||||
<property name="geometry">
|
||||
<rect>
|
||||
<x>10</x>
|
||||
<y>30</y>
|
||||
<width>718</width>
|
||||
<height>56</height>
|
||||
</rect>
|
||||
</property>
|
||||
<layout class="QHBoxLayout" name="horizontalLayout_2">
|
||||
<item>
|
||||
<layout class="QHBoxLayout" name="horizontalLayout">
|
||||
<item>
|
||||
<widget class="QLabel" name="gd_precision_label">
|
||||
<property name="text">
|
||||
<string>Gradient descent precision:</string>
|
||||
</property>
|
||||
</widget>
|
||||
</item>
|
||||
<item>
|
||||
<widget class="QDoubleSpinBox" name="gd_dSpinBox">
|
||||
<property name="decimals">
|
||||
<number>4</number>
|
||||
</property>
|
||||
</widget>
|
||||
</item>
|
||||
</layout>
|
||||
</item>
|
||||
<item>
|
||||
<widget class="QCheckBox" name="use_weights_checkBox">
|
||||
<property name="layoutDirection">
|
||||
<enum>Qt::RightToLeft</enum>
|
||||
</property>
|
||||
<property name="text">
|
||||
<string>use weights</string>
|
||||
</property>
|
||||
</widget>
|
||||
</item>
|
||||
</layout>
|
||||
</widget>
|
||||
</widget>
|
||||
</item>
|
||||
<item>
|
||||
<widget class="QSplitter" name="splitter">
|
||||
<property name="orientation">
|
||||
<enum>Qt::Horizontal</enum>
|
||||
</property>
|
||||
<widget class="QGroupBox" name="Type_groupBox">
|
||||
<property name="enabled">
|
||||
<bool>true</bool>
|
||||
</property>
|
||||
<property name="layoutDirection">
|
||||
<enum>Qt::LeftToRight</enum>
|
||||
</property>
|
||||
<property name="title">
|
||||
<string>Remeshing by type </string>
|
||||
</property>
|
||||
<property name="flat">
|
||||
<bool>false</bool>
|
||||
</property>
|
||||
<property name="checkable">
|
||||
<bool>false</bool>
|
||||
</property>
|
||||
<layout class="QVBoxLayout" name="verticalLayout">
|
||||
<item>
|
||||
<layout class="QGridLayout" name="gridLayout">
|
||||
<item row="6" column="0">
|
||||
<widget class="QCheckBox" name="Area_checkBox">
|
||||
<property name="text">
|
||||
<string>area</string>
|
||||
</property>
|
||||
</widget>
|
||||
</item>
|
||||
<item row="6" column="1">
|
||||
<widget class="QSpinBox" name="Area_spinBox"/>
|
||||
</item>
|
||||
<item row="0" column="0">
|
||||
<widget class="QLabel" name="label">
|
||||
<property name="text">
|
||||
<string>Type</string>
|
||||
</property>
|
||||
<property name="alignment">
|
||||
<set>Qt::AlignCenter</set>
|
||||
</property>
|
||||
</widget>
|
||||
</item>
|
||||
<item row="0" column="1">
|
||||
<widget class="QLabel" name="label_2">
|
||||
<property name="text">
|
||||
<string>Iterations</string>
|
||||
</property>
|
||||
<property name="alignment">
|
||||
<set>Qt::AlignCenter</set>
|
||||
</property>
|
||||
</widget>
|
||||
</item>
|
||||
<item row="5" column="1">
|
||||
<widget class="QSpinBox" name="Angle_spinBox"/>
|
||||
</item>
|
||||
<item row="5" column="0">
|
||||
<widget class="QCheckBox" name="Angle_checkBox">
|
||||
<property name="sizePolicy">
|
||||
<sizepolicy hsizetype="Minimum" vsizetype="Minimum">
|
||||
<horstretch>0</horstretch>
|
||||
<verstretch>0</verstretch>
|
||||
</sizepolicy>
|
||||
</property>
|
||||
<property name="text">
|
||||
<string>angle</string>
|
||||
</property>
|
||||
</widget>
|
||||
</item>
|
||||
</layout>
|
||||
</item>
|
||||
<item>
|
||||
<widget class="QPushButton" name="Apply_button">
|
||||
<property name="text">
|
||||
<string>Apply</string>
|
||||
</property>
|
||||
</widget>
|
||||
</item>
|
||||
</layout>
|
||||
</widget>
|
||||
<widget class="QGroupBox" name="Compatible_groupBox">
|
||||
<property name="title">
|
||||
<string>Compatible remeshing algorithm </string>
|
||||
</property>
|
||||
<layout class="QVBoxLayout" name="verticalLayout_2">
|
||||
<item>
|
||||
<layout class="QGridLayout" name="gridLayout_3">
|
||||
<item row="1" column="0">
|
||||
<widget class="QLabel" name="distance_label">
|
||||
<property name="layoutDirection">
|
||||
<enum>Qt::RightToLeft</enum>
|
||||
</property>
|
||||
<property name="text">
|
||||
<string>Hausdorff distance:</string>
|
||||
</property>
|
||||
<property name="alignment">
|
||||
<set>Qt::AlignRight|Qt::AlignTrailing|Qt::AlignVCenter</set>
|
||||
</property>
|
||||
</widget>
|
||||
</item>
|
||||
<item row="1" column="1">
|
||||
<widget class="QDoubleSpinBox" name="dist_dSpinBox">
|
||||
<property name="decimals">
|
||||
<number>4</number>
|
||||
</property>
|
||||
</widget>
|
||||
</item>
|
||||
<item row="2" column="0">
|
||||
<widget class="QLabel" name="iterations_label">
|
||||
<property name="layoutDirection">
|
||||
<enum>Qt::RightToLeft</enum>
|
||||
</property>
|
||||
<property name="text">
|
||||
<string>max iterations: </string>
|
||||
</property>
|
||||
<property name="alignment">
|
||||
<set>Qt::AlignRight|Qt::AlignTrailing|Qt::AlignVCenter</set>
|
||||
</property>
|
||||
</widget>
|
||||
</item>
|
||||
<item row="2" column="1">
|
||||
<widget class="QSpinBox" name="iterations_spinBox">
|
||||
<property name="sizePolicy">
|
||||
<sizepolicy hsizetype="Minimum" vsizetype="Fixed">
|
||||
<horstretch>0</horstretch>
|
||||
<verstretch>0</verstretch>
|
||||
</sizepolicy>
|
||||
</property>
|
||||
</widget>
|
||||
</item>
|
||||
<item row="0" column="0" colspan="2">
|
||||
<widget class="QLabel" name="label_3">
|
||||
<property name="text">
|
||||
<string>Convergence criteria</string>
|
||||
</property>
|
||||
</widget>
|
||||
</item>
|
||||
</layout>
|
||||
</item>
|
||||
<item>
|
||||
<widget class="QPushButton" name="Run_convergence_button">
|
||||
<property name="text">
|
||||
<string>Run to convergence</string>
|
||||
</property>
|
||||
</widget>
|
||||
</item>
|
||||
</layout>
|
||||
</widget>
|
||||
</widget>
|
||||
</item>
|
||||
</layout>
|
||||
</widget>
|
||||
</item>
|
||||
<item>
|
||||
<widget class="QGroupBox" name="Shape_groupBox">
|
||||
<property name="title">
|
||||
<string>Shape smoothing</string>
|
||||
</property>
|
||||
<layout class="QHBoxLayout" name="horizontalLayout_4">
|
||||
<item>
|
||||
<widget class="QPushButton" name="Curvature_Button">
|
||||
<property name="text">
|
||||
<string>Apply curvature flow algorithm</string>
|
||||
</property>
|
||||
</widget>
|
||||
</item>
|
||||
<item>
|
||||
<layout class="QHBoxLayout" name="horizontalLayout_3">
|
||||
<item>
|
||||
<widget class="QLabel" name="iterations_label_2">
|
||||
<property name="text">
|
||||
<string>Iterations:</string>
|
||||
</property>
|
||||
</widget>
|
||||
</item>
|
||||
<item>
|
||||
<widget class="QSpinBox" name="curv_iterations_spinBox"/>
|
||||
</item>
|
||||
</layout>
|
||||
</item>
|
||||
</layout>
|
||||
</widget>
|
||||
</item>
|
||||
</layout>
|
||||
</widget>
|
||||
</widget>
|
||||
<resources/>
|
||||
<connections/>
|
||||
</ui>
|
||||
|
|
@ -212,6 +212,7 @@ triangulate_facets_plugin \
|
|||
trivial_plugin \
|
||||
vtk_plugin \
|
||||
xyz_plugin \
|
||||
smoothing_plugin \
|
||||
all
|
||||
do
|
||||
if ${MAKE_CMD} -f Makefile help | grep "$target" > /dev/null; then
|
||||
|
|
|
|||
Loading…
Reference in New Issue