View Javadoc

    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.euclidean.threed;
    18
    19  import java.text.MessageFormat;
    20  import java.util.ArrayList;
    21  import java.util.Arrays;
    22  import java.util.Collection;
    23  import java.util.List;
    24  import java.util.function.BiFunction;
    25
    26  import org.apache.commons.geometry.core.partitioning.HyperplaneBoundedRegion;
    27  import org.apache.commons.geometry.core.partitioning.Split;
    28  import org.apache.commons.geometry.core.partitioning.SplitLocation;
    29  import org.apache.commons.geometry.euclidean.internal.EuclideanUtils;
    30  import org.apache.commons.geometry.euclidean.threed.line.Line3D;
    31  import org.apache.commons.geometry.euclidean.threed.line.LineConvexSubset3D;
    32  import org.apache.commons.geometry.euclidean.twod.ConvexArea;
    33  import org.apache.commons.geometry.euclidean.twod.Line;
    34  import org.apache.commons.geometry.euclidean.twod.LineConvexSubset;
    35  import org.apache.commons.geometry.euclidean.twod.Lines;
    36  import org.apache.commons.geometry.euclidean.twod.RegionBSPTree2D;
    37  import org.apache.commons.geometry.euclidean.twod.Vector2D;
    38  import org.apache.commons.geometry.euclidean.twod.path.LinePath;
    39  import org.apache.commons.numbers.core.Precision;
    40
    41  /** Class containing factory methods for constructing {@link Plane} and {@link PlaneSubset} instances.
    42   */
    43  public final class Planes {
    44
    45      /** Utility class; no instantiation. */
    46      private Planes() {
    47      }
    48
    49      /** Build a plane from a point and two (on plane) vectors.
    50       * @param p the provided point (on plane)
    51       * @param u u vector (on plane)
    52       * @param v v vector (on plane)
    53       * @param precision precision context used to compare floating point values
    54       * @return a new plane
    55       * @throws IllegalArgumentException if the norm of the given values is zero, NaN, or infinite.
    56       */
    57      public static EmbeddingPlane fromPointAndPlaneVectors(final Vector3D p, final Vector3D u, final Vector3D v,
    58              final Precision.DoubleEquivalence precision) {
    59          final Vector3D.Unit uNorm = u.normalize();
    60          final Vector3D.Unit vNorm = uNorm.orthogonal(v);
    61          final Vector3D.Unit wNorm = uNorm.cross(vNorm).normalize();
    62          final double originOffset = -p.dot(wNorm);
    63
    64          return new EmbeddingPlane(uNorm, vNorm, wNorm, originOffset, precision);
    65      }
    66
    67      /** Build a plane from a normal.
    68       * Chooses origin as point on plane.
    69       * @param normal normal direction to the plane
    70       * @param precision precision context used to compare floating point values
    71       * @return a new plane
    72       * @throws IllegalArgumentException if the norm of the given values is zero, NaN, or infinite.
    73       */
    74      public static Plane fromNormal(final Vector3D normal, final Precision.DoubleEquivalence precision) {
    75          return fromPointAndNormal(Vector3D.ZERO, normal, precision);
    76      }
    77
    78      /** Build a plane from a point and a normal.
    79       *
    80       * @param p point belonging to the plane
    81       * @param normal normal direction to the plane
    82       * @param precision precision context used to compare floating point values
    83       * @return a new plane
    84       * @throws IllegalArgumentException if the norm of the given values is zero, NaN, or infinite.
    85       */
    86      public static Plane fromPointAndNormal(final Vector3D p, final Vector3D normal,
    87              final Precision.DoubleEquivalence precision) {
    88          final Vector3D.Unit unitNormal = normal.normalize();
    89          final double originOffset = -p.dot(unitNormal);
    90
    91          return new Plane(unitNormal, originOffset, precision);
    92      }
    93
    94      /** Build a plane from three points.
    95       * <p>
    96       * The plane is oriented in the direction of {@code (p2-p1) ^ (p3-p1)}
    97       * </p>
    98       *
    99       * @param p1 first point belonging to the plane
    100      * @param p2 second point belonging to the plane
    101      * @param p3 third point belonging to the plane
    102      * @param precision precision context used to compare floating point values
    103      * @return a new plane
    104      * @throws IllegalArgumentException if the points do not define a unique plane
    105      */
    106     public static Plane fromPoints(final Vector3D p1, final Vector3D p2, final Vector3D p3,
    107             final Precision.DoubleEquivalence precision) {
    108         return fromPoints(Arrays.asList(p1, p2, p3), precision);
    109     }
    110
    111     /** Construct a plane from a collection of points lying on the plane. The plane orientation is
    112      * determined by the overall orientation of the point sequence. For example, if the points wind
    113      * around the z-axis in a counter-clockwise direction, then the plane normal will point up the
    114      * +z axis. If the points wind in the opposite direction, then the plane normal will point down
    115      * the -z axis. The {@code u} vector for the plane is set to the first non-zero vector between
    116      * points in the sequence (ie, the first direction in the path).
    117      *
    118      * @param pts collection of sequenced points lying on the plane
    119      * @param precision precision context used to compare floating point values
    120      * @return a new plane containing the given points
    121      * @throws IllegalArgumentException if the given collection does not contain at least 3 points or the
    122      *      points do not define a unique plane
    123      */
    124     public static Plane fromPoints(final Collection<Vector3D> pts, final Precision.DoubleEquivalence precision) {
    125         return new PlaneBuilder(pts, precision).build();
    126     }
    127
    128     /** Create a new plane subset from a plane and an embedded convex subspace area.
    129      * @param plane embedding plane for the area
    130      * @param area area embedded in the plane
    131      * @return a new convex sub plane instance
    132      */
    133     public static PlaneConvexSubset subsetFromConvexArea(final EmbeddingPlane plane, final ConvexArea area) {
    134         if (area.isFinite()) {
    135             // prefer a vertex-based representation for finite areas
    136             final List<Vector3D> vertices = plane.toSpace(area.getVertices());
    137             return fromConvexPlanarVertices(plane, vertices);
    138         }
    139
    140         return new EmbeddedAreaPlaneConvexSubset(plane, area);
    141     }
    142
    143     /** Create a new convex polygon from the given sequence of vertices. The vertices must define a unique
    144      * plane, meaning that at least 3 unique vertices must be given. The given sequence is assumed to be closed,
    145      * ie that an edge exists between the last vertex and the first.
    146      * @param pts collection of points defining the convex polygon
    147      * @param precision precision context used to compare floating point values
    148      * @return a new convex polygon defined by the given sequence of vertices
    149      * @throws IllegalArgumentException if fewer than 3 vertices are given or the vertices do not define a
    150      *       unique plane
    151      * @see #fromPoints(Collection, Precision.DoubleEquivalence)
    152      */
    153     public static ConvexPolygon3D convexPolygonFromVertices(final Collection<Vector3D> pts,
    154             final Precision.DoubleEquivalence precision) {
    155         final List<Vector3D> vertices = new ArrayList<>(pts.size());
    156         final Plane plane = new PlaneBuilder(pts, precision).buildForConvexPolygon(vertices);
    157
    158         // make sure that the first point is not repeated at the end
    159         final Vector3D firstPt = vertices.get(0);
    160         final Vector3D lastPt = vertices.get(vertices.size() - 1);
    161         if (firstPt.eq(lastPt, precision)) {
    162             vertices.remove(vertices.size() - 1);
    163         }
    164
    165         if (vertices.size() == EuclideanUtils.TRIANGLE_VERTEX_COUNT) {
    166             return new SimpleTriangle3D(plane, vertices.get(0), vertices.get(1), vertices.get(2));
    167         }
    168         return new VertexListConvexPolygon3D(plane, vertices);
    169     }
    170
    171     /** Construct a triangle from three vertices. The triangle plane is oriented such that the points
    172      * are arranged in a counter-clockwise order when looking down the plane normal.
    173      * @param p1 first vertex
    174      * @param p2 second vertex
    175      * @param p3 third vertex
    176      * @param precision precision context used for floating point comparisons
    177      * @return a triangle constructed from the three vertices
    178      * @throws IllegalArgumentException if the points do not define a unique plane
    179      */
    180     public static Triangle3D triangleFromVertices(final Vector3D p1, final Vector3D p2, final Vector3D p3,
    181             final Precision.DoubleEquivalence precision) {
    182         final Plane plane = fromPoints(p1, p2, p3, precision);
    183         return new SimpleTriangle3D(plane, p1, p2, p3);
    184     }
    185
    186     /** Construct a list of {@link Triangle3D} instances from a set of vertices and arrays of face indices.
    187      * For example, the following code constructs a list of triangles forming a square pyramid.
    188      * <pre>
    189      * Precision.DoubleEquivalence precision = Precision.doubleEquivalenceOfEpsilon(1e-10);
    190      *
    191      * Vector3D[] vertices = {
    192      *      Vector3D.ZERO,
    193      *      Vector3D.of(1, 0, 0),
    194      *      Vector3D.of(1, 1, 0),
    195      *      Vector3D.of(0, 1, 0),
    196      *      Vector3D.of(0.5, 0.5, 4)
    197      * };
    198      *
    199      * int[][] faceIndices = {
    200      *      {0, 2, 1},
    201      *      {0, 3, 2},
    202      *      {0, 1, 4},
    203      *      {1, 2, 4},
    204      *      {2, 3, 4},
    205      *      {3, 0, 4}
    206      * };
    207      *
    208      * List&lt;Triangle3D&gt; triangles = Planes.indexedTriangles(vertices, faceIndices, TEST_PRECISION);
    209      * </pre>
    210      * @param vertices vertices available for use in triangle construction
    211      * @param faceIndices array of indices for each triangular face; each entry in the array is an array of
    212      *      3 index values into {@code vertices}, defining the 3 vertices that will be used to construct the
    213      *      triangle
    214      * @param precision precision context used for floating point comparisons
    215      * @return a list of triangles constructed from the set of vertices and face indices
    216      * @throws IllegalArgumentException if any face index array does not contain exactly 3 elements or a set
    217      *      of 3 vertices do not define a plane
    218      * @throws IndexOutOfBoundsException if any index into {@code vertices} is out of bounds
    219      */
    220     public static List<Triangle3D> indexedTriangles(final Vector3D[] vertices, final int[][] faceIndices,
    221             final Precision.DoubleEquivalence precision) {
    222         return indexedTriangles(Arrays.asList(vertices), faceIndices, precision);
    223     }
    224
    225     /** Construct a list of {@link Triangle3D} instances from a set of vertices and arrays of face indices.
    226      * @param vertices vertices available for use in triangle construction
    227      * @param faceIndices array of indices for each triangular face; each entry in the array is an array of
    228      *      3 index values into {@code vertices}, defining the 3 vertices that will be used to construct the
    229      *      triangle
    230      * @param precision precision context used for floating point comparisons
    231      * @return a list of triangles constructed from the set of vertices and face indices
    232      * @throws IllegalArgumentException if any face index array does not contain exactly 3 elements or a set
    233      *      of 3 vertices do not define a plane
    234      * @throws IndexOutOfBoundsException if any index into {@code vertices} is out of bounds
    235      * @see #indexedTriangles(Vector3D[], int[][], Precision.DoubleEquivalence)
    236      */
    237     public static List<Triangle3D> indexedTriangles(final List<? extends Vector3D> vertices, final int[][] faceIndices,
    238             final Precision.DoubleEquivalence precision) {
    239
    240         final int numFaces = faceIndices.length;
    241         final List<Triangle3D> triangles = new ArrayList<>(numFaces);
    242
    243         int[] face;
    244         for (int i = 0; i < numFaces; ++i) {
    245             face = faceIndices[i];
    246             if (face.length != EuclideanUtils.TRIANGLE_VERTEX_COUNT) {
    247                 throw new IllegalArgumentException(MessageFormat.format(
    248                         "Invalid number of vertex indices for face at index {0}: expected {1} but found {2}",
    249                         i, EuclideanUtils.TRIANGLE_VERTEX_COUNT, face.length));
    250             }
    251
    252             triangles.add(triangleFromVertices(
    253                         vertices.get(face[0]),
    254                         vertices.get(face[1]),
    255                         vertices.get(face[2]),
    256                         precision
    257                     ));
    258         }
    259
    260         return triangles;
    261     }
    262
    263     /** Construct a list of {@link ConvexPolygon3D} instances from a set of vertices and arrays of face indices. Each
    264      * face must contain at least 3 vertices but the number of vertices per face does not need to be constant.
    265      * For example, the following code constructs a list of convex polygons forming a square pyramid.
    266      * Note that the first face (the pyramid base) uses a different number of vertices than the other faces.
    267      * <pre>
    268      * Precision.DoubleEquivalence precision = Precision.doubleEquivalenceOfEpsilon(1e-10);
    269      *
    270      * Vector3D[] vertices = {
    271      *      Vector3D.ZERO,
    272      *      Vector3D.of(1, 0, 0),
    273      *      Vector3D.of(1, 1, 0),
    274      *      Vector3D.of(0, 1, 0),
    275      *      Vector3D.of(0.5, 0.5, 4)
    276      * };
    277      *
    278      * int[][] faceIndices = {
    279      *      {0, 3, 2, 1}, // square base
    280      *      {0, 1, 4},
    281      *      {1, 2, 4},
    282      *      {2, 3, 4},
    283      *      {3, 0, 4}
    284      * };
    285      *
    286      * List&lt;ConvexPolygon3D&gt; polygons = Planes.indexedConvexPolygons(vertices, faceIndices, precision);
    287      * </pre>
    288      * @param vertices vertices available for use in convex polygon construction
    289      * @param faceIndices array of indices for each triangular face; each entry in the array is an array of
    290      *      at least 3 index values into {@code vertices}, defining the vertices that will be used to construct the
    291      *      convex polygon
    292      * @param precision precision context used for floating point comparisons
    293      * @return a list of convex polygons constructed from the set of vertices and face indices
    294      * @throws IllegalArgumentException if any face index array does not contain at least 3 elements or a set
    295      *      of vertices do not define a planar convex polygon
    296      * @throws IndexOutOfBoundsException if any index into {@code vertices} is out of bounds
    297      */
    298     public static List<ConvexPolygon3D> indexedConvexPolygons(final Vector3D[] vertices, final int[][] faceIndices,
    299             final Precision.DoubleEquivalence precision) {
    300         return indexedConvexPolygons(Arrays.asList(vertices), faceIndices, precision);
    301     }
    302
    303     /** Construct a list of {@link ConvexPolygon3D} instances from a set of vertices and arrays of face indices. Each
    304      * face must contain at least 3 vertices but the number of vertices per face does not need to be constant.
    305      * @param vertices vertices available for use in convex polygon construction
    306      * @param faceIndices array of indices for each triangular face; each entry in the array is an array of
    307      *      at least 3 index values into {@code vertices}, defining the vertices that will be used to construct the
    308      *      convex polygon
    309      * @param precision precision context used for floating point comparisons
    310      * @return a list of convex polygons constructed from the set of vertices and face indices
    311      * @throws IllegalArgumentException if any face index array does not contain at least 3 elements or a set
    312      *      of vertices do not define a planar convex polygon
    313      * @throws IndexOutOfBoundsException if any index into {@code vertices} is out of bounds
    314      * @see #indexedConvexPolygons(Vector3D[], int[][], Precision.DoubleEquivalence)
    315      */
    316     public static List<ConvexPolygon3D> indexedConvexPolygons(final List<? extends Vector3D> vertices,
    317             final int[][] faceIndices, final Precision.DoubleEquivalence precision) {
    318         final int numFaces = faceIndices.length;
    319         final List<ConvexPolygon3D> polygons = new ArrayList<>(numFaces);
    320         final List<Vector3D> faceVertices = new ArrayList<>();
    321
    322         int[] face;
    323         for (int i = 0; i < numFaces; ++i) {
    324             face = faceIndices[i];
    325             if (face.length < EuclideanUtils.TRIANGLE_VERTEX_COUNT) {
    326                 throw new IllegalArgumentException(MessageFormat.format(
    327                         "Invalid number of vertex indices for face at index {0}: required at least {1} but found {2}",
    328                         i, EuclideanUtils.TRIANGLE_VERTEX_COUNT, face.length));
    329             }
    330
    331             for (final int vertexIndex : face) {
    332                 faceVertices.add(vertices.get(vertexIndex));
    333             }
    334
    335             polygons.add(convexPolygonFromVertices(
    336                         faceVertices,
    337                         precision
    338                     ));
    339
    340             faceVertices.clear();
    341         }
    342
    343         return polygons;
    344     }
    345
    346     /** Get the boundaries of a 3D region created by extruding a polygon defined by a list of vertices. The ends
    347      * ("top" and "bottom") of the extruded 3D region are flat while the sides follow the boundaries of the original
    348      * 2D region.
    349      * @param vertices vertices forming the 2D polygon to extrude
    350      * @param plane plane to extrude the 2D polygon from
    351      * @param extrusionVector vector to extrude the polygon vertices through
    352      * @param precision precision context used to construct the 3D region boundaries
    353      * @return the boundaries of the extruded 3D region
    354      * @throws IllegalStateException if {@code vertices} contains only a single unique vertex
    355      * @throws IllegalArgumentException if regions of non-zero size cannot be produced with the
    356      *      given plane and extrusion vector. This occurs when the extrusion vector has zero length
    357      *      or is orthogonal to the plane normal
    358      * @see LinePath#fromVertexLoop(Collection, Precision.DoubleEquivalence)
    359      * @see #extrude(LinePath, EmbeddingPlane, Vector3D, Precision.DoubleEquivalence)
    360      */
    361     public static List<PlaneConvexSubset> extrudeVertexLoop(final List<Vector2D> vertices,
    362             final EmbeddingPlane plane, final Vector3D extrusionVector, final Precision.DoubleEquivalence precision) {
    363         final LinePath path = LinePath.fromVertexLoop(vertices, precision);
    364         return extrude(path, plane, extrusionVector, precision);
    365     }
    366
    367     /** Get the boundaries of the 3D region created by extruding a 2D line path. The ends ("top" and "bottom") of
    368      * the extruded 3D region are flat while the sides follow the boundaries of the original 2D region. The path is
    369      * converted to a BSP tree before extrusion.
    370      * @param path path to extrude
    371      * @param plane plane to extrude the path from
    372      * @param extrusionVector vector to extrude the polygon points through
    373      * @param precision precision precision context used to construct the 3D region boundaries
    374      * @return the boundaries of the extruded 3D region
    375      * @throws IllegalArgumentException if regions of non-zero size cannot be produced with the
    376      *      given plane and extrusion vector. This occurs when the extrusion vector has zero length
    377      *      or is orthogonal to the plane normal
    378      * @see #extrude(RegionBSPTree2D, EmbeddingPlane, Vector3D, Precision.DoubleEquivalence)
    379      */
    380     public static List<PlaneConvexSubset> extrude(final LinePath path, final EmbeddingPlane plane,
    381             final Vector3D extrusionVector, final Precision.DoubleEquivalence precision) {
    382         return extrude(path.toTree(), plane, extrusionVector, precision);
    383     }
    384
    385     /** Get the boundaries of the 3D region created by extruding a 2D region. The ends ("top" and "bottom") of
    386      * the extruded 3D region are flat while the sides follow the boundaries of the original 2D region.
    387      * @param region region to extrude
    388      * @param plane plane to extrude the region from
    389      * @param extrusionVector vector to extrude the region points through
    390      * @param precision precision precision context used to construct the 3D region boundaries
    391      * @return the boundaries of the extruded 3D region
    392      * @throws IllegalArgumentException if regions of non-zero size cannot be produced with the
    393      *      given plane and extrusion vector. This occurs when the extrusion vector has zero length
    394      *      or is orthogonal to the plane normal
    395      */
    396     public static List<PlaneConvexSubset> extrude(final RegionBSPTree2D region, final EmbeddingPlane plane,
    397             final Vector3D extrusionVector, final Precision.DoubleEquivalence precision) {
    398         return new PlaneRegionExtruder(plane, extrusionVector, precision).extrude(region);
    399     }
    400
    401     /** Get the unique intersection of the plane subset with the given line. Null is
    402      * returned if no unique intersection point exists (ie, the line and plane are
    403      * parallel or coincident) or the line does not intersect the plane subset.
    404      * @param planeSubset plane subset to intersect with
    405      * @param line line to intersect with this plane subset
    406      * @return the unique intersection point between the line and this plane subset
    407      *      or null if no such point exists.
    408      */
    409     static Vector3D intersection(final PlaneSubset planeSubset, final Line3D line) {
    410         final Vector3D pt = planeSubset.getPlane().intersection(line);
    411         return (pt != null && planeSubset.contains(pt)) ? pt : null;
    412     }
    413
    414     /** Get the unique intersection of the plane subset with the given line subset. Null
    415      * is returned if the underlying line and plane do not have a unique intersection
    416      * point (ie, they are parallel or coincident) or the intersection point is unique
    417      * but is not contained in both the line subset and plane subset.
    418      * @param planeSubset plane subset to intersect with
    419      * @param lineSubset line subset to intersect with
    420      * @return the unique intersection point between this plane subset and the argument or
    421      *      null if no such point exists.
    422      */
    423     static Vector3D intersection(final PlaneSubset planeSubset, final LineConvexSubset3D lineSubset) {
    424         final Vector3D pt = intersection(planeSubset, lineSubset.getLine());
    425         return (pt != null && lineSubset.contains(pt)) ? pt : null;
    426     }
    427
    428     /** Validate that the actual plane contains the same points as the expected plane, throwing an exception if not.
    429      * The subspace orientations of embedding planes are not considered.
    430      * @param expected the expected plane
    431      * @param actual the actual plane
    432      * @throws IllegalArgumentException if the actual plane is not equivalent to the expected plane
    433      */
    434     static void validatePlanesEquivalent(final Plane expected, final Plane actual) {
    435         if (!expected.eq(actual, expected.getPrecision())) {
    436             throw new IllegalArgumentException("Arguments do not represent the same plane. Expected " +
    437                     expected + " but was " + actual + ".");
    438         }
    439     }
    440
    441     /** Generic split method that uses performs the split using the subspace region of the plane subset.
    442      * @param splitter splitting hyperplane
    443      * @param subset the plane subset being split
    444      * @param factory function used to create new plane subset instances
    445      * @param <T> Plane subset implementation type
    446      * @return the result of the split operation
    447      */
    448     static <T extends PlaneSubset> Split<T> subspaceSplit(final Plane splitter, final T subset,
    449             final BiFunction<? super EmbeddingPlane, ? super HyperplaneBoundedRegion<Vector2D>, T> factory) {
    450
    451         final EmbeddingPlane thisPlane = subset.getPlane().getEmbedding();
    452
    453         final Line3D intersection = thisPlane.intersection(splitter);
    454         if (intersection == null) {
    455             return getNonIntersectingSplitResult(splitter, subset);
    456         } else {
    457             final EmbeddingPlane embeddingPlane = subset.getPlane().getEmbedding();
    458
    459             // the lines intersect; split the subregion
    460             final Vector3D intersectionOrigin = intersection.getOrigin();
    461             final Vector2D subspaceP1 = embeddingPlane.toSubspace(intersectionOrigin);
    462             final Vector2D subspaceP2 = embeddingPlane.toSubspace(intersectionOrigin.add(intersection.getDirection()));
    463
    464             final Line subspaceSplitter = Lines.fromPoints(subspaceP1, subspaceP2, thisPlane.getPrecision());
    465
    466             final Split<? extends HyperplaneBoundedRegion<Vector2D>> split =
    467                     subset.getEmbedded().getSubspaceRegion().split(subspaceSplitter);
    468             final SplitLocation subspaceSplitLoc = split.getLocation();
    469
    470             if (SplitLocation.MINUS == subspaceSplitLoc) {
    471                 return new Split<>(subset, null);
    472             } else if (SplitLocation.PLUS == subspaceSplitLoc) {
    473                 return new Split<>(null, subset);
    474             }
    475
    476             final T minus = (split.getMinus() != null) ? factory.apply(thisPlane, split.getMinus()) : null;
    477             final T plus = (split.getPlus() != null) ? factory.apply(thisPlane, split.getPlus()) : null;
    478
    479             return new Split<>(minus, plus);
    480         }
    481     }
    482
    483     /** Get a split result for cases where the splitting plane and the plane containing the subset being split
    484      * do not intersect. Callers are responsible for ensuring that the planes involved do not actually intersect.
    485      * @param <T> Plane subset implementation type
    486      * @param splitter plane performing the splitting
    487      * @param subset subset being split
    488      * @return the split result for the non-intersecting split
    489      */
    490     private static <T extends PlaneSubset> Split<T> getNonIntersectingSplitResult(
    491             final Plane splitter, final T subset) {
    492         final Plane plane = subset.getPlane();
    493
    494         final double offset = splitter.offset(plane);
    495         final int comp = plane.getPrecision().compare(offset, 0.0);
    496
    497         if (comp < 0) {
    498             return new Split<>(subset, null);
    499         } else if (comp > 0) {
    500             return new Split<>(null, subset);
    501         } else {
    502             return new Split<>(null, null);
    503         }
    504     }
    505
    506     /** Construct a convex polygon 3D from a plane and a list of vertices lying in the plane. Callers are
    507      * responsible for ensuring that the vertices lie in the plane and define a convex polygon.
    508      * @param plane the plane containing the convex polygon
    509      * @param vertices vertices defining the closed, convex polygon. The must must contain at least 3 unique
    510      *      vertices and should not include the start vertex at the end of the list.
    511      * @return a new convex polygon instance
    512      * @throws IllegalArgumentException if the size of {@code vertices} if less than 3
    513      */
    514     static ConvexPolygon3D fromConvexPlanarVertices(final Plane plane, final List<Vector3D> vertices) {
    515         final int size = vertices.size();
    516
    517         if (size == EuclideanUtils.TRIANGLE_VERTEX_COUNT) {
    518             return new SimpleTriangle3D(plane, vertices.get(0), vertices.get(1), vertices.get(2));
    519         }
    520
    521         return new VertexListConvexPolygon3D(plane, vertices);
    522     }
    523
    524     /** Convert a convex polygon defined by a plane and list of points into a triangle fan.
    525      * @param plane plane containing the convex polygon
    526      * @param vertices vertices defining the convex polygon
    527      * @return a triangle fan representing the same area as the convex polygon
    528      * @throws IllegalArgumentException if fewer than 3 vertices are given
    529      */
    530     static List<Triangle3D> convexPolygonToTriangleFan(final Plane plane, final List<Vector3D> vertices) {
    531         return EuclideanUtils.convexPolygonToTriangleFan(vertices,
    532                 tri -> new SimpleTriangle3D(plane, tri.get(0), tri.get(1), tri.get(2)));
    533     }
    534
    535     /** Internal helper class used to construct planes from sequences of points. Instances can be also be
    536      * configured to collect lists of unique points found during plane construction and validate that the
    537      * defined region is convex.
    538      */
    539     private static final class PlaneBuilder {
    540
    541         /** The point sequence to build a plane for. */
    542         private final Collection<? extends Vector3D> pts;
    543
    544         /** Precision context used for floating point comparisons. */
    545         private final Precision.DoubleEquivalence precision;
    546
    547         /** The start point from the point sequence. */
    548         private Vector3D startPt;
    549
    550         /** The previous point from the point sequence. */
    551         private Vector3D prevPt;
    552
    553         /** The previous vector from the point sequence, preceding from the {@code startPt} to {@code prevPt}. */
    554         private Vector3D prevVector;
    555
    556         /** The computed {@code normal} vector for the plane. */
    557         private Vector3D.Unit normal;
    558
    559         /** The x component of the sum of all cross products from adjacent vectors in the point sequence. */
    560         private double crossSumX;
    561
    562         /** The y component of the sum of all cross products from adjacent vectors in the point sequence. */
    563         private double crossSumY;
    564
    565         /** The z component of the sum of all cross products from adjacent vectors in the point sequence. */
    566         private double crossSumZ;
    567
    568         /** If true, an exception will be thrown if the point sequence is discovered to be non-convex. */
    569         private boolean requireConvex;
    570
    571         /** List that unique vertices discovered in the input sequence will be added to. */
    572         private List<? super Vector3D> uniqueVertexOutput;
    573
    574         /** Construct a new build instance for the given point sequence and precision context.
    575          * @param pts point sequence
    576          * @param precision precision context used to perform floating point comparisons
    577          */
    578         PlaneBuilder(final Collection<? extends Vector3D> pts, final Precision.DoubleEquivalence precision) {
    579             this.pts = pts;
    580             this.precision = precision;
    581         }
    582
    583         /** Build a plane from the configured point sequence.
    584          * @return a plane built from the configured point sequence
    585          * @throws IllegalArgumentException if the points do not define a plane
    586          */
    587         Plane build() {
    588             if (pts.size() < EuclideanUtils.TRIANGLE_VERTEX_COUNT) {
    589                 throw nonPlanar();
    590             }
    591
    592             pts.forEach(this::processPoint);
    593
    594             return createPlane();
    595         }
    596
    597         /** Build a plane from the configured point sequence, validating that the points form a convex region
    598          * and adding all discovered unique points to the given list.
    599          * @param vertexOutput list that unique points discovered in the point sequence will be added to
    600          * @return a plane created from the configured point sequence
    601          * @throws IllegalArgumentException if the points do not define a plane or the {@code requireConvex}
    602          *      flag is true and the points do not define a convex area
    603          */
    604         Plane buildForConvexPolygon(final List<? super Vector3D> vertexOutput) {
    605             this.requireConvex = true;
    606             this.uniqueVertexOutput = vertexOutput;
    607
    608             return build();
    609         }
    610
    611         /** Process a point from the point sequence.
    612          * @param pt
    613          * @throws IllegalArgumentException if the points do not define a plane or the {@code requireConvex}
    614          *      flag is true and the points do not define a convex area
    615          */
    616         private void processPoint(final Vector3D pt) {
    617             if (prevPt == null) {
    618                 startPt = pt;
    619                 prevPt = pt;
    620
    621                 if (uniqueVertexOutput != null) {
    622                     uniqueVertexOutput.add(pt);
    623                 }
    624
    625             } else if (!prevPt.eq(pt, precision)) { // skip duplicate points
    626                 final Vector3D vec = startPt.vectorTo(pt);
    627
    628                 if (prevVector != null) {
    629                     processCrossProduct(prevVector.cross(vec));
    630                 }
    631
    632                 if (uniqueVertexOutput != null) {
    633                     uniqueVertexOutput.add(pt);
    634                 }
    635
    636                 prevPt = pt;
    637                 prevVector = vec;
    638             }
    639         }
    640
    641         /** Process the computed cross product of two vectors from the input point sequence. The vectors
    642          * start at the first point in the sequence and point to adjacent points later in the sequence.
    643          * @param cross the cross product of two vectors from the input point sequence
    644          * @throws IllegalArgumentException if the points do not define a plane or the {@code requireConvex}
    645          *      flag is true and the points do not define a convex area
    646          */
    647         private void processCrossProduct(final Vector3D cross) {
    648             crossSumX += cross.getX();
    649             crossSumY += cross.getY();
    650             crossSumZ += cross.getZ();
    651
    652             final double crossNorm = cross.norm();
    653
    654             if (!precision.eqZero(crossNorm)) {
    655                 // the cross product has non-zero magnitude
    656                 if (normal == null) {
    657                     // save the first non-zero cross product as our normal
    658                     normal = cross.normalize();
    659                 } else {
    660                     final double crossDot = normal.dot(cross) / crossNorm;
    661
    662                     // check non-planar before non-convex since the former is a more general type
    663                     // of issue
    664                     if (!precision.eq(1.0, Math.abs(crossDot))) {
    665                         throw nonPlanar();
    666                     } else if (requireConvex && crossDot < 0) {
    667                         throw nonConvex();
    668                     }
    669                 }
    670             }
    671         }
    672
    673         /** Construct the plane instance using the value gathered during point processing.
    674          * @return the created plane instance
    675          * @throws IllegalArgumentException if the point do not define a plane
    676          */
    677         private Plane createPlane() {
    678             if (normal == null) {
    679                 throw nonPlanar();
    680             }
    681
    682             // flip the normal if needed to match the overall orientation of the points
    683             if (normal.dot(Vector3D.of(crossSumX, crossSumY, crossSumZ)) < 0) {
    684                 normal = normal.negate();
    685             }
    686
    687             // construct the plane
    688             final double originOffset = -startPt.dot(normal);
    689
    690             return new Plane(normal, originOffset, precision);
    691         }
    692
    693         /** Return an exception with a message stating that the points given to this builder do not
    694          * define a plane.
    695          * @return an exception stating that the points do not define a plane
    696          */
    697         private IllegalArgumentException nonPlanar() {
    698             return new IllegalArgumentException("Points do not define a plane: " + pts);
    699         }
    700
    701         /** Return an exception with a message stating that the points given to this builder do not
    702          * define a convex region.
    703          * @return an exception stating that the points do not define a plane
    704          */
    705         private IllegalArgumentException nonConvex() {
    706             return new IllegalArgumentException("Points do not define a convex region: " + pts);
    707         }
    708     }
    709
    710     /** Class designed to create 3D regions by taking a 2D region and extruding from a base plane
    711      * through an extrusion vector. The ends ("top" and "bottom") of the extruded 3D region are flat
    712      * while the sides follow the boundaries of the original 2D region.
    713      */
    714     private static final class PlaneRegionExtruder {
    715         /** Base plane to extrude from. */
    716         private final EmbeddingPlane basePlane;
    717
    718         /** Vector to extrude along; the extruded plane is translated from the base plane by this amount. */
    719         private final Vector3D extrusionVector;
    720
    721         /** True if the extrusion vector points to the plus side of the base plane. */
    722         private final boolean extrudingOnPlusSide;
    723
    724         /** Precision context used to create boundaries. */
    725         private final Precision.DoubleEquivalence precision;
    726
    727         /** Construct a new instance that performs extrusions from {@code basePlane} along {@code extrusionVector}.
    728          * @param basePlane base plane to extrude from
    729          * @param extrusionVector vector to extrude along
    730          * @param precision precision context used to construct boundaries
    731          * @throws IllegalArgumentException if the given extrusion vector and plane produce regions
    732          *      of zero size
    733          */
    734         PlaneRegionExtruder(final EmbeddingPlane basePlane, final Vector3D extrusionVector,
    735                 final Precision.DoubleEquivalence precision) {
    736
    737             this.basePlane = basePlane;
    738
    739             // Extruded plane; this forms the end of the 3D region opposite the base plane.
    740             final EmbeddingPlane extrudedPlane = basePlane.translate(extrusionVector);
    741
    742             if (basePlane.contains(extrudedPlane)) {
    743                 throw new IllegalArgumentException(
    744                         "Extrusion vector produces regions of zero size: extrusionVector= " +
    745                                 extrusionVector + ",  plane= " + basePlane);
    746             }
    747
    748             this.extrusionVector = extrusionVector;
    749             this.extrudingOnPlusSide = basePlane.getNormal().dot(extrusionVector) > 0;
    750
    751             this.precision = precision;
    752         }
    753
    754         /** Extrude the given 2D BSP tree using the configured base plane and extrusion vector.
    755          * @param subspaceRegion region to extrude
    756          * @return the boundaries of the extruded region
    757          */
    758         public List<PlaneConvexSubset> extrude(final RegionBSPTree2D subspaceRegion) {
    759             final List<PlaneConvexSubset> extrudedBoundaries = new ArrayList<>();
    760
    761             // add the boundaries
    762             addEnds(subspaceRegion, extrudedBoundaries);
    763             addSides(subspaceRegion, extrudedBoundaries);
    764
    765             return extrudedBoundaries;
    766         }
    767
    768         /** Add the end ("top" and "bottom") of the extruded subspace region to the result list.
    769          * @param subspaceRegion subspace region being extruded.
    770          * @param result list to add the boundary results to
    771          */
    772         private void addEnds(final RegionBSPTree2D subspaceRegion, final List<? super PlaneConvexSubset> result) {
    773             // add the base boundaries
    774             final List<ConvexArea> baseAreas = subspaceRegion.toConvex();
    775
    776             final List<PlaneConvexSubset> baseList = new ArrayList<>(baseAreas.size());
    777             final List<PlaneConvexSubset> extrudedList = new ArrayList<>(baseAreas.size());
    778
    779             final AffineTransformMatrix3D extrudeTransform = AffineTransformMatrix3D.createTranslation(extrusionVector);
    780
    781             PlaneConvexSubset base;
    782             for (final ConvexArea area : baseAreas) {
    783                 base = subsetFromConvexArea(basePlane, area);
    784                 if (extrudingOnPlusSide) {
    785                     base = base.reverse();
    786                 }
    787
    788                 baseList.add(base);
    789                 extrudedList.add(base.transform(extrudeTransform).reverse());
    790             }
    791
    792             result.addAll(baseList);
    793             result.addAll(extrudedList);
    794         }
    795
    796         /** Add the side boundaries of the extruded region to the result list.
    797          * @param subspaceRegion subspace region being extruded.
    798          * @param result list to add the boundary results to
    799          */
    800         private void addSides(final RegionBSPTree2D subspaceRegion, final List<? super PlaneConvexSubset> result) {
    801             Vector2D subStartPt;
    802             Vector2D subEndPt;
    803
    804             PlaneConvexSubset boundary;
    805             for (final LinePath path : subspaceRegion.getBoundaryPaths()) {
    806                 for (final LineConvexSubset lineSubset : path.getElements()) {
    807                     subStartPt = lineSubset.getStartPoint();
    808                     subEndPt = lineSubset.getEndPoint();
    809
    810                     boundary = (subStartPt != null && subEndPt != null) ?
    811                             extrudeSideFinite(basePlane.toSpace(subStartPt), basePlane.toSpace(subEndPt)) :
    812                             extrudeSideInfinite(lineSubset);
    813
    814                     result.add(boundary);
    815                 }
    816             }
    817         }
    818
    819         /** Extrude a single, finite boundary forming one of the sides of the extruded region.
    820          * @param startPt start point of the boundary
    821          * @param endPt end point of the boundary
    822          * @return the extruded region side boundary
    823          */
    824         private ConvexPolygon3D extrudeSideFinite(final Vector3D startPt, final Vector3D endPt) {
    825             final Vector3D extrudedStartPt = startPt.add(extrusionVector);
    826             final Vector3D extrudedEndPt = endPt.add(extrusionVector);
    827
    828             final List<Vector3D> vertices = extrudingOnPlusSide ?
    829                     Arrays.asList(startPt, endPt, extrudedEndPt, extrudedStartPt) :
    830                     Arrays.asList(startPt, extrudedStartPt, extrudedEndPt, endPt);
    831
    832             return convexPolygonFromVertices(vertices, precision);
    833         }
    834
    835         /** Extrude a single, infinite boundary forming one of the sides of the extruded region.
    836          * @param lineSubset line subset to extrude
    837          * @return the extruded region side boundary
    838          */
    839         private PlaneConvexSubset extrudeSideInfinite(final LineConvexSubset lineSubset) {
    840             final Vector2D subLinePt = lineSubset.getLine().getOrigin();
    841             final Vector2D subLineDir = lineSubset.getLine().getDirection();
    842
    843             final Vector3D linePt = basePlane.toSpace(subLinePt);
    844             final Vector3D lineDir = linePt.vectorTo(basePlane.toSpace(subLinePt.add(subLineDir)));
    845
    846             final EmbeddingPlane sidePlane;
    847             if (extrudingOnPlusSide) {
    848                 sidePlane = fromPointAndPlaneVectors(linePt, lineDir, extrusionVector, precision);
    849             } else {
    850                 sidePlane = fromPointAndPlaneVectors(linePt, extrusionVector, lineDir, precision);
    851             }
    852
    853             final Vector2D sideLineOrigin = sidePlane.toSubspace(linePt);
    854             final Vector2D sideLineDir = sideLineOrigin.vectorTo(sidePlane.toSubspace(linePt.add(lineDir)));
    855
    856             final Vector2D extrudedSideLineOrigin = sidePlane.toSubspace(linePt.add(extrusionVector));
    857
    858             final Vector2D sideExtrusionDir = sidePlane.toSubspace(sidePlane.getOrigin().add(extrusionVector))
    859                     .normalize();
    860
    861             // construct a list of lines forming the bounds of the extruded subspace region
    862             final List<Line> lines = new ArrayList<>();
    863
    864             // add the top and bottom lines (original and extruded)
    865             if (extrudingOnPlusSide) {
    866                 lines.add(Lines.fromPointAndDirection(sideLineOrigin, sideLineDir, precision));
    867                 lines.add(Lines.fromPointAndDirection(extrudedSideLineOrigin, sideLineDir.negate(), precision));
    868             } else {
    869                 lines.add(Lines.fromPointAndDirection(sideLineOrigin, sideLineDir.negate(), precision));
    870                 lines.add(Lines.fromPointAndDirection(extrudedSideLineOrigin, sideLineDir, precision));
    871             }
    872
    873             // if we have a point on the original line, then connect the two
    874             final Vector2D startPt = lineSubset.getStartPoint();
    875             final Vector2D endPt = lineSubset.getEndPoint();
    876             if (startPt != null) {
    877                 lines.add(Lines.fromPointAndDirection(
    878                         sidePlane.toSubspace(basePlane.toSpace(startPt)),
    879                         extrudingOnPlusSide ? sideExtrusionDir.negate() : sideExtrusionDir,
    880                         precision));
    881             } else if (endPt != null) {
    882                 lines.add(Lines.fromPointAndDirection(
    883                         sidePlane.toSubspace(basePlane.toSpace(endPt)),
    884                         extrudingOnPlusSide ? sideExtrusionDir : sideExtrusionDir.negate(),
    885                         precision));
    886             }
    887
    888             return subsetFromConvexArea(sidePlane, ConvexArea.fromBounds(lines));
    889         }
    890     }
    891 }