AbstractConvexHyperplaneBoundedRegion.java

  1. /*
  2.  * Licensed to the Apache Software Foundation (ASF) under one or more
  3.  * contributor license agreements.  See the NOTICE file distributed with
  4.  * this work for additional information regarding copyright ownership.
  5.  * The ASF licenses this file to You under the Apache License, Version 2.0
  6.  * (the "License"); you may not use this file except in compliance with
  7.  * the License.  You may obtain a copy of the License at
  8.  *
  9.  *      http://www.apache.org/licenses/LICENSE-2.0
  10.  *
  11.  * Unless required by applicable law or agreed to in writing, software
  12.  * distributed under the License is distributed on an "AS IS" BASIS,
  13.  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
  14.  * See the License for the specific language governing permissions and
  15.  * limitations under the License.
  16.  */
  17. package org.apache.commons.geometry.core.partitioning;

  18. import java.util.ArrayList;
  19. import java.util.Collections;
  20. import java.util.Iterator;
  21. import java.util.List;
  22. import java.util.function.Function;

  23. import org.apache.commons.geometry.core.Point;
  24. import org.apache.commons.geometry.core.RegionLocation;
  25. import org.apache.commons.geometry.core.Transform;

  26. /** Base class for convex hyperplane-bounded regions. This class provides generic implementations of many
  27.  * algorithms related to convex regions.
  28.  * @param <P> Point implementation type
  29.  * @param <S> Hyperplane convex subset implementation type
  30.  */
  31. public abstract class AbstractConvexHyperplaneBoundedRegion<P extends Point<P>, S extends HyperplaneConvexSubset<P>>
  32.     implements HyperplaneBoundedRegion<P> {
  33.     /** List of boundaries for the region. */
  34.     private final List<S> boundaries;

  35.     /** Simple constructor. Callers are responsible for ensuring that the given list of boundaries
  36.      * define a convex region. No validation is performed.
  37.      * @param boundaries the boundaries of the convex region
  38.      */
  39.     protected AbstractConvexHyperplaneBoundedRegion(final List<S> boundaries) {
  40.         this.boundaries = Collections.unmodifiableList(boundaries);
  41.     }

  42.     /** Get the boundaries of the convex region. The exact ordering of the boundaries
  43.      * is not guaranteed.
  44.      * @return the boundaries of the convex region
  45.      */
  46.     public List<S> getBoundaries() {
  47.         return boundaries;
  48.     }

  49.     /** {@inheritDoc} */
  50.     @Override
  51.     public boolean isFull() {
  52.         // no boundaries => no outside
  53.         return boundaries.isEmpty();
  54.     }

  55.     /** {@inheritDoc}
  56.      *
  57.      * <p>This method always returns false.</p>
  58.      */
  59.     @Override
  60.     public boolean isEmpty() {
  61.         return false;
  62.     }

  63.     /** {@inheritDoc} */
  64.     @Override
  65.     public double getBoundarySize() {
  66.         double sum = 0.0;
  67.         for (final S boundary : boundaries) {
  68.             sum += boundary.getSize();
  69.         }

  70.         return sum;
  71.     }

  72.     /** {@inheritDoc} */
  73.     @Override
  74.     public RegionLocation classify(final P pt) {
  75.         boolean isOn = false;

  76.         HyperplaneLocation loc;
  77.         for (final S boundary : boundaries) {
  78.             loc = boundary.getHyperplane().classify(pt);

  79.             if (loc == HyperplaneLocation.PLUS) {
  80.                 return RegionLocation.OUTSIDE;
  81.             } else if (loc == HyperplaneLocation.ON) {
  82.                 isOn = true;
  83.             }
  84.         }

  85.         return isOn ? RegionLocation.BOUNDARY : RegionLocation.INSIDE;
  86.     }

  87.     /** {@inheritDoc} */
  88.     @Override
  89.     public P project(final P pt) {

  90.         P projected;
  91.         double dist;

  92.         P closestPt = null;
  93.         double closestDist = Double.POSITIVE_INFINITY;

  94.         for (final S boundary : boundaries) {
  95.             projected = boundary.closest(pt);
  96.             dist = pt.distance(projected);

  97.             if (projected != null && (closestPt == null || dist < closestDist)) {
  98.                 closestPt = projected;
  99.                 closestDist = dist;
  100.             }
  101.         }

  102.         return closestPt;
  103.     }

  104.     /** Trim the given hyperplane subset to the portion contained inside this instance.
  105.      * @param sub hyperplane subset to trim. Null is returned if the subset does not intersect the instance.
  106.      * @return portion of the argument that lies entirely inside the region represented by
  107.      *      this instance, or null if it does not intersect.
  108.      */
  109.     public HyperplaneConvexSubset<P> trim(final HyperplaneConvexSubset<P> sub) {
  110.         HyperplaneConvexSubset<P> remaining = sub;
  111.         for (final S boundary : boundaries) {
  112.             remaining = remaining.split(boundary.getHyperplane()).getMinus();
  113.             if (remaining == null) {
  114.                 break;
  115.             }
  116.         }

  117.         return remaining;
  118.     }

  119.     /** {@inheritDoc} */
  120.     @Override
  121.     public String toString() {
  122.         final StringBuilder sb = new StringBuilder();
  123.         sb.append(this.getClass().getSimpleName())
  124.             .append("[boundaries= ")
  125.             .append(boundaries)
  126.             .append(']');

  127.         return sb.toString();
  128.     }

  129.     /** Generic, internal transform method. Subclasses should use this to implement their own transform methods.
  130.      * @param transform the transform to apply to the instance
  131.      * @param thisInstance a reference to the current instance; this is passed as
  132.      *      an argument in order to allow it to be a generic type
  133.      * @param boundaryType the type used for the boundary hyperplane subsets
  134.      * @param factory function used to create new convex region instances
  135.      * @param <R> Region implementation type
  136.      * @return the result of the transform operation
  137.      */
  138.     protected <R extends AbstractConvexHyperplaneBoundedRegion<P, S>> R transformInternal(
  139.             final Transform<P> transform, final R thisInstance, final Class<S> boundaryType,
  140.             final Function<? super List<S>, R> factory) {

  141.         if (isFull()) {
  142.             return thisInstance;
  143.         }

  144.         final List<S> origBoundaries = getBoundaries();

  145.         final int size = origBoundaries.size();
  146.         final List<S> tBoundaries = new ArrayList<>(size);

  147.         // determine if the hyperplanes should be reversed
  148.         final S boundary = origBoundaries.get(0);
  149.         HyperplaneConvexSubset<P> tBoundary = boundary.transform(transform);

  150.         final boolean reverseDirection = swapsInsideOutside(transform);

  151.         // transform all of the segments
  152.         if (reverseDirection) {
  153.             tBoundary = tBoundary.reverse();
  154.         }
  155.         tBoundaries.add(boundaryType.cast(tBoundary));

  156.         for (int i = 1; i < origBoundaries.size(); ++i) {
  157.             tBoundary = origBoundaries.get(i).transform(transform);

  158.             if (reverseDirection) {
  159.                 tBoundary = tBoundary.reverse();
  160.             }

  161.             tBoundaries.add(boundaryType.cast(tBoundary));
  162.         }

  163.         return factory.apply(tBoundaries);
  164.     }

  165.     /** Return true if the given transform swaps the inside and outside of
  166.      * the region.
  167.      *
  168.      * <p>The default behavior of this method is to return true if the transform
  169.      * does not preserve spatial orientation (ie, {@link Transform#preservesOrientation()}
  170.      * is false). Subclasses may need to override this method to implement the correct
  171.      * behavior for their space and dimension.</p>
  172.      * @param transform transform to check
  173.      * @return true if the given transform swaps the interior and exterior of
  174.      *      the region
  175.      */
  176.     protected boolean swapsInsideOutside(final Transform<P> transform) {
  177.         return !transform.preservesOrientation();
  178.     }

  179.     /** Generic, internal split method. Subclasses should call this from their {@link #split(Hyperplane)} methods.
  180.      * @param splitter splitting hyperplane
  181.      * @param thisInstance a reference to the current instance; this is passed as
  182.      *      an argument in order to allow it to be a generic type
  183.      * @param boundaryType the type used for the boundary hyperplane subsets
  184.      * @param factory function used to create new convex region instances
  185.      * @param <R> Region implementation type
  186.      * @return the result of the split operation
  187.      */
  188.     protected <R extends AbstractConvexHyperplaneBoundedRegion<P, S>> Split<R> splitInternal(
  189.             final Hyperplane<P> splitter, final R thisInstance, final Class<S> boundaryType,
  190.             final Function<List<S>, R> factory) {

  191.         return isFull() ?
  192.                 splitInternalFull(splitter, boundaryType, factory) :
  193.                 splitInternalNonFull(splitter, thisInstance, boundaryType, factory);
  194.     }

  195.     /** Internal split method for use with full regions, i.e. regions that cover the entire space.
  196.      * @param splitter splitting hyperplane
  197.      * @param boundaryType the type used for the boundary hyperplane subsets
  198.      * @param factory function used to create new convex region instances
  199.      * @param <R> Region implementation type
  200.      * @return the result of the split operation
  201.      */
  202.     private <R extends AbstractConvexHyperplaneBoundedRegion<P, S>> Split<R> splitInternalFull(
  203.             final Hyperplane<P> splitter, final Class<S> boundaryType, final Function<? super List<S>, R> factory) {

  204.         final R minus = factory.apply(Collections.singletonList(boundaryType.cast(splitter.span())));
  205.         final R plus = factory.apply(Collections.singletonList(boundaryType.cast(splitter.reverse().span())));

  206.         return new Split<>(minus, plus);
  207.     }

  208.     /** Internal split method for use with non-full regions, i.e. regions that do not cover the entire space.
  209.      * @param splitter splitting hyperplane
  210.      * @param thisInstance a reference to the current instance; this is passed as
  211.      *      an argument in order to allow it to be a generic type
  212.      * @param boundaryType the type used for the boundary hyperplane subsets
  213.      * @param factory function used to create new convex region instances
  214.      * @param <R> Region implementation type
  215.      * @return the result of the split operation
  216.      */
  217.     private <R extends AbstractConvexHyperplaneBoundedRegion<P, S>> Split<R> splitInternalNonFull(
  218.             final Hyperplane<P> splitter, final R thisInstance, final Class<S> boundaryType,
  219.             final Function<? super List<S>, R> factory) {

  220.         final HyperplaneConvexSubset<P> trimmedSplitter = trim(splitter.span());

  221.         if (trimmedSplitter == null) {
  222.             // The splitter lies entirely outside of the region; we need
  223.             // to determine whether we lie on the plus or minus side of the splitter.

  224.             final SplitLocation regionLoc = determineRegionPlusMinusLocation(splitter);
  225.             return regionLoc == SplitLocation.MINUS ?
  226.                     new Split<>(thisInstance, null) :
  227.                     new Split<>(null, thisInstance);
  228.         }

  229.         // the splitter passes through the region; split the other region boundaries
  230.         // by the splitter
  231.         final ArrayList<S> minusBoundaries = new ArrayList<>();
  232.         final ArrayList<S> plusBoundaries = new ArrayList<>();

  233.         splitBoundaries(splitter, boundaryType, minusBoundaries, plusBoundaries);

  234.         // if the splitter was trimmed by the region boundaries, double-check that the split boundaries
  235.         // actually lie on both sides of the splitter; this is another case where floating point errors
  236.         // can cause a discrepancy between the results of splitting the splitter by the boundaries and
  237.         // splitting the boundaries by the splitter
  238.         if (!trimmedSplitter.isFull()) {
  239.             if (minusBoundaries.isEmpty()) {
  240.                 if (plusBoundaries.isEmpty()) {
  241.                     return new Split<>(null, null);
  242.                 }
  243.                 return new Split<>(null, thisInstance);
  244.             } else if (plusBoundaries.isEmpty()) {
  245.                 return new Split<>(thisInstance, null);
  246.             }
  247.         }

  248.         // we have a consistent region split; create the new plus and minus regions
  249.         minusBoundaries.add(boundaryType.cast(trimmedSplitter));
  250.         plusBoundaries.add(boundaryType.cast(trimmedSplitter.reverse()));

  251.         minusBoundaries.trimToSize();
  252.         plusBoundaries.trimToSize();

  253.         return new Split<>(factory.apply(minusBoundaries), factory.apply(plusBoundaries));
  254.     }

  255.     /** Determine whether the region lies on the plus or minus side of the given splitter. It is assumed
  256.      * that (1) the region is not full, and (2) the given splitter does not pass through the region.
  257.      *
  258.      * <p>In theory, this is a very simple operation: one need only test a single region boundary
  259.      * to see if it lies on the plus or minus side of the splitter. In practice, however, accumulated
  260.      * floating point errors can cause discrepancies between the splitting operations, causing
  261.      * boundaries to be classified as lying on both sides of the splitter when they should only lie on one.
  262.      * Therefore, this method examines as many boundaries as needed in order to determine the best response.
  263.      * The algorithm proceeds as follows:
  264.      * <ol>
  265.      *  <li>If any boundary lies completely on the minus or plus side of the splitter, then
  266.      *      {@link SplitLocation#MINUS MINUS} or {@link SplitLocation#PLUS PLUS} is returned, respectively.</li>
  267.      *  <li>If any boundary is coincident with the splitter ({@link SplitLocation#NEITHER NEITHER}), then
  268.      *      {@link SplitLocation#MINUS MINUS} is returned if the boundary hyperplane has the same orientation
  269.      *      as the splitter, otherwise {@link SplitLocation#PLUS PLUS}.</li>
  270.      *  <li>If no boundaries match the above conditions, then the sizes of the split boundaries are compared. If
  271.      *      the sum of the sizes of the boundaries on the minus side is greater than the sum of the sizes of
  272.      *      the boundaries on the plus size, then {@link SplitLocation#MINUS MINUS} is returned. Otherwise,
  273.      *      {@link SplitLocation#PLUS PLUS} is returned.
  274.      * </ol>
  275.      * @param splitter splitter to classify the region against; the splitter is assumed to lie
  276.      *      completely outside of the region
  277.      * @return {@link SplitLocation#MINUS} if the region lies on the minus side of the splitter and
  278.      *      {@link SplitLocation#PLUS} if the region lies on the plus side of the splitter
  279.      */
  280.     private SplitLocation determineRegionPlusMinusLocation(final Hyperplane<P> splitter) {
  281.         double minusSize = 0;
  282.         double plusSize = 0;

  283.         Split<? extends HyperplaneConvexSubset<P>> split;
  284.         SplitLocation loc;

  285.         for (final S boundary : boundaries) {
  286.             split = boundary.split(splitter);
  287.             loc = split.getLocation();

  288.             if (loc == SplitLocation.MINUS || loc == SplitLocation.PLUS) {
  289.                 return loc;
  290.             } else if (loc == SplitLocation.NEITHER) {
  291.                 return splitter.similarOrientation(boundary.getHyperplane()) ?
  292.                         SplitLocation.MINUS :
  293.                         SplitLocation.PLUS;
  294.             } else {
  295.                 minusSize += split.getMinus().getSize();
  296.                 plusSize += split.getPlus().getSize();
  297.             }
  298.         }

  299.         return minusSize > plusSize ? SplitLocation.MINUS : SplitLocation.PLUS;
  300.     }

  301.     /** Split the boundaries of the region by the given hyperplane, adding the split parts into the
  302.      * corresponding lists.
  303.      * @param splitter splitting hyperplane
  304.      * @param boundaryType the type used for the boundary hyperplane subsets
  305.      * @param minusBoundaries list that will contain the portions of the boundaries on the minus side
  306.      *      of the splitting hyperplane
  307.      * @param plusBoundaries list that will contain the portions of the boundaries on the plus side of
  308.      *      the splitting hyperplane
  309.      */
  310.     private void splitBoundaries(final Hyperplane<P> splitter, final Class<S> boundaryType,
  311.             final List<S> minusBoundaries, final List<S> plusBoundaries) {

  312.         Split<? extends HyperplaneConvexSubset<P>> split;
  313.         HyperplaneConvexSubset<P> minusBoundary;
  314.         HyperplaneConvexSubset<P> plusBoundary;

  315.         for (final S boundary : boundaries) {
  316.             split = boundary.split(splitter);

  317.             minusBoundary = split.getMinus();
  318.             plusBoundary = split.getPlus();

  319.             if (minusBoundary != null) {
  320.                 minusBoundaries.add(boundaryType.cast(minusBoundary));
  321.             }

  322.             if (plusBoundary != null) {
  323.                 plusBoundaries.add(boundaryType.cast(plusBoundary));
  324.             }
  325.         }
  326.     }

  327.     /** Internal class encapsulating the logic for building convex region boundaries from collections of hyperplanes.
  328.      * @param <P> Point implementation type
  329.      * @param <S> Hyperplane convex subset implementation type
  330.      */
  331.     protected static class ConvexRegionBoundaryBuilder<P extends Point<P>, S extends HyperplaneConvexSubset<P>> {

  332.         /** Hyperplane convex subset implementation type. */
  333.         private final Class<S> subsetType;

  334.         /** Construct a new instance for building convex region boundaries with the given hyperplane
  335.          * convex subset implementation type.
  336.          * @param subsetType Hyperplane convex subset implementation type
  337.          */
  338.         public ConvexRegionBoundaryBuilder(final Class<S> subsetType) {
  339.             this.subsetType = subsetType;
  340.         }

  341.         /** Compute a list of hyperplane convex subsets representing the boundaries of the convex region
  342.          * bounded by the given collection of hyperplanes.
  343.          * @param bounds hyperplanes defining the convex region
  344.          * @return a list of hyperplane convex subsets representing the boundaries of the convex region
  345.          * @throws IllegalArgumentException if the given hyperplanes do not form a convex region
  346.          */
  347.         public List<S> build(final Iterable<? extends Hyperplane<P>> bounds) {

  348.             final List<S> boundaries = new ArrayList<>();

  349.             // cut each hyperplane by every other hyperplane in order to get the region boundaries
  350.             int boundIdx = -1;
  351.             HyperplaneConvexSubset<P> boundary;

  352.             for (final Hyperplane<P> currentBound : bounds) {
  353.                 ++boundIdx;

  354.                 boundary = splitBound(currentBound, bounds, boundIdx);
  355.                 if (boundary != null) {
  356.                     boundaries.add(subsetType.cast(boundary));
  357.                 }
  358.             }

  359.             if (boundIdx > 0 && boundaries.isEmpty()) {
  360.                 // nothing was added
  361.                 throw nonConvexException(bounds);
  362.             }

  363.             return boundaries;
  364.         }

  365.         /** Split the given bounding hyperplane by all of the other hyperplanes in the given collection, returning the
  366.          * remaining hyperplane subset.
  367.          * @param currentBound the bound to split; this value is assumed to have come from {@code bounds}
  368.          * @param bounds collection of bounds to use to split {@code currentBound}
  369.          * @param currentBoundIdx the index of {@code currentBound} in {@code bounds}
  370.          * @return the part of {@code currentBound}'s hyperplane subset that lies on the minus side of all of the
  371.          *      splitting hyperplanes
  372.          * @throws IllegalArgumentException if the hyperplanes do not form a convex region
  373.          */
  374.         private HyperplaneConvexSubset<P> splitBound(final Hyperplane<P> currentBound,
  375.                 final Iterable<? extends Hyperplane<P>> bounds, final int currentBoundIdx) {

  376.             HyperplaneConvexSubset<P> boundary = currentBound.span();

  377.             final Iterator<? extends Hyperplane<P>> boundsIt = bounds.iterator();

  378.             Hyperplane<P> splitter;
  379.             int splitterIdx = -1;

  380.             while (boundsIt.hasNext() && boundary != null) {
  381.                 splitter = boundsIt.next();
  382.                 ++splitterIdx;

  383.                 if (currentBound == splitter) {
  384.                     // do not split the bound with itself

  385.                     if (currentBoundIdx > splitterIdx) {
  386.                         // this hyperplane is duplicated in the list; skip all but the
  387.                         // first insertion of its hyperplane subset
  388.                         return null;
  389.                     }
  390.                 } else {
  391.                     // split the boundary
  392.                     final Split<? extends HyperplaneConvexSubset<P>> split = boundary.split(splitter);

  393.                     if (split.getLocation() != SplitLocation.NEITHER) {
  394.                         // retain the minus portion of the split
  395.                         boundary = split.getMinus();
  396.                     } else if (!currentBound.similarOrientation(splitter)) {
  397.                         // two or more splitters are coincident and have opposite
  398.                         // orientations, meaning that no area is on the minus side
  399.                         // of both
  400.                         throw nonConvexException(bounds);
  401.                     } else if (currentBoundIdx > splitterIdx) {
  402.                         // two or more hyperplanes are equivalent; only use the boundary
  403.                         // from the first one and return null for this one
  404.                         return null;
  405.                     }
  406.                 }
  407.             }

  408.             return boundary;
  409.         }

  410.         /** Return an exception indicating that the given collection of hyperplanes do not produce a convex region.
  411.          * @param bounds collection of hyperplanes
  412.          * @return an exception indicating that the given collection of hyperplanes do not produce a convex region
  413.          */
  414.         private IllegalArgumentException nonConvexException(final Iterable<? extends Hyperplane<P>> bounds) {
  415.             return new IllegalArgumentException("Bounding hyperplanes do not produce a convex region: " + bounds);
  416.         }
  417.     }
  418. }