1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
|
/* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT.
* See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details.
* Author(s): Clément Maria
*
* Copyright (C) 2014 Inria
*
* Modification(s):
* - YYYY/MM Author: Description of the modification
*/
#ifndef DISTANCE_FUNCTIONS_H_
#define DISTANCE_FUNCTIONS_H_
#include <gudhi/Debug_utils.h>
#include <gudhi/Miniball.hpp>
#include <boost/range/metafunctions.hpp>
#include <boost/range/size.hpp>
#include <CGAL/Epeck_d.h>
#include <cmath> // for std::sqrt
#include <type_traits> // for std::decay
#include <iterator> // for std::begin, std::end
#include <utility>
namespace Gudhi {
/** @file
* @brief Global distance functions
*/
/** @brief Compute the Euclidean distance between two Points given by a range of coordinates. The points are assumed to
* have the same dimension. */
class Euclidean_distance {
public:
// boost::range_value is not SFINAE-friendly so we cannot use it in the return type
template< typename Point >
typename std::iterator_traits<typename boost::range_iterator<Point>::type>::value_type
operator()(const Point& p1, const Point& p2) const {
auto it1 = std::begin(p1);
auto it2 = std::begin(p2);
typedef typename boost::range_value<Point>::type NT;
NT dist = 0;
for (; it1 != std::end(p1); ++it1, ++it2) {
GUDHI_CHECK(it2 != std::end(p2), "inconsistent point dimensions");
NT tmp = *it1 - *it2;
dist += tmp*tmp;
}
GUDHI_CHECK(it2 == std::end(p2), "inconsistent point dimensions");
using std::sqrt;
return sqrt(dist);
}
template< typename T >
T operator() (const std::pair< T, T >& f, const std::pair< T, T >& s) const {
T dx = f.first - s.first;
T dy = f.second - s.second;
using std::sqrt;
return sqrt(dx*dx + dy*dy);
}
};
/** @brief Compute the radius of the minimal enclosing ball between Points given by a range of coordinates.
* The points are assumed to have the same dimension. */
class Minimal_enclosing_ball_radius {
public:
/** \brief Enclosing ball radius from two points using CGAL.
*
* @param[in] point_1
* @param[in] point_2
* @return Enclosing ball radius for the two points.
* \tparam Point must be a Kernel::Point_d from CGAL.
*
*/
template< typename Kernel = CGAL::Epeck_d<CGAL::Dynamic_dimension_tag>,
typename Point= typename Kernel::Point_d>
double operator()(const Point& point_1, const Point& point_2) const {
Kernel kernel_;
return std::sqrt(CGAL::to_double(kernel_.squared_distance_d_object()(point_1, point_2))) / 2.;
}
/** \brief Minimal_enclosing_ball_radius from two points.
*
* @param[in] point_1 First point.
* @param[in] point_2 second point.
* @return The minimal enclosing ball radius for the two points (aka. Euclidean distance / 2.).
*
* \tparam Point must be a range of Cartesian coordinates.
*
*/
template< typename Point >
typename std::iterator_traits<typename boost::range_iterator<Point>::type>::value_type
operator()(const Point& point_1, const Point& point_2) const {
std::clog << "#" << *point_1.begin() << "##" << *point_2.begin() << std::endl;
return Euclidean_distance()(point_1, point_2) / 2.;
}
/** \brief Enclosing ball radius from a point cloud using CGAL.
*
* @param[in] point_cloud The points.
* @return Enclosing ball radius for the points.
* \tparam Point_cloud must be a range of Kernel::Point_d points from CGAL.
*
*/
template< typename Kernel = CGAL::Epeck_d<CGAL::Dynamic_dimension_tag>,
typename Point= typename Kernel::Point_d,
typename Point_cloud = std::vector<Point>>
double operator()(const Point_cloud& point_cloud) const {
Kernel kernel_;
return std::sqrt(CGAL::to_double(kernel_.compute_squared_radius_d_object()(point_cloud.begin(), point_cloud.end())));
}
/** \brief Minimal_enclosing_ball_radius from a point cloud.
*
* @param[in] point_cloud The points.
* @return The minimal enclosing ball radius for the points.
*
* \tparam Point_cloud must be a range of points with Cartesian coordinates.
* Point_cloud is a range over a range of Coordinate.
*
*/
template< typename Point_cloud,
typename Point_iterator = typename boost::range_const_iterator<Point_cloud>::type,
typename Point = typename std::iterator_traits<Point_iterator>::value_type,
typename Coordinate_iterator = typename boost::range_const_iterator<Point>::type,
typename Coordinate = typename std::iterator_traits<Coordinate_iterator>::value_type>
Coordinate
operator()(const Point_cloud& point_cloud) const {
using Min_sphere = Miniball::Miniball<Miniball::CoordAccessor<Point_iterator, Coordinate_iterator>>;
Min_sphere ms(boost::size(*point_cloud.begin()), point_cloud.begin(), point_cloud.end());
#ifdef DEBUG_TRACES
std::clog << "Minimal_enclosing_ball_radius = " << std::sqrt(ms.squared_radius()) << " | nb points = "
<< boost::size(point_cloud) << " | dimension = "
<< boost::size(*point_cloud.begin()) << std::endl;
#endif // DEBUG_TRACES
return std::sqrt(ms.squared_radius());
}
};
} // namespace Gudhi
#endif // DISTANCE_FUNCTIONS_H_
|
