CGAL Cartesian grid

Question

In my code, I organize objects into a regular Cartesian grid (such as 10x10). Often given a point, I need to test whether the point intersects grid and if so, which bins contain the point. I already have my own implementation but I don't like to hassle with precision issues.

So, does CGAL has a 2D regular Cartesian grid?

Answer 1

You can use CGAL::points_on_square_grid_2 to generate the grid points. CGAL kernels provide Kernel::CompareXY_2 functors, which you can use to figure out the exact location of your query point on the grid. For example you can sort your grid points and then use std::lower_bound followed by CGAL::orientation or CGAL::collinear on the appropriate elements of your range. You could also build an arrangement, but this would be an overkill.

Here is a sample code.

#include <CGAL/Exact_predicates_exact_constructions_kernel.h>
#include <CGAL/point_generators_2.h>
#include <CGAL/random_selection.h>
#include <CGAL/Polygon_2_algorithms.h>

using namespace CGAL;
using K= Exact_predicates_exact_constructions_kernel;
using Point =K::Point_2;
using Creator = Creator_uniform_2<double, Point>;
using Grid = std::vector<Point>;

const int gridSide = 3;

void locate_point (Point p, Grid grid);

int main ()
{
  Grid points;
  points_on_square_grid_2(gridSide * gridSide, gridSide * gridSide, std::back_inserter(points), Creator());

  std::sort(points.begin(), points.end(), K::Less_xy_2());

  std::cout << "Grid points:\n";
  for (auto& p:points)
    std::cout << p << '\n';

  std::cout << "\ncorner points:\n";
  Grid cornerPoints{points[0], points[gridSide - 1], points[gridSide * gridSide - 1],
                    points[gridSide * (gridSide - 1)]};
  for (auto& p:cornerPoints)
    std::cout << p << '\n';

  std::cout << '\n';

  Point p1{-8, -8};
  Point p2{-10, 3};
  Point p3{-9, -8};
  Point p4{0, 4};
  Point p5{1, 5};

  locate_point(p1, points);
  locate_point(p2, points);
  locate_point(p3, points);
  locate_point(p4, points);
  locate_point(p5, points);

}

void locate_point (Point p, Grid grid)
{
  if (grid.empty())
    {
      std::cout << "Point " << p << " not in grid";
      return;
    }

  // check if point is in grid
  Grid cornerPoints{grid[0], grid[gridSide - 1], grid[gridSide * gridSide - 1], grid[gridSide * (gridSide - 1)]};
  auto point_is = CGAL::bounded_side_2(cornerPoints.begin(), cornerPoints.end(), p);
  switch (point_is)
    {
  case CGAL::ON_UNBOUNDED_SIDE:
    std::cout << "Point " << p << " not in grid\n";
      return;
  case CGAL::ON_BOUNDARY:
    std::cout << "Point " << p << " on grid boundary\n";
      return;
  case CGAL::ON_BOUNDED_SIDE:
    std::cout << "Point " << p << " is in grid\n";
    }

  auto f = std::lower_bound(grid.begin(), grid.end(), p, K::Less_xy_2());
  auto g = std::find_if(f, grid.end(), [&p] (const Point& gridpoint)
  { return K::Less_y_2()(p, gridpoint); });

  if (CGAL::collinear(p, *g, *(g - 1)))
    {
      std::cout << "Point " << p << " on grid side between points " << *(g - 1) << " and " << *g << '\n';
      return;
    }

  std::cout << "Point " << p << " in bin whose upper right point is " << *g << '\n';
  return;

}

Output:

Grid points:
-9 -9
-9 0
-9 9
0 -9
0 0
0 9
9 -9
9 0
9 9

corner points:
-9 -9
-9 9
9 9
9 -9

Point -8 -8 is in grid
Point -8 -8 in bin whose upper right point is 0 0
Point -10 3 not in grid
Point -9 -8 on grid boundary
Point 0 4 is in grid
Point 0 4 on grid side between points 0 0 and 0 9
Point 1 5 is in grid
Point 1 5 in bin whose upper right point is 9 9