predicates_helpers.hpp

Go to the documentation of this file.

#pragma once
 
#include "implementation/duality.hpp"

 
 
namespace pgl {
 
namespace detail {

template <class ANumber, class BNumber>
using Exact1DNumber = std::conditional_t<
    std::is_floating_point_v<ANumber> || std::is_floating_point_v<BNumber>,
    double, Rational<BigInt>>;

 
template <class T>
struct is_area_cut_target : std::false_type {};
 
template <class PointType, class Label>
struct is_area_cut_target<Rectangle<PointType, Label>> : std::true_type {};
 
template <class PointType, class Label>
struct is_area_cut_target<Triangle<PointType, Label>> : std::true_type {};
 
template <class T>
inline constexpr bool is_area_cut_target_v = is_area_cut_target<std::remove_cvref_t<T>>::value;

 
template <class Polygon, class OtherSegment>
constexpr bool polygonBoundaryContainsSegment(const Polygon& polygon, const OtherSegment& other) {
    const auto boundary = polygon.edges();
    for (const auto& edge : boundary) {
        if (edge.contains(other)) {
            return true;
        }
    }
    return false;
}

template <class RectangleType, class FirstPoint, class SecondPoint>
constexpr bool lineIntersectsRectangle(const RectangleType& rectangle, const FirstPoint& first, const SecondPoint& second) {
    // The line's two points appear in all four signs, so filtering them once
    // converts each of their coordinates a single time instead of four.
    using Coordinate = sign_coordinate_t<typename FirstPoint::NumberType,
                                         typename SecondPoint::NumberType,
                                         typename RectangleType::NumberType>;
    const auto a = filtered<Coordinate>(first);
    const auto b = filtered<Coordinate>(second);
    bool has_positive = false;
    bool has_negative = false;
    const auto vertices = rectangle.vertices();
    for (const auto& vertex : vertices) {
        const auto side = orientationSignOf(a, b, filtered<Coordinate>(vertex)).value();
        if (side == std::partial_ordering::equivalent) {
            return true;
        }
        has_positive = has_positive || side == std::partial_ordering::greater;
        has_negative = has_negative || side == std::partial_ordering::less;
    }
    return has_positive && has_negative;
}

template <class RectangleType, class FirstPoint, class SecondPoint>
constexpr bool segmentIntersectsRectangle(const RectangleType& rectangle, const FirstPoint& first, const SecondPoint& second) {
    const auto& low = rectangle.min();
    const auto& high = rectangle.max();
 
    const bool rising_x = first.x() < second.x();
    if ((rising_x ? first.x() : second.x()) > high.x() ||
        (rising_x ? second.x() : first.x()) < low.x()) {
        return false;
    }
    const bool rising_y = first.y() < second.y();
    if ((rising_y ? first.y() : second.y()) > high.y() ||
        (rising_y ? second.y() : first.y()) < low.y()) {
        return false;
    }
    return lineIntersectsRectangle(rectangle, first, second);
}

template <class RectangleType, class FirstPoint, class SecondPoint>
constexpr bool lineIntersectsRectangleInterior(const RectangleType& rectangle, const FirstPoint& first, const SecondPoint& second) {
    bool has_positive = false;
    bool has_negative = false;
    const auto vertices = rectangle.vertices();
    for (const auto& vertex : vertices) {
        const auto side = orientationSign(first, second, vertex);
        has_positive = has_positive || side == std::partial_ordering::greater;
        has_negative = has_negative || side == std::partial_ordering::less;
    }
    return has_positive && has_negative;
}

template <class RectangleType, class FirstPoint, class SecondPoint>
constexpr bool segmentIntersectsRectangleInteriorExact(const RectangleType& rectangle, const FirstPoint& first, const SecondPoint& second) {
    using Coordinate = promoted_number_t<std::common_type_t<
        std::remove_cvref_t<decltype(rectangle.min().x())>,
        std::remove_cvref_t<decltype(first.x())>,
        std::remove_cvref_t<decltype(second.x())>>>;
    using SegmentPoint = Point<Coordinate>;
    using TestSegment = Segment<SegmentPoint>;
 
    const TestSegment segment{SegmentPoint(first), SegmentPoint(second)};
    if (segment.isDegenerate()) {
        return false;
    }
 
    if (rectangle.interiorContains(first) || rectangle.interiorContains(second)) {
        return true;
    }
 
    if (!lineIntersectsRectangleInterior(rectangle, first, second)) {
        return false;
    }
 
    const auto edges = rectangle.edges();
    std::array<bool, 4> edge_hits{};
    int distinct_boundary_contacts = 0;
    for (std::size_t i = 0; i < edges.size(); ++i) {
        edge_hits[i] = segment.intersects(edges[i]);
        distinct_boundary_contacts += edge_hits[i] ? 1 : 0;
    }
 
    if (distinct_boundary_contacts < 2) {
        return false;
    }
 
    const auto vertices = rectangle.vertices();
    if (segment.contains(vertices[0]) && edge_hits[0] && edge_hits[3]) {
        --distinct_boundary_contacts;
    }
    if (segment.contains(vertices[1]) && edge_hits[0] && edge_hits[1]) {
        --distinct_boundary_contacts;
    }
    if (segment.contains(vertices[2]) && edge_hits[1] && edge_hits[2]) {
        --distinct_boundary_contacts;
    }
    if (segment.contains(vertices[3]) && edge_hits[2] && edge_hits[3]) {
        --distinct_boundary_contacts;
    }
 
    return distinct_boundary_contacts >= 2;
}

template <class Remover, MonotoneChainConcept Chain>
constexpr bool separatesChain(const Remover& remover, const Chain& chain) {
    const std::size_t n = chain.size();
    if (n < 2) {
        // Removing anything from at most one point cannot disconnect it.
        return false;
    }
    bool active = false;     // a component is in progress
    bool touched = false;    // ... and it touches an extreme vertex
    for (std::size_t i = 0; i + 1 < n; ++i) {
        const Segment<typename Chain::PointType> edge(chain[i], chain[i + 1]);
        if (edge.intersects(remover)) {
            const bool connected = active && remover.contains(chain[i]);
            if (!connected) {
                if (active && !touched) {
                    return true;
                }
                touched = false;
            }
            active = true;
            if (i == 0 && remover.contains(chain[0])) {
                touched = true;
            }
            if (i + 2 == n && remover.contains(chain[n - 1])) {
                touched = true;
            }
        } else {
            if (active && !touched) {
                return true;
            }
            active = false;
            touched = false;
        }
    }
    return active && !touched;
}


template <class Number>
using region_exact_number_t =
    std::conditional_t<std::is_floating_point_v<Number>, double, Rational<BigInt>>;

template <class Region, HalfplaneConcept OtherHalfplane>
constexpr bool regionInsideHalfplane(const Region& region, const OtherHalfplane& halfplane) {
    std::remove_cvref_t<Region> copy(region);
    return !copy.insert(halfplane);
}

template <class Region, HalfplaneConcept OtherHalfplane>
constexpr bool regionInsideHalfplaneInterior(const Region& region, const OtherHalfplane& halfplane) {
    return regionInsideHalfplane(region, halfplane) && !region.intersects(halfplane.asLine());
}

template <class Region, PointConcept FirstPoint, PointConcept SecondPoint>
constexpr bool regionStrictlyOnBothSides(const Region& region, const FirstPoint& first, const SecondPoint& second) {
    const Halfplane<std::remove_cvref_t<FirstPoint>> left(first, second);
    return !regionInsideHalfplane(region, left) && !regionInsideHalfplane(region, left.opposite());
}

template <class Region>
constexpr auto degenerateRegionCarrier(const Region& region) {
    using E = region_exact_number_t<typename Region::NumberType>;
    using EPoint = Point<E, typename Region::PointType::LabelType>;
    using Carrier = std::variant<EPoint, Segment<EPoint>, Ray<EPoint>, Line<EPoint>>;
    if (region.isBounded()) {
        // A point or a segment: the convex hull of the (partly coincident)
        // implicit vertices.
        const auto verts = region.template vertices<E>();
        EPoint lo = verts[0];
        EPoint hi = verts[0];
        for (const auto& v : verts) {
            if (v < lo) {
                lo = v;
            }
            if (hi < v) {
                hi = v;
            }
        }
        if (lo == hi) {
            return Carrier(lo);
        }
        return Carrier(Segment<EPoint>(lo, hi));
    }
    // An unbounded degenerate region is a ray or a line, and then some
    // boundary edge equals the whole region.
    for (std::size_t i = 0; i < region.size(); ++i) {
        const auto e = region.template edge<E>(i);
        if (const auto* ray = std::get_if<Ray<EPoint>>(&e)) {
            return Carrier(*ray);
        }
        if (const auto* line = std::get_if<Line<EPoint>>(&e)) {
            return Carrier(*line);
        }
    }
    // Unreachable for a degenerate region; keep a deterministic fallback.
    return Carrier(EPoint(region[0].source()));
}

template <class Region, RectangleConcept BoundingRectangle>
constexpr std::remove_cvref_t<Region> regionClippedToBox(const Region& region, const BoundingRectangle& bounds) {
    using Result = std::remove_cvref_t<Region>;
    using RegionPoint = typename Result::PointType;
    using N = typename Result::NumberType;
    using RegionHalfplane = typename Result::HalfplaneType;
    const N margin(2);
    // A fractional bound is rounded away from the box before the margin is
    // added, so the corners are whole numbers whatever N is. The box only has
    // to contain @p bounds, and rounding outward keeps that; what it buys is
    // the depth of everything downstream. Clipping against a corner that is
    // itself a deep fraction — a disk's bounding box under rational
    // coordinates carries twelve-digit numerators over eight-digit
    // denominators — hands every clipped vertex those denominators, and the
    // degree-four predicates over them then run on numbers with hundreds of
    // digits. An integral corner leaves the clipped vertices as shallow as the
    // region's own. For an integral N this is what the toward-zero cast
    // already did.
    const auto roundOutward = [](const auto& value, bool up) -> N {
        if constexpr (pgl::is_Rational_v<std::remove_cvref_t<decltype(value)>>) {
            // Both parts are wanted, so reduce once and read them off that copy:
            // numerator() and denominator() each run their own gcd otherwise.
            const auto reduced = value.simplified();
            const auto numerator = reduced.numerator();
            const auto denominator = reduced.denominator();  // always positive
            auto quotient = numerator / denominator;        // truncates toward zero
            const bool exact = quotient * denominator == numerator;
            if (!exact && (numerator < 0) == !up) {
                quotient = up ? quotient + 1 : quotient - 1;
            }
            return N(quotient);
        } else {
            return N(value);
        }
    };
    const RegionPoint lo(roundOutward(bounds.min().x(), false) - margin,
                         roundOutward(bounds.min().y(), false) - margin);
    const RegionPoint hi(roundOutward(bounds.max().x(), true) + margin,
                         roundOutward(bounds.max().y(), true) + margin);
    const RegionPoint lohi(lo.x(), hi.y());
    const RegionPoint hilo(hi.x(), lo.y());
    Result result(region);
    result.insert(RegionHalfplane(lo, hilo));
    result.insert(RegionHalfplane(hilo, hi));
    result.insert(RegionHalfplane(hi, lohi));
    result.insert(RegionHalfplane(lohi, lo));
    return result;
}

template <class HoledRegion, class EdgePredicate>
constexpr bool everyHoledRegionEdge(const HoledRegion& region, EdgePredicate&& predicate) {
    for (const auto& edge : region.outer().edgesView()) {
        if (!predicate(edge)) {
            return false;
        }
    }
    for (const auto& hole : region.holes()) {
        for (const auto& edge : hole.edgesView()) {
            if (!predicate(edge)) {
                return false;
            }
        }
    }
    return true;
}
 
}  // namespace detail
 
}  // namespace pgl