001 /*
002 * Licensed to the Apache Software Foundation (ASF) under one or more
003 * contributor license agreements. See the NOTICE file distributed with
004 * this work for additional information regarding copyright ownership.
005 * The ASF licenses this file to You under the Apache License, Version 2.0
006 * (the "License"); you may not use this file except in compliance with
007 * the License. You may obtain a copy of the License at
008 *
009 * http://www.apache.org/licenses/LICENSE-2.0
010 *
011 * Unless required by applicable law or agreed to in writing, software
012 * distributed under the License is distributed on an "AS IS" BASIS,
013 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
014 * See the License for the specific language governing permissions and
015 * limitations under the License.
016 */
017 package org.apache.commons.math3.geometry.euclidean.threed;
018
019 import java.util.ArrayList;
020
021 import org.apache.commons.math3.geometry.euclidean.twod.Euclidean2D;
022 import org.apache.commons.math3.geometry.euclidean.twod.PolygonsSet;
023 import org.apache.commons.math3.geometry.euclidean.twod.Vector2D;
024 import org.apache.commons.math3.geometry.partitioning.AbstractSubHyperplane;
025 import org.apache.commons.math3.geometry.partitioning.BSPTree;
026 import org.apache.commons.math3.geometry.partitioning.BSPTreeVisitor;
027 import org.apache.commons.math3.geometry.partitioning.BoundaryAttribute;
028 import org.apache.commons.math3.geometry.partitioning.RegionFactory;
029 import org.apache.commons.math3.geometry.partitioning.SubHyperplane;
030 import org.apache.commons.math3.util.FastMath;
031
032 /** Extractor for {@link PolygonsSet polyhedrons sets} outlines.
033 * <p>This class extracts the 2D outlines from {{@link PolygonsSet
034 * polyhedrons sets} in a specified projection plane.</p>
035 * @version $Id: OutlineExtractor.java 1416643 2012-12-03 19:37:14Z tn $
036 * @since 3.0
037 */
038 public class OutlineExtractor {
039
040 /** Abscissa axis of the projection plane. */
041 private Vector3D u;
042
043 /** Ordinate axis of the projection plane. */
044 private Vector3D v;
045
046 /** Normal of the projection plane (viewing direction). */
047 private Vector3D w;
048
049 /** Build an extractor for a specific projection plane.
050 * @param u abscissa axis of the projection point
051 * @param v ordinate axis of the projection point
052 */
053 public OutlineExtractor(final Vector3D u, final Vector3D v) {
054 this.u = u;
055 this.v = v;
056 w = Vector3D.crossProduct(u, v);
057 }
058
059 /** Extract the outline of a polyhedrons set.
060 * @param polyhedronsSet polyhedrons set whose outline must be extracted
061 * @return an outline, as an array of loops.
062 */
063 public Vector2D[][] getOutline(final PolyhedronsSet polyhedronsSet) {
064
065 // project all boundary facets into one polygons set
066 final BoundaryProjector projector = new BoundaryProjector();
067 polyhedronsSet.getTree(true).visit(projector);
068 final PolygonsSet projected = projector.getProjected();
069
070 // Remove the spurious intermediate vertices from the outline
071 final Vector2D[][] outline = projected.getVertices();
072 for (int i = 0; i < outline.length; ++i) {
073 final Vector2D[] rawLoop = outline[i];
074 int end = rawLoop.length;
075 int j = 0;
076 while (j < end) {
077 if (pointIsBetween(rawLoop, end, j)) {
078 // the point should be removed
079 for (int k = j; k < (end - 1); ++k) {
080 rawLoop[k] = rawLoop[k + 1];
081 }
082 --end;
083 } else {
084 // the point remains in the loop
085 ++j;
086 }
087 }
088 if (end != rawLoop.length) {
089 // resize the array
090 outline[i] = new Vector2D[end];
091 System.arraycopy(rawLoop, 0, outline[i], 0, end);
092 }
093 }
094
095 return outline;
096
097 }
098
099 /** Check if a point is geometrically between its neighbour in an array.
100 * <p>The neighbours 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 neighbours
106 */
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 /** Simple constructor.
128 */
129 public BoundaryProjector() {
130 projected = new PolygonsSet(new BSPTree<Euclidean2D>(Boolean.FALSE));
131 }
132
133 /** {@inheritDoc} */
134 public Order visitOrder(final BSPTree<Euclidean3D> node) {
135 return Order.MINUS_SUB_PLUS;
136 }
137
138 /** {@inheritDoc} */
139 public void visitInternalNode(final BSPTree<Euclidean3D> node) {
140 @SuppressWarnings("unchecked")
141 final BoundaryAttribute<Euclidean3D> attribute =
142 (BoundaryAttribute<Euclidean3D>) node.getAttribute();
143 if (attribute.getPlusOutside() != null) {
144 addContribution(attribute.getPlusOutside(), false);
145 }
146 if (attribute.getPlusInside() != null) {
147 addContribution(attribute.getPlusInside(), true);
148 }
149 }
150
151 /** {@inheritDoc} */
152 public void visitLeafNode(final BSPTree<Euclidean3D> node) {
153 }
154
155 /** Add he contribution of a boundary facet.
156 * @param facet boundary facet
157 * @param reversed if true, the facet has the inside on its plus side
158 */
159 private void addContribution(final SubHyperplane<Euclidean3D> facet, final boolean reversed) {
160
161 // extract the vertices of the facet
162 @SuppressWarnings("unchecked")
163 final AbstractSubHyperplane<Euclidean3D, Euclidean2D> absFacet =
164 (AbstractSubHyperplane<Euclidean3D, Euclidean2D>) facet;
165 final Plane plane = (Plane) facet.getHyperplane();
166
167 final double scal = plane.getNormal().dotProduct(w);
168 if (FastMath.abs(scal) > 1.0e-3) {
169 Vector2D[][] vertices =
170 ((PolygonsSet) absFacet.getRemainingRegion()).getVertices();
171
172 if ((scal < 0) ^ reversed) {
173 // the facet is seen from the inside,
174 // we need to invert its boundary orientation
175 final Vector2D[][] newVertices = new Vector2D[vertices.length][];
176 for (int i = 0; i < vertices.length; ++i) {
177 final Vector2D[] loop = vertices[i];
178 final Vector2D[] newLoop = new Vector2D[loop.length];
179 if (loop[0] == null) {
180 newLoop[0] = null;
181 for (int j = 1; j < loop.length; ++j) {
182 newLoop[j] = loop[loop.length - j];
183 }
184 } else {
185 for (int j = 0; j < loop.length; ++j) {
186 newLoop[j] = loop[loop.length - (j + 1)];
187 }
188 }
189 newVertices[i] = newLoop;
190 }
191
192 // use the reverted vertices
193 vertices = newVertices;
194
195 }
196
197 // compute the projection of the facet in the outline plane
198 final ArrayList<SubHyperplane<Euclidean2D>> edges = new ArrayList<SubHyperplane<Euclidean2D>>();
199 for (Vector2D[] loop : vertices) {
200 final boolean closed = loop[0] != null;
201 int previous = closed ? (loop.length - 1) : 1;
202 Vector3D previous3D = plane.toSpace(loop[previous]);
203 int current = (previous + 1) % loop.length;
204 Vector2D pPoint = new Vector2D(previous3D.dotProduct(u),
205 previous3D.dotProduct(v));
206 while (current < loop.length) {
207
208 final Vector3D current3D = plane.toSpace(loop[current]);
209 final Vector2D cPoint = new Vector2D(current3D.dotProduct(u),
210 current3D.dotProduct(v));
211 final org.apache.commons.math3.geometry.euclidean.twod.Line line =
212 new org.apache.commons.math3.geometry.euclidean.twod.Line(pPoint, cPoint);
213 SubHyperplane<Euclidean2D> edge = line.wholeHyperplane();
214
215 if (closed || (previous != 1)) {
216 // the previous point is a real vertex
217 // it defines one bounding point of the edge
218 final double angle = line.getAngle() + 0.5 * FastMath.PI;
219 final org.apache.commons.math3.geometry.euclidean.twod.Line l =
220 new org.apache.commons.math3.geometry.euclidean.twod.Line(pPoint, angle);
221 edge = edge.split(l).getPlus();
222 }
223
224 if (closed || (current != (loop.length - 1))) {
225 // the current point is a real vertex
226 // it defines one bounding point of the edge
227 final double angle = line.getAngle() + 0.5 * FastMath.PI;
228 final org.apache.commons.math3.geometry.euclidean.twod.Line l =
229 new org.apache.commons.math3.geometry.euclidean.twod.Line(cPoint, angle);
230 edge = edge.split(l).getMinus();
231 }
232
233 edges.add(edge);
234
235 previous = current++;
236 previous3D = current3D;
237 pPoint = cPoint;
238
239 }
240 }
241 final PolygonsSet projectedFacet = new PolygonsSet(edges);
242
243 // add the contribution of the facet to the global outline
244 projected = (PolygonsSet) new RegionFactory<Euclidean2D>().union(projected, projectedFacet);
245
246 }
247 }
248
249 /** Get the projection of the polyhedrons set on the plane.
250 * @return projection of the polyhedrons set on the plane
251 */
252 public PolygonsSet getProjected() {
253 return projected;
254 }
255
256 }
257
258 }