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.util.ArrayList;
    20  import java.util.List;
    21  import java.util.stream.Stream;
    22  import java.util.stream.StreamSupport;
    23
    24  import org.apache.commons.geometry.core.partitioning.Hyperplane;
    25  import org.apache.commons.geometry.core.partitioning.HyperplaneConvexSubset;
    26  import org.apache.commons.geometry.core.partitioning.HyperplaneSubset;
    27  import org.apache.commons.geometry.core.partitioning.Split;
    28  import org.apache.commons.geometry.core.partitioning.bsp.AbstractBSPTree;
    29  import org.apache.commons.geometry.core.partitioning.bsp.AbstractPartitionedRegionBuilder;
    30  import org.apache.commons.geometry.core.partitioning.bsp.AbstractRegionBSPTree;
    31  import org.apache.commons.geometry.core.partitioning.bsp.BSPTreeVisitor;
    32  import org.apache.commons.geometry.core.partitioning.bsp.RegionCutBoundary;
    33  import org.apache.commons.geometry.euclidean.threed.line.Line3D;
    34  import org.apache.commons.geometry.euclidean.threed.line.LineConvexSubset3D;
    35  import org.apache.commons.geometry.euclidean.threed.line.LinecastPoint3D;
    36  import org.apache.commons.numbers.core.Precision;
    37
    38  /** Binary space partitioning (BSP) tree representing a region in three dimensional
    39   * Euclidean space.
    40   */
    41  public final class RegionBSPTree3D extends AbstractRegionBSPTree<Vector3D, RegionBSPTree3D.RegionNode3D>
    42      implements BoundarySource3D {
    43
    44      /** Create a new, empty region. */
    45      public RegionBSPTree3D() {
    46          this(false);
    47      }
    48
    49      /** Create a new region. If {@code full} is true, then the region will
    50       * represent the entire 3D space. Otherwise, it will be empty.
    51       * @param full whether or not the region should contain the entire
    52       *      3D space or be empty
    53       */
    54      public RegionBSPTree3D(final boolean full) {
    55          super(full);
    56      }
    57
    58      /** Return a deep copy of this instance.
    59       * @return a deep copy of this instance.
    60       * @see #copy(org.apache.commons.geometry.core.partitioning.bsp.BSPTree)
    61       */
    62      public RegionBSPTree3D copy() {
    63          final RegionBSPTree3D result = RegionBSPTree3D.empty();
    64          result.copy(this);
    65
    66          return result;
    67      }
    68
    69      /** {@inheritDoc} */
    70      @Override
    71      public Iterable<PlaneConvexSubset> boundaries() {
    72          return createBoundaryIterable(PlaneConvexSubset.class::cast);
    73      }
    74
    75      /** {@inheritDoc} */
    76      @Override
    77      public Stream<PlaneConvexSubset> boundaryStream() {
    78          return StreamSupport.stream(boundaries().spliterator(), false);
    79      }
    80
    81      /** {@inheritDoc} */
    82      @Override
    83      public List<PlaneConvexSubset> getBoundaries() {
    84          return createBoundaryList(PlaneConvexSubset.class::cast);
    85      }
    86
    87      /** Return a list of {@link ConvexVolume}s representing the same region
    88       * as this instance. One convex volume is returned for each interior leaf
    89       * node in the tree.
    90       * @return a list of convex volumes representing the same region as this
    91       *      instance
    92       */
    93      public List<ConvexVolume> toConvex() {
    94          final List<ConvexVolume> result = new ArrayList<>();
    95
    96          toConvexRecursive(getRoot(), ConvexVolume.full(), result);
    97
    98          return result;
    99      }
    100
    101     /** Recursive method to compute the convex volumes of all inside leaf nodes in the subtree rooted at the given
    102      * node. The computed convex volumes are added to the given list.
    103      * @param node root of the subtree to compute the convex volumes for
    104      * @param nodeVolume the convex volume for the current node; this will be split by the node's cut hyperplane to
    105      *      form the convex volumes for any child nodes
    106      * @param result list containing the results of the computation
    107      */
    108     private void toConvexRecursive(final RegionNode3D node, final ConvexVolume nodeVolume,
    109             final List<? super ConvexVolume> result) {
    110
    111         if (node.isLeaf()) {
    112             // base case; only add to the result list if the node is inside
    113             if (node.isInside()) {
    114                 result.add(nodeVolume);
    115             }
    116         } else {
    117             // recurse
    118             final Split<ConvexVolume> split = nodeVolume.split(node.getCutHyperplane());
    119
    120             toConvexRecursive(node.getMinus(), split.getMinus(), result);
    121             toConvexRecursive(node.getPlus(), split.getPlus(), result);
    122         }
    123     }
    124
    125     /** {@inheritDoc} */
    126     @Override
    127     public Split<RegionBSPTree3D> split(final Hyperplane<Vector3D> splitter) {
    128         return split(splitter, RegionBSPTree3D.empty(), RegionBSPTree3D.empty());
    129     }
    130
    131     /** {@inheritDoc} */
    132     @Override
    133     public Vector3D project(final Vector3D pt) {
    134         // use our custom projector so that we can disambiguate points that are
    135         // actually equidistant from the target point
    136         final BoundaryProjector3D projector = new BoundaryProjector3D(pt);
    137         accept(projector);
    138
    139         return projector.getProjected();
    140     }
    141
    142     /** Return the current instance.
    143      */
    144     @Override
    145     public RegionBSPTree3D toTree() {
    146         return this;
    147     }
    148
    149     /** {@inheritDoc} */
    150     @Override
    151     public List<LinecastPoint3D> linecast(final LineConvexSubset3D subset) {
    152         final LinecastVisitor visitor = new LinecastVisitor(subset, false);
    153         accept(visitor);
    154
    155         return visitor.getResults();
    156     }
    157
    158     /** {@inheritDoc} */
    159     @Override
    160     public LinecastPoint3D linecastFirst(final LineConvexSubset3D subset) {
    161         final LinecastVisitor visitor = new LinecastVisitor(subset, true);
    162         accept(visitor);
    163
    164         return visitor.getFirstResult();
    165     }
    166
    167     /** {@inheritDoc} */
    168     @Override
    169     protected RegionSizeProperties<Vector3D> computeRegionSizeProperties() {
    170         // handle simple cases
    171         if (isFull()) {
    172             return new RegionSizeProperties<>(Double.POSITIVE_INFINITY, null);
    173         } else if (isEmpty()) {
    174             return new RegionSizeProperties<>(0, null);
    175         }
    176
    177         final RegionSizePropertiesVisitor visitor = new RegionSizePropertiesVisitor();
    178         accept(visitor);
    179
    180         return visitor.getRegionSizeProperties();
    181     }
    182
    183     /** {@inheritDoc} */
    184     @Override
    185     protected RegionNode3D createNode() {
    186         return new RegionNode3D(this);
    187     }
    188
    189     /** Return a new instance containing all of 3D space.
    190      * @return a new instance containing all of 3D space.
    191      */
    192     public static RegionBSPTree3D full() {
    193         return new RegionBSPTree3D(true);
    194     }
    195
    196     /** Return a new, empty instance. The represented region is completely empty.
    197      * @return a new, empty instance.
    198      */
    199     public static RegionBSPTree3D empty() {
    200         return new RegionBSPTree3D(false);
    201     }
    202
    203     /** Construct a new tree from the given boundaries. If no boundaries
    204      * are present, the returned tree is empty.
    205      * @param boundaries boundaries to construct the tree from
    206      * @return a new tree instance constructed from the given boundaries
    207      * @see #from(Iterable, boolean)
    208      */
    209     public static RegionBSPTree3D from(final Iterable<? extends PlaneConvexSubset> boundaries) {
    210         return from(boundaries, false);
    211     }
    212
    213     /** Construct a new tree from the given boundaries. If {@code full} is true, then
    214      * the initial tree before boundary insertion contains the entire space. Otherwise,
    215      * it is empty.
    216      * @param boundaries boundaries to construct the tree from
    217      * @param full if true, the initial tree will contain the entire space
    218      * @return a new tree instance constructed from the given boundaries
    219      */
    220     public static RegionBSPTree3D from(final Iterable<? extends PlaneConvexSubset> boundaries, final boolean full) {
    221         final RegionBSPTree3D tree = new RegionBSPTree3D(full);
    222         tree.insert(boundaries);
    223
    224         return tree;
    225     }
    226
    227     /** Create a new {@link PartitionedRegionBuilder3D} instance which can be used to build balanced
    228      * BSP trees from region boundaries.
    229      * @return a new {@link PartitionedRegionBuilder3D} instance
    230      */
    231     public static PartitionedRegionBuilder3D partitionedRegionBuilder() {
    232         return new PartitionedRegionBuilder3D();
    233     }
    234
    235     /** BSP tree node for three dimensional Euclidean space.
    236      */
    237     public static final class RegionNode3D extends AbstractRegionBSPTree.AbstractRegionNode<Vector3D, RegionNode3D> {
    238         /** Simple constructor.
    239          * @param tree the owning tree instance
    240          */
    241         RegionNode3D(final AbstractBSPTree<Vector3D, RegionNode3D> tree) {
    242             super(tree);
    243         }
    244
    245         /** Get the region represented by this node. The returned region contains
    246          * the entire area contained in this node, regardless of the attributes of
    247          * any child nodes.
    248          * @return the region represented by this node
    249          */
    250         public ConvexVolume getNodeRegion() {
    251             ConvexVolume volume = ConvexVolume.full();
    252
    253             RegionNode3D child = this;
    254             RegionNode3D parent;
    255
    256             while ((parent = child.getParent()) != null) {
    257                 final Split<ConvexVolume> split = volume.split(parent.getCutHyperplane());
    258
    259                 volume = child.isMinus() ? split.getMinus() : split.getPlus();
    260
    261                 child = parent;
    262             }
    263
    264             return volume;
    265         }
    266
    267         /** {@inheritDoc} */
    268         @Override
    269         protected RegionNode3D getSelf() {
    270             return this;
    271         }
    272     }
    273
    274     /** Class used to build regions in Euclidean 3D space by inserting boundaries into a BSP
    275      * tree containing "partitions", i.e. structural cuts where both sides of the cut have the same region location.
    276      * When partitions are chosen that effectively divide the region boundaries at each partition level, the
    277      * constructed tree is shallower and more balanced than one constructed from the region boundaries alone,
    278      * resulting in improved performance. For example, consider a mesh approximation of a sphere. The region is
    279      * convex so each boundary has all of the other boundaries on its minus side; the plus sides are all empty.
    280      * When these boundaries are inserted directly into a tree, the tree degenerates into a simple linked list of
    281      * nodes with a height directly proportional to the number of boundaries. This means that many operations on the
    282      * tree, such as inside/outside testing of points, involve iterating through each and every region boundary. In
    283      * contrast, if a partition is first inserted that passes through the sphere center, the first BSP tree node
    284      * contains region nodes on its plus <em>and</em> minus sides, cutting the height of the tree in half. Operations
    285      * such as inside/outside testing are then able to skip half of the tree nodes with a single test on the
    286      * root node, resulting in drastically improved performance. Insertion of additional partitions (using a grid
    287      * layout, for example) can produce even shallower trees, although there is a point unique to each boundary set at
    288      * which the addition of more partitions begins to decrease instead of increase performance.
    289      *
    290      * <h2>Usage</h2>
    291      * <p>Usage of this class consists of two phases: (1) <em>partition insertion</em> and (2) <em>boundary
    292      * insertion</em>. Instances begin in the <em>partition insertion</em> phase. Here, partitions can be inserted
    293      * into the empty tree using {@link PartitionedRegionBuilder3D#insertPartition(PlaneConvexSubset) insertPartition}
    294      * or similar methods. The {@link org.apache.commons.geometry.core.partitioning.bsp.RegionCutRule#INHERIT INHERIT}
    295      * cut rule is used internally to insert the cut so the represented region remains empty even as partitions are
    296      * inserted.
    297      * </p>
    298      *
    299      * <p>The instance moves into the <em>boundary insertion</em> phase when the caller inserts the first region
    300      * boundary, using {@link PartitionedRegionBuilder3D#insertBoundary(PlaneConvexSubset) insertBoundary} or
    301      * similar methods. Attempting to insert a partition after this point results in an {@code IllegalStateException}.
    302      * This ensures that partitioning cuts are always located higher up the tree than boundary cuts.</p>
    303      *
    304      * <p>After all boundaries are inserted, the {@link PartitionedRegionBuilder3D#build() build} method is used
    305      * to perform final processing and return the computed tree.</p>
    306      */
    307     public static final class PartitionedRegionBuilder3D
    308         extends AbstractPartitionedRegionBuilder<Vector3D, RegionNode3D> {
    309
    310         /** Construct a new builder instance.
    311          */
    312         private PartitionedRegionBuilder3D() {
    313             super(RegionBSPTree3D.empty());
    314         }
    315
    316         /** Insert a partition plane.
    317          * @param partition partition to insert
    318          * @return this instance
    319          * @throws IllegalStateException if a boundary has previously been inserted
    320          */
    321         public PartitionedRegionBuilder3D insertPartition(final Plane partition) {
    322             return insertPartition(partition.span());
    323         }
    324
    325         /** Insert a plane convex subset as a partition.
    326          * @param partition partition to insert
    327          * @return this instance
    328          * @throws IllegalStateException if a boundary has previously been inserted
    329          */
    330         public PartitionedRegionBuilder3D insertPartition(final PlaneConvexSubset partition) {
    331             insertPartitionInternal(partition);
    332
    333             return this;
    334         }
    335
    336         /** Insert a set of three axis aligned planes intersecting at the given point as partitions.
    337          * The planes all contain the {@code center} point and have the normals {@code +x}, {@code +y},
    338          * and {@code +z} in that order. If inserted into an empty tree, this will partition the space
    339          * into 8 sections.
    340          * @param center center point for the partitions; all 3 inserted planes intersect at this point
    341          * @param precision precision context used to construct the planes
    342          * @return this instance
    343          * @throws IllegalStateException if a boundary has previously been inserted
    344          */
    345         public PartitionedRegionBuilder3D insertAxisAlignedPartitions(final Vector3D center,
    346                 final Precision.DoubleEquivalence precision) {
    347
    348             insertPartition(Planes.fromPointAndNormal(center, Vector3D.Unit.PLUS_X, precision));
    349             insertPartition(Planes.fromPointAndNormal(center, Vector3D.Unit.PLUS_Y, precision));
    350             insertPartition(Planes.fromPointAndNormal(center, Vector3D.Unit.PLUS_Z, precision));
    351
    352             return this;
    353         }
    354
    355         /** Insert a 3D grid of partitions. The partitions are constructed recursively: at each level a set of
    356          * three axis-aligned partitioning planes are inserted using
    357          * {@link #insertAxisAlignedPartitions(Vector3D, Precision.DoubleEquivalence) insertAxisAlignedPartitions}.
    358          * The algorithm then recurses using bounding boxes from the min point to the center and from the center
    359          * point to the max. Note that this means no partitions are ever inserted directly on the boundaries of
    360          * the given bounding box. This is intentional and done to allow this method to be called directly with the
    361          * bounding box from a set of boundaries to be inserted without unnecessarily adding partitions that will
    362          * never have region boundaries on both sides.
    363          * @param bounds bounding box for the grid
    364          * @param level recursion level for the grid; each level subdivides each grid cube into 8 sections, making the
    365          *      total number of grid cubes equal to {@code 8 ^ level}
    366          * @param precision precision context used to construct the partition planes
    367          * @return this instance
    368          * @throws IllegalStateException if a boundary has previously been inserted
    369          */
    370         public PartitionedRegionBuilder3D insertAxisAlignedGrid(final Bounds3D bounds, final int level,
    371                 final Precision.DoubleEquivalence precision) {
    372
    373             insertAxisAlignedGridRecursive(bounds.getMin(), bounds.getMax(), level, precision);
    374
    375             return this;
    376         }
    377
    378         /** Recursively insert axis-aligned grid partitions.
    379          * @param min min point for the grid cube to partition
    380          * @param max max point for the grid cube to partition
    381          * @param level current recursion level
    382          * @param precision precision context used to construct the partition planes
    383          */
    384         private void insertAxisAlignedGridRecursive(final Vector3D min, final Vector3D max, final int level,
    385                 final Precision.DoubleEquivalence precision) {
    386             if (level > 0) {
    387                 final Vector3D center = min.lerp(max, 0.5);
    388
    389                 insertAxisAlignedPartitions(center, precision);
    390
    391                 final int nextLevel = level - 1;
    392                 insertAxisAlignedGridRecursive(min, center, nextLevel, precision);
    393                 insertAxisAlignedGridRecursive(center, max, nextLevel, precision);
    394             }
    395         }
    396
    397         /** Insert a region boundary.
    398          * @param boundary region boundary to insert
    399          * @return this instance
    400          */
    401         public PartitionedRegionBuilder3D insertBoundary(final PlaneConvexSubset boundary) {
    402             insertBoundaryInternal(boundary);
    403
    404             return this;
    405         }
    406
    407         /** Insert a collection of region boundaries.
    408          * @param boundaries boundaries to insert
    409          * @return this instance
    410          */
    411         public PartitionedRegionBuilder3D insertBoundaries(final Iterable<? extends PlaneConvexSubset> boundaries) {
    412             for (final PlaneConvexSubset boundary : boundaries) {
    413                 insertBoundaryInternal(boundary);
    414             }
    415
    416             return this;
    417         }
    418
    419         /** Insert all boundaries from the given source.
    420          * @param boundarySrc source of boundaries to insert
    421          * @return this instance
    422          */
    423         public PartitionedRegionBuilder3D insertBoundaries(final BoundarySource3D boundarySrc) {
    424             try (Stream<PlaneConvexSubset> stream = boundarySrc.boundaryStream()) {
    425                 stream.forEach(this::insertBoundaryInternal);
    426             }
    427
    428             return this;
    429         }
    430
    431         /** Build and return the region BSP tree.
    432          * @return the region BSP tree
    433          */
    434         public RegionBSPTree3D build() {
    435             return (RegionBSPTree3D) buildInternal();
    436         }
    437     }
    438
    439     /** Class used to project points onto the 3D region boundary.
    440      */
    441     private static final class BoundaryProjector3D extends BoundaryProjector<Vector3D, RegionNode3D> {
    442         /** Simple constructor.
    443          * @param point the point to project onto the region's boundary
    444          */
    445         private BoundaryProjector3D(final Vector3D point) {
    446             super(point);
    447         }
    448
    449         /** {@inheritDoc} */
    450         @Override
    451         protected Vector3D disambiguateClosestPoint(final Vector3D target, final Vector3D a, final Vector3D b) {
    452             // return the point with the smallest coordinate values
    453             final int cmp = Vector3D.COORDINATE_ASCENDING_ORDER.compare(a, b);
    454             return cmp < 0 ? a : b;
    455         }
    456     }
    457
    458     /** Visitor for computing geometric properties for 3D BSP tree instances.
    459      *  The volume of the region is computed using the equation
    460      *  <code>V = (1/3)*&Sigma;<sub>F</sub>[(C<sub>F</sub>&sdot;N<sub>F</sub>)*area(F)]</code>,
    461      *  where <code>F</code> represents each face in the region, <code>C<sub>F</sub></code>
    462      *  represents the centroid of the face, and <code>N<sub>F</sub></code> represents the
    463      *  normal of the face. (More details can be found in the article
    464      *  <a href="https://en.wikipedia.org/wiki/Polyhedron#Volume">here</a>.)
    465      *  This essentially splits up the region into pyramids with a 2D face forming
    466      *  the base of each pyramid. The centroid is computed in a similar way. The centroid
    467      *  of each pyramid is calculated using the fact that it is located 3/4 of the way along the
    468      *  line from the apex to the base. The region centroid then becomes the volume-weighted
    469      *  average of these pyramid centers.
    470      *  @see <a href="https://en.wikipedia.org/wiki/Polyhedron#Volume">Polyhedron#Volume</a>
    471      */
    472     private static final class RegionSizePropertiesVisitor implements BSPTreeVisitor<Vector3D, RegionNode3D> {
    473
    474         /** Accumulator for boundary volume contributions. */
    475         private double volumeSum;
    476
    477         /** Centroid contribution x coordinate accumulator. */
    478         private double sumX;
    479
    480         /** Centroid contribution y coordinate accumulator. */
    481         private double sumY;
    482
    483         /** Centroid contribution z coordinate accumulator. */
    484         private double sumZ;
    485
    486         /** {@inheritDoc} */
    487         @Override
    488         public Result visit(final RegionNode3D node) {
    489             if (node.isInternal()) {
    490                 final RegionCutBoundary<Vector3D> boundary = node.getCutBoundary();
    491
    492                 for (final HyperplaneConvexSubset<Vector3D> outsideFacing : boundary.getOutsideFacing()) {
    493                     addBoundaryContribution(outsideFacing, false);
    494                 }
    495
    496                 for (final HyperplaneConvexSubset<Vector3D> insideFacing : boundary.getInsideFacing()) {
    497                     addBoundaryContribution(insideFacing, true);
    498                 }
    499             }
    500
    501             return Result.CONTINUE;
    502         }
    503
    504         /** Return the computed size properties for the visited region.
    505          * @return the computed size properties for the visited region.
    506          */
    507         public RegionSizeProperties<Vector3D> getRegionSizeProperties() {
    508             double size = Double.POSITIVE_INFINITY;
    509             Vector3D centroid = null;
    510
    511             // we only have a finite size if the volume sum is finite and positive
    512             // (negative indicates a finite outside surrounded by an infinite inside)
    513             if (Double.isFinite(volumeSum) && volumeSum > 0.0) {
    514                 // apply the 1/3 pyramid volume scaling factor
    515                 size = volumeSum / 3.0;
    516
    517                 // Since the volume we used when adding together the boundary contributions
    518                 // was 3x the actual pyramid size, we'll multiply by 1/4 here instead
    519                 // of 3/4 to adjust for the actual centroid position in each pyramid.
    520                 final double centroidScale = 1.0 / (4 * size);
    521                 centroid =  Vector3D.of(
    522                         sumX * centroidScale,
    523                         sumY * centroidScale,
    524                         sumZ * centroidScale);
    525             }
    526
    527             return new RegionSizeProperties<>(size, centroid);
    528         }
    529
    530         /** Add the contribution of the given node cut boundary. If {@code reverse} is true,
    531          * the volume of the contribution is reversed before being added to the total.
    532          * @param boundary node cut boundary
    533          * @param reverse if true, the boundary contribution is reversed before being added to the total.
    534          */
    535         private void addBoundaryContribution(final HyperplaneSubset<Vector3D> boundary, final boolean reverse) {
    536             final PlaneSubset boundarySubset = (PlaneSubset) boundary;
    537
    538             final Plane boundaryPlane = boundarySubset.getPlane();
    539             final double boundaryArea = boundarySubset.getSize();
    540             final Vector3D boundaryCentroid = boundarySubset.getCentroid();
    541
    542             if (Double.isInfinite(boundaryArea)) {
    543                 volumeSum = Double.POSITIVE_INFINITY;
    544             } else if (boundaryCentroid != null) {
    545                 // the volume here is actually 3x the actual pyramid volume; we'll apply
    546                 // the final scaling all at once at the end
    547                 double scaledVolume = boundaryArea * boundaryCentroid.dot(boundaryPlane.getNormal());
    548                 if (reverse) {
    549                     scaledVolume = -scaledVolume;
    550                 }
    551
    552                 volumeSum += scaledVolume;
    553
    554                 sumX += scaledVolume * boundaryCentroid.getX();
    555                 sumY += scaledVolume * boundaryCentroid.getY();
    556                 sumZ += scaledVolume * boundaryCentroid.getZ();
    557             }
    558         }
    559     }
    560
    561     /** BSP tree visitor that performs a linecast operation against the boundaries of the visited tree.
    562      */
    563     private static final class LinecastVisitor implements BSPTreeVisitor<Vector3D, RegionNode3D> {
    564
    565         /** The line subset to intersect with the boundaries of the BSP tree. */
    566         private final LineConvexSubset3D linecastSubset;
    567
    568         /** If true, the visitor will stop visiting the tree once the first linecast
    569          * point is determined.
    570          */
    571         private final boolean firstOnly;
    572
    573         /** The minimum abscissa found during the search. */
    574         private double minAbscissa = Double.POSITIVE_INFINITY;
    575
    576         /** List of results from the linecast operation. */
    577         private final List<LinecastPoint3D> results = new ArrayList<>();
    578
    579         /** Create a new instance with the given intersecting line convex subset.
    580          * @param linecastSubset line subset to intersect with the BSP tree region boundary
    581          * @param firstOnly if true, the visitor will stop visiting the tree once the first
    582          *      linecast point is determined
    583          */
    584         LinecastVisitor(final LineConvexSubset3D linecastSubset, final boolean firstOnly) {
    585             this.linecastSubset = linecastSubset;
    586             this.firstOnly = firstOnly;
    587         }
    588
    589         /** Get the first {@link org.apache.commons.geometry.euclidean.twod.LinecastPoint2D}
    590          * resulting from the linecast operation.
    591          * @return the first linecast result point
    592          */
    593         public LinecastPoint3D getFirstResult() {
    594             final List<LinecastPoint3D> sortedResults = getResults();
    595
    596             return sortedResults.isEmpty() ?
    597                     null :
    598                     sortedResults.get(0);
    599         }
    600
    601         /** Get a list containing the results of the linecast operation. The list is
    602          * sorted and filtered.
    603          * @return list of sorted and filtered results from the linecast operation
    604          */
    605         public List<LinecastPoint3D> getResults() {
    606             LinecastPoint3D.sortAndFilter(results);
    607
    608             return results;
    609         }
    610
    611         /** {@inheritDoc} */
    612         @Override
    613         public Order visitOrder(final RegionNode3D internalNode) {
    614             final Plane cut = (Plane) internalNode.getCutHyperplane();
    615             final Line3D line = linecastSubset.getLine();
    616
    617             final boolean plusIsNear = line.getDirection().dot(cut.getNormal()) < 0;
    618
    619             return plusIsNear ?
    620                     Order.PLUS_NODE_MINUS :
    621                     Order.MINUS_NODE_PLUS;
    622         }
    623
    624         /** {@inheritDoc} */
    625         @Override
    626         public Result visit(final RegionNode3D node) {
    627             if (node.isInternal()) {
    628                 // check if the line subset intersects the node cut hyperplane
    629                 final Line3D line = linecastSubset.getLine();
    630                 final Vector3D pt = ((Plane) node.getCutHyperplane()).intersection(line);
    631
    632                 if (pt != null) {
    633                     if (firstOnly && !results.isEmpty() &&
    634                             line.getPrecision().compare(minAbscissa, line.abscissa(pt)) < 0) {
    635                         // we have results and we are now sure that no other intersection points will be
    636                         // found that are closer or at the same position on the intersecting line.
    637                         return Result.TERMINATE;
    638                     } else if (linecastSubset.contains(pt)) {
    639                         // we've potentially found a new linecast point; add it to the list of potential
    640                         // results
    641                         final LinecastPoint3D potentialResult = computeLinecastPoint(pt, node);
    642                         if (potentialResult != null) {
    643                             results.add(potentialResult);
    644
    645                             // update the min abscissa
    646                             minAbscissa = Math.min(minAbscissa, potentialResult.getAbscissa());
    647                         }
    648                     }
    649                 }
    650             }
    651
    652             return Result.CONTINUE;
    653         }
    654
    655         /** Compute the linecast point for the given intersection point and tree node, returning null
    656          * if the point does not actually lie on the region boundary.
    657          * @param pt intersection point
    658          * @param node node containing the cut that the linecast line intersected with
    659          * @return a new linecast point instance or null if the intersection point does not lie
    660          *      on the region boundary
    661          */
    662         private LinecastPoint3D computeLinecastPoint(final Vector3D pt, final RegionNode3D node) {
    663             final Plane cut = (Plane) node.getCutHyperplane();
    664             final RegionCutBoundary<Vector3D> boundary = node.getCutBoundary();
    665
    666             boolean onBoundary = false;
    667             boolean negateNormal = false;
    668
    669             if (boundary.containsInsideFacing(pt)) {
    670                 // on inside-facing boundary
    671                 onBoundary = true;
    672                 negateNormal = true;
    673             } else  if (boundary.containsOutsideFacing(pt)) {
    674                 // on outside-facing boundary
    675                 onBoundary = true;
    676             }
    677
    678             if (onBoundary) {
    679                 Vector3D normal = cut.getNormal();
    680                 if (negateNormal) {
    681                     normal = normal.negate();
    682                 }
    683
    684                 return new LinecastPoint3D(pt, normal, linecastSubset.getLine());
    685             }
    686
    687             return null;
    688         }
    689     }
    690 }