CGAL 5.4 - 2D Generalized Barycentric Coordinates
CGAL::Barycentric_coordinates::Harmonic_coordinates_2< VertexRange, DiscretizedDomain, GeomTraits, PointMap > Class Template Reference

#include <CGAL/Barycentric_coordinates_2/Harmonic_coordinates_2.h>

Definition

template<typename VertexRange, typename DiscretizedDomain, typename GeomTraits, typename PointMap = CGAL::Identity_property_map<typename GeomTraits::Point_2>>
class CGAL::Barycentric_coordinates::Harmonic_coordinates_2< VertexRange, DiscretizedDomain, GeomTraits, PointMap >

2D harmonic coordinates.

This class implements 2D harmonic coordinate functions ( [2], [7] ), which can be evaluated at any point inside a simple polygon.

Harmonic coordinates are well-defined and non-negative in the closure of any simple polygon, however they cannot be computed analytically and hence they are approximated. The classical way to approximate these coordinates is by discretizing over the space of piecewise linear functions with respect to a partition of the polygon's interior domain, e.g. a triangulation.

Once computed at the vertices of the discretized domain, the coordinate functions can be evaluated at any point inside a polygon by locating the finite element that contains a query point and linearly interpolating within this element. See more details in the user manual here.

Template Parameters
VertexRangea model of ConstRange whose iterator type is RandomAccessIterator
DiscretizedDomaina model of DiscretizedDomain_2. For the moment, we only support domains whose partition's finite elements are triangles.
GeomTraitsa model of BarycentricTraits_2
PointMapa model of ReadablePropertyMap whose key type is VertexRange::value_type and value type is Point_2. The default is CGAL::Identity_property_map.
Examples:
Barycentric_coordinates_2/harmonic_coordinates.cpp, and Barycentric_coordinates_2/shape_deformation.cpp.

Types

typedef GeomTraits::FT FT
 Number type.
 
typedef GeomTraits::Point_2 Point_2
 Point type.
 

Initialization

 Harmonic_coordinates_2 (const VertexRange &polygon, const DiscretizedDomain &domain, const GeomTraits traits=GeomTraits(), const PointMap point_map=PointMap())
 initializes all internal data structures. More...
 

Access

template<typename OutIterator >
OutIterator operator() (const Point_2 &query, OutIterator c_begin)
 evaluates 2D harmonic coordinates. More...
 
template<typename OutIterator >
OutIterator operator() (const std::size_t query_index, OutIterator c_begin)
 returns 2D harmonic coordinates at one domain vertex. More...
 
template<typename OutIterator >
OutIterator operator() (OutIterator c_begin)
 returns 2D harmonic coordinates at all domain vertices. More...
 

Computation

void compute ()
 computes 2D harmonic coordinates at the vertices of the input domain.
 

Step by Step Computation

void setup ()
 computes all necessary harmonic data. More...
 
void factorize ()
 factorizes the matrix A. More...
 
void solve ()
 solves the linear system Ax = b. More...
 

Memory Management

void clear ()
 clears all internal data structures.
 
void release_memory ()
 releases all memory that is used internally.
 

Constructor & Destructor Documentation

◆ Harmonic_coordinates_2()

template<typename VertexRange , typename DiscretizedDomain , typename GeomTraits , typename PointMap = CGAL::Identity_property_map<typename GeomTraits::Point_2>>
CGAL::Barycentric_coordinates::Harmonic_coordinates_2< VertexRange, DiscretizedDomain, GeomTraits, PointMap >::Harmonic_coordinates_2 ( const VertexRange &  polygon,
const DiscretizedDomain &  domain,
const GeomTraits  traits = GeomTraits(),
const PointMap  point_map = PointMap() 
)

initializes all internal data structures.

This class implements the behavior of harmonic coordinates for 2D query points.

Parameters
polygonan instance of VertexRange with the vertices of a simple polygon
domainan instance of DiscretizedDomain with a partition of the interior part of a simple polygon
traitsa traits class with geometric objects, predicates, and constructions; the default initialization is provided
point_mapan instance of PointMap that maps a vertex from polygon to Point_2; the default initialization is provided
Precondition
polygon.size() >= 3
polygon is simple

Member Function Documentation

◆ factorize()

template<typename VertexRange , typename DiscretizedDomain , typename GeomTraits , typename PointMap = CGAL::Identity_property_map<typename GeomTraits::Point_2>>
void CGAL::Barycentric_coordinates::Harmonic_coordinates_2< VertexRange, DiscretizedDomain, GeomTraits, PointMap >::factorize ( )

factorizes the matrix A.

Precondition
setup() is called

◆ operator()() [1/3]

template<typename VertexRange , typename DiscretizedDomain , typename GeomTraits , typename PointMap = CGAL::Identity_property_map<typename GeomTraits::Point_2>>
template<typename OutIterator >
OutIterator CGAL::Barycentric_coordinates::Harmonic_coordinates_2< VertexRange, DiscretizedDomain, GeomTraits, PointMap >::operator() ( const Point_2 query,
OutIterator  c_begin 
)

evaluates 2D harmonic coordinates.

This function fills c_begin with harmonic coordinates evaluated at the query point with respect to the vertices of the input polygon. Evaluation is performed by locating the finite element in the input domain that contains query and then linearly interpolating harmonic coordinates within this element.

If query does not belong to the input domain or the located element has more than 3 vertices, all coordinates are set to zero.

The number of returned coordinates equals to the number of polygon vertices.

After the coordinates \(b_i\) with \(i = 1\dots n\) are computed, where \(n\) is the number of polygon vertices, the query point \(q\) can be obtained as \(q = \sum_{i = 1}^{n}b_ip_i\), where \(p_i\) are the polygon vertices.

Template Parameters
OutIteratora model of OutputIterator that accepts values of type FT
Parameters
querya query point
c_beginthe beginning of the destination range with the computed coordinates
Returns
an output iterator to the element in the destination range, one past the last coordinate stored

◆ operator()() [2/3]

template<typename VertexRange , typename DiscretizedDomain , typename GeomTraits , typename PointMap = CGAL::Identity_property_map<typename GeomTraits::Point_2>>
template<typename OutIterator >
OutIterator CGAL::Barycentric_coordinates::Harmonic_coordinates_2< VertexRange, DiscretizedDomain, GeomTraits, PointMap >::operator() ( const std::size_t  query_index,
OutIterator  c_begin 
)

returns 2D harmonic coordinates at one domain vertex.

This function fills c_begin with harmonic coordinates computed at the vertex of the input domain with the index query_index.

The number of returned coordinates equals to the number of polygon vertices.

After the coordinates \(b_i\) with \(i = 1\dots n\) are computed, where \(n\) is the number of polygon vertices, the partition vertex \(q\) with the index query_index can be obtained as \(q = \sum_{i = 1}^{n}b_ip_i\), where \(p_i\) are the polygon vertices.

Template Parameters
OutIteratora model of OutputIterator that accepts values of type FT
Parameters
query_indexa domain's vertex index
c_beginthe beginning of the destination range with the computed coordinates
Returns
an output iterator to the element in the destination range, one past the last coordinate stored
Precondition
query_index >= 0 && query_index < domain.number_of_vertices()

◆ operator()() [3/3]

template<typename VertexRange , typename DiscretizedDomain , typename GeomTraits , typename PointMap = CGAL::Identity_property_map<typename GeomTraits::Point_2>>
template<typename OutIterator >
OutIterator CGAL::Barycentric_coordinates::Harmonic_coordinates_2< VertexRange, DiscretizedDomain, GeomTraits, PointMap >::operator() ( OutIterator  c_begin)

returns 2D harmonic coordinates at all domain vertices.

This function fills c_begin with harmonic coordinates computed at the vertices of the input domain.

The number of returned coordinates equals to the number of input domain vertices.

Template Parameters
OutIteratora model of OutputIterator that accepts values of type std::vector<FT>
Parameters
c_beginthe beginning of the destination range with the computed coordinates
Returns
an output iterator to the element in the destination range, one past the last coordinate set stored

◆ setup()

template<typename VertexRange , typename DiscretizedDomain , typename GeomTraits , typename PointMap = CGAL::Identity_property_map<typename GeomTraits::Point_2>>
void CGAL::Barycentric_coordinates::Harmonic_coordinates_2< VertexRange, DiscretizedDomain, GeomTraits, PointMap >::setup ( )

computes all necessary harmonic data.

This function fills in the left side matrix A and the right side vector b of the linear system Ax = b. The matrix A is a sparse symmetric positive definite matrix. The solution vector x of this system gives the harmonic coordinates at the vertices of the input domain.

◆ solve()

template<typename VertexRange , typename DiscretizedDomain , typename GeomTraits , typename PointMap = CGAL::Identity_property_map<typename GeomTraits::Point_2>>
void CGAL::Barycentric_coordinates::Harmonic_coordinates_2< VertexRange, DiscretizedDomain, GeomTraits, PointMap >::solve ( )

solves the linear system Ax = b.

Precondition
factorize() is called