Search code examples
computational-geometrycgal

Check whether a given point is on the (un)bounded side/boundary of a circle (CGAL)


I'm using Kernel::Circle_2 with CGAL::Arr_circle_segment_traits_2.

Given a point (of the nested type Point_2 of this set of traits), I would like to check whether it's on the bounded side, unbounded side or on the boundary of a given circle.

There is the function called bounded_side of the class Circle_2, yet it supports only points of Kernel::Point_2. When I'm using CGAL::to_double() to convert the point to this class, I lose accuracy.

Is there another way to perform this check? I store the information in a 2D_Arrangement, if that helps.


Solution

  • You can use the code following. Note that the coordinate of the 2D points are of type Sqrt_extension.

        #include <CGAL/Exact_predicates_exact_constructions_kernel.h>
        #include <CGAL/Arr_circle_segment_traits_2.h>
    
        typedef CGAL::Exact_predicates_exact_constructions_kernel K;
        typedef CGAL::Arr_circle_segment_traits_2<K> Traits;
        typedef Traits::Point_2::CoordNT Sqrt_extension;
    
        CGAL::Bounded_side
        incircle(const typename K::Circle_2& circle,
                 const typename Traits::Point_2& p)
        {
          const K::Point_2& center = circle.center();
          K::FT sq_rad = circle.squared_radius();
    
          switch(CGAL::compare(  CGAL::square(p.x()-center.x())-sq_rad,
                                -CGAL::square(p.y()-center.y()) ) )
          {
            case CGAL::LARGER:
              return CGAL::ON_UNBOUNDED_SIDE;
            case CGAL::SMALLER:
              return CGAL::ON_BOUNDED_SIDE;
            case CGAL::EQUAL:
              break;
          }
          return CGAL::ON_BOUNDARY;
        }
    
        int main()
        {
    
          K::Circle_2 circle(K::Point_2(0,0), 2);
    
          Traits::Point_2 out(Sqrt_extension(1,2,3) , Sqrt_extension(4,5,6));
          CGAL_assertion( incircle(circle, out) == CGAL::ON_UNBOUNDED_SIDE );
    
          Traits::Point_2 in(Sqrt_extension(1) , Sqrt_extension(0));
          CGAL_assertion( incircle(circle, in) == CGAL::ON_BOUNDED_SIDE );
    
          Traits::Point_2 bnd(Sqrt_extension(0,1,2) , Sqrt_extension(0));
          CGAL_assertion( incircle(circle, bnd) == CGAL::ON_BOUNDARY );
    
          return 0;
        }