/*
* Licensed to the Apache Software Foundation (ASF) under one or more
* contributor license agreements. See the NOTICE file distributed with
* this work for additional information regarding copyright ownership.
* The ASF licenses this file to You under the Apache License, Version 2.0
* (the "License"); you may not use this file except in compliance with
* the License. You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*/
package org.apache.commons.math3.geometry.enclosing;
import java.util.ArrayList;
import java.util.List;
import org.apache.commons.math3.exception.MathInternalError;
import org.apache.commons.math3.geometry.Point;
import org.apache.commons.math3.geometry.Space;
Class implementing Emo Welzl algorithm to find the smallest enclosing ball in linear time.
The class implements the algorithm described in paper Smallest
Enclosing Disks (Balls and Ellipsoids) by Emo Welzl, Lecture Notes in Computer Science
555 (1991) 359-370. The pivoting improvement published in the paper Fast and
Robust Smallest Enclosing Balls, by Bernd Gärtner and further modified in
paper
Efficient Computation of Smallest Enclosing Balls in Three Dimensions by Linus Källberg
to avoid performing local copies of data have been included.
Type parameters: - <S> – Space type.
- <P> – Point type.
Since: 3.3
/** Class implementing Emo Welzl algorithm to find the smallest enclosing ball in linear time.
* <p>
* The class implements the algorithm described in paper <a
* href="http://www.inf.ethz.ch/personal/emo/PublFiles/SmallEnclDisk_LNCS555_91.pdf">Smallest
* Enclosing Disks (Balls and Ellipsoids)</a> by Emo Welzl, Lecture Notes in Computer Science
* 555 (1991) 359-370. The pivoting improvement published in the paper <a
* href="http://www.inf.ethz.ch/personal/gaertner/texts/own_work/esa99_final.pdf">Fast and
* Robust Smallest Enclosing Balls</a>, by Bernd Gärtner and further modified in
* paper <a
* href=http://www.idt.mdh.se/kurser/ct3340/ht12/MINICONFERENCE/FinalPapers/ircse12_submission_30.pdf">
* Efficient Computation of Smallest Enclosing Balls in Three Dimensions</a> by Linus Källberg
* to avoid performing local copies of data have been included.
* </p>
* @param <S> Space type.
* @param <P> Point type.
* @since 3.3
*/
public class WelzlEncloser<S extends Space, P extends Point<S>> implements Encloser<S, P> {
Tolerance below which points are consider to be identical. /** Tolerance below which points are consider to be identical. */
private final double tolerance;
Generator for balls on support. /** Generator for balls on support. */
private final SupportBallGenerator<S, P> generator;
Simple constructor.
Params: - tolerance – below which points are consider to be identical
- generator – generator for balls on support
/** Simple constructor.
* @param tolerance below which points are consider to be identical
* @param generator generator for balls on support
*/
public WelzlEncloser(final double tolerance, final SupportBallGenerator<S, P> generator) {
this.tolerance = tolerance;
this.generator = generator;
}
{@inheritDoc} /** {@inheritDoc} */
public EnclosingBall<S, P> enclose(final Iterable<P> points) {
if (points == null || !points.iterator().hasNext()) {
// return an empty ball
return generator.ballOnSupport(new ArrayList<P>());
}
// Emo Welzl algorithm with Bernd Gärtner and Linus Källberg improvements
return pivotingBall(points);
}
Compute enclosing ball using Gärtner's pivoting heuristic.
Params: - points – points to be enclosed
Returns: enclosing ball
/** Compute enclosing ball using Gärtner's pivoting heuristic.
* @param points points to be enclosed
* @return enclosing ball
*/
private EnclosingBall<S, P> pivotingBall(final Iterable<P> points) {
final P first = points.iterator().next();
final List<P> extreme = new ArrayList<P>(first.getSpace().getDimension() + 1);
final List<P> support = new ArrayList<P>(first.getSpace().getDimension() + 1);
// start with only first point selected as a candidate support
extreme.add(first);
EnclosingBall<S, P> ball = moveToFrontBall(extreme, extreme.size(), support);
while (true) {
// select the point farthest to current ball
final P farthest = selectFarthest(points, ball);
if (ball.contains(farthest, tolerance)) {
// we have found a ball containing all points
return ball;
}
// recurse search, restricted to the small subset containing support and farthest point
support.clear();
support.add(farthest);
EnclosingBall<S, P> savedBall = ball;
ball = moveToFrontBall(extreme, extreme.size(), support);
if (ball.getRadius() < savedBall.getRadius()) {
// this should never happen
throw new MathInternalError();
}
// it was an interesting point, move it to the front
// according to Gärtner's heuristic
extreme.add(0, farthest);
// prune the least interesting points
extreme.subList(ball.getSupportSize(), extreme.size()).clear();
}
}
Compute enclosing ball using Welzl's move to front heuristic.
Params: - extreme – subset of extreme points
- nbExtreme – number of extreme points to consider
- support – points that must belong to the ball support
Returns: enclosing ball, for the extreme subset only
/** Compute enclosing ball using Welzl's move to front heuristic.
* @param extreme subset of extreme points
* @param nbExtreme number of extreme points to consider
* @param support points that must belong to the ball support
* @return enclosing ball, for the extreme subset only
*/
private EnclosingBall<S, P> moveToFrontBall(final List<P> extreme, final int nbExtreme,
final List<P> support) {
// create a new ball on the prescribed support
EnclosingBall<S, P> ball = generator.ballOnSupport(support);
if (ball.getSupportSize() <= ball.getCenter().getSpace().getDimension()) {
for (int i = 0; i < nbExtreme; ++i) {
final P pi = extreme.get(i);
if (!ball.contains(pi, tolerance)) {
// we have found an outside point,
// enlarge the ball by adding it to the support
support.add(pi);
ball = moveToFrontBall(extreme, i, support);
support.remove(support.size() - 1);
// it was an interesting point, move it to the front
// according to Welzl's heuristic
for (int j = i; j > 0; --j) {
extreme.set(j, extreme.get(j - 1));
}
extreme.set(0, pi);
}
}
}
return ball;
}
Select the point farthest to the current ball.
Params: - points – points to be enclosed
- ball – current ball
Returns: farthest point
/** Select the point farthest to the current ball.
* @param points points to be enclosed
* @param ball current ball
* @return farthest point
*/
public P selectFarthest(final Iterable<P> points, final EnclosingBall<S, P> ball) {
final P center = ball.getCenter();
P farthest = null;
double dMax = -1.0;
for (final P point : points) {
final double d = point.distance(center);
if (d > dMax) {
farthest = point;
dMax = d;
}
}
return farthest;
}
}