File: internal_functions_on_circle_3.h

package info (click to toggle)
cgal 4.0-5
  • links: PTS
  • area: main
  • in suites: wheezy
  • size: 65,068 kB
  • sloc: cpp: 500,870; ansic: 102,544; sh: 321; python: 92; makefile: 75; xml: 2
file content (196 lines) | stat: -rw-r--r-- 8,212 bytes parent folder | download
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
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
// Copyright (c) 2008  INRIA Sophia-Antipolis (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: svn+ssh://scm.gforge.inria.fr/svn/cgal/branches/next/Circular_kernel_3/include/CGAL/Circular_kernel_3/internal_functions_on_circle_3.h $
// $Id: internal_functions_on_circle_3.h 67117 2012-01-13 18:14:48Z lrineau $
//
// Author(s) : Monique Teillaud, Sylvain Pion, Pedro Machado, 
//             Sebastien Loriot, Julien Hazebrouck, Damien Leroy

// Partially supported by the IST Programme of the EU as a 
// STREP (FET Open) Project under Contract No  IST-006413 
// (ACS -- Algorithms for Complex Shapes)

#ifndef CGAL_SPHERICAL_KERNEL_PREDICATES_ON_CIRCLE_3_H
#define CGAL_SPHERICAL_KERNEL_PREDICATES_ON_CIRCLE_3_H

#include <CGAL/Circle_type.h>
#include <CGAL/Circular_kernel_3/internal_functions_on_plane_3.h>
#include <CGAL/Circular_kernel_3/internal_functions_on_circular_arc_point_3.h>
#include <CGAL/Root_of_traits.h>

namespace CGAL {
  namespace SphericalFunctors {
            
    template < class SK >
    typename SK::Circle_3
    construct_circle_3(const typename SK::Polynomials_for_circle_3 &eq)
    {
      typedef typename SK::Point_3 Point_3;
      typedef typename SK::Plane_3 Plane_3;
      typedef typename SK::Circle_3 Circle_3;
      typedef typename SK::Sphere_3 Sphere_3;
      typedef typename SK::FT FT;
      Sphere_3 s = construct_sphere_3<SK>(eq.first);
      Plane_3 p = construct_plane_3<SK>(eq.second);
      const FT d2 = CGAL::square(p.a()*s.center().x() + 
                                 p.b()*s.center().y() + 
                                 p.c()*s.center().z() + p.d()) /
       (CGAL::square(p.a()) + CGAL::square(p.b()) + CGAL::square(p.c()));
      // We do not accept circles with radius 0 (should we?)
      CGAL_kernel_precondition(d2 <  s.squared_radius());
      // d2 < s.squared_radius()
      Point_3 center = p.projection(s.center());
      return Circle_3(center,s.squared_radius() - d2,p);
    }


    template < class SK >
    CGAL::Circle_type
    classify_circle_3(const typename SK::Circle_3& circle,const typename SK::Sphere_3& sphere)
    {
      typedef typename SK::Algebraic_kernel::Root_for_spheres_2_3   Root_for_spheres_2_3;
      typedef typename SK::Circular_arc_point_3                     Circular_arc_point_3;
      
      CGAL_kernel_precondition(SK().has_on_3_object()(sphere,circle));
      
      //if circle is a great circle, it can only be a bipolar or a threaded.
      if (circle.center()==sphere.center()){
        if (circle.supporting_plane().orthogonal_vector().z()==0)
          return CGAL::BIPOLAR;
        return CGAL::THREADED;
      }
      
      typename SK::Root_of_2 radius=CGAL::make_sqrt(sphere.squared_radius());
      Circular_arc_point_3 north_pole( Root_for_spheres_2_3(sphere.center().x(),sphere.center().y(),sphere.center().z()+radius) );
      Circular_arc_point_3 south_pole( Root_for_spheres_2_3(sphere.center().x(),sphere.center().y(),sphere.center().z()-radius) );


      const typename SK::Sphere_3& supp_sphere=circle.diametral_sphere();
      typename SK::Bounded_side_3 bounded_side=SK().bounded_side_3_object();
      

      CGAL::Bounded_side side_of_north=bounded_side(supp_sphere,north_pole);
      CGAL::Bounded_side side_of_south=bounded_side(supp_sphere,south_pole);
      
      if (side_of_north==ON_BOUNDARY || side_of_south==ON_BOUNDARY)
        return CGAL::POLAR;

      if (side_of_north==ON_BOUNDED_SIDE || side_of_south==ON_BOUNDED_SIDE)
        return CGAL::THREADED;      
      
      CGAL_kernel_precondition(side_of_north==ON_UNBOUNDED_SIDE && side_of_south==ON_UNBOUNDED_SIDE);
      
      return CGAL::NORMAL;
    }


    template < class SK >
    inline typename SK::FT
    extremal_points_z_coordinate(const typename SK::Circle_3& circle,const typename SK::Sphere_3& sphere)
    {
      CGAL_kernel_precondition(SK().has_on_3_object()(sphere,circle));
      CGAL_kernel_precondition(classify_circle_3<SK>(circle,sphere)==CGAL::NORMAL);
      
      const typename SK::Point_3& circle_center=circle.center();
      const typename SK::Point_3& sphere_center=sphere.center();

      return
      typename SK::FT(2) * (circle_center-sphere_center).z() * sphere.squared_radius() 
               / ( SK().compute_squared_distance_3_object()(circle_center,sphere_center) + sphere.squared_radius()-circle.squared_radius() )
               + sphere_center.z();
    }
    
    template < class SK, class Output_iterator >
    Output_iterator theta_extremal_points(const typename SK::Circle_3& circle,const typename SK::Sphere_3& sphere,Output_iterator out_it){
      CGAL_kernel_precondition(classify_circle_3<SK>(circle,sphere)==NORMAL);
      CGAL_kernel_precondition(SK().has_on_3_object()(sphere,circle));
      
      typename SK::FT z_coord=extremal_points_z_coordinate<SK>(circle,sphere);
      
      typename SK::Plane_3 plane(0,0,1,-z_coord);
      std::vector<CGAL::Object> inters;
      
      intersect_3<SK>(circle,plane,std::back_inserter(inters));      
      CGAL_kernel_precondition(inters.size()==2);
      const std::pair<typename SK::Circular_arc_point_3,unsigned>* pt[2]={NULL,NULL};
      pt[0]=object_cast<std::pair<typename SK::Circular_arc_point_3,unsigned> >(&inters[0]);
      pt[1]=object_cast<std::pair<typename SK::Circular_arc_point_3,unsigned> >(&inters[1]);
      CGAL_kernel_precondition(pt[0]!=NULL);
      CGAL_kernel_precondition(pt[1]!=NULL);
      
      if ( compare_theta_of_pts<SK>(pt[0]->first,pt[1]->first,sphere) == SMALLER){
        *out_it++=pt[0]->first;
        *out_it++=pt[1]->first;
      }
      else{
        *out_it++=pt[1]->first;
        *out_it++=pt[0]->first;
      }
        
      return out_it;
    }
    
    template < class SK >
    typename SK::Circular_arc_point_3 theta_extremal_point(const typename SK::Circle_3& circle,const typename SK::Sphere_3& sphere,bool is_smallest){
      typename SK::Circular_arc_point_3 pts[2];
      theta_extremal_points(circle,sphere,pts);
      if (is_smallest)
        return pts[0];
      return pts[1];
    }
    
    template < class SK,class Output_iterator>
    Output_iterator make_circle_theta_monotone(const typename SK::Circle_3& circle,const typename SK::Sphere_3& sphere,Output_iterator out_it){
      CGAL::Circle_type type=classify_circle_3<SK>(circle,sphere);
      switch (type){
        case THREADED:
        {
          *out_it++=typename SK::Circular_arc_3(circle);
          break;
        }  
        case POLAR:{
          typename SK::Vector_3 ortho=circle.center()-sphere.center();
          CGAL_kernel_precondition(ortho.z()!=0);
          bool is_north_pole=ortho.z()>0;
          typename SK::Root_of_2 radius = (is_north_pole?1:-1)* make_sqrt(sphere.squared_radius());
          typename SK::Circular_arc_point_3 source_target(
            typename SK::Algebraic_kernel::Root_for_spheres_2_3(
              sphere.center().x(),
              sphere.center().y(),
              sphere.center().z()+radius
            )
          );
          *out_it++=typename SK::Circular_arc_3(circle,source_target);
          break;
        }
        case NORMAL:{
          typename SK::Circular_arc_point_3 ints[2];
          theta_extremal_points(circle,sphere,ints);
          *out_it++=typename SK::Circular_arc_3(circle,ints[0],ints[1]);
          *out_it++=typename SK::Circular_arc_3(circle,ints[1],ints[0]);
        }
        break;
        case BIPOLAR:
          CGAL_kernel_precondition(!"This function does not accept bipolar circle as input.");
      }
      return out_it;
    }
    

    
  }//SphericalFunctors
}//CGAL

#endif //CGAL_SPHERICAL_KERNEL_PREDICATES_ON_CIRCLE_3_H