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.math3.geometry.euclidean.threed; 18 19 import java.util.ArrayList; 20 21 import org.apache.commons.math3.geometry.Point; 22 import org.apache.commons.math3.geometry.euclidean.twod.Euclidean2D; 23 import org.apache.commons.math3.geometry.euclidean.twod.PolygonsSet; 24 import org.apache.commons.math3.geometry.euclidean.twod.Vector2D; 25 import org.apache.commons.math3.geometry.partitioning.AbstractSubHyperplane; 26 import org.apache.commons.math3.geometry.partitioning.BSPTree; 27 import org.apache.commons.math3.geometry.partitioning.BSPTreeVisitor; 28 import org.apache.commons.math3.geometry.partitioning.BoundaryAttribute; 29 import org.apache.commons.math3.geometry.partitioning.RegionFactory; 30 import org.apache.commons.math3.geometry.partitioning.SubHyperplane; 31 import org.apache.commons.math3.util.FastMath; 32 33 /** Extractor for {@link PolygonsSet polyhedrons sets} outlines. 34 * <p>This class extracts the 2D outlines from {{@link PolygonsSet 35 * polyhedrons sets} in a specified projection plane.</p> 36 * @since 3.0 37 */ 38 public class OutlineExtractor { 39 40 /** Abscissa axis of the projection plane. */ 41 private Vector3D u; 42 43 /** Ordinate axis of the projection plane. */ 44 private Vector3D v; 45 46 /** Normal of the projection plane (viewing direction). */ 47 private Vector3D w; 48 49 /** Build an extractor for a specific projection plane. 50 * @param u abscissa axis of the projection point 51 * @param v ordinate axis of the projection point 52 */ OutlineExtractor(final Vector3D u, final Vector3D v)53 public OutlineExtractor(final Vector3D u, final Vector3D v) { 54 this.u = u; 55 this.v = v; 56 w = Vector3D.crossProduct(u, v); 57 } 58 59 /** Extract the outline of a polyhedrons set. 60 * @param polyhedronsSet polyhedrons set whose outline must be extracted 61 * @return an outline, as an array of loops. 62 */ getOutline(final PolyhedronsSet polyhedronsSet)63 public Vector2D[][] getOutline(final PolyhedronsSet polyhedronsSet) { 64 65 // project all boundary facets into one polygons set 66 final BoundaryProjector projector = new BoundaryProjector(polyhedronsSet.getTolerance()); 67 polyhedronsSet.getTree(true).visit(projector); 68 final PolygonsSet projected = projector.getProjected(); 69 70 // Remove the spurious intermediate vertices from the outline 71 final Vector2D[][] outline = projected.getVertices(); 72 for (int i = 0; i < outline.length; ++i) { 73 final Vector2D[] rawLoop = outline[i]; 74 int end = rawLoop.length; 75 int j = 0; 76 while (j < end) { 77 if (pointIsBetween(rawLoop, end, j)) { 78 // the point should be removed 79 for (int k = j; k < (end - 1); ++k) { 80 rawLoop[k] = rawLoop[k + 1]; 81 } 82 --end; 83 } else { 84 // the point remains in the loop 85 ++j; 86 } 87 } 88 if (end != rawLoop.length) { 89 // resize the array 90 outline[i] = new Vector2D[end]; 91 System.arraycopy(rawLoop, 0, outline[i], 0, end); 92 } 93 } 94 95 return outline; 96 97 } 98 99 /** Check if a point is geometrically between its neighbor in an array. 100 * <p>The neighbors are computed considering the array is a loop 101 * (i.e. point at index (n-1) is before point at index 0)</p> 102 * @param loop points array 103 * @param n number of points to consider in the array 104 * @param i index of the point to check (must be between 0 and n-1) 105 * @return true if the point is exactly between its neighbors 106 */ pointIsBetween(final Vector2D[] loop, final int n, final int i)107 private boolean pointIsBetween(final Vector2D[] loop, final int n, final int i) { 108 final Vector2D previous = loop[(i + n - 1) % n]; 109 final Vector2D current = loop[i]; 110 final Vector2D next = loop[(i + 1) % n]; 111 final double dx1 = current.getX() - previous.getX(); 112 final double dy1 = current.getY() - previous.getY(); 113 final double dx2 = next.getX() - current.getX(); 114 final double dy2 = next.getY() - current.getY(); 115 final double cross = dx1 * dy2 - dx2 * dy1; 116 final double dot = dx1 * dx2 + dy1 * dy2; 117 final double d1d2 = FastMath.sqrt((dx1 * dx1 + dy1 * dy1) * (dx2 * dx2 + dy2 * dy2)); 118 return (FastMath.abs(cross) <= (1.0e-6 * d1d2)) && (dot >= 0.0); 119 } 120 121 /** Visitor projecting the boundary facets on a plane. */ 122 private class BoundaryProjector implements BSPTreeVisitor<Euclidean3D> { 123 124 /** Projection of the polyhedrons set on the plane. */ 125 private PolygonsSet projected; 126 127 /** Tolerance below which points are considered identical. */ 128 private final double tolerance; 129 130 /** Simple constructor. 131 * @param tolerance tolerance below which points are considered identical 132 */ BoundaryProjector(final double tolerance)133 BoundaryProjector(final double tolerance) { 134 this.projected = new PolygonsSet(new BSPTree<Euclidean2D>(Boolean.FALSE), tolerance); 135 this.tolerance = tolerance; 136 } 137 138 /** {@inheritDoc} */ visitOrder(final BSPTree<Euclidean3D> node)139 public Order visitOrder(final BSPTree<Euclidean3D> node) { 140 return Order.MINUS_SUB_PLUS; 141 } 142 143 /** {@inheritDoc} */ visitInternalNode(final BSPTree<Euclidean3D> node)144 public void visitInternalNode(final BSPTree<Euclidean3D> node) { 145 @SuppressWarnings("unchecked") 146 final BoundaryAttribute<Euclidean3D> attribute = 147 (BoundaryAttribute<Euclidean3D>) node.getAttribute(); 148 if (attribute.getPlusOutside() != null) { 149 addContribution(attribute.getPlusOutside(), false); 150 } 151 if (attribute.getPlusInside() != null) { 152 addContribution(attribute.getPlusInside(), true); 153 } 154 } 155 156 /** {@inheritDoc} */ visitLeafNode(final BSPTree<Euclidean3D> node)157 public void visitLeafNode(final BSPTree<Euclidean3D> node) { 158 } 159 160 /** Add he contribution of a boundary facet. 161 * @param facet boundary facet 162 * @param reversed if true, the facet has the inside on its plus side 163 */ addContribution(final SubHyperplane<Euclidean3D> facet, final boolean reversed)164 private void addContribution(final SubHyperplane<Euclidean3D> facet, final boolean reversed) { 165 166 // extract the vertices of the facet 167 @SuppressWarnings("unchecked") 168 final AbstractSubHyperplane<Euclidean3D, Euclidean2D> absFacet = 169 (AbstractSubHyperplane<Euclidean3D, Euclidean2D>) facet; 170 final Plane plane = (Plane) facet.getHyperplane(); 171 172 final double scal = plane.getNormal().dotProduct(w); 173 if (FastMath.abs(scal) > 1.0e-3) { 174 Vector2D[][] vertices = 175 ((PolygonsSet) absFacet.getRemainingRegion()).getVertices(); 176 177 if ((scal < 0) ^ reversed) { 178 // the facet is seen from the inside, 179 // we need to invert its boundary orientation 180 final Vector2D[][] newVertices = new Vector2D[vertices.length][]; 181 for (int i = 0; i < vertices.length; ++i) { 182 final Vector2D[] loop = vertices[i]; 183 final Vector2D[] newLoop = new Vector2D[loop.length]; 184 if (loop[0] == null) { 185 newLoop[0] = null; 186 for (int j = 1; j < loop.length; ++j) { 187 newLoop[j] = loop[loop.length - j]; 188 } 189 } else { 190 for (int j = 0; j < loop.length; ++j) { 191 newLoop[j] = loop[loop.length - (j + 1)]; 192 } 193 } 194 newVertices[i] = newLoop; 195 } 196 197 // use the reverted vertices 198 vertices = newVertices; 199 200 } 201 202 // compute the projection of the facet in the outline plane 203 final ArrayList<SubHyperplane<Euclidean2D>> edges = new ArrayList<SubHyperplane<Euclidean2D>>(); 204 for (Vector2D[] loop : vertices) { 205 final boolean closed = loop[0] != null; 206 int previous = closed ? (loop.length - 1) : 1; 207 Vector3D previous3D = plane.toSpace((Point<Euclidean2D>) loop[previous]); 208 int current = (previous + 1) % loop.length; 209 Vector2D pPoint = new Vector2D(previous3D.dotProduct(u), 210 previous3D.dotProduct(v)); 211 while (current < loop.length) { 212 213 final Vector3D current3D = plane.toSpace((Point<Euclidean2D>) loop[current]); 214 final Vector2D cPoint = new Vector2D(current3D.dotProduct(u), 215 current3D.dotProduct(v)); 216 final org.apache.commons.math3.geometry.euclidean.twod.Line line = 217 new org.apache.commons.math3.geometry.euclidean.twod.Line(pPoint, cPoint, tolerance); 218 SubHyperplane<Euclidean2D> edge = line.wholeHyperplane(); 219 220 if (closed || (previous != 1)) { 221 // the previous point is a real vertex 222 // it defines one bounding point of the edge 223 final double angle = line.getAngle() + 0.5 * FastMath.PI; 224 final org.apache.commons.math3.geometry.euclidean.twod.Line l = 225 new org.apache.commons.math3.geometry.euclidean.twod.Line(pPoint, angle, tolerance); 226 edge = edge.split(l).getPlus(); 227 } 228 229 if (closed || (current != (loop.length - 1))) { 230 // the current point is a real vertex 231 // it defines one bounding point of the edge 232 final double angle = line.getAngle() + 0.5 * FastMath.PI; 233 final org.apache.commons.math3.geometry.euclidean.twod.Line l = 234 new org.apache.commons.math3.geometry.euclidean.twod.Line(cPoint, angle, tolerance); 235 edge = edge.split(l).getMinus(); 236 } 237 238 edges.add(edge); 239 240 previous = current++; 241 previous3D = current3D; 242 pPoint = cPoint; 243 244 } 245 } 246 final PolygonsSet projectedFacet = new PolygonsSet(edges, tolerance); 247 248 // add the contribution of the facet to the global outline 249 projected = (PolygonsSet) new RegionFactory<Euclidean2D>().union(projected, projectedFacet); 250 251 } 252 } 253 254 /** Get the projection of the polyhedrons set on the plane. 255 * @return projection of the polyhedrons set on the plane 256 */ getProjected()257 public PolygonsSet getProjected() { 258 return projected; 259 } 260 261 } 262 263 } 264