Commit c0ea7924 authored by Benedikt Zoennchen's avatar Benedikt Zoennchen
Browse files

implementation of a new merging method which is supported by an R-Tree data structure.

parent bec78e97
package org.vadere.gui.topographycreator.control;
import org.apache.commons.lang3.tuple.Pair;
import org.vadere.gui.topographycreator.model.IDrawPanelModel;
import org.vadere.simulator.projects.migration.GeometryCleaner;
import org.vadere.state.attributes.scenario.AttributesObstacle;
import org.vadere.state.scenario.Obstacle;
import org.vadere.meshing.WeilerAtherton;
import org.vadere.util.geometry.shapes.VPolygon;
import org.vadere.util.geometry.shapes.VRectangle;
import org.vadere.util.logging.Logger;
import java.awt.event.ActionEvent;
import java.util.ArrayList;
......@@ -22,6 +24,7 @@ import javax.swing.undo.UndoableEditSupport;
public class ActionMergeObstacles extends TopographyAction {
private final UndoableEditSupport undoSupport;
private static Logger logger = Logger.getLogger(ActionMergeObstacles.class);
public ActionMergeObstacles(String name, ImageIcon icon, IDrawPanelModel panelModel,
UndoableEditSupport undoSupport) {
......@@ -41,19 +44,35 @@ public class ActionMergeObstacles extends TopographyAction {
.map(shape -> ((VPolygon)shape))
.collect(Collectors.toList());
WeilerAtherton weilerAtherton = new WeilerAtherton(polygons);
List<VPolygon> mergedPolygons = weilerAtherton.cup();
GeometryCleaner topographyCleaner = new GeometryCleaner(
new VRectangle(getScenarioPanelModel().getTopographyBound()), polygons, 0.01);
logger.info("start merging, this can require some time.");
Pair<VPolygon, List<VPolygon>> mergedPolygons = topographyCleaner.clean();
//polygons = WeilerAtherton.magnet(polygons, new VRectangle(getScenarioPanelModel().getBounds()));
logger.info("merging process finisehd.");
/*WeilerAtherton weilerAtherton = new WeilerAtherton(polygons);
List<VPolygon> mergedPolygons = weilerAtherton.cup();*/
// remove polygon obstacles
getScenarioPanelModel().removeObstacleIf(obstacle ->
obstacle.getShape() instanceof VPolygon || obstacle.getShape() instanceof VRectangle);
// add merged obstacles
mergedPolygons
mergedPolygons.getRight()
.stream()
.map(polygon -> new Obstacle(new AttributesObstacle(-1, polygon)))
.forEach(obstacle -> getScenarioPanelModel().addShape(obstacle));
//WeilerAtherton weilerAtherton = new WeilerAtherton(Arrays.asList(mergedPolygons.getLeft(), new VPolygon(getScenarioPanelModel().getBounds())));
//weilerAtherton.cap().get();
//getScenarioPanelModel().addShape(new Obstacle(new AttributesObstacle(-1, weilerAtherton.cap().get())));
List<Obstacle> after = new ArrayList<>(getScenarioPanelModel().getObstacles());
......
......@@ -81,9 +81,9 @@ public class IncrementalTriangulation<V extends IVertex, E extends IHalfEdge, F
private Predicate<E> illegalPredicate;
private static Logger log = Logger.getLogger(IncrementalTriangulation.class);
static {
/*static {
ITriConnectivity.log.setDebug();
}
}*/
/**
* Construct a triangulation using an empty mesh.
......@@ -398,7 +398,7 @@ public class IncrementalTriangulation<V extends IVertex, E extends IHalfEdge, F
}
@Override
public E insertVertex(@NotNull V vertex, @NotNull F face) {
public E insertVertex(@NotNull final V vertex, @NotNull final F face) {
if(!initialized) {
init();
}
......
......@@ -17,7 +17,6 @@ import java.util.Collection;
import java.util.Collections;
import java.util.HashSet;
import java.util.List;
import java.util.Optional;
import java.util.Set;
import java.util.stream.Collectors;
import java.util.stream.Stream;
......
......@@ -2316,7 +2316,7 @@ public interface IMesh<V extends IVertex, E extends IHalfEdge, F extends IFace>
*
* @return a string representing the mesh
*/
default String toPythonTriangulation(@NotNull final Function<V, Double> evalPoint) {
default String toPythonTriangulation(@Nullable final Function<V, Double> evalPoint) {
garbageCollection();
StringBuilder builder = new StringBuilder();
List<V> vertices = getVertices();
......@@ -2341,12 +2341,14 @@ public interface IMesh<V extends IVertex, E extends IHalfEdge, F extends IFace>
builder.append("]\n");
// [z1, z2, ...] z = value
builder.append("z = [");
for(V v : vertices) {
builder.append(evalPoint.apply(v) + ",");
if(evalPoint != null) {
builder.append("z = [");
for(V v : vertices) {
builder.append(evalPoint.apply(v) + ",");
}
builder.delete(builder.length()-1, builder.length());
builder.append("]\n");
}
builder.delete(builder.length()-1, builder.length());
builder.append("]\n");
// [[vId1, vId2, vId3], ...]
List<F> faces = getFaces();
......
......@@ -1020,6 +1020,23 @@ public interface IPolyConnectivity<V extends IVertex, E extends IHalfEdge, F ext
return GeometryUtils.isRightOf(v1.getX(), v1.getY(), v2.getX(), v2.getY(), x1, y1);
}
/**
* Returns true if the point (x1, y1) is left of the half-edge in O(1). The half-edge is directed
* and ends in its point.
*
* Does not change the connectivity.
*
* @param x1 the x-coordinate of the point
* @param y1 the y-coordinate of the point
* @param edge the half-edge
* @return true if the point (x1, y1) is left of the half-edge, false otherwise
*/
default boolean isLeftOfRobust(final double x1, final double y1, @NotNull final E edge) {
V v1 = getMesh().getVertex(getMesh().getPrev(edge));
V v2 = getMesh().getVertex(edge);
return GeometryUtils.isLeftOf(v1.getX(), v1.getY(), v2.getX(), v2.getY(), x1, y1);
}
/**
* Returns true if the point (x1, y1) is left of the half-edge in O(1). The half-edge is directed
* and ends in its point.
......
......@@ -2612,7 +2612,7 @@ public interface ITriConnectivity<V extends IVertex, E extends IHalfEdge, F exte
// find the entering edge
for(E e : getMesh().getEdgeIt(startFace)) {
// line intersection
if(intersects(q, pDirection, e) && isLeftOf(p.getX(), p.getY(), e)) {
if(intersects(q, pDirection, e) && isLeftOfRobust(p.getX(), p.getY(), e)) {
inEdge = e;
break;
}
......
......@@ -23,7 +23,6 @@ import java.util.Collections;
import java.util.HashMap;
import java.util.HashSet;
import java.util.LinkedList;
import java.util.List;
import java.util.Map;
import java.util.Optional;
import java.util.Set;
......@@ -237,7 +236,7 @@ public class GenConstrainedDelaunayTriangulator<V extends IVertex, E extends IHa
VLine vEdge = triangulation.getMesh().toLine(edge);
// to be save
// to be save, TODO inconsistent geometry check which may lead to a deadlock (isLeftOfRobust is "robust" while intersectLineSegment is not")
if(GeometryUtils.intersectLineSegment(vEdge.x1, vEdge.y1, vEdge.x2, vEdge.y2, v1.getX(), v1.getY(), v2.getX(), v2.getY())) {
E next = triangulation.getMesh().getNext(edge);
E prev = triangulation.getMesh().getPrev(edge);
......
......@@ -10,7 +10,6 @@ import org.vadere.meshing.mesh.inter.IIncrementalTriangulation;
import org.vadere.meshing.mesh.inter.IVertex;
import org.vadere.meshing.mesh.triangulation.triangulator.inter.IRefiner;
import org.vadere.meshing.utils.io.tex.TexGraphGenerator;
import org.vadere.util.geometry.Geometry;
import org.vadere.util.geometry.GeometryUtils;
import org.vadere.util.logging.Logger;
import org.vadere.util.math.IDistanceFunction;
......
......@@ -119,7 +119,7 @@ public class TestBoyerWatson {
IMesh<PVertex, PHalfEdge, PFace> mesh = delaunayTriangulation.getMesh();
PHalfEdge edge = mesh.getEdge(face);
// triangulations are always ccw ordered!
// triangulations are always ccwRobust ordered!
assertTrue(GeometryUtils.isCCW(mesh.getVertex(edge), mesh.getVertex(mesh.getNext(edge)), mesh.getVertex(mesh.getPrev(edge))));
delaunayTriangulation.splitTriangle(face, centerPoint, false);
......
package org.vadere.simulator.projects.migration;
import org.apache.commons.lang3.tuple.Pair;
import org.jetbrains.annotations.NotNull;
import org.vadere.meshing.mesh.gen.IncrementalTriangulation;
import org.vadere.meshing.mesh.gen.PFace;
import org.vadere.meshing.mesh.gen.PHalfEdge;
import org.vadere.meshing.mesh.gen.PMesh;
import org.vadere.meshing.mesh.gen.PVertex;
import org.vadere.meshing.mesh.triangulation.triangulator.gen.GenConstrainedDelaunayTriangulator;
import org.vadere.util.geometry.GeometryUtils;
import org.vadere.util.geometry.PlanarGraphGenerator;
import org.vadere.util.geometry.shapes.IPoint;
import org.vadere.util.geometry.shapes.VLine;
import org.vadere.util.geometry.shapes.VPoint;
import org.vadere.util.geometry.shapes.VPolygon;
import org.vadere.util.geometry.shapes.VRectangle;
import org.vadere.util.logging.Logger;
import org.vadere.util.math.IDistanceFunction;
import java.util.ArrayList;
import java.util.Collection;
import java.util.List;
import java.util.function.Predicate;
import java.util.stream.Collectors;
/**
* The {@link GeometryCleaner} eliminates most of problematic obstacle definitions.
* It eliminates all intersecting or co-linear lines i.e. if two obstacles intersect
* or align (with a tolerance) they will be merged together.
* For a user defined tolerance {@link GeometryCleaner#tol} points are identical if
* they are closer than {@link GeometryCleaner#tol}.
*
* @author Benedikt Zoennchen
*/
public class GeometryCleaner {
private final static Logger logger = Logger.getLogger(GeometryCleaner.class);
private final VRectangle boundingBox;
private final VRectangle bound;
private final List<VPolygon> polygons;
private final double tol;
public GeometryCleaner(@NotNull final VRectangle bound, @NotNull final Collection<VPolygon> polygons, final double tol) {
this.bound = bound;
this.polygons = new ArrayList<>(polygons);
this.tol = tol;
List<VPoint> points = polygons.stream().flatMap(poly -> poly.getPath().stream()).collect(Collectors.toList());
points.addAll(bound.getPath());
this.boundingBox = GeometryUtils.boundRelative(points);
}
/**
* This algorithm merges points which are very close (dependent on {@link GeometryCleaner#tol}.
*/
private List<VPolygon> magnet() {
List<VPolygon> clone = new ArrayList<>(polygons.size());
IncrementalTriangulation<PVertex, PHalfEdge, PFace> dt = new IncrementalTriangulation<>(new PMesh(), boundingBox);
for(int i = 0; i < polygons.size(); i++) {
List<VPoint> points = polygons.get(i).getPoints();
List<IPoint> clonePoints = new ArrayList<>(points.size());
for(int j = 0; j < points.size(); j++) {
IPoint point = points.get(j);
PHalfEdge edge = dt.insert(point);
IPoint insertedPoint = dt.getMesh().getPoint(edge);
PVertex v = dt.getMesh().getVertex(edge);
PVertex umin = null;
double len = Double.MAX_VALUE;
for(PVertex u : dt.getMesh().getAdjacentVertexIt(v)) {
if(len > u.distance(v)) {
umin = u;
len = u.distance(v);
}
}
if(!insertedPoint.equals(point)) {
point = insertedPoint;
}
if(v.distanceSq(umin) < tol * tol) {
point = new VPoint(umin.getX(), umin.getY());
}
clonePoints.add(point);
}
clone.add(GeometryUtils.toPolygon(clonePoints));
}
return clone;
}
/**
* Removes all intersecting and co-linear lines by splitting them into multiple line or
* by removing identical lines.
*
* @param polygons polygons of the topography (obstacle shapes)
* @param bound the bounding of the topography
*
* @return a list of non-intersecting lines
*/
private Collection<VLine> removeIntersectingLines(@NotNull final List<VPolygon> polygons, @NotNull final VRectangle bound) {
List<VLine> lines = polygons.stream().flatMap(p -> p.getLinePath().stream()).collect(Collectors.toList());
for(VLine line : bound.getLines()) {
lines.add(line);
}
logger.debug("resolve line intersection of " + lines.size() + " lines.");
PlanarGraphGenerator planarGraphGenerator = new PlanarGraphGenerator(lines, 0.01);
Collection<VLine> allLines = planarGraphGenerator.generate();
return allLines;
}
/**
* Computes all merged {@link VPolygon}s and the segment-bounding {@link VPolygon} based on a list of non-intersecting lines
* and the list of (non-merged) {@link VPolygon}. Note that this method does not support polygons with holes.
*
* @param lines non-intersecting lines
* @param polygons original (intersecting) polygons
*
* @return the segment-bounding polygon and a list of merged polygons
*/
private Pair<VPolygon, List<VPolygon>> mergePolygons(@NotNull final Collection<VLine> lines, @NotNull final List<VPolygon> polygons) {
// 1. compute the contrained Delaunay triangulation (all non-intersecting lines are constrained)
IncrementalTriangulation<PVertex, PHalfEdge, PFace> dt = new IncrementalTriangulation<>(new PMesh(), boundingBox);
GenConstrainedDelaunayTriangulator<PVertex, PHalfEdge, PFace> cdt = new GenConstrainedDelaunayTriangulator<>(dt, lines, false);
cdt.generate(false);
var triangulation = cdt.getTriangulation();
// 2. definition of a distance function which defines points inside an obstacle (or outside the topography bound)
IDistanceFunction distanceFunction = IDistanceFunction.create(bound, polygons);
// 3. compute the segment-bounding simple polygon by using the distance function
Predicate<PFace> removePredicate = face -> distanceFunction.apply(triangulation.getMesh().toTriangle(face).midPoint()) > 0;
triangulation.shrinkBorder(removePredicate, true);
VPolygon boundingPolygon = GeometryUtils.toPolygon(triangulation.getMesh().getPoints(triangulation.getMesh().getBorder()));
// 4. compute all holes. Each hole will be the union of intersecting polygons. A hole is generated by removing
List<PFace> faces = triangulation.getMesh().getFaces();
for(PFace face : faces) {
if(!triangulation.getMesh().isBorder(face) && !triangulation.getMesh().isDestroyed(face) && !triangulation.getMesh().isHole(face)) {
triangulation.createHole(face, f -> distanceFunction.apply(triangulation.getMesh().toTriangle(f).midPoint()) > 0, true);
}
}
// 5. convert all holes (merged polygons) back to VPolygon(s)
List<VPolygon> result = new ArrayList<>();
for (PFace hole : triangulation.getMesh().getHoles()) {
VPolygon poly = GeometryUtils.toPolygon(triangulation.getMesh().streamVertices(hole).map(v -> triangulation.getMesh().toPoint(v)).collect(Collectors.toList()));
result.add(poly);
}
//String str = cdt.getMesh().toPythonTriangulation(null);
//System.out.println(str);
//System.out.println(TexGraphGenerator.toTikz(triangulation.getMesh(),1.0f, true));
return Pair.of(boundingPolygon, result);
}
/**
* Cleans the geometry of the topography by:
* 1. merging very close points
* 2. merging very close or intersecting polygons
* 2.1 elimination of intersecting lines
* 2.2 merging of intersecting polygons
*
* @return the segment-bounding polygon and a list of merged polygons
*/
public Pair<VPolygon, List<VPolygon>> clean() {
logger.debug("start magnet algorithm");
List<VPolygon> magnetResult = magnet();
logger.debug("start planar graph generation");
Collection<VLine> lines = removeIntersectingLines(magnetResult, bound);
logger.debug("start polygon merging");
Pair<VPolygon, List<VPolygon>> mergedPolygons = mergePolygons(lines, magnetResult);
return mergedPolygons;
}
}
......@@ -44,7 +44,6 @@
<lwjgl.version>3.1.2</lwjgl.version>
</properties>
<!-- end generated by lwjgl -->
<dependencies>
<!-- GeoLib: JTS -->
<dependency>
......@@ -53,6 +52,13 @@
<version>1.13</version>
</dependency>
<!-- https://mvnrepository.com/artifact/com.github.davidmoten/rtree -->
<!-- R-Tree Lib -->
<dependency>
<groupId>com.github.davidmoten</groupId>
<artifactId>rtree</artifactId>
<version>0.8.7</version>
</dependency>
<!-- generated by lwjgl -->
<dependency>
<groupId>org.lwjgl</groupId>
......
......@@ -109,6 +109,10 @@ public class GeometryUtils {
return Optional.empty();
}
public static VPoint intersectionPoint(@NotNull final VLine line1, @NotNull final VLine line2) {
return intersectionPoint(line1.x1, line1.y1, line1.x2, line1.y2, line2.x1, line2.y1, line2.x2, line2.y2);
}
/**
* Tests if there is a intersection between the {@link VRectangle} boundary and the line-segment defined by ((x1, y1), (x2, y2)).
*
......@@ -1811,7 +1815,7 @@ public class GeometryUtils {
* @param by y-coordinate of b
* @return the projection of a onto b
*/
public static IPoint projectOnto(double ax, double ay, double bx, double by) {
public static VPoint projectOnto(double ax, double ay, double bx, double by) {
assert bx * bx + by * by > GeometryUtils.DOUBLE_EPS;
double blen = Math.sqrt(bx * bx + by * by);
double bxn = bx / blen;
......@@ -1819,7 +1823,7 @@ public class GeometryUtils {
// scalar product
double alpha = ax * bxn + ay * byn;
IPoint a1 = new VPoint(bxn * alpha, byn * alpha);
VPoint a1 = new VPoint(bxn * alpha, byn * alpha);
return a1;
}
......@@ -1835,7 +1839,7 @@ public class GeometryUtils {
*
* @return he projection of a onto the line (p,q)
*/
public static IPoint projectOntoLine(double ax, double ay, double px, double py, double qx, double qy) {
public static VPoint projectOntoLine(double ax, double ay, double px, double py, double qx, double qy) {
double bx = qx - px;
double by = qy - py;
double apx = ax - px;
......
......@@ -41,7 +41,7 @@ public class GrahamScan {
// sort by y first, x secondly to get an extreme point
Arrays.sort(points, new VPointCoordinateComparator());
// now sort by ccw relative to first extreme point
// now sort by ccwRobust relative to first extreme point
Arrays.sort(points, 1, numberOfPoints, new VPointPolarComparator(points[0]));
// first extreme point will be in the hull
......@@ -178,21 +178,21 @@ public class GrahamScan {
double dx2 = o2.x - point.x;
double dy2 = o2.y - point.y;
// o1 above point and o2 below point => o1 before o2 ccw
// o1 above point and o2 below point => o1 before o2 ccwRobust
if (dy1 >= 0 && dy2 < 0) {
return -1;
}
// o1 is below point and o2 above point => o2 before o1 ccw
// o1 is below point and o2 above point => o2 before o1 ccwRobust
else if (dy2 >= 0 && dy1 < 0) {
return 1;
}
// horizontal line between o1 and o2, so look at dx
else if (dy1 == 0 && dy2 == 0) {
// o1 right of point and o2 left of point => o1 before o2 ccw
// o1 right of point and o2 left of point => o1 before o2 ccwRobust
if (dx1 >= 0 && dx2 < 0) {
return -1;
}
// o2 right of point and o2 left of point => o2 before o1 ccw
// o2 right of point and o2 left of point => o2 before o1 ccwRobust
else if (dx2 >= 0 && dx1 < 0) {
return 1;
}
......
package org.vadere.util.geometry;
import com.github.davidmoten.rtree.RTree;
import com.github.davidmoten.rtree.geometry.Rectangle;
import com.github.davidmoten.rtree.geometry.internal.RectangleDouble;
import org.jetbrains.annotations.NotNull;
import org.vadere.util.geometry.shapes.VLine;
import org.vadere.util.geometry.shapes.VPoint;
import java.util.ArrayList;
import java.util.Collection;
import java.util.List;
/**
* This class generates a planar graph based on a list of possibly intersecting
* and co-linear lines. By line we mean line-segment. There are multiple cases for line intersection:
*
* 1. two lines intersect at one point which is not an endpoint of the line (general case)
* 2. two lines are identical (the have the same end point)
* 3. two lines share infinitely many points but are not identical
* 4. two lines share exactly one endpoint
*
* @author Benedikt Zoennchen
*/
public class PlanarGraphGenerator {
private RTree<String, VLine> lineRTree;
private RTree<String, VLine> unresolvedLines;
private ArrayList<VLine> allLines;
private double tol;
private double ccwTol;
/**
* Default constructor.
*
* @param lines collection of intersecting lines
* @param tol a distance tolerance at which two points are seen as identical
*/
public PlanarGraphGenerator(@NotNull final Collection<VLine> lines, final double tol) {
this.allLines = new ArrayList<>();
this.tol = tol;
this.ccwTol = tol;
this.allLines.addAll(lines);
}
private void reset() {
this.lineRTree = RTree.create();
this.unresolvedLines = RTree.create();
for(VLine line : allLines) {
lineRTree = lineRTree.add(line.toString(), line);
unresolvedLines = unresolvedLines.add(line.toString(), line);
}
}
/**
* Computes and returns non-intersecting lines.
* The computation is supported by a R-Tree data structure such that
* we avoid to test all pairs of lines.
*
* @return non-intersecting lines.
*/
public Collection<VLine> generate() {
reset();
boolean hasChanged;
do {
hasChanged = false;
ArrayList<VLine> allLines = new ArrayList<>();
unresolvedLines.entries().map(e -> e.geometry()).forEach(line -> allLines.add(line));
//int size = unresolvedLines.size();
//System.out.println(size);
for(VLine line1 : allLines) {
Rectangle mbr = line1.mbr();
RectangleDouble bufferedMbr = RectangleDouble.create(mbr.x1() - tol, mbr.y1() - tol, mbr.x2() + tol, mbr.y2() + tol);
ArrayList<VLine> otherLines = new ArrayList<>();
unresolvedLines.search(bufferedMbr).map(e -> e.geometry()).forEach(l -> otherLines.add(l));
//lineRTree.entries().map(e -> e.geometry()).forEach(line -> otherLines.add(line));
for(VLine line2 : otherLines) {
if(line1 != line2) {
deleteLine(line1);
deleteLine(line2);
List<VLine> result = new ArrayList<>();
hasChanged = adjustLines(line1, line2, result);
for(VLine line : result) {
addLine(line);
}
if(hasChanged) {
break;
}
}
}
if(hasChanged) {
break;
} else {
unresolvedLines = unresolvedLines.delete(line1.toString(), line1);