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 * https://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
18 /*
19 * This is not the original file distributed by the Apache Software Foundation
20 * It has been modified by the Hipparchus project
21 */
22 package org.hipparchus.geometry.euclidean.oned;
23
24 import java.util.ArrayList;
25 import java.util.Collection;
26 import java.util.Iterator;
27 import java.util.List;
28 import java.util.NoSuchElementException;
29
30 import org.hipparchus.geometry.Point;
31 import org.hipparchus.geometry.partitioning.AbstractRegion;
32 import org.hipparchus.geometry.partitioning.BSPTree;
33 import org.hipparchus.geometry.partitioning.BoundaryProjection;
34 import org.hipparchus.geometry.partitioning.SubHyperplane;
35 import org.hipparchus.util.Precision;
36
37 /** This class represents a 1D region: a set of intervals.
38 */
39 public class IntervalsSet extends AbstractRegion<Euclidean1D, Euclidean1D> implements Iterable<double[]> {
40
41 /** Build an intervals set representing the whole real line.
42 * @param tolerance tolerance below which points are considered identical.
43 */
44 public IntervalsSet(final double tolerance) {
45 super(tolerance);
46 }
47
48 /** Build an intervals set corresponding to a single interval.
49 * @param lower lower bound of the interval, must be lesser or equal
50 * to {@code upper} (may be {@code Double.NEGATIVE_INFINITY})
51 * @param upper upper bound of the interval, must be greater or equal
52 * to {@code lower} (may be {@code Double.POSITIVE_INFINITY})
53 * @param tolerance tolerance below which points are considered identical.
54 */
55 public IntervalsSet(final double lower, final double upper, final double tolerance) {
56 super(buildTree(lower, upper, tolerance), tolerance);
57 }
58
59 /** Build an intervals set from an inside/outside BSP tree.
60 * <p>The leaf nodes of the BSP tree <em>must</em> have a
61 * {@code Boolean} attribute representing the inside status of
62 * the corresponding cell (true for inside cells, false for outside
63 * cells). In order to avoid building too many small objects, it is
64 * recommended to use the predefined constants
65 * {@code Boolean.TRUE} and {@code Boolean.FALSE}</p>
66 * @param tree inside/outside BSP tree representing the intervals set
67 * @param tolerance tolerance below which points are considered identical.
68 */
69 public IntervalsSet(final BSPTree<Euclidean1D> tree, final double tolerance) {
70 super(tree, tolerance);
71 }
72
73 /** Build an intervals set from a Boundary REPresentation (B-rep).
74 * <p>The boundary is provided as a collection of {@link
75 * SubHyperplane sub-hyperplanes}. Each sub-hyperplane has the
76 * interior part of the region on its minus side and the exterior on
77 * its plus side.</p>
78 * <p>The boundary elements can be in any order, and can form
79 * several non-connected sets (like for example polygons with holes
80 * or a set of disjoints polyhedrons considered as a whole). In
81 * fact, the elements do not even need to be connected together
82 * (their topological connections are not used here). However, if the
83 * boundary does not really separate an inside open from an outside
84 * open (open having here its topological meaning), then subsequent
85 * calls to the {@link
86 * org.hipparchus.geometry.partitioning.Region#checkPoint(org.hipparchus.geometry.Point)
87 * checkPoint} method will not be meaningful anymore.</p>
88 * <p>If the boundary is empty, the region will represent the whole
89 * space.</p>
90 * @param boundary collection of boundary elements
91 * @param tolerance tolerance below which points are considered identical.
92 */
93 public IntervalsSet(final Collection<SubHyperplane<Euclidean1D>> boundary,
94 final double tolerance) {
95 super(boundary, tolerance);
96 }
97
98 /** Build an inside/outside tree representing a single interval.
99 * @param lower lower bound of the interval, must be lesser or equal
100 * to {@code upper} (may be {@code Double.NEGATIVE_INFINITY})
101 * @param upper upper bound of the interval, must be greater or equal
102 * to {@code lower} (may be {@code Double.POSITIVE_INFINITY})
103 * @param tolerance tolerance below which points are considered identical.
104 * @return the built tree
105 */
106 private static BSPTree<Euclidean1D> buildTree(final double lower, final double upper,
107 final double tolerance) {
108 if (Double.isInfinite(lower) && (lower < 0)) {
109 if (Double.isInfinite(upper) && (upper > 0)) {
110 // the tree must cover the whole real line
111 return new BSPTree<Euclidean1D>(Boolean.TRUE);
112 }
113 // the tree must be open on the negative infinity side
114 final SubHyperplane<Euclidean1D> upperCut =
115 new OrientedPoint(new Vector1D(upper), true, tolerance).wholeHyperplane();
116 return new BSPTree<Euclidean1D>(upperCut,
117 new BSPTree<Euclidean1D>(Boolean.FALSE),
118 new BSPTree<Euclidean1D>(Boolean.TRUE),
119 null);
120 }
121 final SubHyperplane<Euclidean1D> lowerCut =
122 new OrientedPoint(new Vector1D(lower), false, tolerance).wholeHyperplane();
123 if (Double.isInfinite(upper) && (upper > 0)) {
124 // the tree must be open on the positive infinity side
125 return new BSPTree<Euclidean1D>(lowerCut,
126 new BSPTree<Euclidean1D>(Boolean.FALSE),
127 new BSPTree<Euclidean1D>(Boolean.TRUE),
128 null);
129 }
130
131 // the tree must be bounded on the two sides
132 final SubHyperplane<Euclidean1D> upperCut =
133 new OrientedPoint(new Vector1D(upper), true, tolerance).wholeHyperplane();
134 return new BSPTree<Euclidean1D>(lowerCut,
135 new BSPTree<Euclidean1D>(Boolean.FALSE),
136 new BSPTree<Euclidean1D>(upperCut,
137 new BSPTree<Euclidean1D>(Boolean.FALSE),
138 new BSPTree<Euclidean1D>(Boolean.TRUE),
139 null),
140 null);
141
142 }
143
144 /** {@inheritDoc} */
145 @Override
146 public IntervalsSet buildNew(final BSPTree<Euclidean1D> tree) {
147 return new IntervalsSet(tree, getTolerance());
148 }
149
150 /** {@inheritDoc} */
151 @Override
152 protected void computeGeometricalProperties() {
153 if (getTree(false).getCut() == null) {
154 setBarycenter((Point<Euclidean1D>) Vector1D.NaN);
155 setSize(((Boolean) getTree(false).getAttribute()) ? Double.POSITIVE_INFINITY : 0);
156 } else {
157 double size = 0.0;
158 double sum = 0.0;
159 for (final Interval interval : asList()) {
160 size += interval.getSize();
161 sum += interval.getSize() * interval.getBarycenter();
162 }
163 setSize(size);
164 if (Double.isInfinite(size)) {
165 setBarycenter((Point<Euclidean1D>) Vector1D.NaN);
166 } else if (size >= Precision.SAFE_MIN) {
167 setBarycenter((Point<Euclidean1D>) new Vector1D(sum / size));
168 } else {
169 setBarycenter((Point<Euclidean1D>) ((OrientedPoint) getTree(false).getCut().getHyperplane()).getLocation());
170 }
171 }
172 }
173
174 /** Get the lowest value belonging to the instance.
175 * @return lowest value belonging to the instance
176 * ({@code Double.NEGATIVE_INFINITY} if the instance doesn't
177 * have any low bound, {@code Double.POSITIVE_INFINITY} if the
178 * instance is empty)
179 */
180 public double getInf() {
181 BSPTree<Euclidean1D> node = getTree(false);
182 double inf = Double.POSITIVE_INFINITY;
183 while (node.getCut() != null) {
184 final OrientedPoint op = (OrientedPoint) node.getCut().getHyperplane();
185 inf = op.getLocation().getX();
186 node = op.isDirect() ? node.getMinus() : node.getPlus();
187 }
188 return ((Boolean) node.getAttribute()) ? Double.NEGATIVE_INFINITY : inf;
189 }
190
191 /** Get the highest value belonging to the instance.
192 * @return highest value belonging to the instance
193 * ({@code Double.POSITIVE_INFINITY} if the instance doesn't
194 * have any high bound, {@code Double.NEGATIVE_INFINITY} if the
195 * instance is empty)
196 */
197 public double getSup() {
198 BSPTree<Euclidean1D> node = getTree(false);
199 double sup = Double.NEGATIVE_INFINITY;
200 while (node.getCut() != null) {
201 final OrientedPoint op = (OrientedPoint) node.getCut().getHyperplane();
202 sup = op.getLocation().getX();
203 node = op.isDirect() ? node.getPlus() : node.getMinus();
204 }
205 return ((Boolean) node.getAttribute()) ? Double.POSITIVE_INFINITY : sup;
206 }
207
208 /** {@inheritDoc}
209 */
210 @Override
211 public BoundaryProjection<Euclidean1D> projectToBoundary(final Point<Euclidean1D> point) {
212
213 // get position of test point
214 final double x = ((Vector1D) point).getX();
215
216 double previous = Double.NEGATIVE_INFINITY;
217 for (final double[] a : this) {
218 if (x < a[0]) {
219 // the test point lies between the previous and the current intervals
220 // offset will be positive
221 final double previousOffset = x - previous;
222 final double currentOffset = a[0] - x;
223 if (previousOffset < currentOffset) {
224 return new BoundaryProjection<Euclidean1D>(point, finiteOrNullPoint(previous), previousOffset);
225 } else {
226 return new BoundaryProjection<Euclidean1D>(point, finiteOrNullPoint(a[0]), currentOffset);
227 }
228 } else if (x <= a[1]) {
229 // the test point lies within the current interval
230 // offset will be negative
231 final double offset0 = a[0] - x;
232 final double offset1 = x - a[1];
233 if (offset0 < offset1) {
234 return new BoundaryProjection<Euclidean1D>(point, finiteOrNullPoint(a[1]), offset1);
235 } else {
236 return new BoundaryProjection<Euclidean1D>(point, finiteOrNullPoint(a[0]), offset0);
237 }
238 }
239 previous = a[1];
240 }
241
242 // the test point if past the last sub-interval
243 return new BoundaryProjection<Euclidean1D>(point, finiteOrNullPoint(previous), x - previous);
244
245 }
246
247 /** Build a finite point.
248 * @param x abscissa of the point
249 * @return a new point for finite abscissa, null otherwise
250 */
251 private Vector1D finiteOrNullPoint(final double x) {
252 return Double.isInfinite(x) ? null : new Vector1D(x);
253 }
254
255 /** Build an ordered list of intervals representing the instance.
256 * <p>This method builds this intervals set as an ordered list of
257 * {@link Interval Interval} elements. If the intervals set has no
258 * lower limit, the first interval will have its low bound equal to
259 * {@code Double.NEGATIVE_INFINITY}. If the intervals set has
260 * no upper limit, the last interval will have its upper bound equal
261 * to {@code Double.POSITIVE_INFINITY}. An empty tree will
262 * build an empty list while a tree representing the whole real line
263 * will build a one element list with both bounds being
264 * infinite.</p>
265 * @return a new ordered list containing {@link Interval Interval}
266 * elements
267 */
268 public List<Interval> asList() {
269 final List<Interval> list = new ArrayList<>();
270 for (final double[] a : this) {
271 list.add(new Interval(a[0], a[1]));
272 }
273 return list;
274 }
275
276 /** Get the first leaf node of a tree.
277 * @param root tree root
278 * @return first leaf node
279 */
280 private BSPTree<Euclidean1D> getFirstLeaf(final BSPTree<Euclidean1D> root) {
281
282 if (root.getCut() == null) {
283 return root;
284 }
285
286 // find the smallest internal node
287 BSPTree<Euclidean1D> smallest = null;
288 for (BSPTree<Euclidean1D> n = root; n != null; n = previousInternalNode(n)) {
289 smallest = n;
290 }
291
292 return leafBefore(smallest);
293
294 }
295
296 /** Get the node corresponding to the first interval boundary.
297 * @return smallest internal node,
298 * or null if there are no internal nodes (i.e. the set is either empty or covers the real line)
299 */
300 private BSPTree<Euclidean1D> getFirstIntervalBoundary() {
301
302 // start search at the tree root
303 BSPTree<Euclidean1D> node = getTree(false);
304 if (node.getCut() == null) {
305 return null;
306 }
307
308 // walk tree until we find the smallest internal node
309 node = getFirstLeaf(node).getParent();
310
311 // walk tree until we find an interval boundary
312 while (node != null && !(isIntervalStart(node) || isIntervalEnd(node))) {
313 node = nextInternalNode(node);
314 }
315
316 return node;
317
318 }
319
320 /** Check if an internal node corresponds to the start abscissa of an interval.
321 * @param node internal node to check
322 * @return true if the node corresponds to the start abscissa of an interval
323 */
324 private boolean isIntervalStart(final BSPTree<Euclidean1D> node) {
325
326 if ((Boolean) leafBefore(node).getAttribute()) {
327 // it has an inside cell before it, it may end an interval but not start it
328 return false;
329 }
330
331 if (!(Boolean) leafAfter(node).getAttribute()) {
332 // it has an outside cell after it, it is a dummy cut away from real intervals
333 return false;
334 }
335
336 // the cell has an outside before and an inside after it
337 // it is the start of an interval
338 return true;
339
340 }
341
342 /** Check if an internal node corresponds to the end abscissa of an interval.
343 * @param node internal node to check
344 * @return true if the node corresponds to the end abscissa of an interval
345 */
346 private boolean isIntervalEnd(final BSPTree<Euclidean1D> node) {
347
348 if (!(Boolean) leafBefore(node).getAttribute()) {
349 // it has an outside cell before it, it may start an interval but not end it
350 return false;
351 }
352
353 if ((Boolean) leafAfter(node).getAttribute()) {
354 // it has an inside cell after it, it is a dummy cut in the middle of an interval
355 return false;
356 }
357
358 // the cell has an inside before and an outside after it
359 // it is the end of an interval
360 return true;
361
362 }
363
364 /** Get the next internal node.
365 * @param node current internal node
366 * @return next internal node in ascending order, or null
367 * if this is the last internal node
368 */
369 private BSPTree<Euclidean1D> nextInternalNode(BSPTree<Euclidean1D> node) {
370
371 if (childAfter(node).getCut() != null) {
372 // the next node is in the sub-tree
373 return leafAfter(node).getParent();
374 }
375
376 // there is nothing left deeper in the tree, we backtrack
377 while (isAfterParent(node)) {
378 node = node.getParent();
379 }
380 return node.getParent();
381
382 }
383
384 /** Get the previous internal node.
385 * @param node current internal node
386 * @return previous internal node in ascending order, or null
387 * if this is the first internal node
388 */
389 private BSPTree<Euclidean1D> previousInternalNode(BSPTree<Euclidean1D> node) {
390
391 if (childBefore(node).getCut() != null) {
392 // the next node is in the sub-tree
393 return leafBefore(node).getParent();
394 }
395
396 // there is nothing left deeper in the tree, we backtrack
397 while (isBeforeParent(node)) {
398 node = node.getParent();
399 }
400 return node.getParent();
401
402 }
403
404 /** Find the leaf node just before an internal node.
405 * @param node internal node at which the sub-tree starts
406 * @return leaf node just before the internal node
407 */
408 private BSPTree<Euclidean1D> leafBefore(BSPTree<Euclidean1D> node) {
409
410 node = childBefore(node);
411 while (node.getCut() != null) {
412 node = childAfter(node);
413 }
414
415 return node;
416
417 }
418
419 /** Find the leaf node just after an internal node.
420 * @param node internal node at which the sub-tree starts
421 * @return leaf node just after the internal node
422 */
423 private BSPTree<Euclidean1D> leafAfter(BSPTree<Euclidean1D> node) {
424
425 node = childAfter(node);
426 while (node.getCut() != null) {
427 node = childBefore(node);
428 }
429
430 return node;
431
432 }
433
434 /** Check if a node is the child before its parent in ascending order.
435 * @param node child node considered
436 * @return true is the node has a parent end is before it in ascending order
437 */
438 private boolean isBeforeParent(final BSPTree<Euclidean1D> node) {
439 final BSPTree<Euclidean1D> parent = node.getParent();
440 if (parent == null) {
441 return false;
442 } else {
443 return node == childBefore(parent);
444 }
445 }
446
447 /** Check if a node is the child after its parent in ascending order.
448 * @param node child node considered
449 * @return true is the node has a parent end is after it in ascending order
450 */
451 private boolean isAfterParent(final BSPTree<Euclidean1D> node) {
452 final BSPTree<Euclidean1D> parent = node.getParent();
453 if (parent == null) {
454 return false;
455 } else {
456 return node == childAfter(parent);
457 }
458 }
459
460 /** Find the child node just before an internal node.
461 * @param node internal node at which the sub-tree starts
462 * @return child node just before the internal node
463 */
464 private BSPTree<Euclidean1D> childBefore(BSPTree<Euclidean1D> node) {
465 if (isDirect(node)) {
466 // smaller abscissas are on minus side, larger abscissas are on plus side
467 return node.getMinus();
468 } else {
469 // smaller abscissas are on plus side, larger abscissas are on minus side
470 return node.getPlus();
471 }
472 }
473
474 /** Find the child node just after an internal node.
475 * @param node internal node at which the sub-tree starts
476 * @return child node just after the internal node
477 */
478 private BSPTree<Euclidean1D> childAfter(BSPTree<Euclidean1D> node) {
479 if (isDirect(node)) {
480 // smaller abscissas are on minus side, larger abscissas are on plus side
481 return node.getPlus();
482 } else {
483 // smaller abscissas are on plus side, larger abscissas are on minus side
484 return node.getMinus();
485 }
486 }
487
488 /** Check if an internal node has a direct oriented point.
489 * @param node internal node to check
490 * @return true if the oriented point is direct
491 */
492 private boolean isDirect(final BSPTree<Euclidean1D> node) {
493 return ((OrientedPoint) node.getCut().getHyperplane()).isDirect();
494 }
495
496 /** Get the abscissa of an internal node.
497 * @param node internal node to check
498 * @return abscissa
499 */
500 private double getAngle(final BSPTree<Euclidean1D> node) {
501 return ((OrientedPoint) node.getCut().getHyperplane()).getLocation().getX();
502 }
503
504 /** {@inheritDoc}
505 * <p>
506 * The iterator returns the limit values of sub-intervals in ascending order.
507 * </p>
508 * <p>
509 * The iterator does <em>not</em> support the optional {@code remove} operation.
510 * </p>
511 */
512 @Override
513 public Iterator<double[]> iterator() {
514 return new SubIntervalsIterator();
515 }
516
517 /** Local iterator for sub-intervals. */
518 private class SubIntervalsIterator implements Iterator<double[]> {
519
520 /** Current node. */
521 private BSPTree<Euclidean1D> current;
522
523 /** Sub-interval no yet returned. */
524 private double[] pending;
525
526 /** Simple constructor.
527 */
528 SubIntervalsIterator() {
529
530 current = getFirstIntervalBoundary();
531
532 if (current == null) {
533 // all the leaf tree nodes share the same inside/outside status
534 if ((Boolean) getFirstLeaf(getTree(false)).getAttribute()) {
535 // it is an inside node, it represents the full real line
536 pending = new double[] {
537 Double.NEGATIVE_INFINITY, Double.POSITIVE_INFINITY
538 };
539 } else {
540 pending = null;
541 }
542 } else if (isIntervalEnd(current)) {
543 // the first boundary is an interval end,
544 // so the first interval starts at infinity
545 pending = new double[] {
546 Double.NEGATIVE_INFINITY, getAngle(current)
547 };
548 } else {
549 selectPending();
550 }
551 }
552
553 /** Walk the tree to select the pending sub-interval.
554 */
555 private void selectPending() {
556
557 // look for the start of the interval
558 BSPTree<Euclidean1D> start = current;
559 while (start != null && !isIntervalStart(start)) {
560 start = nextInternalNode(start);
561 }
562
563 if (start == null) {
564 // we have exhausted the iterator
565 current = null;
566 pending = null;
567 return;
568 }
569
570 // look for the end of the interval
571 BSPTree<Euclidean1D> end = start;
572 while (end != null && !isIntervalEnd(end)) {
573 end = nextInternalNode(end);
574 }
575
576 if (end != null) {
577
578 // we have identified the interval
579 pending = new double[] {
580 getAngle(start), getAngle(end)
581 };
582
583 // prepare search for next interval
584 current = end;
585
586 } else {
587
588 // the final interval is open toward infinity
589 pending = new double[] {
590 getAngle(start), Double.POSITIVE_INFINITY
591 };
592
593 // there won't be any other intervals
594 current = null;
595
596 }
597
598 }
599
600 /** {@inheritDoc} */
601 @Override
602 public boolean hasNext() {
603 return pending != null;
604 }
605
606 /** {@inheritDoc} */
607 @Override
608 public double[] next() {
609 if (pending == null) {
610 throw new NoSuchElementException();
611 }
612 final double[] next = pending;
613 selectPending();
614 return next;
615 }
616
617 /** {@inheritDoc} */
618 @Override
619 public void remove() {
620 throw new UnsupportedOperationException();
621 }
622
623 }
624
625 }