File difference.cpp
File List > algorithm > difference.cpp
Go to the documentation of this file
// Copyright (c) 2012-2013, IGN France.
// Copyright (c) 2012-2022, Oslandia.
// SPDX-License-Identifier: LGPL-2.0-or-later
#include "SFCGAL/algorithm/difference.h"
#include "SFCGAL/Exception.h"
#include "SFCGAL/Polygon.h"
#include "SFCGAL/TriangulatedSurface.h"
#include "SFCGAL/algorithm/differencePrimitives.h"
#include "SFCGAL/algorithm/isValid.h"
#include "SFCGAL/detail/GeometrySet.h"
#include "SFCGAL/triangulate/triangulatePolygon.h"
#include <CGAL/Boolean_set_operations_2.h>
#include <CGAL/Exact_predicates_exact_constructions_kernel.h>
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/box_intersection_d.h>
//
// Intersection kernel
using namespace SFCGAL::detail;
namespace SFCGAL::algorithm {
template <int Dim>
struct CollisionMapper {
using PrimitiveHandleSet = std::vector<PrimitiveHandle<Dim> *>;
using Map = std::map<PrimitiveHandle<Dim> *, PrimitiveHandleSet>;
CollisionMapper(Map &map) : _map(map){};
void
operator()(const typename PrimitiveBox<Dim>::Type &a,
const typename PrimitiveBox<Dim>::Type &b)
{
_map[a.handle()].push_back(b.handle());
}
private:
Map &_map;
};
template <typename OutputIteratorType>
auto
difference(const Point_2 &primitive, const PrimitiveHandle<2> &pb,
OutputIteratorType out) -> OutputIteratorType
{
switch (pb.handle.which()) {
case PrimitivePoint:
difference(primitive, *pb.as<Point_2>(), out);
break;
case PrimitiveSegment:
difference(primitive, *pb.as<Segment_2>(), out);
break;
case PrimitiveSurface:
difference(primitive, *pb.as<PolygonWH_2>(), out);
break;
}
return out;
}
template <typename OutputIteratorType>
auto
difference(const Segment_2 &primitive, const PrimitiveHandle<2> &pb,
OutputIteratorType out) -> OutputIteratorType
{
switch (pb.handle.which()) {
case PrimitivePoint:
*out++ = primitive;
break;
case PrimitiveSegment:
difference(primitive, *pb.as<Segment_2>(), out);
break;
case PrimitiveSurface:
difference(primitive, *pb.as<PolygonWH_2>(), out);
break;
}
return out;
}
template <typename OutputIteratorType>
auto
difference(const PolygonWH_2 &primitive, const PrimitiveHandle<2> &pb,
OutputIteratorType out) -> OutputIteratorType
{
switch (pb.handle.which()) {
case PrimitivePoint:
*out++ = primitive;
break;
case PrimitiveSegment:
*out++ = primitive;
break;
case PrimitiveSurface:
difference(primitive, *pb.as<PolygonWH_2>(), out);
break;
}
return out;
}
template <typename OutputIteratorType>
auto
difference(const Point_3 &primitive, const PrimitiveHandle<3> &pb,
OutputIteratorType out) -> OutputIteratorType
{
switch (pb.handle.which()) {
case PrimitivePoint:
difference(primitive, *pb.as<Point_3>(), out);
break;
case PrimitiveSegment:
difference(primitive, *pb.as<Segment_3>(), out);
break;
case PrimitiveSurface:
difference(primitive, *pb.as<Triangle_3>(), out);
break;
case PrimitiveVolume:
difference(primitive, *pb.as<MarkedPolyhedron>(), out);
break;
}
return out;
}
template <typename OutputIteratorType>
auto
difference(const Segment_3 &primitive, const PrimitiveHandle<3> &pb,
OutputIteratorType out) -> OutputIteratorType
{
switch (pb.handle.which()) {
case PrimitivePoint:
*out++ = primitive;
break;
case PrimitiveSegment:
difference(primitive, *pb.as<Segment_3>(), out);
break;
case PrimitiveSurface:
difference(primitive, *pb.as<Triangle_3>(), out);
break;
case PrimitiveVolume:
difference(primitive, *pb.as<MarkedPolyhedron>(), out);
break;
}
return out;
}
template <typename OutputIteratorType>
auto
difference(const Triangle_3 &primitive, const PrimitiveHandle<3> &pb,
OutputIteratorType out) -> OutputIteratorType
{
switch (pb.handle.which()) {
case PrimitivePoint:
*out++ = primitive;
break;
case PrimitiveSegment:
*out++ = primitive;
break;
case PrimitiveSurface:
difference(primitive, *pb.as<Triangle_3>(), out);
break;
case PrimitiveVolume:
difference(primitive, *pb.as<MarkedPolyhedron>(), out);
break;
}
return out;
}
template <typename OutputIteratorType>
auto
difference(const MarkedPolyhedron &primitive, const PrimitiveHandle<3> &pb,
OutputIteratorType out) -> OutputIteratorType
{
switch (pb.handle.which()) {
case PrimitivePoint:
*out++ = primitive;
break;
case PrimitiveSegment:
*out++ = primitive;
break;
case PrimitiveSurface:
*out++ = primitive;
break;
case PrimitiveVolume:
difference(primitive, *pb.as<MarkedPolyhedron>(), out);
break;
}
return out;
}
template <typename Primitive, typename PrimitiveHandleConstIterator>
auto
difference(const Primitive &primitive, PrimitiveHandleConstIterator begin,
PrimitiveHandleConstIterator end) -> std::vector<Primitive>
{
std::vector<Primitive> primitives;
primitives.push_back(primitive);
for (PrimitiveHandleConstIterator b = begin; b != end; ++b) {
std::vector<Primitive> new_primitives;
for (auto a = primitives.begin(); a != primitives.end(); ++a) {
difference(*a, *(*b), std::back_inserter(new_primitives));
}
primitives.swap(new_primitives);
}
return primitives;
}
// just performs the type switch for the primitive to substract from
//
void
appendDifference(const PrimitiveHandle<2> &pa,
CollisionMapper<2>::PrimitiveHandleSet::const_iterator begin,
CollisionMapper<2>::PrimitiveHandleSet::const_iterator end,
GeometrySet<2> &output)
{
switch (pa.handle.which()) {
case PrimitivePoint: {
std::vector<Point_2> res = difference(*pa.as<Point_2>(), begin, end);
output.addPoints(res.begin(), res.end());
return;
}
case PrimitiveSegment: {
std::vector<Segment_2> res = difference(*pa.as<Segment_2>(), begin, end);
output.addSegments(res.begin(), res.end());
return;
}
case PrimitiveSurface: {
std::vector<PolygonWH_2> res =
difference(*pa.as<PolygonWH_2>(), begin, end);
output.addSurfaces(res.begin(), res.end());
return;
}
}
}
void
appendDifference(const PrimitiveHandle<3> &pa,
CollisionMapper<3>::PrimitiveHandleSet::const_iterator begin,
CollisionMapper<3>::PrimitiveHandleSet::const_iterator end,
GeometrySet<3> &output)
{
switch (pa.handle.which()) {
case PrimitivePoint: {
std::vector<Point_3> res = difference(*pa.as<Point_3>(), begin, end);
output.addPoints(res.begin(), res.end());
return;
}
case PrimitiveSegment: {
std::vector<Segment_3> res = difference(*pa.as<Segment_3>(), begin, end);
output.addSegments(res.begin(), res.end());
break;
}
case PrimitiveSurface: {
std::vector<Triangle_3> res = difference(*pa.as<Triangle_3>(), begin, end);
output.addSurfaces(res.begin(), res.end());
break;
}
case PrimitiveVolume: {
std::vector<MarkedPolyhedron> res =
difference(*pa.as<MarkedPolyhedron>(), begin, end);
output.addVolumes(res.begin(), res.end());
break;
}
}
}
void
post_difference(const GeometrySet<2> &input, GeometrySet<2> &output)
{
//
// reverse orientation of polygons if needed
for (const auto &it : input.surfaces()) {
const PolygonWH_2 &p = it.primitive();
Polygon_2 outer = p.outer_boundary();
if (outer.orientation() == CGAL::CLOCKWISE) {
outer.reverse_orientation();
}
std::list<CGAL::Polygon_2<Kernel>> rings;
for (auto hit = p.holes_begin(); hit != p.holes_end(); ++hit) {
rings.push_back(*hit);
if (hit->orientation() == CGAL::COUNTERCLOCKWISE) {
rings.back().reverse_orientation();
}
}
output.surfaces().emplace_back(
PolygonWH_2(outer, rings.begin(), rings.end()));
}
output.points() = input.points();
output.segments() = input.segments();
output.volumes() = input.volumes();
}
void
post_difference(const GeometrySet<3> &input, GeometrySet<3> &output)
{
// nothing special to do
output = input;
}
template <int Dim>
void
difference(const GeometrySet<Dim> &a, const GeometrySet<Dim> &b,
GeometrySet<Dim> &output)
{
using BoxCollection = typename SFCGAL::detail::BoxCollection<Dim>::Type;
typename SFCGAL::detail::HandleCollection<Dim>::Type ahandles;
typename SFCGAL::detail::HandleCollection<Dim>::Type bhandles;
BoxCollection aboxes;
BoxCollection bboxes;
a.computeBoundingBoxes(ahandles, aboxes);
b.computeBoundingBoxes(bhandles, bboxes);
// here we use box_intersection_d to build the list of operations
// that actually need to be performed
GeometrySet<Dim> temp;
GeometrySet<Dim> temp2;
typename CollisionMapper<Dim>::Map map;
CollisionMapper<Dim> const cb(map);
CGAL::box_intersection_d(aboxes.begin(), aboxes.end(), bboxes.begin(),
bboxes.end(), cb);
// now we have in cb a map of operations to perform
// we can put in the result right away all the primitives
// that are not keys in this map
{
auto ait = aboxes.begin();
const auto end = aboxes.end();
for (; ait != end; ++ait) {
if (map.find(ait->handle()) == map.end()) {
temp.addPrimitive(*ait->handle());
}
}
}
// then we delegate the operations according to type
{
auto cbit = map.begin();
const auto end = map.end();
for (; cbit != end; ++cbit) {
appendDifference(*cbit->first, cbit->second.begin(), cbit->second.end(),
temp);
}
}
post_difference(temp, temp2);
output.merge(temp2);
}
template void
difference<2>(const GeometrySet<2> &a, const GeometrySet<2> &b,
GeometrySet<2> &);
template void
difference<3>(const GeometrySet<3> &a, const GeometrySet<3> &b,
GeometrySet<3> &);
auto
difference(const Geometry &ga, const Geometry &gb, NoValidityCheck /*unused*/)
-> std::unique_ptr<Geometry>
{
GeometrySet<2> const gsa(ga);
GeometrySet<2> const gsb(gb);
GeometrySet<2> output;
algorithm::difference(gsa, gsb, output);
GeometrySet<2> filtered;
output.filterCovered(filtered);
return filtered.recompose();
}
auto
difference(const Geometry &ga, const Geometry &gb) -> std::unique_ptr<Geometry>
{
SFCGAL_ASSERT_GEOMETRY_VALIDITY_2D(ga);
SFCGAL_ASSERT_GEOMETRY_VALIDITY_2D(gb);
return difference(ga, gb, NoValidityCheck());
}
auto
difference3D(const Geometry &ga, const Geometry &gb, NoValidityCheck /*unused*/)
-> std::unique_ptr<Geometry>
{
GeometrySet<3> const gsa(ga);
GeometrySet<3> const gsb(gb);
GeometrySet<3> output;
algorithm::difference(gsa, gsb, output);
GeometrySet<3> filtered;
output.filterCovered(filtered);
return filtered.recompose();
}
auto
difference3D(const Geometry &ga, const Geometry &gb)
-> std::unique_ptr<Geometry>
{
SFCGAL_ASSERT_GEOMETRY_VALIDITY_3D(ga);
SFCGAL_ASSERT_GEOMETRY_VALIDITY_3D(gb);
return difference3D(ga, gb, NoValidityCheck());
}
} // namespace SFCGAL::algorithm