songsenand
/
cgal
mirror of https://github.com/CGAL/cgal
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

Merge pull request #6378 from GilesBathgate/Nef_3-performance_sphere_circle-GilesBathgate 

Nef_3: Sphere circle constructors.
This commit is contained in:
Sébastien Loriot 2022-03-29 17:45:52 +02:00
parent b34bc3f8f7 1c4b191c71
commit 99f0598f3e
9 changed files with 228 additions and 164 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

21
Cartesian_kernel/include/CGAL/Cartesian/function_objects.h
View File

@ -3031,6 +3031,18 @@ namespace CartesianKernelFunctors {
return construct_vector(vx, vy, vz);
}
Vector_3
operator()( Origin, const Point_3& q, const Point_3& r ) const
{
// Cross product oq * or
FT vx = q.y()*r.z() - r.y()*q.z();
FT vy = q.z()*r.x() - r.z()*q.x();
FT vz = q.x()*r.y() - r.x()*q.y();
typename K::Construct_vector_3 construct_vector;
return construct_vector(vx, vy, vz);
}
};
template <typename K>
@ -4307,6 +4319,15 @@ namespace CartesianKernelFunctors {
w.x(), w.y(), w.z());
}
result_type
operator()( Origin, const Point_3& u,
const Point_3& v, const Point_3& w) const
{
return orientationC3(u.x(), u.y(), u.z(),
v.x(), v.y(), v.z(),
w.x(), w.y(), w.z());
}
result_type
operator()( const Tetrahedron_3& t) const
{

8
Filtered_kernel/include/CGAL/Epic_converter.h
View File

@ -49,7 +49,15 @@ class Epic_converter {
typedef typename IK::FT IK_FT;
public:
std::pair<Origin, bool> operator()(Origin o) const
{
return std::make_pair(o, true);
}
std::pair<Null_vector, bool> operator()(Null_vector n) const
{
return std::make_pair(n, true);
}
std::pair<double,bool> operator()(const typename IK::FT n) const
{

4
Filtered_kernel/include/CGAL/Lazy.h
View File

@ -799,7 +799,7 @@ struct Approx_converter
template < typename T >
decltype(auto)
operator()(const T&t) const
{ return t.approx(); }
{ return approx(t); }
const Null_vector&
operator()(const Null_vector& n) const
@ -824,7 +824,7 @@ struct Exact_converter
template < typename T >
decltype(auto)
operator()(const T&t) const
{ return t.exact(); }
{ return exact(t); }
const Null_vector&
operator()(const Null_vector& n) const

160
Filtered_kernel/include/CGAL/Static_filtered_predicate.h
View File

@ -7,7 +7,6 @@
// $Id$
// SPDX-License-Identifier: LGPL-3.0-or-later OR LicenseRef-Commercial
//
//
// Author(s) : Andreas Fabri, Laurent Rineau
#ifndef CGAL_STATIC_FILTERED_PREDICATE_H
@ -29,9 +28,7 @@ public:
result_type operator()(const A1& a1) const
{
CGAL::Epic_converter<AK> convert;
typedef typename Kernel_traits<A1>::type EK;
typedef typename Type_mapper<A1,EK,Exact_predicates_inexact_constructions_kernel>::type T1;
std::pair<T1,bool> aa1 = convert(a1.approx());
auto aa1 = convert(approx(a1));
if(! aa1.second){
return fp(a1);
}
@ -39,83 +36,36 @@ public:
return epicp(aa1.first);
}
template <typename A1>
result_type operator()(const A1& a1, const Null_vector& v) const
{
CGAL::Epic_converter<AK> convert;
typedef typename Kernel_traits<A1>::type EK;
typedef typename Type_mapper<A1,EK,Exact_predicates_inexact_constructions_kernel>::type T1;
std::pair<T1,bool> aa1 = convert(a1.approx());
if(! aa1.second){
return fp(a1, v);
}
return epicp(aa1.first, v);
}
template <typename A1, typename A2>
result_type operator()(const A1& a1, const A2& a2) const
{
CGAL::Epic_converter<AK> convert;
typedef typename Kernel_traits<A1>::type EK;
typedef typename Type_mapper<A1,EK, Exact_predicates_inexact_constructions_kernel>::type T1;
std::pair<T1,bool> aa1 = convert(approx(a1));
auto aa1 = convert(approx(a1));
if(! aa1.second){
return fp(a1, a2);
}
typedef typename Type_mapper<A2,EK, Exact_predicates_inexact_constructions_kernel>::type T2;
std::pair<T2,bool> aa2 = convert(approx(a2));
auto aa2 = convert(approx(a2));
if(! aa2.second){
return fp(a1, a2);
}
return epicp(aa1.first, aa2.first);
}
// We need these two specializations as in general we determine
// the kernel for the template argument A1, and this does not work for Bbox_2 and Bbox_3
template <typename A2>
result_type operator()(const Bbox_2& bb, const A2& a2) const
{
CGAL::Epic_converter<AK> convert;
typedef typename Kernel_traits<A2>::type EK;
typedef typename Type_mapper<A2,EK, Exact_predicates_inexact_constructions_kernel>::type T2;
std::pair<T2,bool> aa2 = convert(approx(a2));
if(! aa2.second){
return fp(bb, a2);
}
return epicp(bb, aa2.first);
}
template <typename A2>
result_type operator()(const Bbox_3& bb, const A2& a2) const
{
CGAL::Epic_converter<AK> convert;
typedef typename Kernel_traits<A2>::type EK;
typedef typename Type_mapper<A2,EK, Exact_predicates_inexact_constructions_kernel>::type T2;
std::pair<T2,bool> aa2 = convert(approx(a2));
if(! aa2.second){
return fp(bb, a2);
}
return epicp(bb, aa2.first);
}
template <typename A1, typename A2, typename A3>
result_type operator()(const A1& a1, const A2& a2, const A3& a3) const
{
CGAL::Epic_converter<AK> convert;
typedef typename Kernel_traits<A1>::type EK;
typedef typename Type_mapper<A1,EK, Exact_predicates_inexact_constructions_kernel>::type T1;
std::pair<T1,bool> aa1 = convert(a1.approx());
auto aa1 = convert(approx(a1));
if(! aa1.second){
return fp(a1, a2, a3);
}
typedef typename Type_mapper<A2, EK, Exact_predicates_inexact_constructions_kernel>::type T2;
std::pair<T2,bool> aa2 = convert(a2.approx());
auto aa2 = convert(approx(a2));
if(! aa2.second){
return fp(a1, a2, a3);
}
typedef typename Type_mapper<A3,EK, Exact_predicates_inexact_constructions_kernel>::type T3;
std::pair<T3,bool> aa3 = convert(a3.approx());
auto aa3 = convert(approx(a3));
if(! aa3.second){
return fp(a1, a2, a3);
}
@ -127,24 +77,22 @@ public:
result_type operator()(const A1& a1, const A2& a2, const A3& a3, const A4& a4) const
{
CGAL::Epic_converter<AK> convert;
typedef typename Kernel_traits<A1>::type EK;
typedef typename Type_mapper<A1,EK,Exact_predicates_inexact_constructions_kernel>::type T1;
std::pair<T1,bool> aa1 = convert(a1.approx());
auto aa1 = convert(approx(a1));
if(! aa1.second){
return fp(a1, a2, a3, a4);
}
typedef typename Type_mapper<A2,EK,Exact_predicates_inexact_constructions_kernel>::type T2;
std::pair<T2,bool> aa2 = convert(a2.approx());
auto aa2 = convert(approx(a2));
if(! aa2.second){
return fp(a1, a2, a3, a4);
}
typedef typename Type_mapper<A3,EK,Exact_predicates_inexact_constructions_kernel>::type T3;
std::pair<T3,bool> aa3 = convert(a3.approx());
auto aa3 = convert(approx(a3));
if(! aa3.second){
return fp(a1, a2, a3, a4);
}
typedef typename Type_mapper<A4,EK,Exact_predicates_inexact_constructions_kernel>::type T4;
std::pair<T4,bool> aa4 = convert(a4.approx());
auto aa4 = convert(approx(a4));
if(! aa4.second){
return fp(a1, a2, a3, a4);
}
@ -155,29 +103,23 @@ public:
result_type operator()(const A1& a1, const A2& a2, const A3& a3, const A4& a4, const A5& a5) const
{
CGAL::Epic_converter<AK> convert;
typedef typename Kernel_traits<A1>::type EK;
typedef typename Type_mapper<A1,EK,Exact_predicates_inexact_constructions_kernel>::type T1;
std::pair<T1,bool> aa1 = convert(a1.approx());
auto aa1 = convert(approx(a1));
if(! aa1.second){
return fp(a1, a2, a3, a4, a5);
}
typedef typename Type_mapper<A2,EK,Exact_predicates_inexact_constructions_kernel>::type T2;
std::pair<T2,bool> aa2 = convert(a2.approx());
auto aa2 = convert(approx(a2));
if(! aa2.second){
return fp(a1, a2, a3, a4, a5);
}
typedef typename Type_mapper<A3,EK,Exact_predicates_inexact_constructions_kernel>::type T3;
std::pair<T3,bool> aa3 = convert(a3.approx());
auto aa3 = convert(approx(a3));
if(! aa3.second){
return fp(a1, a2, a3, a4, a5);
}
typedef typename Type_mapper<A4,EK,Exact_predicates_inexact_constructions_kernel>::type T4;
std::pair<T4,bool> aa4 = convert(a4.approx());
auto aa4 = convert(approx(a4));
if(! aa4.second){
return fp(a1, a2, a3, a4, a5);
}
typedef typename Type_mapper<A5,EK,Exact_predicates_inexact_constructions_kernel>::type T5;
std::pair<T5,bool> aa5 = convert(a5.approx());
auto aa5 = convert(approx(a5));
if(! aa5.second){
return fp(a1, a2, a3, a4, a5);
}
@ -188,34 +130,27 @@ public:
result_type operator()(const A1& a1, const A2& a2, const A3& a3, const A4& a4, const A5& a5, const A6& a6) const
{
CGAL::Epic_converter<AK> convert;
typedef typename Kernel_traits<A1>::type EK;
typedef typename Type_mapper<A1,EK,Exact_predicates_inexact_constructions_kernel>::type T1;
std::pair<T1,bool> aa1 = convert(a1.approx());
auto aa1 = convert(approx(a1));
if(! aa1.second){
return fp(a1, a2, a3, a4, a5, a6);
}
typedef typename Type_mapper<A2,EK,Exact_predicates_inexact_constructions_kernel>::type T2;
std::pair<T2,bool> aa2 = convert(a2.approx());
auto aa2 = convert(approx(a2));
if(! aa2.second){
return fp(a1, a2, a3, a4, a5, a6);
}
typedef typename Type_mapper<A3,EK,Exact_predicates_inexact_constructions_kernel>::type T3;
std::pair<T3,bool> aa3 = convert(a3.approx());
auto aa3 = convert(approx(a3));
if(! aa3.second){
return fp(a1, a2, a3, a4, a5, a6);
}
typedef typename Type_mapper<A4,EK,Exact_predicates_inexact_constructions_kernel>::type T4;
std::pair<T4,bool> aa4 = convert(a4.approx());
auto aa4 = convert(approx(a4));
if(! aa4.second){
return fp(a1, a2, a3, a4, a5, a6);
}
typedef typename Type_mapper<A5,EK,Exact_predicates_inexact_constructions_kernel>::type T5;
std::pair<T5,bool> aa5 = convert(a5.approx());
auto aa5 = convert(approx(a5));
if(! aa5.second){
return fp(a1, a2, a3, a4, a5, a6);
}
typedef typename Type_mapper<A6,EK,Exact_predicates_inexact_constructions_kernel>::type T6;
std::pair<T6,bool> aa6 = convert(a6.approx());
auto aa6 = convert(approx(a6));
if(! aa6.second){
return fp(a1, a2, a3, a4, a5, a6);
}
@ -226,39 +161,31 @@ public:
result_type operator()(const A1& a1, const A2& a2, const A3& a3, const A4& a4, const A5& a5, const A6& a6, const A6& a7) const
{
CGAL::Epic_converter<AK> convert;
typedef typename Kernel_traits<A1>::type EK;
typedef typename Type_mapper<A1,EK,Exact_predicates_inexact_constructions_kernel>::type T1;
std::pair<T1,bool> aa1 = convert(a1.approx());
auto aa1 = convert(approx(a1));
if(! aa1.second){
return fp(a1, a2, a3, a4, a5, a6, a7);
}
typedef typename Type_mapper<A2,EK,Exact_predicates_inexact_constructions_kernel>::type T2;
std::pair<T2,bool> aa2 = convert(a2.approx());
auto aa2 = convert(approx(a2));
if(! aa2.second){
return fp(a1, a2, a3, a4, a5, a6, a7);
}
typedef typename Type_mapper<A3,EK,Exact_predicates_inexact_constructions_kernel>::type T3;
std::pair<T3,bool> aa3 = convert(a3.approx());
auto aa3 = convert(approx(a3));
if(! aa3.second){
return fp(a1, a2, a3, a4, a5, a6, a7);
}
typedef typename Type_mapper<A4,EK,Exact_predicates_inexact_constructions_kernel>::type T4;
std::pair<T4,bool> aa4 = convert(a4.approx());
auto aa4 = convert(approx(a4));
if(! aa4.second){
return fp(a1, a2, a3, a4, a5, a6, a7);
}
typedef typename Type_mapper<A5,EK,Exact_predicates_inexact_constructions_kernel>::type T5;
std::pair<T5,bool> aa5 = convert(a5.approx());
auto aa5 = convert(approx(a5));
if(! aa5.second){
return fp(a1, a2, a3, a4, a5, a6, a7);
}
typedef typename Type_mapper<A6,EK,Exact_predicates_inexact_constructions_kernel>::type T6;
std::pair<T6,bool> aa6 = convert(a6.approx());
auto aa6 = convert(approx(a6));
if(! aa6.second){
return fp(a1, a2, a3, a4, a5, a6, a7);
}
typedef typename Type_mapper<A7,EK,Exact_predicates_inexact_constructions_kernel>::type T7;
std::pair<T7,bool> aa7 = convert(a7.approx());
auto aa7 = convert(approx(a7));
if(! aa7.second){
return fp(a1, a2, a3, a4, a5, a6, a7);
}
@ -270,44 +197,35 @@ public:
result_type operator()(const A1& a1, const A2& a2, const A3& a3, const A4& a4, const A5& a5, const A6& a6, const A7& a7, const A8& a8) const
{
CGAL::Epic_converter<AK> convert;
typedef typename Kernel_traits<A1>::type EK;
typedef typename Type_mapper<A1,EK,Exact_predicates_inexact_constructions_kernel>::type T1;
std::pair<T1,bool> aa1 = convert(a1.approx());
auto aa1 = convert(approx(a1));
if(! aa1.second){
return fp(a1, a2, a3, a4, a5, a6, a7, a8);
}
typedef typename Type_mapper<A2,EK,Exact_predicates_inexact_constructions_kernel>::type T2;
std::pair<T2,bool> aa2 = convert(a2.approx());
auto aa2 = convert(approx(a2));
if(! aa2.second){
return fp(a1, a2, a3, a4, a5, a6, a7, a8);
}
typedef typename Type_mapper<A3,EK,Exact_predicates_inexact_constructions_kernel>::type T3;
std::pair<T3,bool> aa3 = convert(a3.approx());
auto aa3 = convert(approx(a3));
if(! aa3.second){
return fp(a1, a2, a3, a4, a5, a6, a7, a8);
}
typedef typename Type_mapper<A4,EK,Exact_predicates_inexact_constructions_kernel>::type T4;
std::pair<T4,bool> aa4 = convert(a4.approx());
auto aa4 = convert(approx(a4));
if(! aa4.second){
return fp(a1, a2, a3, a4, a5, a6, a7, a8);
}
typedef typename Type_mapper<A5,EK,Exact_predicates_inexact_constructions_kernel>::type T5;
std::pair<T5,bool> aa5 = convert(a5.approx());
auto aa5 = convert(approx(a5));
if(! aa5.second){
return fp(a1, a2, a3, a4, a5, a6, a7, a8);
}
typedef typename Type_mapper<A6,EK,Exact_predicates_inexact_constructions_kernel>::type T6;
std::pair<T6,bool> aa6 = convert(a6.approx());
auto aa6 = convert(approx(a6));
if(! aa6.second){
return fp(a1, a2, a3, a5, a5, a6, a7, a8);
}
typedef typename Type_mapper<A7,EK,Exact_predicates_inexact_constructions_kernel>::type T7;
std::pair<T7,bool> aa7 = convert(a7.approx());
auto aa7 = convert(approx(a7));
if(! aa7.second){
return fp(a1, a2, a3, a5, a5, a6, a7, a8);
}
typedef typename Type_mapper<A8,EK,Exact_predicates_inexact_constructions_kernel>::type T8;
std::pair<T8,bool> aa8 = convert(a8.approx());
auto aa8 = convert(approx(a8));
if(! aa8.second){
return fp(a1, a2, a3, a4, a5, a6, a7, a8);
}

108
Kernel_23/test/Kernel_23/origin_3.cpp Normal file
View File

@ -0,0 +1,108 @@
//#define CGAL_PROFILE
#include <CGAL/Exact_predicates_exact_constructions_kernel.h>
#include <CGAL/Timer.h>
#include <vector>
#include <iostream>
typedef CGAL::Exact_predicates_exact_constructions_kernel K;
typedef K::Point_3 Point_3;
typedef K::Vector_3 Vector_3;
typedef CGAL::Timer Timer;
int main(int argc, char* argv[] )
{
std::ifstream ifs((argc>1)? argv[1]:CGAL::data_file_path("points_3/cube.xyz"));
std::vector<Point_3> points;
Point_3 p;
while(ifs >> p){
points.push_back(p);
}
std::cout << "Orientation_3" << std::endl;
Timer t;
const std::size_t N = points.size()-3;
const K::Orientation_3 orientation = K().orientation_3_object();
int positive = 0;
t.start();
{
std::cout << "overload with 4 points" << std::endl;
for(std::size_t k = 0; k < 100; ++k)
for(std::size_t i = 0; i < N; ++i){
Point_3 o(CGAL::ORIGIN);
if(orientation(o, points[i], points[i+1], points[i+2]) == CGAL::POSITIVE){
++positive;
}
}
}
t.stop();
std::cout << t.time() << " sec." << std::endl;
t.reset();
t.start();
{
std::cout << "overload with origin and 3 points" << std::endl;
for (std::size_t k = 0; k < 100; ++k)
for(std::size_t i = 0; i < N; ++i){
if(orientation(CGAL::ORIGIN, points[i], points[i+1], points[i+2]) == CGAL::POSITIVE){
--positive;
}
}
}
t.stop();
if(positive != 0){
std::cout << "Not the same results for Orientation_3"<< std::endl;
assert(false);
}
std::cout << t.time() << " sec." << std::endl;
std::cout << "Construct_orthogonal_vector_3" << std::endl;
const K::Construct_orthogonal_vector_3 construct_orthogonal_vector = K().construct_orthogonal_vector_3_object();
double sumx1 = 0, sumx2 = 0;
t.start();
{
std::cout << "overload with 3 points" << std::endl;
for(std::size_t k = 0; k < 100; ++k)
for(std::size_t i = 0; i < N; ++i){
Point_3 o(CGAL::ORIGIN);
Vector_3 v = construct_orthogonal_vector(o, points[i], points[i+1]);
sumx1 += CGAL::to_double(v.approx().x());
}
}
t.stop();
std::cout << t.time() << " sec." << std::endl;
t.reset();
t.start();
{
std::cout << "overload with origin and 2 points" << std::endl;
for (std::size_t k = 0; k < 100; ++k)
for(std::size_t i = 0; i < N; ++i){
Vector_3 v = construct_orthogonal_vector(CGAL::ORIGIN, points[i], points[i+1]);
sumx2 += CGAL::to_double(v.approx().x());
}
}
t.stop();
if(sumx1 != sumx2){
std::cout << "Not the same results for Construct_orthogonal_vector" << std::endl;
assert(false);
}
std::cout << t.time() << " sec." << std::endl;
return 0;
}

6
Nef_3/include/CGAL/Nef_3/SNC_constructor.h
View File

@ -2037,7 +2037,7 @@ class SNC_constructor<SNC_indexed_items, SNC_structure_>
D.link_as_isolated_vertex(v2, f2);
D.link_as_loop(l,f1);
D.link_as_loop(l->twin(),f2);
l->circle() = Sphere_circle(faces_p->plane());
l->circle() = Sphere_circle(CGAL::ORIGIN,faces_p->plane());
l->twin()->circle() = l->circle().opposite();
f2->mark() = mf2;
l->mark() = l->twin()->mark() = ml;
@ -2053,7 +2053,7 @@ class SNC_constructor<SNC_indexed_items, SNC_structure_>
SHalfedge_handle se1;
SHalfedge_handle se2;
SFace_handle sf;
Sphere_circle c(f->plane());
Sphere_circle c(CGAL::ORIGIN,f->plane());
SHalfedge_handle next_edge;
SHalfedge_around_svertex_const_circulator ec(E.out_edges(e)), ee(ec);
@ -2104,7 +2104,7 @@ class SNC_constructor<SNC_indexed_items, SNC_structure_>
se1 = D.new_shalfedge_pair(ec2->twin(), en->twin(), -1, 1);
CGAL_NEF_TRACEN("new edge pair " << ec2->twin()->source()->vector() <<
" -> " << en->twin()->source()->vector());
se1->circle() = Sphere_circle(faces_p->plane());
se1->circle() = Sphere_circle(CGAL::ORIGIN,faces_p->plane());
se1->twin()->circle() = se1->circle().opposite();
se1->mark() = se1->twin()->mark() = BOP(mark_of_right_sface[ec2], faces_p->mark(), inv);

7
Nef_3/include/CGAL/Nef_3/polygon_mesh_to_nef_3.h
View File

@ -171,7 +171,6 @@ template <class PolygonMesh, class SNC_structure, class FaceIndexMap, class Half
void polygon_mesh_to_nef_3(PolygonMesh& P, SNC_structure& S, FaceIndexMap fimap, HalfedgeIndexMap himap)
{
typedef typename boost::property_map<PolygonMesh, vertex_point_t>::type PMap;
typedef typename SNC_structure::Plane_3 Plane;
typedef typename SNC_structure::Vector_3 Vector_3;
typedef typename boost::graph_traits<PolygonMesh>::vertex_descriptor vertex_descriptor;
typedef typename boost::graph_traits<PolygonMesh>::halfedge_descriptor halfedge_descriptor;
@ -249,8 +248,7 @@ void polygon_mesh_to_nef_3(PolygonMesh& P, SNC_structure& S, FaceIndexMap fimap,
with_border = true;
else {
std::size_t i = get(fimap,face(pe_prev,P));
Plane ss_plane( CGAL::ORIGIN, normals[i]);
Sphere_circle ss_circle(ss_plane);
Sphere_circle ss_circle(CGAL::ORIGIN, normals[i]);
CGAL_assertion_code(if(num_edges[i] > 3) {
CGAL_assertion(ss_circle.has_on(sp));
CGAL_assertion(ss_circle.has_on(sv_prev->point()));
@ -282,8 +280,7 @@ void polygon_mesh_to_nef_3(PolygonMesh& P, SNC_structure& S, FaceIndexMap fimap,
e = sv_prev->out_sedge();
} else {
std::size_t i = get(fimap,face(pe_prev,P));
Plane ss_plane( CGAL::ORIGIN, normals[i]);
Sphere_circle ss_circle(ss_plane);
Sphere_circle ss_circle(CGAL::ORIGIN, normals[i]);
CGAL_assertion_code(if(num_edges[i] > 3) {
CGAL_assertion(ss_circle.has_on(sp_0));

39
Nef_S2/include/CGAL/Nef_S2/Sphere_circle.h
View File

@ -17,6 +17,7 @@
#include <CGAL/basic.h>
#include <CGAL/Nef_S2/Sphere_point.h>
namespace CGAL {
@ -47,6 +48,7 @@ typedef std::pair< Sphere_segment<R>,Sphere_segment<R> >
typedef typename R_::Plane_3 Plane_3;
typedef typename R_::Line_3 Line_3;
typedef typename R_::Point_3 Point_3;
typedef typename R_::Vector_3 Vector_3;
typedef Sphere_circle<R_> Self;
typedef typename R_::Plane_3 Base;
@ -61,15 +63,15 @@ $q$ are not antipodal on $S_2$, then this circle is unique and oriented
such that a walk along |\Mvar| meets $p$ just before the shorter segment
between $p$ and $q$. If $p$ and $q$ are antipodal of each other then we
create any great circle that contains $p$ and $q$.}*/
{ Point_3 p1(0,0,0), p4 = CGAL::ORIGIN + ((Base*) this)->orthogonal_vector();
{
if ( p != q.antipode() ) {
if (R_().orientation_3_object()(p1,Point_3(p),
Point_3(q), p4) != CGAL::POSITIVE )
*this = Self(opposite());
Point_3 po = CGAL::ORIGIN + Base::orthogonal_vector();
if (R_().orientation_3_object()(CGAL::ORIGIN, Point_3(p),
Point_3(q), po) != CGAL::POSITIVE )
*this = opposite();
} else {
/* previous method was: *this = Self(Plane_3(p1,q-p));
but p, q don't belong to he plane ((0,0,0), q-p) */
/* previous method was: *this = Self(Plane_3((0,0,0),q-p));
but p, q don't belong to the plane ((0,0,0), q-p) */
if(!Line_3(p,q).has_on(Point_3(1,0,0)))
*this = Self(Plane_3(p,q,Point_3(1,0,0)));
else
@ -78,16 +80,22 @@ create any great circle that contains $p$ and $q$.}*/
}
}
Sphere_circle(const Plane_3& h) : Base(h)
Sphere_circle(const Plane_3& h) : Base(h)
/*{\Mcreate creates the circle of $S_2$ corresponding to the plane
|h|. If |h| does not contain the origin, then |\Mvar| becomes the
circle parallel to |h| containing the origin.}*/
{
if(h.d() != 0) *this = Plane_3(h.a(),h.b(),h.c(),RT(0));
if(!is_zero(h.d()))
*this = Plane_3(h.a(),h.b(),h.c(),RT(0));
}
Sphere_circle(const RT& x, const RT& y, const RT& z): Base(x,y,z,0) {}
Sphere_circle(const Origin& o, const Vector_3& v) : Base(o,v) {}
Sphere_circle(const Origin&, const Plane_3& h) : Base(h.a(),h.b(),h.c(),RT(0))
{/* Even if |h| does not contain the origin, the circle will contain the origin
and be parallel to |h| */}
Sphere_circle(const RT& x, const RT& y, const RT& z) : Base(x,y,z,0) {}
/*{\Mcreate creates the circle orthogonal to the vector $(x,y,z)$.}*/
Sphere_circle(Sphere_circle<R> c, const Sphere_point<R>& p)
@ -102,10 +110,9 @@ Sphere_circle(Sphere_circle<R> c, const Sphere_point<R>& p)
}
/*{\Moperations 4 2}*/
Sphere_circle<R> opposite() const
/*{\Mop returns the opposite of |\Mvar|.}*/
{ return Base::opposite(); }
{ return Sphere_circle<R>(Base::opposite(),Assume_d_equal_0{}); }
bool has_on(const Sphere_point<R>& p) const
/*{\Mop returns true iff |\Mvar| contains |p|.}*/
@ -117,12 +124,12 @@ Plane_3 plane() const { return Base(*this); }
Plane_3 plane_through(const Point_3& p) const
/*{\Mop returns the plane parallel to |\Mvar| that
contains point |p|.}*/
{ return Plane_3(p,((Base*) this)->orthogonal_direction()); }
{ return Plane_3(p,Base::orthogonal_vector()); }
Sphere_point<R> orthogonal_pole() const
/*{\Mop returns the point that is the pole of the
hemisphere left of |\Mvar|.}*/
{ return CGAL::ORIGIN+((Base*) this)->orthogonal_vector(); }
{ return CGAL::ORIGIN+Base::orthogonal_vector(); }
Sphere_segment_pair split_at(const Sphere_point<R>& p) const;
/*{\Mop returns the pair of circle segments that is the result
@ -134,6 +141,10 @@ of splitting |\Mvar| at the $x$-$y$-coordinate plane if |\Mvar|
is not part of it. Otherwise |\Mvar| is split at the
$x$-$z$-coordinate plane.}*/
private:
struct Assume_d_equal_0{};
Sphere_circle(const Plane_3& h, Assume_d_equal_0) : Base(h){}
}; // Sphere_circle<R>
/*{\Mtext\headerline{Global functions}}*/

39
Nef_S2/include/CGAL/Nef_S2/Sphere_segment.h
View File

@ -37,14 +37,15 @@ Sphere_segment_rep() { ps_ = pt_ = Point(); c_ = Circle(); }
Sphere_segment_rep(const Point& p1, const Point& p2,
bool shorter_arc=true) :
ps_(p1), pt_(p2), c_(Plane_3(p1,p2,Point_3(CGAL::ORIGIN)))
ps_(p1), pt_(p2),
c_(CGAL::ORIGIN,R_().construct_orthogonal_vector_3_object()(CGAL::ORIGIN,p1,p2))
{ // warning stays as reminder that one gets an arbitrary plane equation
// in this degenerate case
CGAL_warning(p1 != p2.antipode());
CGAL_assertion(p1 != p2.antipode());
if ( p1 == p2 ) {
Plane_3 h(Point_3(CGAL::ORIGIN),(p1-CGAL::ORIGIN));
c_ = Sphere_circle<R_>(Plane_3(Point_3(CGAL::ORIGIN),h.base1()));
Plane_3 h(CGAL::ORIGIN,p1-CGAL::ORIGIN);
c_ = Sphere_circle<R_>(CGAL::ORIGIN,h.base1());
}
if (!shorter_arc) c_ = c_.opposite();
CGAL_exactness_assertion(c_.has_on(p1) && c_.has_on(p2));
@ -59,7 +60,7 @@ Sphere_segment_rep(const Circle& c1,
{ CGAL_assertion(!equal_as_sets(c1,c2));
ps_ = intersection(c1,c2);
pt_ = ps_.antipode();
if ( R_::orientation(Point_3(CGAL::ORIGIN),ps_,pt_,
if ( R_().orientation_3_object()(CGAL::ORIGIN,ps_,pt_,
CGAL::ORIGIN + c_.orthogonal_vector()) !=
CGAL::POSITIVE ) std::swap(ps_,pt_);
}
@ -176,7 +177,7 @@ void split_halfcircle(Sphere_segment<R>& s1,
/*{\Mop splits a halfcircle into two equally sized segments.
\precond |\Mvar| is a halfcircle.}*/
{ CGAL_assertion( is_halfcircle() );
Plane_3 h(Point_3(CGAL::ORIGIN),(target()-CGAL::ORIGIN));
Plane_3 h(CGAL::ORIGIN,(target()-CGAL::ORIGIN));
Sphere_point<R> p =
CGAL::intersection(sphere_circle(),Sphere_circle<R>(h));
if ( !has_on_after_intersection(p) ) p = p.antipode();
@ -187,17 +188,17 @@ void split_halfcircle(Sphere_segment<R>& s1,
bool is_short() const
/*{\Mop a segment is short iff it is shorter than a halfcircle.}*/
{
return R().orientation_3_object()(Point_3(CGAL::ORIGIN),
Point_3(source()),
Point_3(target()),
return R().orientation_3_object()(CGAL::ORIGIN,
source(),
target(),
orthogonal_pole())
== CGAL::POSITIVE; }
bool is_long() const
/*{\Mop a segment is long iff it is longer than a halfcircle.}*/
{ return R().orientation_3_object()(Point_3(CGAL::ORIGIN),
Point_3(source()),
Point_3(target()),
{ return R().orientation_3_object()(CGAL::ORIGIN,
source(),
target(),
orthogonal_pole())
== CGAL::NEGATIVE; }
@ -229,17 +230,17 @@ Point_3 orthogonal_pole() const
{ return CGAL::ORIGIN + sphere_circle().orthogonal_vector(); }
CGAL::Orientation source_orientation(const CGAL::Sphere_point<R>& p) const
{ return orientation(Point_3(CGAL::ORIGIN),
orthogonal_pole(),
source(),
p);
{ return R().orientation_3_object()(CGAL::ORIGIN,
orthogonal_pole(),
source(),
p);
}
CGAL::Orientation target_orientation(const CGAL::Sphere_point<R>& p) const
{ return orientation(Point_3(CGAL::ORIGIN),
target(),
orthogonal_pole(),
p);
{ return R().orientation_3_object()(CGAL::ORIGIN,
target(),
orthogonal_pole(),
p);
}
};