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.twod;
018    
019    import java.util.ArrayList;
020    import java.util.Collection;
021    import java.util.List;
022    
023    import org.apache.commons.math3.exception.MathInternalError;
024    import org.apache.commons.math3.geometry.euclidean.oned.Euclidean1D;
025    import org.apache.commons.math3.geometry.euclidean.oned.Interval;
026    import org.apache.commons.math3.geometry.euclidean.oned.IntervalsSet;
027    import org.apache.commons.math3.geometry.euclidean.oned.Vector1D;
028    import org.apache.commons.math3.geometry.partitioning.AbstractRegion;
029    import org.apache.commons.math3.geometry.partitioning.AbstractSubHyperplane;
030    import org.apache.commons.math3.geometry.partitioning.BSPTree;
031    import org.apache.commons.math3.geometry.partitioning.BSPTreeVisitor;
032    import org.apache.commons.math3.geometry.partitioning.BoundaryAttribute;
033    import org.apache.commons.math3.geometry.partitioning.Side;
034    import org.apache.commons.math3.geometry.partitioning.SubHyperplane;
035    import org.apache.commons.math3.geometry.partitioning.utilities.AVLTree;
036    import org.apache.commons.math3.geometry.partitioning.utilities.OrderedTuple;
037    import org.apache.commons.math3.util.FastMath;
038    
039    /** This class represents a 2D region: a set of polygons.
040     * @version $Id: PolygonsSet.java 1422195 2012-12-15 06:45:18Z psteitz $
041     * @since 3.0
042     */
043    public class PolygonsSet extends AbstractRegion<Euclidean2D, Euclidean1D> {
044    
045        /** Vertices organized as boundary loops. */
046        private Vector2D[][] vertices;
047    
048        /** Build a polygons set representing the whole real line.
049         */
050        public PolygonsSet() {
051            super();
052        }
053    
054        /** Build a polygons set from a BSP tree.
055         * <p>The leaf nodes of the BSP tree <em>must</em> have a
056         * {@code Boolean} attribute representing the inside status of
057         * the corresponding cell (true for inside cells, false for outside
058         * cells). In order to avoid building too many small objects, it is
059         * recommended to use the predefined constants
060         * {@code Boolean.TRUE} and {@code Boolean.FALSE}</p>
061         * @param tree inside/outside BSP tree representing the region
062         */
063        public PolygonsSet(final BSPTree<Euclidean2D> tree) {
064            super(tree);
065        }
066    
067        /** Build a polygons set from a Boundary REPresentation (B-rep).
068         * <p>The boundary is provided as a collection of {@link
069         * SubHyperplane sub-hyperplanes}. Each sub-hyperplane has the
070         * interior part of the region on its minus side and the exterior on
071         * its plus side.</p>
072         * <p>The boundary elements can be in any order, and can form
073         * several non-connected sets (like for example polygons with holes
074         * or a set of disjoint polyhedrons considered as a whole). In
075         * fact, the elements do not even need to be connected together
076         * (their topological connections are not used here). However, if the
077         * boundary does not really separate an inside open from an outside
078         * open (open having here its topological meaning), then subsequent
079         * calls to the {@link
080         * org.apache.commons.math3.geometry.partitioning.Region#checkPoint(org.apache.commons.math3.geometry.Vector)
081         * checkPoint} method will not be meaningful anymore.</p>
082         * <p>If the boundary is empty, the region will represent the whole
083         * space.</p>
084         * @param boundary collection of boundary elements, as a
085         * collection of {@link SubHyperplane SubHyperplane} objects
086         */
087        public PolygonsSet(final Collection<SubHyperplane<Euclidean2D>> boundary) {
088            super(boundary);
089        }
090    
091        /** Build a parallellepipedic box.
092         * @param xMin low bound along the x direction
093         * @param xMax high bound along the x direction
094         * @param yMin low bound along the y direction
095         * @param yMax high bound along the y direction
096         */
097        public PolygonsSet(final double xMin, final double xMax,
098                           final double yMin, final double yMax) {
099            super(boxBoundary(xMin, xMax, yMin, yMax));
100        }
101    
102        /** Build a polygon from a simple list of vertices.
103         * <p>The boundary is provided as a list of points considering to
104         * represent the vertices of a simple loop. The interior part of the
105         * region is on the left side of this path and the exterior is on its
106         * right side.</p>
107         * <p>This constructor does not handle polygons with a boundary
108         * forming several disconnected paths (such as polygons with holes).</p>
109         * <p>For cases where this simple constructor applies, it is expected to
110         * be numerically more robust than the {@link #PolygonsSet(Collection) general
111         * constructor} using {@link SubHyperplane subhyperplanes}.</p>
112         * <p>If the list is empty, the region will represent the whole
113         * space.</p>
114         * <p>
115         * Polygons with thin pikes or dents are inherently difficult to handle because
116         * they involve lines with almost opposite directions at some vertices. Polygons
117         * whose vertices come from some physical measurement with noise are also
118         * difficult because an edge that should be straight may be broken in lots of
119         * different pieces with almost equal directions. In both cases, computing the
120         * lines intersections is not numerically robust due to the almost 0 or almost
121         * &pi; angle. Such cases need to carefully adjust the {@code hyperplaneThickness}
122         * parameter. A too small value would often lead to completely wrong polygons
123         * with large area wrongly identified as inside or outside. Large values are
124         * often much safer. As a rule of thumb, a value slightly below the size of the
125         * most accurate detail needed is a good value for the {@code hyperplaneThickness}
126         * parameter.
127         * </p>
128         * @param hyperplaneThickness tolerance below which points are considered to
129         * belong to the hyperplane (which is therefore more a slab)
130         * @param vertices vertices of the simple loop boundary
131         * @since 3.1
132         */
133        public PolygonsSet(final double hyperplaneThickness, final Vector2D ... vertices) {
134            super(verticesToTree(hyperplaneThickness, vertices));
135        }
136    
137        /** Create a list of hyperplanes representing the boundary of a box.
138         * @param xMin low bound along the x direction
139         * @param xMax high bound along the x direction
140         * @param yMin low bound along the y direction
141         * @param yMax high bound along the y direction
142         * @return boundary of the box
143         */
144        private static Line[] boxBoundary(final double xMin, final double xMax,
145                                          final double yMin, final double yMax) {
146            final Vector2D minMin = new Vector2D(xMin, yMin);
147            final Vector2D minMax = new Vector2D(xMin, yMax);
148            final Vector2D maxMin = new Vector2D(xMax, yMin);
149            final Vector2D maxMax = new Vector2D(xMax, yMax);
150            return new Line[] {
151                new Line(minMin, maxMin),
152                new Line(maxMin, maxMax),
153                new Line(maxMax, minMax),
154                new Line(minMax, minMin)
155            };
156        }
157    
158        /** Build the BSP tree of a polygons set from a simple list of vertices.
159         * <p>The boundary is provided as a list of points considering to
160         * represent the vertices of a simple loop. The interior part of the
161         * region is on the left side of this path and the exterior is on its
162         * right side.</p>
163         * <p>This constructor does not handle polygons with a boundary
164         * forming several disconnected paths (such as polygons with holes).</p>
165         * <p>For cases where this simple constructor applies, it is expected to
166         * be numerically more robust than the {@link #PolygonsSet(Collection) general
167         * constructor} using {@link SubHyperplane subhyperplanes}.</p>
168         * @param hyperplaneThickness tolerance below which points are consider to
169         * belong to the hyperplane (which is therefore more a slab)
170         * @param vertices vertices of the simple loop boundary
171         * @return the BSP tree of the input vertices
172         */
173        private static BSPTree<Euclidean2D> verticesToTree(final double hyperplaneThickness,
174                                                           final Vector2D ... vertices) {
175    
176            final int n = vertices.length;
177            if (n == 0) {
178                // the tree represents the whole space
179                return new BSPTree<Euclidean2D>(Boolean.TRUE);
180            }
181    
182            // build the vertices
183            final Vertex[] vArray = new Vertex[n];
184            for (int i = 0; i < n; ++i) {
185                vArray[i] = new Vertex(vertices[i]);
186            }
187    
188            // build the edges
189            List<Edge> edges = new ArrayList<Edge>();
190            for (int i = 0; i < n; ++i) {
191    
192                // get the endpoints of the edge
193                final Vertex start = vArray[i];
194                final Vertex end   = vArray[(i + 1) % n];
195    
196                // get the line supporting the edge, taking care not to recreate it
197                // if it was already created earlier due to another edge being aligned
198                // with the current one
199                Line line = start.sharedLineWith(end);
200                if (line == null) {
201                    line = new Line(start.getLocation(), end.getLocation());
202                }
203    
204                // create the edge and store it
205                edges.add(new Edge(start, end, line));
206    
207                // check if another vertex also happens to be on this line
208                for (final Vertex vertex : vArray) {
209                    if (vertex != start && vertex != end &&
210                        FastMath.abs(line.getOffset(vertex.getLocation())) <= hyperplaneThickness) {
211                        vertex.bindWith(line);
212                    }
213                }
214    
215            }
216    
217            // build the tree top-down
218            final BSPTree<Euclidean2D> tree = new BSPTree<Euclidean2D>();
219            insertEdges(hyperplaneThickness, tree, edges);
220    
221            return tree;
222    
223        }
224    
225        /** Recursively build a tree by inserting cut sub-hyperplanes.
226         * @param hyperplaneThickness tolerance below which points are consider to
227         * belong to the hyperplane (which is therefore more a slab)
228         * @param node current tree node (it is a leaf node at the beginning
229         * of the call)
230         * @param edges list of edges to insert in the cell defined by this node
231         * (excluding edges not belonging to the cell defined by this node)
232         */
233        private static void insertEdges(final double hyperplaneThickness,
234                                        final BSPTree<Euclidean2D> node,
235                                        final List<Edge> edges) {
236    
237            // find an edge with an hyperplane that can be inserted in the node
238            int index = 0;
239            Edge inserted =null;
240            while (inserted == null && index < edges.size()) {
241                inserted = edges.get(index++);
242                if (inserted.getNode() == null) {
243                    if (node.insertCut(inserted.getLine())) {
244                        inserted.setNode(node);
245                    } else {
246                        inserted = null;
247                    }
248                } else {
249                    inserted = null;
250                }
251            }
252    
253            if (inserted == null) {
254                // no suitable edge was found, the node remains a leaf node
255                // we need to set its inside/outside boolean indicator
256                final BSPTree<Euclidean2D> parent = node.getParent();
257                if (parent == null || node == parent.getMinus()) {
258                    node.setAttribute(Boolean.TRUE);
259                } else {
260                    node.setAttribute(Boolean.FALSE);
261                }
262                return;
263            }
264    
265            // we have split the node by inserted an edge as a cut sub-hyperplane
266            // distribute the remaining edges in the two sub-trees
267            final List<Edge> plusList  = new ArrayList<Edge>();
268            final List<Edge> minusList = new ArrayList<Edge>();
269            for (final Edge edge : edges) {
270                if (edge != inserted) {
271                    final double startOffset = inserted.getLine().getOffset(edge.getStart().getLocation());
272                    final double endOffset   = inserted.getLine().getOffset(edge.getEnd().getLocation());
273                    Side startSide = (FastMath.abs(startOffset) <= hyperplaneThickness) ?
274                                     Side.HYPER : ((startOffset < 0) ? Side.MINUS : Side.PLUS);
275                    Side endSide   = (FastMath.abs(endOffset) <= hyperplaneThickness) ?
276                                     Side.HYPER : ((endOffset < 0) ? Side.MINUS : Side.PLUS);
277                    switch (startSide) {
278                        case PLUS:
279                            if (endSide == Side.MINUS) {
280                                // we need to insert a split point on the hyperplane
281                                final Vertex splitPoint = edge.split(inserted.getLine());
282                                minusList.add(splitPoint.getOutgoing());
283                                plusList.add(splitPoint.getIncoming());
284                            } else {
285                                plusList.add(edge);
286                            }
287                            break;
288                        case MINUS:
289                            if (endSide == Side.PLUS) {
290                                // we need to insert a split point on the hyperplane
291                                final Vertex splitPoint = edge.split(inserted.getLine());
292                                minusList.add(splitPoint.getIncoming());
293                                plusList.add(splitPoint.getOutgoing());
294                            } else {
295                                minusList.add(edge);
296                            }
297                            break;
298                        default:
299                            if (endSide == Side.PLUS) {
300                                plusList.add(edge);
301                            } else if (endSide == Side.MINUS) {
302                                minusList.add(edge);
303                            }
304                            break;
305                    }
306                }
307            }
308    
309            // recurse through lower levels
310            if (!plusList.isEmpty()) {
311                insertEdges(hyperplaneThickness, node.getPlus(),  plusList);
312            } else {
313                node.getPlus().setAttribute(Boolean.FALSE);
314            }
315            if (!minusList.isEmpty()) {
316                insertEdges(hyperplaneThickness, node.getMinus(), minusList);
317            } else {
318                node.getMinus().setAttribute(Boolean.TRUE);
319            }
320    
321        }
322    
323        /** Internal class for holding vertices while they are processed to build a BSP tree. */
324        private static class Vertex {
325    
326            /** Vertex location. */
327            private final Vector2D location;
328    
329            /** Incoming edge. */
330            private Edge incoming;
331    
332            /** Outgoing edge. */
333            private Edge outgoing;
334    
335            /** Lines bound with this vertex. */
336            private final List<Line> lines;
337    
338            /** Build a non-processed vertex not owned by any node yet.
339             * @param location vertex location
340             */
341            public Vertex(final Vector2D location) {
342                this.location = location;
343                this.incoming = null;
344                this.outgoing = null;
345                this.lines    = new ArrayList<Line>();
346            }
347    
348            /** Get Vertex location.
349             * @return vertex location
350             */
351            public Vector2D getLocation() {
352                return location;
353            }
354    
355            /** Bind a line considered to contain this vertex.
356             * @param line line to bind with this vertex
357             */
358            public void bindWith(final Line line) {
359                lines.add(line);
360            }
361    
362            /** Get the common line bound with both the instance and another vertex, if any.
363             * <p>
364             * When two vertices are both bound to the same line, this means they are
365             * already handled by node associated with this line, so there is no need
366             * to create a cut hyperplane for them.
367             * </p>
368             * @param vertex other vertex to check instance against
369             * @return line bound with both the instance and another vertex, or null if the
370             * two vertices do not share a line yet
371             */
372            public Line sharedLineWith(final Vertex vertex) {
373                for (final Line line1 : lines) {
374                    for (final Line line2 : vertex.lines) {
375                        if (line1 == line2) {
376                            return line1;
377                        }
378                    }
379                }
380                return null;
381            }
382    
383            /** Set incoming edge.
384             * <p>
385             * The line supporting the incoming edge is automatically bound
386             * with the instance.
387             * </p>
388             * @param incoming incoming edge
389             */
390            public void setIncoming(final Edge incoming) {
391                this.incoming = incoming;
392                bindWith(incoming.getLine());
393            }
394    
395            /** Get incoming edge.
396             * @return incoming edge
397             */
398            public Edge getIncoming() {
399                return incoming;
400            }
401    
402            /** Set outgoing edge.
403             * <p>
404             * The line supporting the outgoing edge is automatically bound
405             * with the instance.
406             * </p>
407             * @param outgoing outgoing edge
408             */
409            public void setOutgoing(final Edge outgoing) {
410                this.outgoing = outgoing;
411                bindWith(outgoing.getLine());
412            }
413    
414            /** Get outgoing edge.
415             * @return outgoing edge
416             */
417            public Edge getOutgoing() {
418                return outgoing;
419            }
420    
421        }
422    
423        /** Internal class for holding edges while they are processed to build a BSP tree. */
424        private static class Edge {
425    
426            /** Start vertex. */
427            private final Vertex start;
428    
429            /** End vertex. */
430            private final Vertex end;
431    
432            /** Line supporting the edge. */
433            private final Line line;
434    
435            /** Node whose cut hyperplane contains this edge. */
436            private BSPTree<Euclidean2D> node;
437    
438            /** Build an edge not contained in any node yet.
439             * @param start start vertex
440             * @param end end vertex
441             * @param line line supporting the edge
442             */
443            public Edge(final Vertex start, final Vertex end, final Line line) {
444    
445                this.start = start;
446                this.end   = end;
447                this.line  = line;
448                this.node  = null;
449    
450                // connect the vertices back to the edge
451                start.setOutgoing(this);
452                end.setIncoming(this);
453    
454            }
455    
456            /** Get start vertex.
457             * @return start vertex
458             */
459            public Vertex getStart() {
460                return start;
461            }
462    
463            /** Get end vertex.
464             * @return end vertex
465             */
466            public Vertex getEnd() {
467                return end;
468            }
469    
470            /** Get the line supporting this edge.
471             * @return line supporting this edge
472             */
473            public Line getLine() {
474                return line;
475            }
476    
477            /** Set the node whose cut hyperplane contains this edge.
478             * @param node node whose cut hyperplane contains this edge
479             */
480            public void setNode(final BSPTree<Euclidean2D> node) {
481                this.node = node;
482            }
483    
484            /** Get the node whose cut hyperplane contains this edge.
485             * @return node whose cut hyperplane contains this edge
486             * (null if edge has not yet been inserted into the BSP tree)
487             */
488            public BSPTree<Euclidean2D> getNode() {
489                return node;
490            }
491    
492            /** Split the edge.
493             * <p>
494             * Once split, this edge is not referenced anymore by the vertices,
495             * it is replaced by the two half-edges and an intermediate splitting
496             * vertex is introduced to connect these two halves.
497             * </p>
498             * @param splitLine line splitting the edge in two halves
499             * @return split vertex (its incoming and outgoing edges are the two halves)
500             */
501            public Vertex split(final Line splitLine) {
502                final Vertex splitVertex = new Vertex(line.intersection(splitLine));
503                splitVertex.bindWith(splitLine);
504                final Edge startHalf = new Edge(start, splitVertex, line);
505                final Edge endHalf   = new Edge(splitVertex, end, line);
506                startHalf.node = node;
507                endHalf.node   = node;
508                return splitVertex;
509            }
510    
511        }
512    
513        /** {@inheritDoc} */
514        @Override
515        public PolygonsSet buildNew(final BSPTree<Euclidean2D> tree) {
516            return new PolygonsSet(tree);
517        }
518    
519        /** {@inheritDoc} */
520        @Override
521        protected void computeGeometricalProperties() {
522    
523            final Vector2D[][] v = getVertices();
524    
525            if (v.length == 0) {
526                final BSPTree<Euclidean2D> tree = getTree(false);
527                if (tree.getCut() == null && (Boolean) tree.getAttribute()) {
528                    // the instance covers the whole space
529                    setSize(Double.POSITIVE_INFINITY);
530                    setBarycenter(Vector2D.NaN);
531                } else {
532                    setSize(0);
533                    setBarycenter(new Vector2D(0, 0));
534                }
535            } else if (v[0][0] == null) {
536                // there is at least one open-loop: the polygon is infinite
537                setSize(Double.POSITIVE_INFINITY);
538                setBarycenter(Vector2D.NaN);
539            } else {
540                // all loops are closed, we compute some integrals around the shape
541    
542                double sum  = 0;
543                double sumX = 0;
544                double sumY = 0;
545    
546                for (Vector2D[] loop : v) {
547                    double x1 = loop[loop.length - 1].getX();
548                    double y1 = loop[loop.length - 1].getY();
549                    for (final Vector2D point : loop) {
550                        final double x0 = x1;
551                        final double y0 = y1;
552                        x1 = point.getX();
553                        y1 = point.getY();
554                        final double factor = x0 * y1 - y0 * x1;
555                        sum  += factor;
556                        sumX += factor * (x0 + x1);
557                        sumY += factor * (y0 + y1);
558                    }
559                }
560    
561                if (sum < 0) {
562                    // the polygon as a finite outside surrounded by an infinite inside
563                    setSize(Double.POSITIVE_INFINITY);
564                    setBarycenter(Vector2D.NaN);
565                } else {
566                    setSize(sum / 2);
567                    setBarycenter(new Vector2D(sumX / (3 * sum), sumY / (3 * sum)));
568                }
569    
570            }
571    
572        }
573    
574        /** Get the vertices of the polygon.
575         * <p>The polygon boundary can be represented as an array of loops,
576         * each loop being itself an array of vertices.</p>
577         * <p>In order to identify open loops which start and end by
578         * infinite edges, the open loops arrays start with a null point. In
579         * this case, the first non null point and the last point of the
580         * array do not represent real vertices, they are dummy points
581         * intended only to get the direction of the first and last edge. An
582         * open loop consisting of a single infinite line will therefore be
583         * represented by a three elements array with one null point
584         * followed by two dummy points. The open loops are always the first
585         * ones in the loops array.</p>
586         * <p>If the polygon has no boundary at all, a zero length loop
587         * array will be returned.</p>
588         * <p>All line segments in the various loops have the inside of the
589         * region on their left side and the outside on their right side
590         * when moving in the underlying line direction. This means that
591         * closed loops surrounding finite areas obey the direct
592         * trigonometric orientation.</p>
593         * @return vertices of the polygon, organized as oriented boundary
594         * loops with the open loops first (the returned value is guaranteed
595         * to be non-null)
596         */
597        public Vector2D[][] getVertices() {
598            if (vertices == null) {
599                if (getTree(false).getCut() == null) {
600                    vertices = new Vector2D[0][];
601                } else {
602    
603                    // sort the segments according to their start point
604                    final SegmentsBuilder visitor = new SegmentsBuilder();
605                    getTree(true).visit(visitor);
606                    final AVLTree<ComparableSegment> sorted = visitor.getSorted();
607    
608                    // identify the loops, starting from the open ones
609                    // (their start segments are naturally at the sorted set beginning)
610                    final ArrayList<List<ComparableSegment>> loops = new ArrayList<List<ComparableSegment>>();
611                    while (!sorted.isEmpty()) {
612                        final AVLTree<ComparableSegment>.Node node = sorted.getSmallest();
613                        final List<ComparableSegment> loop = followLoop(node, sorted);
614                        if (loop != null) {
615                            loops.add(loop);
616                        }
617                    }
618    
619                    // tranform the loops in an array of arrays of points
620                    vertices = new Vector2D[loops.size()][];
621                    int i = 0;
622    
623                    for (final List<ComparableSegment> loop : loops) {
624                        if (loop.size() < 2) {
625                            // single infinite line
626                            final Line line = loop.get(0).getLine();
627                            vertices[i++] = new Vector2D[] {
628                                null,
629                                line.toSpace(new Vector1D(-Float.MAX_VALUE)),
630                                line.toSpace(new Vector1D(+Float.MAX_VALUE))
631                            };
632                        } else if (loop.get(0).getStart() == null) {
633                            // open loop with at least one real point
634                            final Vector2D[] array = new Vector2D[loop.size() + 2];
635                            int j = 0;
636                            for (Segment segment : loop) {
637    
638                                if (j == 0) {
639                                    // null point and first dummy point
640                                    double x = segment.getLine().toSubSpace(segment.getEnd()).getX();
641                                    x -= FastMath.max(1.0, FastMath.abs(x / 2));
642                                    array[j++] = null;
643                                    array[j++] = segment.getLine().toSpace(new Vector1D(x));
644                                }
645    
646                                if (j < (array.length - 1)) {
647                                    // current point
648                                    array[j++] = segment.getEnd();
649                                }
650    
651                                if (j == (array.length - 1)) {
652                                    // last dummy point
653                                    double x = segment.getLine().toSubSpace(segment.getStart()).getX();
654                                    x += FastMath.max(1.0, FastMath.abs(x / 2));
655                                    array[j++] = segment.getLine().toSpace(new Vector1D(x));
656                                }
657    
658                            }
659                            vertices[i++] = array;
660                        } else {
661                            final Vector2D[] array = new Vector2D[loop.size()];
662                            int j = 0;
663                            for (Segment segment : loop) {
664                                array[j++] = segment.getStart();
665                            }
666                            vertices[i++] = array;
667                        }
668                    }
669    
670                }
671            }
672    
673            return vertices.clone();
674    
675        }
676    
677        /** Follow a boundary loop.
678         * @param node node containing the segment starting the loop
679         * @param sorted set of segments belonging to the boundary, sorted by
680         * start points (contains {@code node})
681         * @return a list of connected sub-hyperplanes starting at
682         * {@code node}
683         */
684        private List<ComparableSegment> followLoop(final AVLTree<ComparableSegment>.Node node,
685                                                   final AVLTree<ComparableSegment> sorted) {
686    
687            final ArrayList<ComparableSegment> loop = new ArrayList<ComparableSegment>();
688            ComparableSegment segment = node.getElement();
689            loop.add(segment);
690            final Vector2D globalStart = segment.getStart();
691            Vector2D end = segment.getEnd();
692            node.delete();
693    
694            // is this an open or a closed loop ?
695            final boolean open = segment.getStart() == null;
696    
697            while ((end != null) && (open || (globalStart.distance(end) > 1.0e-10))) {
698    
699                // search the sub-hyperplane starting where the previous one ended
700                AVLTree<ComparableSegment>.Node selectedNode = null;
701                ComparableSegment       selectedSegment  = null;
702                double                  selectedDistance = Double.POSITIVE_INFINITY;
703                final ComparableSegment lowerLeft        = new ComparableSegment(end, -1.0e-10, -1.0e-10);
704                final ComparableSegment upperRight       = new ComparableSegment(end, +1.0e-10, +1.0e-10);
705                for (AVLTree<ComparableSegment>.Node n = sorted.getNotSmaller(lowerLeft);
706                     (n != null) && (n.getElement().compareTo(upperRight) <= 0);
707                     n = n.getNext()) {
708                    segment = n.getElement();
709                    final double distance = end.distance(segment.getStart());
710                    if (distance < selectedDistance) {
711                        selectedNode     = n;
712                        selectedSegment  = segment;
713                        selectedDistance = distance;
714                    }
715                }
716    
717                if (selectedDistance > 1.0e-10) {
718                    // this is a degenerated loop, it probably comes from a very
719                    // tiny region with some segments smaller than the threshold, we
720                    // simply ignore it
721                    return null;
722                }
723    
724                end = selectedSegment.getEnd();
725                loop.add(selectedSegment);
726                selectedNode.delete();
727    
728            }
729    
730            if ((loop.size() == 2) && !open) {
731                // this is a degenerated infinitely thin loop, we simply ignore it
732                return null;
733            }
734    
735            if ((end == null) && !open) {
736                throw new MathInternalError();
737            }
738    
739            return loop;
740    
741        }
742    
743        /** Private extension of Segment allowing comparison. */
744        private static class ComparableSegment extends Segment implements Comparable<ComparableSegment> {
745    
746            /** Sorting key. */
747            private OrderedTuple sortingKey;
748    
749            /** Build a segment.
750             * @param start start point of the segment
751             * @param end end point of the segment
752             * @param line line containing the segment
753             */
754            public ComparableSegment(final Vector2D start, final Vector2D end, final Line line) {
755                super(start, end, line);
756                sortingKey = (start == null) ?
757                             new OrderedTuple(Double.NEGATIVE_INFINITY, Double.NEGATIVE_INFINITY) :
758                             new OrderedTuple(start.getX(), start.getY());
759            }
760    
761            /** Build a dummy segment.
762             * <p>
763             * The object built is not a real segment, only the sorting key is used to
764             * allow searching in the neighborhood of a point. This is an horrible hack ...
765             * </p>
766             * @param start start point of the segment
767             * @param dx abscissa offset from the start point
768             * @param dy ordinate offset from the start point
769             */
770            public ComparableSegment(final Vector2D start, final double dx, final double dy) {
771                super(null, null, null);
772                sortingKey = new OrderedTuple(start.getX() + dx, start.getY() + dy);
773            }
774    
775            /** {@inheritDoc} */
776            public int compareTo(final ComparableSegment o) {
777                return sortingKey.compareTo(o.sortingKey);
778            }
779    
780            /** {@inheritDoc} */
781            @Override
782            public boolean equals(final Object other) {
783                if (this == other) {
784                    return true;
785                } else if (other instanceof ComparableSegment) {
786                    return compareTo((ComparableSegment) other) == 0;
787                } else {
788                    return false;
789                }
790            }
791    
792            /** {@inheritDoc} */
793            @Override
794            public int hashCode() {
795                return getStart().hashCode() ^ getEnd().hashCode() ^
796                       getLine().hashCode() ^ sortingKey.hashCode();
797            }
798    
799        }
800    
801        /** Visitor building segments. */
802        private static class SegmentsBuilder implements BSPTreeVisitor<Euclidean2D> {
803    
804            /** Sorted segments. */
805            private AVLTree<ComparableSegment> sorted;
806    
807            /** Simple constructor. */
808            public SegmentsBuilder() {
809                sorted = new AVLTree<ComparableSegment>();
810            }
811    
812            /** {@inheritDoc} */
813            public Order visitOrder(final BSPTree<Euclidean2D> node) {
814                return Order.MINUS_SUB_PLUS;
815            }
816    
817            /** {@inheritDoc} */
818            public void visitInternalNode(final BSPTree<Euclidean2D> node) {
819                @SuppressWarnings("unchecked")
820                final BoundaryAttribute<Euclidean2D> attribute = (BoundaryAttribute<Euclidean2D>) node.getAttribute();
821                if (attribute.getPlusOutside() != null) {
822                    addContribution(attribute.getPlusOutside(), false);
823                }
824                if (attribute.getPlusInside() != null) {
825                    addContribution(attribute.getPlusInside(), true);
826                }
827            }
828    
829            /** {@inheritDoc} */
830            public void visitLeafNode(final BSPTree<Euclidean2D> node) {
831            }
832    
833            /** Add he contribution of a boundary facet.
834             * @param sub boundary facet
835             * @param reversed if true, the facet has the inside on its plus side
836             */
837            private void addContribution(final SubHyperplane<Euclidean2D> sub, final boolean reversed) {
838                @SuppressWarnings("unchecked")
839                final AbstractSubHyperplane<Euclidean2D, Euclidean1D> absSub =
840                    (AbstractSubHyperplane<Euclidean2D, Euclidean1D>) sub;
841                final Line line      = (Line) sub.getHyperplane();
842                final List<Interval> intervals = ((IntervalsSet) absSub.getRemainingRegion()).asList();
843                for (final Interval i : intervals) {
844                    final Vector2D start = Double.isInfinite(i.getInf()) ?
845                                          null : (Vector2D) line.toSpace(new Vector1D(i.getInf()));
846                    final Vector2D end   = Double.isInfinite(i.getSup()) ?
847                                          null : (Vector2D) line.toSpace(new Vector1D(i.getSup()));
848                    if (reversed) {
849                        sorted.insert(new ComparableSegment(end, start, line.getReverse()));
850                    } else {
851                        sorted.insert(new ComparableSegment(start, end, line));
852                    }
853                }
854            }
855    
856            /** Get the sorted segments.
857             * @return sorted segments
858             */
859            public AVLTree<ComparableSegment> getSorted() {
860                return sorted;
861            }
862    
863        }
864    
865    }