Collectives™ on Stack Overflow

Find centralized, trusted content and collaborate around the technologies you use most.

Learn more about Collectives

Teams

Q&A for work

Connect and share knowledge within a single location that is structured and easy to search.

Learn more about Teams

I'm using CGAL (v. 4.14-2) to compute the visibility region (polygon) within a simple polygon from one of its vertices, using the Simple_polygon_visibility_2 class. On a particular combination of polygon/vertex, I'm getting the following assertion error:

terminate called after throwing an instance of 'CGAL::Assertion_exception'
  what():  CGAL ERROR: assertion violation!
Expr: k+1<vertices.size()
File: /usr/include/CGAL/Simple_polygon_visibility_2.h
Line: 678

This occurs in the method scan_edges of the class Simple_polygon_visibility_2, it seems like there should be some intersection with the given segment/ray, but none was found.

When I run the same code with the class Triangular_expansion_visibility_2, it seems to work, giving the output:

Regularized visibility region of q has 8 edges:
[7968 492 -> 7938 428]
[7884 408 -> 7968 492]
[8040 428 -> 7884 408]
[8090.99 428 -> 8040 428]
[8105.55 458.865 -> 8090.99 428]
[8090 456 -> 8105.55 458.865]
[7968 446 -> 8090 456]
[7938 428 -> 7968 446]

Here's a minimal working example:

#include <CGAL/Arr_naive_point_location.h>
#include <CGAL/Arr_segment_traits_2.h>
#include <CGAL/Arrangement_2.h>
#include <CGAL/Exact_predicates_exact_constructions_kernel.h>
#include <CGAL/Polygon_2.h>
#include <CGAL/Simple_polygon_visibility_2.h>
#include <CGAL/Triangular_expansion_visibility_2.h>
#include <fstream>
#include <iostream>
#include <string>
#include <list>
#include <vector>
typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel;
typedef Kernel::Point_2 Point_2;
typedef Kernel::Segment_2 Segment_2;
typedef CGAL::Polygon_2<Kernel, std::list<Point_2>> Polygon_2;
typedef CGAL::Arr_segment_traits_2<Kernel> Traits_2;
typedef CGAL::Arrangement_2<Traits_2> Arrangement_2;
typedef Arrangement_2::Face_handle Face_handle;
typedef Arrangement_2::Vertex_handle ArrVertex_handle;
typedef Arrangement_2::Edge_const_iterator Edge_const_iterator;
typedef Arrangement_2::Ccb_halfedge_circulator Ccb_halfedge_circulator;
using namespace std;
vector<Point_2> readInput() {
  vector<Point_2> Points;
  ifstream File;
  File.open("polygon");
  string Line;
  long int Idx, XCoord, YCoord;
  while (getline(File, Line)) {
    istringstream Iss(Line);
    Iss >> Idx >> XCoord >> YCoord;
    Points.push_back({XCoord, YCoord});
  return Points;
int main() {
  // create environment
  std::vector<Segment_2> segments;
  Arrangement_2 env;
  Segment_2 CurrSegment;
  Polygon_2 AuxPoly;
  auto Points = readInput();
  auto QueryPt = Points[0];
  ArrVertex_handle CurrPtHandle =
      env.insert_in_face_interior(QueryPt, env.unbounded_face());
  AuxPoly.push_back(Points[0]);
  ArrVertex_handle NextPtHandle = CurrPtHandle;
  ArrVertex_handle QueryPtHandle = CurrPtHandle;
  for (auto i = 1; i < Points.size(); ++i) {
    NextPtHandle = env.insert_in_face_interior(Points[i], env.unbounded_face());
    CurrSegment = Segment_2(CurrPtHandle->point(), NextPtHandle->point());
    env.insert_at_vertices(CurrSegment, CurrPtHandle, NextPtHandle);
    CurrPtHandle = NextPtHandle;
    AuxPoly.push_back(Points[i]);
  assert(AuxPoly.is_simple());
  CurrSegment = Segment_2(CurrPtHandle->point(), QueryPtHandle->point());
  auto HalfEHandle =
      env.insert_at_vertices(CurrSegment, CurrPtHandle, QueryPtHandle);
  if (HalfEHandle->face()->is_unbounded()) {
    // there are exactly two incident halfedges in the query point
    auto FirstHalfE = HalfEHandle->target()->incident_halfedges();
    HalfEHandle = (FirstHalfE != HalfEHandle) ? FirstHalfE : next(FirstHalfE, 1);
  // compute non regularized visibility area
  // Define visibiliy object type that computes regularized visibility area
  typedef CGAL::Simple_polygon_visibility_2<Arrangement_2, CGAL::Tag_true> RSPV;
  // typedef CGAL::Triangular_expansion_visibility_2<Arrangement_2, CGAL::Tag_true>  RSPV;
  Arrangement_2 regular_output;
  RSPV regular_visibility(env);
  regular_visibility.compute_visibility(QueryPtHandle->point(), HalfEHandle, regular_output);
  std::cout << "Regularized visibility region of q has "
            << regular_output.number_of_edges() << " edges:" << std::endl;
  for (Edge_const_iterator eit = regular_output.edges_begin();
       eit != regular_output.edges_end(); ++eit)
    std::cout << "[" << eit->source()->point() << " -> "
              << eit->target()->point() << "]" << std::endl;
  return 0;

an auxiliary file arr_print.h to print the arrangement:

#ifndef _PRINT_ARR_H_
#define _PRINT_ARR_H_
#include <iostream>
//-----------------------------------------------------------------------------
// Print all neighboring vertices to a given arrangement vertex.
template<class Arrangement>
void print_neighboring_vertices (typename Arrangement::Vertex_const_handle v)
  if (v->is_isolated())
    std::cout << "The vertex (" << v->point() << ") is isolated" << std::endl;
    return;
  typename Arrangement::Halfedge_around_vertex_const_circulator  first, curr;
  typename Arrangement::Vertex_const_handle                      u;
  std::cout << "The neighbors of the vertex (" << v->point() << ") are:";
  first = curr = v->incident_halfedges();
    // Note that the current halfedge is (u -> v):
    u = curr->source();
    std::cout << " (" << u->point() << ")";
    ++curr;
  } while (curr != first);
  std::cout << std::endl;
  return;
//-----------------------------------------------------------------------------
// Print all vertices (points) and edges (curves) along a connected component
// boundary.
template<class Arrangement>
void print_ccb (typename Arrangement::Ccb_halfedge_const_circulator circ)
  typename Arrangement::Ccb_halfedge_const_circulator  curr = circ;
  typename Arrangement::Halfedge_const_handle          he;
  std::cout << "(" << curr->source()->point() << ")";
    he = curr;
    std::cout << "   [" << he->curve() << "]   "
              << "(" << he->target()->point() << ")";
    ++curr;
  } while (curr != circ);
  std::cout << std::endl;
  return;
//-----------------------------------------------------------------------------
// Print the boundary description of an arrangement face.
template<class Arrangement>
void print_face (typename Arrangement::Face_const_handle f)
  // Print the outer boundary.
  if (f->is_unbounded())
    std::cout << "Unbounded face. " << std::endl;
    std::cout << "Outer boundary: ";
    print_ccb<Arrangement> (f->outer_ccb());
  // Print the boundary of each of the holes.
  typename Arrangement::Hole_const_iterator  hole;
  int                                         index = 1;
  for (hole = f->holes_begin(); hole != f->holes_end(); ++hole, ++index)
    std::cout << "    Hole #" << index << ": ";
    print_ccb<Arrangement> (*hole);
  // Print the isolated vertices.
  typename Arrangement::Isolated_vertex_const_iterator  iv;
  for (iv = f->isolated_vertices_begin(), index = 1;
       iv != f->isolated_vertices_end(); ++iv, ++index)
    std::cout << "    Isolated vertex #" << index << ": "
              << "(" << iv->point() << ")" << std::endl;
  return;
//-----------------------------------------------------------------------------
// Print the given arrangement.
template<class Arrangement>
void print_arrangement (const Arrangement& arr)
  // CGAL_precondition (arr.is_valid());
  // Print the arrangement vertices.
  typename Arrangement::Vertex_const_iterator  vit;
  std::cout << arr.number_of_vertices() << " vertices:" << std::endl;
  for (vit = arr.vertices_begin(); vit != arr.vertices_end(); ++vit)
    std::cout << "(" << vit->point() << ")";
    if (vit->is_isolated())
      std::cout << " - Isolated." << std::endl;
      std::cout << " - degree " << vit->degree() << std::endl;
  // Print the arrangement edges.
  typename Arrangement::Edge_const_iterator    eit;
  std::cout << arr.number_of_edges() << " edges:" << std::endl;
  for (eit = arr.edges_begin(); eit != arr.edges_end(); ++eit)
    std::cout << "[" << eit->curve() << "]" << std::endl;
  // Print the arrangement faces.
  typename Arrangement::Face_const_iterator    fit;
  std::cout << arr.number_of_faces() << " faces:" << std::endl;
  for (fit = arr.faces_begin(); fit != arr.faces_end(); ++fit)
    print_face<Arrangement> (fit);
  return;
#endif

The input polygon:

input polygon

The input file has in each line the polygon vertex id, followed by its x- and-coordinates. The query point is the one with id 1716.

Could anyone help me with this issue?

Thank you all in advance.

Your output shows you used rational numbers for points - for example 1318831/163 etc, but I don't see that in your code – HEKTO Nov 13, 2019 at 19:48 I believe it's a bug in the Simple_polygon_visibility_2 class. If you shift your query point (7938,428) vertically a litlle - for example, to (7938,429) - then all will work fine. Conclusion - this implementation can't handle vertices with the same y-coordinate, and you have three of them with y = 428. I think you should shrink your data file (6000 points is too many) and submit this bug to the CGAL team. – HEKTO Nov 15, 2019 at 19:20

Thanks for contributing an answer to Stack Overflow!

  • Please be sure to answer the question. Provide details and share your research!

But avoid

  • Asking for help, clarification, or responding to other answers.
  • Making statements based on opinion; back them up with references or personal experience.

To learn more, see our tips on writing great answers.