mirror of https://github.com/CGAL/cgal
Merge pull request #1947 from janetournois/Mesh_3-add_lipschitz_sizing-example-jtournois
Mesh_3 - add an example for `Lipschitz_sizing`
This commit is contained in:
Sebastien Loriot
GitHub
commit
cce6932af9
|
|
@ -809,6 +809,21 @@ a 3D medical image.
|
|||
View of a 3D mesh produced from a 3D image with different size for different organs. Code from subsection \ref Mesh_3_subsubsection_examples_3d_image_variable generates this file.
|
||||
\cgalFigureEnd
|
||||
|
||||
\subsubsection Mesh_3LipschitzSizingField Lipschitz Sizing Field
|
||||
|
||||
\anchor Mesh_3_subsubsection_examples_3d_image_lip_sizing
|
||||
|
||||
The following example shows how to use another custom sizing function, that is
|
||||
`k`-Lipschitz. For each subdomain, the user provides the parameter `k`,
|
||||
a minimal size and maximal size for cells.
|
||||
|
||||
\cgalExample{Mesh_3/mesh_polyhedral_domain_with_lipschitz_sizing.cpp}
|
||||
|
||||
\cgalFigureBegin{figurefandisk_lip_sizing,fandisk_lip.jpg}
|
||||
View of a 3D mesh produced with a Lipschitz sizing field,
|
||||
for a single subdomain. Code from subsection \ref Mesh_3_subsubsection_examples_3d_image_lip_sizing generates this mesh.
|
||||
\cgalFigureEnd
|
||||
|
||||
\subsection Mesh_3MeshingDomainswithSharpFeatures Meshing Domains with Sharp Features
|
||||
|
||||
\anchor Mesh_3_subsection_example_polyhedral_with_edges
|
||||
|
|
|
|||
|
|
@ -17,6 +17,7 @@
|
|||
\example Mesh_3/mesh_polyhedral_domain.cpp
|
||||
\example Mesh_3/remesh_polyhedral_surface.cpp
|
||||
\example Mesh_3/mesh_polyhedral_domain_with_features.cpp
|
||||
\example Mesh_3/mesh_polyhedral_domain_with_lipschitz_sizing.cpp
|
||||
\example Mesh_3/mesh_two_implicit_spheres_with_balls.cpp
|
||||
\example Mesh_3/mesh_3D_gray_image.cpp
|
||||
\example Mesh_3/mesh_3D_gray_vtk_image.cpp
|
||||
|
|
|
|||
Binary file not shown.
|
After Width: | Height: | Size: 109 KiB |
|
|
@ -75,6 +75,8 @@ if ( CGAL_FOUND )
|
|||
create_single_source_cgal_program( "mesh_polyhedral_domain.cpp" )
|
||||
create_single_source_cgal_program( "remesh_polyhedral_surface.cpp" )
|
||||
create_single_source_cgal_program( "mesh_polyhedral_domain_with_features.cpp" )
|
||||
create_single_source_cgal_program( "mesh_polyhedral_domain_with_lipschitz_sizing.cpp" )
|
||||
|
||||
if( WITH_CGAL_ImageIO )
|
||||
if( VTK_FOUND AND "${VTK_VERSION_MAJOR}" GREATER "5" )
|
||||
add_executable ( mesh_3D_gray_vtk_image mesh_3D_gray_vtk_image.cpp )
|
||||
|
|
|
|||
|
|
@ -0,0 +1,80 @@
|
|||
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
|
||||
|
||||
#include <CGAL/Mesh_triangulation_3.h>
|
||||
#include <CGAL/Mesh_complex_3_in_triangulation_3.h>
|
||||
#include <CGAL/Mesh_criteria_3.h>
|
||||
|
||||
#include <CGAL/Polyhedral_mesh_domain_with_features_3.h>
|
||||
#include <CGAL/make_mesh_3.h>
|
||||
|
||||
#include <CGAL/Mesh_3/experimental/Lipschitz_sizing_polyhedron.h>
|
||||
|
||||
//Kernel
|
||||
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
|
||||
typedef K::FT FT;
|
||||
|
||||
// Domain
|
||||
typedef CGAL::Mesh_polyhedron_3<K>::type Polyhedron;
|
||||
typedef CGAL::Polyhedral_mesh_domain_with_features_3<K> Mesh_domain;
|
||||
|
||||
#ifdef CGAL_CONCURRENT_MESH_3
|
||||
typedef CGAL::Parallel_tag Concurrency_tag;
|
||||
#else
|
||||
typedef CGAL::Sequential_tag Concurrency_tag;
|
||||
#endif
|
||||
|
||||
// Triangulation
|
||||
typedef CGAL::Mesh_triangulation_3<Mesh_domain,CGAL::Default,Concurrency_tag>::type Tr;
|
||||
typedef CGAL::Mesh_complex_3_in_triangulation_3<
|
||||
Tr,Mesh_domain::Corner_index,Mesh_domain::Curve_segment_index> C3t3;
|
||||
|
||||
// Criteria
|
||||
typedef CGAL::Mesh_criteria_3<Tr> Mesh_criteria;
|
||||
|
||||
// Sizing field
|
||||
typedef CGAL::Mesh_3::Lipschitz_sizing<K, Mesh_domain, Mesh_domain::AABB_tree> Lip_sizing;
|
||||
|
||||
// To avoid verbose function and named parameters call
|
||||
using namespace CGAL::parameters;
|
||||
|
||||
int main(int argc, char*argv[])
|
||||
{
|
||||
const char* fname = (argc>1) ? argv[1] : "data/fandisk.off";
|
||||
std::ifstream input(fname);
|
||||
Polyhedron polyhedron;
|
||||
input >> polyhedron;
|
||||
if (input.fail()){
|
||||
std::cerr << "Error: Cannot read file " << fname << std::endl;
|
||||
return EXIT_FAILURE;
|
||||
}
|
||||
// Create domain
|
||||
Mesh_domain domain(polyhedron);
|
||||
|
||||
// Get sharp features
|
||||
domain.detect_features();
|
||||
|
||||
// Create Lipschitz sizing field
|
||||
Lip_sizing lip_sizing(domain, &domain.aabb_tree());
|
||||
|
||||
FT min_size = 0.02;
|
||||
lip_sizing.add_parameters_for_subdomain(1, //subdomain id
|
||||
0.3, //k
|
||||
min_size,//min_size
|
||||
0.5); //max_size
|
||||
|
||||
// Mesh criteria
|
||||
Mesh_criteria criteria(edge_size = min_size,
|
||||
facet_angle = 25,
|
||||
facet_size = min_size,
|
||||
facet_distance = 0.005,
|
||||
cell_radius_edge_ratio = 3,
|
||||
cell_size = lip_sizing);
|
||||
// Mesh generation
|
||||
C3t3 c3t3 = CGAL::make_mesh_3<C3t3>(domain, criteria);
|
||||
|
||||
// Output
|
||||
std::ofstream medit_file("out.mesh");
|
||||
c3t3.output_to_medit(medit_file);
|
||||
|
||||
return EXIT_SUCCESS;
|
||||
}
|
||||
|
|
@ -0,0 +1,469 @@
|
|||
// Copyright (c) 2016 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) : Jane Tournois
|
||||
//
|
||||
|
||||
#ifndef CGAL_LIPSCHITZ_SIZING_H
|
||||
#define CGAL_LIPSCHITZ_SIZING_H
|
||||
|
||||
#include <CGAL/license/Mesh_3.h>
|
||||
|
||||
#include <CGAL/AABB_tree.h>
|
||||
#include <CGAL/AABB_traits.h>
|
||||
#include <CGAL/AABB_triangle_primitive.h>
|
||||
|
||||
#include <CGAL/Mesh_3/experimental/AABB_filtered_projection_traits.h>
|
||||
#include <CGAL/Mesh_3/experimental/Get_facet_patch_id.h>
|
||||
|
||||
#include <CGAL/Mesh_3/experimental/Lipschitz_sizing_parameters.h>
|
||||
|
||||
#include <CGAL/Default.h>
|
||||
#include <CGAL/array.h>
|
||||
#include <CGAL/Bbox_3.h>
|
||||
|
||||
#include <boost/shared_ptr.hpp>
|
||||
|
||||
#include <list>
|
||||
#include <sstream>
|
||||
#include <string>
|
||||
#include <fstream>
|
||||
#include <vector>
|
||||
|
||||
namespace CGAL
|
||||
{
|
||||
namespace Mesh_3
|
||||
{
|
||||
|
||||
template <class Kernel, class MeshDomain
|
||||
, typename AABBTreeTemplate = CGAL::Default
|
||||
#ifdef CGAL_MESH_3_EXPERIMENTAL_USE_PATCHES_IDS
|
||||
, typename Get_facet_patch_id_ = CGAL::Default
|
||||
, typename Patches_ids_ = CGAL::Default
|
||||
#endif
|
||||
>
|
||||
class Lipschitz_sizing
|
||||
{
|
||||
public:
|
||||
typedef Kernel K;
|
||||
typedef typename Kernel::FT FT;
|
||||
typedef typename Kernel::Triangle_3 Triangle;
|
||||
typedef typename Kernel::Point_3 Point_3;
|
||||
|
||||
typedef typename std::list<Triangle>::iterator Tr_iterator;
|
||||
typedef CGAL::AABB_triangle_primitive<K, Tr_iterator> Primitive;
|
||||
typedef CGAL::AABB_traits<K, Primitive> AABB_tr_traits;
|
||||
typedef CGAL::AABB_tree<AABB_tr_traits> AABB_tree;
|
||||
|
||||
typedef typename CGAL::Default::Get<AABBTreeTemplate, AABB_tree>::type Tree;
|
||||
|
||||
typedef typename MeshDomain::Index Index;
|
||||
typedef typename MeshDomain::Subdomain_index Subdomain_index;
|
||||
typedef typename MeshDomain::Surface_patch_index Surface_patch_index;
|
||||
|
||||
typedef CGAL::Lipschitz_sizing_parameters<MeshDomain, FT> Parameters;
|
||||
|
||||
#ifdef CGAL_MESH_3_EXPERIMENTAL_USE_PATCHES_IDS
|
||||
typedef typename CGAL::Default::Get<Patches_ids_,
|
||||
typename MeshDomain::Surface_patch_index_set>::type Patches_ids;
|
||||
typedef std::vector<Patches_ids> Patches_ids_map;
|
||||
|
||||
typedef typename CGAL::Default::Get<
|
||||
Get_facet_patch_id_,
|
||||
CGAL::Mesh_3::Get_facet_patch_id<typename Tree::Primitive>
|
||||
>::type Get_facet_patch_id;
|
||||
|
||||
typedef CGAL::Mesh_3::Filtered_projection_traits<
|
||||
typename Tree::AABB_traits, Get_facet_patch_id> AABB_filtered_traits;
|
||||
|
||||
private:
|
||||
typedef CGAL::Search_traits_3<Kernel> KdTreeTraits;
|
||||
typedef CGAL::Orthogonal_k_neighbor_search<KdTreeTraits> Neighbor_search;
|
||||
typedef typename Neighbor_search::Tree Kd_tree;
|
||||
#endif
|
||||
|
||||
private:
|
||||
//only one of these aabb_trees is needed
|
||||
const Tree* m_ptree;
|
||||
boost::shared_ptr<Tree> m_own_ptree;
|
||||
|
||||
const MeshDomain& m_domain;
|
||||
Parameters m_params;
|
||||
|
||||
const CGAL::cpp11::array<double, 3>& m_vxyz;
|
||||
const CGAL::Bbox_3& m_bbox;
|
||||
const bool m_domain_is_a_box;
|
||||
|
||||
#ifdef CGAL_MESH_3_EXPERIMENTAL_USE_PATCHES_IDS
|
||||
//help to accelerate aabb_tree queries in m_ptree
|
||||
boost::shared_ptr<Kd_tree> m_kd_tree;
|
||||
|
||||
Get_facet_patch_id m_get_facet_patch_id;
|
||||
const Patches_ids_map& patches_ids_map;
|
||||
#endif
|
||||
|
||||
public:
|
||||
Lipschitz_sizing(const MeshDomain& domain)
|
||||
: m_ptree(NULL)
|
||||
, m_own_ptree()
|
||||
, m_domain(domain)
|
||||
, m_params(domain)
|
||||
{
|
||||
}
|
||||
|
||||
Lipschitz_sizing(const MeshDomain& domain
|
||||
, const Tree* ptree
|
||||
, const CGAL::cpp11::array<double, 3>& vxyz
|
||||
, const CGAL::Bbox_3& bbox
|
||||
, const bool domain_is_a_box
|
||||
#ifdef CGAL_MESH_3_EXPERIMENTAL_USE_PATCHES_IDS
|
||||
, const Patches_ids_map& patches_ids_map
|
||||
#endif
|
||||
)
|
||||
: m_ptree(ptree)
|
||||
, m_own_ptree()
|
||||
, m_domain(domain)
|
||||
, m_params(domain)
|
||||
, m_vxyz(vxyz)
|
||||
, m_bbox(bbox)
|
||||
, m_domain_is_a_box(domain_is_a_box)
|
||||
#ifdef CGAL_MESH_3_EXPERIMENTAL_USE_PATCHES_IDS
|
||||
, m_get_facet_patch_id()
|
||||
, patches_ids_map(patches_ids_map)
|
||||
#endif
|
||||
{
|
||||
}
|
||||
|
||||
FT operator()(const Point_3& p, const int dim, const Index& index) const
|
||||
{
|
||||
CGAL_assertion(!m_params.empty());
|
||||
#ifdef CGAL_MESH_3_LIPSCHITZ_SIZING_VERBOSE
|
||||
std::cout << "D = " << dim << "\t";
|
||||
#endif
|
||||
if (dim == 3)
|
||||
{
|
||||
return size_in_subdomain(p, m_domain.subdomain_index(index));
|
||||
}
|
||||
else if (dim == 2)
|
||||
{
|
||||
#ifdef CGAL_MESH_3_EXPERIMENTAL_USE_PATCHES_IDS
|
||||
Surface_patch_index sp_index = m_domain.surface_patch_index(index);
|
||||
|
||||
if(!is_on_cube_boundary(sp_index)
|
||||
&& !is_on_cube_boundary(p))
|
||||
{
|
||||
#ifdef CGAL_MESH_3_LIPSCHITZ_SIZING_VERBOSE
|
||||
std::cout << " \n";
|
||||
#endif
|
||||
FT size_max;
|
||||
m_params.get_parameters(sp_index, size_max);
|
||||
return size_max;
|
||||
}
|
||||
else
|
||||
{
|
||||
#ifdef CGAL_MESH_3_LIPSCHITZ_SIZING_VERBOSE
|
||||
std::cout << "(on cube) ";
|
||||
#endif
|
||||
const std::pair<Subdomain_index, Subdomain_index>& index
|
||||
= m_params.incident_subdomains(sp_index);
|
||||
|
||||
#ifdef CGAL_MESH_3_IMAGE_EXAMPLE
|
||||
if (index.first == INT_MIN)
|
||||
return size_in_subdomain(p, index.second);
|
||||
else
|
||||
return size_in_subdomain(p, index.first);
|
||||
#else //==POLYHEDRAL EXAMPLE
|
||||
if (!is_in_domain(index.first))
|
||||
return size_in_subdomain(p, index.second);
|
||||
else
|
||||
return size_in_subdomain(p, index.first);
|
||||
#endif //CGAL_MESH_3_IMAGE_EXAMPLE
|
||||
}
|
||||
|
||||
#else //CGAL_MESH_3_EXPERIMENTAL_USE_PATCHES_IDS
|
||||
CGAL_assertion(false);
|
||||
return 0.;
|
||||
#endif //CGAL_MESH_3_EXPERIMENTAL_USE_PATCHES_IDS
|
||||
}
|
||||
|
||||
#ifndef CGAL_MESH_3_IMAGE_EXAMPLE
|
||||
else if (dim == 1)
|
||||
{
|
||||
#ifdef CGAL_MESH_3_EXPERIMENTAL_USE_PATCHES_IDS
|
||||
const typename MeshDomain::Curve_segment_index& curve_id =
|
||||
m_domain.curve_segment_index(index);
|
||||
const Patches_ids& ids = patches_ids_map[curve_id];
|
||||
|
||||
if (m_domain_is_a_box && ids.size() == 2)
|
||||
{
|
||||
//we are on an edge of the box
|
||||
//same code as when dim == 2
|
||||
Surface_patch_index spi = *(ids.begin());
|
||||
const std::pair<Subdomain_index, Subdomain_index>& subdomains
|
||||
= m_params.incident_subdomains(spi);
|
||||
if (!is_in_domain(subdomains.first))
|
||||
return size_in_subdomain(p, subdomains.second);
|
||||
else
|
||||
return size_in_subdomain(p, subdomains.first);
|
||||
}
|
||||
return min_size_in_incident_subdomains(ids);
|
||||
#else
|
||||
CGAL_assertion(false);//should not be used for dimension 1
|
||||
return 0.;
|
||||
#endif
|
||||
}
|
||||
else if (dim == 0)
|
||||
{
|
||||
#ifdef CGAL_MESH_3_EXPERIMENTAL_USE_PATCHES_IDS
|
||||
const Patches_ids& ids =
|
||||
(m_domain.corners_incidences_map().find(p)->second);
|
||||
|
||||
if (m_domain_is_a_box && ids.size() == 3)
|
||||
{
|
||||
//we are on a corner of the box
|
||||
//same code as when dim == 2
|
||||
Surface_patch_index spi = *(ids.begin());
|
||||
const std::pair<Subdomain_index, Subdomain_index>& subdomains
|
||||
= m_params.incident_subdomains(spi);
|
||||
if (!is_in_domain(subdomains.first))
|
||||
return size_in_subdomain(p, subdomains.second);
|
||||
else
|
||||
return size_in_subdomain(p, subdomains.first);
|
||||
}
|
||||
|
||||
return min_size_in_incident_subdomains(ids);;
|
||||
#else
|
||||
CGAL_assertion(false);//should not be used for dimension 0
|
||||
return 0;
|
||||
#endif
|
||||
}
|
||||
#endif
|
||||
|
||||
CGAL_assertion(false);
|
||||
return 0.;
|
||||
}
|
||||
|
||||
private:
|
||||
#ifdef CGAL_MESH_3_EXPERIMENTAL_USE_PATCHES_IDS
|
||||
std::vector<Subdomain_index> incident_subdomains(const Patches_ids& ids) const
|
||||
{
|
||||
std::vector<Subdomain_index> vec;
|
||||
BOOST_FOREACH(Surface_patch_index spi, ids)
|
||||
{
|
||||
const std::pair<Subdomain_index, Subdomain_index>& subdomains
|
||||
= m_params.incident_subdomains(spi);
|
||||
|
||||
if (is_in_domain(subdomains.first))
|
||||
vec.push_back(subdomains.first);
|
||||
|
||||
if (is_in_domain(subdomains.second))
|
||||
vec.push_back(subdomains.second);
|
||||
}
|
||||
return vec;
|
||||
}
|
||||
|
||||
FT min_size_in_incident_subdomains(const Patches_ids& ids) const
|
||||
{
|
||||
FT size = static_cast<FT>((std::numeric_limits<double>::max)());
|
||||
BOOST_FOREACH(Surface_patch_index spi, ids)
|
||||
{
|
||||
const std::pair<Subdomain_index, Subdomain_index>& subdomains
|
||||
= m_params.incident_subdomains(spi);
|
||||
|
||||
FT k, size_min, size_max;
|
||||
if (is_in_domain(subdomains.first))
|
||||
{
|
||||
m_params.get_parameters(subdomains.first, k, size_min, size_max);
|
||||
size = (std::min)(size, size_min);
|
||||
}
|
||||
if (is_in_domain(subdomains.second))
|
||||
{
|
||||
m_params.get_parameters(subdomains.second, k, size_min, size_max);
|
||||
size = (std::min)(size, size_min);
|
||||
}
|
||||
}
|
||||
return size;
|
||||
}
|
||||
#endif
|
||||
|
||||
public:
|
||||
template <typename C3T3>
|
||||
void init_aabb_tree_from_c3t3(const C3T3* p_c3t3)
|
||||
{
|
||||
static std::list<Triangle> triangles;
|
||||
for (typename C3T3::Facets_in_complex_iterator
|
||||
fit = p_c3t3->facets_in_complex_begin();
|
||||
fit != p_c3t3->facets_in_complex_end();
|
||||
++fit)
|
||||
{
|
||||
if (!is_on_cube_boundary(*fit))
|
||||
triangles.push_back(p_c3t3->triangulation().triangle(*fit));
|
||||
}
|
||||
|
||||
m_own_ptree.reset(new Tree(triangles.begin(), triangles.end()));
|
||||
m_own_ptree->build();
|
||||
m_own_ptree->accelerate_distance_queries();
|
||||
}
|
||||
|
||||
private:
|
||||
template <typename Facet>
|
||||
bool is_on_cube_boundary(const Facet& f) const
|
||||
{
|
||||
return is_on_cube_boundary(f.first->surface_patch_index(f.second));
|
||||
}
|
||||
|
||||
bool is_on_cube_boundary(const Surface_patch_index& sp_index) const
|
||||
{
|
||||
const std::pair<Subdomain_index, Subdomain_index>& index
|
||||
= m_params.incident_subdomains(sp_index);
|
||||
|
||||
#ifdef CGAL_MESH_3_IMAGE_EXAMPLE
|
||||
return (index.first == INT_MIN || index.second == INT_MIN);
|
||||
#else //POLYHEDRAL EXAMPLE
|
||||
|
||||
if (m_domain_is_a_box)
|
||||
return !is_in_domain(index.first) || !is_in_domain(index.second);
|
||||
else
|
||||
return false;
|
||||
|
||||
#endif //CGAL_MESH_3_IMAGE_EXAMPLE
|
||||
}
|
||||
|
||||
bool is_on_cube_boundary(const Point_3& p) const
|
||||
{
|
||||
#ifndef CGAL_MESH_3_IMAGE_EXAMPLE//POLYHEDRAL EXAMPLE
|
||||
|
||||
if (m_domain_is_a_box)
|
||||
//checks that p is in the outer 'shell' of voxels
|
||||
return p.x() < m_bbox.xmin() + m_vxyz[0]
|
||||
|| p.x() > m_bbox.xmax() - m_vxyz[0]
|
||||
|| p.y() < m_bbox.ymin() + m_vxyz[1]
|
||||
|| p.y() > m_bbox.ymax() - m_vxyz[1]
|
||||
|| p.z() < m_bbox.zmin() + m_vxyz[2]
|
||||
|| p.z() > m_bbox.zmax() - m_vxyz[2];
|
||||
else
|
||||
return false;
|
||||
|
||||
#else //IMAGE EXAMPLE
|
||||
CGAL_USE(p);
|
||||
return false;
|
||||
#endif //IMAGE_EXAMPLE
|
||||
}
|
||||
|
||||
bool is_in_domain(const Subdomain_index& index) const
|
||||
{
|
||||
#ifdef CGAL_MESH_3_IMAGE_EXAMPLE
|
||||
return (index != 0 && index != INT_MIN);
|
||||
#else //POLYHEDRAL EXAMPLE
|
||||
return (index != 0);
|
||||
#endif
|
||||
}
|
||||
|
||||
FT size_in_subdomain(const Point_3& p, const Subdomain_index& index) const
|
||||
{
|
||||
FT k, size_min, size_max;
|
||||
m_params.get_parameters(index, k, size_min, size_max);
|
||||
|
||||
FT sqdist = 0.;
|
||||
if(m_ptree == NULL)
|
||||
{
|
||||
sqdist = m_own_ptree->squared_distance(p);
|
||||
}
|
||||
else
|
||||
{
|
||||
Point_3 closest = compute_closest_point(p);
|
||||
sqdist = CGAL::squared_distance(p, closest);
|
||||
}
|
||||
|
||||
FT size = k * CGAL::sqrt(sqdist) + size_min;
|
||||
return (std::min)(size, size_max);
|
||||
}
|
||||
|
||||
#ifdef CGAL_MESH_3_EXPERIMENTAL_USE_PATCHES_IDS
|
||||
void kd_tree()
|
||||
{
|
||||
typedef typename MeshDomain::Polyhedron Polyhedron;
|
||||
if(m_kd_tree.get() == 0) {
|
||||
m_kd_tree.reset(new Kd_tree);
|
||||
BOOST_FOREACH(std::size_t poly_id, m_domain.inside_polyhedra()) {
|
||||
const Polyhedron& poly = m_domain.polyhedra()[poly_id];
|
||||
BOOST_FOREACH(typename Polyhedron::Vertex_handle v, vertices(poly))
|
||||
{
|
||||
m_kd_tree->insert(v->point());
|
||||
}
|
||||
}
|
||||
BOOST_FOREACH(std::size_t poly_id, m_domain.boundary_polyhedra()) {
|
||||
const Polyhedron& poly = m_domain.polyhedra()[poly_id];
|
||||
BOOST_FOREACH(typename Polyhedron::Vertex_handle v, vertices(poly))
|
||||
{
|
||||
if(!is_on_cube_boundary(v->point()))
|
||||
m_kd_tree->insert(v->point());
|
||||
}
|
||||
}
|
||||
m_kd_tree->build();
|
||||
}
|
||||
}
|
||||
#endif
|
||||
|
||||
Point_3 compute_closest_point(const Point_3& p) const
|
||||
{
|
||||
#ifndef CGAL_MESH_3_IMAGE_EXAMPLE //POLYHEDRAL_EXAMPLE
|
||||
|
||||
#ifdef CGAL_MESH_3_EXPERIMENTAL_USE_PATCHES_IDS
|
||||
const std::vector<Surface_patch_index>& boundary_ids =
|
||||
m_domain.boundary_patches();
|
||||
|
||||
CGAL_STATIC_THREAD_LOCAL_VARIABLE_4(AABB_filtered_traits,
|
||||
projection_traits,
|
||||
boundary_ids.begin(),
|
||||
boundary_ids.end(),
|
||||
m_ptree->traits(),
|
||||
m_get_facet_patch_id);
|
||||
kd_tree();//build it if needed
|
||||
Neighbor_search search(*m_kd_tree, p, 1);
|
||||
projection_traits.reset(search.begin()->first);
|
||||
|
||||
m_ptree->traversal(p, projection_traits);
|
||||
CGAL_assertion(projection_traits.found());
|
||||
return projection_traits.closest_point();
|
||||
|
||||
#else //CGAL_MESH_3_EXPERIMENTAL_USE_PATCHES_IDS
|
||||
return m_ptree->closest_point(p);
|
||||
#endif //CGAL_MESH_3_EXPERIMENTAL_USE_PATCHES_IDS
|
||||
|
||||
#else
|
||||
CGAL_assertion(false);
|
||||
return CGAL::ORIGIN;//not used
|
||||
#endif
|
||||
}
|
||||
|
||||
public:
|
||||
void add_parameters_for_subdomain(const Subdomain_index& id
|
||||
, const FT& k
|
||||
, const FT& size_min
|
||||
, const FT& size_max)
|
||||
{
|
||||
m_params.add_subdomain(id, k, size_min, size_max);
|
||||
}
|
||||
|
||||
};
|
||||
|
||||
}//namespace Mesh_3
|
||||
}//namespace CGAL
|
||||
|
||||
#endif // CGAL_LIPSCHITZ_SIZING_H
|
||||
|
|
@ -0,0 +1,162 @@
|
|||
// Copyright (c) 2016 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) : Jane Tournois
|
||||
//
|
||||
|
||||
#ifndef CGAL_LIP_SIZING_PARAMETERS_H
|
||||
#define CGAL_LIP_SIZING_PARAMETERS_H
|
||||
|
||||
#include <CGAL/license/Mesh_3.h>
|
||||
|
||||
#include <map>
|
||||
#include <limits>
|
||||
|
||||
namespace CGAL
|
||||
{
|
||||
template <typename MeshDomain, typename FT>
|
||||
class Lipschitz_sizing_parameters
|
||||
{
|
||||
struct SubdomainParam
|
||||
{
|
||||
FT m_k;
|
||||
FT m_size_min;//min size in subdomain
|
||||
FT m_size_max;//max size in subdomain
|
||||
|
||||
public:
|
||||
SubdomainParam()
|
||||
: m_k(0.)
|
||||
, m_size_min(0.)
|
||||
, m_size_max(0.)
|
||||
{}
|
||||
SubdomainParam(const FT& k
|
||||
, const FT& size_min
|
||||
, const FT& size_max)
|
||||
: m_k(k)
|
||||
, m_size_min(size_min)
|
||||
, m_size_max(size_max)
|
||||
{}
|
||||
|
||||
const FT& k() const { return m_k; }
|
||||
const FT& size_min() const { return m_size_min; }
|
||||
const FT& size_max() const { return m_size_max; }
|
||||
};
|
||||
|
||||
private:
|
||||
typedef typename MeshDomain::Subdomain_index Subdomain_index;
|
||||
typedef typename MeshDomain::Surface_patch_index Surface_patch_index;
|
||||
|
||||
typedef std::map<Subdomain_index, SubdomainParam> Parameters_map;
|
||||
Parameters_map m_parameters;
|
||||
SubdomainParam m_default_params;
|
||||
const MeshDomain* p_domain_;
|
||||
|
||||
public:
|
||||
Lipschitz_sizing_parameters(const MeshDomain& domain)
|
||||
: p_domain_(&domain)
|
||||
{}
|
||||
|
||||
bool empty() const { return m_parameters.empty(); }
|
||||
std::size_t size() const { return m_parameters.size(); }
|
||||
|
||||
void add_subdomain(const Subdomain_index& index
|
||||
, const FT& k
|
||||
, const FT& size_min
|
||||
, const FT& size_max)
|
||||
{
|
||||
typename Parameters_map::iterator it = m_parameters.find(index);
|
||||
if (it != m_parameters.end())
|
||||
std::cout << "Warning : sizing parameters for subdomain " << index
|
||||
<< "will be overwritten." << std::endl;
|
||||
|
||||
if (index == -1)
|
||||
m_default_params
|
||||
= SubdomainParam(k, size_min, size_max);
|
||||
else
|
||||
m_parameters[index]
|
||||
= SubdomainParam(k, size_min, size_max);
|
||||
}
|
||||
|
||||
void get_parameters(const Subdomain_index& index
|
||||
, FT& k
|
||||
, FT& size_min
|
||||
, FT& size_max) const
|
||||
{
|
||||
typename Parameters_map::const_iterator it
|
||||
= m_parameters.find(index);
|
||||
|
||||
SubdomainParam p = (it != m_parameters.end())
|
||||
? (*it).second
|
||||
: m_default_params;
|
||||
k = p.k();
|
||||
size_min = p.size_min();
|
||||
size_max = p.size_max();
|
||||
}
|
||||
|
||||
#ifdef CGAL_MESH_3_LIPSCHITZ_SIZING_EXPERIMENTAL
|
||||
void get_parameters(const Surface_patch_index& spi
|
||||
, FT& size_max) const
|
||||
{
|
||||
const std::pair<Subdomain_index, Subdomain_index>& index
|
||||
= incident_subdomains(spi);
|
||||
|
||||
typename Parameters_map::const_iterator it1
|
||||
= m_parameters.find(index.first);
|
||||
typename Parameters_map::const_iterator it2
|
||||
= m_parameters.find(index.second);
|
||||
|
||||
SubdomainParam p1
|
||||
= (it1 != m_parameters.end()) ? (*it1).second : m_default_params;
|
||||
SubdomainParam p2
|
||||
= (it2 != m_parameters.end()) ? (*it2).second : m_default_params;
|
||||
|
||||
bool boundary1 = (index.first == 0); //boundary of the domain, inside the cube
|
||||
bool boundary2 = (index.second == 0); //boundary of the domain, inside the cube
|
||||
|
||||
#ifdef CGAL_MESH_3_IMAGE_EXAMPLE
|
||||
boundary1 = boundary1 || (index.first == INT_MIN); //boundary of the cube
|
||||
boundary2 = boundary2 || (index.second == INT_MIN);//boundary of the cube
|
||||
#endif
|
||||
CGAL_assertion(!boundary1 || !boundary2);
|
||||
|
||||
if (!boundary1)
|
||||
{
|
||||
if (!boundary2)
|
||||
size_max = (std::min)(p1.size_min(), p2.size_min());
|
||||
else
|
||||
size_max = p1.size_min();
|
||||
}
|
||||
else
|
||||
size_max = p2.size_min();
|
||||
}
|
||||
|
||||
const std::pair<Subdomain_index, Subdomain_index> &
|
||||
incident_subdomains(const Surface_patch_index& index) const
|
||||
{
|
||||
//#ifdef CGAL_MESH_3_IMAGE_EXAMPLE
|
||||
// return index;
|
||||
//#else //POLYHEDRAL_EXAMPLE
|
||||
return p_domain_->incident_subdomains_indices(index);
|
||||
//#endif
|
||||
}
|
||||
#endif
|
||||
|
||||
};
|
||||
}//namespace CGAL
|
||||
|
||||
#endif //CGAL_LIP_SIZING_PARAMETERS_H
|
||||
|
|
@ -0,0 +1,186 @@
|
|||
// Copyright (c) 2016 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) : Jane Tournois
|
||||
//
|
||||
|
||||
#ifndef CGAL_LIPSCHITZ_SIZING_POLYHEDRON_H
|
||||
#define CGAL_LIPSCHITZ_SIZING_POLYHEDRON_H
|
||||
|
||||
#include <CGAL/license/Mesh_3.h>
|
||||
|
||||
#include <CGAL/AABB_tree.h>
|
||||
#include <CGAL/AABB_traits.h>
|
||||
#include <CGAL/AABB_triangle_primitive.h>
|
||||
|
||||
#include <CGAL/Mesh_3/experimental/Lipschitz_sizing_parameters.h>
|
||||
|
||||
#include <CGAL/Default.h>
|
||||
#include <CGAL/array.h>
|
||||
#include <CGAL/Bbox_3.h>
|
||||
|
||||
#include <boost/shared_ptr.hpp>
|
||||
|
||||
#include <list>
|
||||
#include <sstream>
|
||||
#include <string>
|
||||
#include <fstream>
|
||||
#include <vector>
|
||||
|
||||
namespace CGAL
|
||||
{
|
||||
namespace Mesh_3
|
||||
{
|
||||
|
||||
template <class Kernel, class MeshDomain
|
||||
, typename AABBTreeTemplate = CGAL::Default
|
||||
>
|
||||
class Lipschitz_sizing
|
||||
{
|
||||
public:
|
||||
typedef Kernel K;
|
||||
typedef typename Kernel::FT FT;
|
||||
typedef typename Kernel::Triangle_3 Triangle;
|
||||
typedef typename Kernel::Point_3 Point_3;
|
||||
|
||||
typedef typename std::list<Triangle>::iterator Tr_iterator;
|
||||
typedef CGAL::AABB_triangle_primitive<K, Tr_iterator> Primitive;
|
||||
typedef CGAL::AABB_traits<K, Primitive> AABB_tr_traits;
|
||||
typedef CGAL::AABB_tree<AABB_tr_traits> AABB_tree;
|
||||
|
||||
typedef typename CGAL::Default::Get<AABBTreeTemplate, AABB_tree>::type Tree;
|
||||
|
||||
typedef typename MeshDomain::Index Index;
|
||||
typedef typename MeshDomain::Subdomain_index Subdomain_index;
|
||||
typedef typename MeshDomain::Surface_patch_index Surface_patch_index;
|
||||
|
||||
typedef CGAL::Lipschitz_sizing_parameters<MeshDomain, FT> Parameters;
|
||||
|
||||
|
||||
private:
|
||||
const Tree* m_ptree;
|
||||
boost::shared_ptr<Tree> m_own_ptree;
|
||||
|
||||
const MeshDomain& m_domain;
|
||||
Parameters m_params;
|
||||
|
||||
public:
|
||||
Lipschitz_sizing(const MeshDomain& domain)
|
||||
: m_ptree(NULL)
|
||||
, m_own_ptree()
|
||||
, m_domain(domain)
|
||||
, m_params(domain)
|
||||
{
|
||||
}
|
||||
|
||||
Lipschitz_sizing(const MeshDomain& domain
|
||||
, const Tree* ptree
|
||||
)
|
||||
: m_ptree(ptree)
|
||||
, m_own_ptree()
|
||||
, m_domain(domain)
|
||||
, m_params(domain)
|
||||
{
|
||||
}
|
||||
|
||||
FT operator()(const Point_3& p, const int dim, const Index& index) const
|
||||
{
|
||||
CGAL_assertion(!m_params.empty());
|
||||
if (dim == 3)
|
||||
{
|
||||
return size_in_subdomain(p, m_domain.subdomain_index(index));
|
||||
}
|
||||
else if (dim == 2)
|
||||
{
|
||||
CGAL_assertion(false);//should not be used for dimension 2
|
||||
return 0.;
|
||||
}
|
||||
else if (dim == 1)
|
||||
{
|
||||
CGAL_assertion(false);//should not be used for dimension 1
|
||||
return 0.;
|
||||
}
|
||||
else if (dim == 0)
|
||||
{
|
||||
CGAL_assertion(false);//should not be used for dimension 0
|
||||
return 0;
|
||||
}
|
||||
CGAL_assertion(false);
|
||||
return 0.;
|
||||
}
|
||||
|
||||
public:
|
||||
template <typename C3T3>
|
||||
void init_aabb_tree_from_c3t3(const C3T3* p_c3t3)
|
||||
{
|
||||
static std::list<Triangle> triangles;
|
||||
for (typename C3T3::Facets_in_complex_iterator
|
||||
fit = p_c3t3->facets_in_complex_begin();
|
||||
fit != p_c3t3->facets_in_complex_end();
|
||||
++fit)
|
||||
{
|
||||
triangles.push_back(p_c3t3->triangulation().triangle(*fit));
|
||||
}
|
||||
|
||||
m_own_ptree.reset(new Tree(triangles.begin(), triangles.end()));
|
||||
m_own_ptree->build();
|
||||
m_own_ptree->accelerate_distance_queries();
|
||||
}
|
||||
|
||||
private:
|
||||
FT size_in_subdomain(const Point_3& p, const Subdomain_index& index) const
|
||||
{
|
||||
FT k, size_min, size_max;
|
||||
m_params.get_parameters(index, k, size_min, size_max);
|
||||
|
||||
FT sqdist = 0.;
|
||||
if(m_ptree == NULL)
|
||||
{
|
||||
sqdist = m_own_ptree->squared_distance(p);
|
||||
}
|
||||
else
|
||||
{
|
||||
Point_3 closest = compute_closest_point(p);
|
||||
sqdist = CGAL::squared_distance(p, closest);
|
||||
}
|
||||
|
||||
FT size = k * CGAL::sqrt(sqdist) + size_min;
|
||||
return (std::min)(size, size_max);
|
||||
}
|
||||
|
||||
private:
|
||||
Point_3 compute_closest_point(const Point_3& p) const
|
||||
{
|
||||
return m_ptree->closest_point(p);
|
||||
}
|
||||
|
||||
public:
|
||||
void add_parameters_for_subdomain(const Subdomain_index& id
|
||||
, const FT& k
|
||||
, const FT& size_min
|
||||
, const FT& size_max)
|
||||
{
|
||||
m_params.add_subdomain(id, k, size_min, size_max);
|
||||
}
|
||||
|
||||
};
|
||||
|
||||
}//namespace Mesh_3
|
||||
}//namespace CGAL
|
||||
|
||||
#endif // CGAL_LIPSCHITZ_SIZING_POLYHEDRON_H
|
||||
Loading…
Reference in New Issue