1 | // License: GPL. For details, see LICENSE file.
|
---|
2 | package org.openstreetmap.josm.tools;
|
---|
3 |
|
---|
4 | import java.awt.geom.Area;
|
---|
5 | import java.awt.geom.Line2D;
|
---|
6 | import java.awt.geom.Path2D;
|
---|
7 | import java.awt.geom.PathIterator;
|
---|
8 | import java.awt.geom.Rectangle2D;
|
---|
9 | import java.math.BigDecimal;
|
---|
10 | import java.math.MathContext;
|
---|
11 | import java.util.ArrayList;
|
---|
12 | import java.util.Collection;
|
---|
13 | import java.util.Collections;
|
---|
14 | import java.util.Comparator;
|
---|
15 | import java.util.Iterator;
|
---|
16 | import java.util.LinkedHashSet;
|
---|
17 | import java.util.List;
|
---|
18 | import java.util.Set;
|
---|
19 | import java.util.TreeSet;
|
---|
20 | import java.util.function.Predicate;
|
---|
21 | import java.util.stream.Collectors;
|
---|
22 |
|
---|
23 | import org.openstreetmap.josm.command.AddCommand;
|
---|
24 | import org.openstreetmap.josm.command.ChangeCommand;
|
---|
25 | import org.openstreetmap.josm.command.Command;
|
---|
26 | import org.openstreetmap.josm.data.coor.EastNorth;
|
---|
27 | import org.openstreetmap.josm.data.coor.ILatLon;
|
---|
28 | import org.openstreetmap.josm.data.osm.BBox;
|
---|
29 | import org.openstreetmap.josm.data.osm.DataSet;
|
---|
30 | import org.openstreetmap.josm.data.osm.INode;
|
---|
31 | import org.openstreetmap.josm.data.osm.IPrimitive;
|
---|
32 | import org.openstreetmap.josm.data.osm.IRelation;
|
---|
33 | import org.openstreetmap.josm.data.osm.IWay;
|
---|
34 | import org.openstreetmap.josm.data.osm.MultipolygonBuilder;
|
---|
35 | import org.openstreetmap.josm.data.osm.MultipolygonBuilder.JoinedPolygon;
|
---|
36 | import org.openstreetmap.josm.data.osm.Node;
|
---|
37 | import org.openstreetmap.josm.data.osm.NodePositionComparator;
|
---|
38 | import org.openstreetmap.josm.data.osm.OsmPrimitive;
|
---|
39 | import org.openstreetmap.josm.data.osm.Relation;
|
---|
40 | import org.openstreetmap.josm.data.osm.Way;
|
---|
41 | import org.openstreetmap.josm.data.osm.WaySegment;
|
---|
42 | import org.openstreetmap.josm.data.osm.visitor.paint.relations.Multipolygon;
|
---|
43 | import org.openstreetmap.josm.data.osm.visitor.paint.relations.Multipolygon.PolyData;
|
---|
44 | import org.openstreetmap.josm.data.osm.visitor.paint.relations.MultipolygonCache;
|
---|
45 | import org.openstreetmap.josm.data.projection.Projection;
|
---|
46 | import org.openstreetmap.josm.data.projection.ProjectionRegistry;
|
---|
47 | import org.openstreetmap.josm.data.projection.Projections;
|
---|
48 |
|
---|
49 | /**
|
---|
50 | * Some tools for geometry related tasks.
|
---|
51 | *
|
---|
52 | * @author viesturs
|
---|
53 | */
|
---|
54 | public final class Geometry {
|
---|
55 |
|
---|
56 | private Geometry() {
|
---|
57 | // Hide default constructor for utils classes
|
---|
58 | }
|
---|
59 |
|
---|
60 | /**
|
---|
61 | * The result types for a {@link Geometry#polygonIntersection(Area, Area)} test
|
---|
62 | */
|
---|
63 | public enum PolygonIntersection {
|
---|
64 | /**
|
---|
65 | * The first polygon is inside the second one
|
---|
66 | */
|
---|
67 | FIRST_INSIDE_SECOND,
|
---|
68 | /**
|
---|
69 | * The second one is inside the first
|
---|
70 | */
|
---|
71 | SECOND_INSIDE_FIRST,
|
---|
72 | /**
|
---|
73 | * The polygons do not overlap
|
---|
74 | */
|
---|
75 | OUTSIDE,
|
---|
76 | /**
|
---|
77 | * The polygon borders cross each other
|
---|
78 | */
|
---|
79 | CROSSING
|
---|
80 | }
|
---|
81 |
|
---|
82 | /** threshold value for size of intersection area given in east/north space */
|
---|
83 | public static final double INTERSECTION_EPS_EAST_NORTH = 1e-4;
|
---|
84 |
|
---|
85 | /**
|
---|
86 | * Will find all intersection and add nodes there for list of given ways.
|
---|
87 | * Handles self-intersections too.
|
---|
88 | * And makes commands to add the intersection points to ways.
|
---|
89 | *
|
---|
90 | * Prerequisite: no two nodes have the same coordinates.
|
---|
91 | *
|
---|
92 | * @param ways a list of ways to test
|
---|
93 | * @param test if true, do not build list of Commands, just return nodes
|
---|
94 | * @param cmds list of commands, typically empty when handed to this method.
|
---|
95 | * Will be filled with commands that add intersection nodes to
|
---|
96 | * the ways.
|
---|
97 | * @return list of new nodes, if test is true the list might not contain all intersections
|
---|
98 | */
|
---|
99 | public static Set<Node> addIntersections(List<Way> ways, boolean test, List<Command> cmds) {
|
---|
100 |
|
---|
101 | int n = ways.size();
|
---|
102 | @SuppressWarnings("unchecked")
|
---|
103 | List<Node>[] newNodes = new ArrayList[n];
|
---|
104 | BBox[] wayBounds = new BBox[n];
|
---|
105 | boolean[] changedWays = new boolean[n];
|
---|
106 |
|
---|
107 | Set<Node> intersectionNodes = new LinkedHashSet<>();
|
---|
108 |
|
---|
109 | //copy node arrays for local usage.
|
---|
110 | for (int pos = 0; pos < n; pos++) {
|
---|
111 | newNodes[pos] = new ArrayList<>(ways.get(pos).getNodes());
|
---|
112 | wayBounds[pos] = ways.get(pos).getBBox();
|
---|
113 | changedWays[pos] = false;
|
---|
114 | }
|
---|
115 |
|
---|
116 | DataSet dataset = ways.get(0).getDataSet();
|
---|
117 |
|
---|
118 | //iterate over all way pairs and introduce the intersections
|
---|
119 | Comparator<Node> coordsComparator = new NodePositionComparator();
|
---|
120 | for (int seg1Way = 0; seg1Way < n; seg1Way++) {
|
---|
121 | for (int seg2Way = seg1Way; seg2Way < n; seg2Way++) {
|
---|
122 |
|
---|
123 | //do not waste time on bounds that do not intersect
|
---|
124 | if (!wayBounds[seg1Way].intersects(wayBounds[seg2Way])) {
|
---|
125 | continue;
|
---|
126 | }
|
---|
127 |
|
---|
128 | List<Node> way1Nodes = newNodes[seg1Way];
|
---|
129 | List<Node> way2Nodes = newNodes[seg2Way];
|
---|
130 |
|
---|
131 | //iterate over primary segmemt
|
---|
132 | for (int seg1Pos = 0; seg1Pos + 1 < way1Nodes.size(); seg1Pos++) {
|
---|
133 |
|
---|
134 | //iterate over secondary segment
|
---|
135 | int seg2Start = seg1Way != seg2Way ? 0 : seg1Pos + 2; //skip the adjacent segment
|
---|
136 |
|
---|
137 | for (int seg2Pos = seg2Start; seg2Pos + 1 < way2Nodes.size(); seg2Pos++) {
|
---|
138 |
|
---|
139 | //need to get them again every time, because other segments may be changed
|
---|
140 | Node seg1Node1 = way1Nodes.get(seg1Pos);
|
---|
141 | Node seg1Node2 = way1Nodes.get(seg1Pos + 1);
|
---|
142 | Node seg2Node1 = way2Nodes.get(seg2Pos);
|
---|
143 | Node seg2Node2 = way2Nodes.get(seg2Pos + 1);
|
---|
144 |
|
---|
145 | int commonCount = 0;
|
---|
146 | //test if we have common nodes to add.
|
---|
147 | if (seg1Node1 == seg2Node1 || seg1Node1 == seg2Node2) {
|
---|
148 | commonCount++;
|
---|
149 |
|
---|
150 | if (seg1Way == seg2Way &&
|
---|
151 | seg1Pos == 0 &&
|
---|
152 | seg2Pos == way2Nodes.size() -2) {
|
---|
153 | //do not add - this is first and last segment of the same way.
|
---|
154 | } else {
|
---|
155 | intersectionNodes.add(seg1Node1);
|
---|
156 | }
|
---|
157 | }
|
---|
158 |
|
---|
159 | if (seg1Node2 == seg2Node1 || seg1Node2 == seg2Node2) {
|
---|
160 | commonCount++;
|
---|
161 |
|
---|
162 | intersectionNodes.add(seg1Node2);
|
---|
163 | }
|
---|
164 |
|
---|
165 | //no common nodes - find intersection
|
---|
166 | if (commonCount == 0) {
|
---|
167 | EastNorth intersection = getSegmentSegmentIntersection(
|
---|
168 | seg1Node1.getEastNorth(), seg1Node2.getEastNorth(),
|
---|
169 | seg2Node1.getEastNorth(), seg2Node2.getEastNorth());
|
---|
170 |
|
---|
171 | if (intersection != null) {
|
---|
172 | Node newNode = new Node(ProjectionRegistry.getProjection().eastNorth2latlon(intersection));
|
---|
173 | Node intNode = newNode;
|
---|
174 | boolean insertInSeg1 = false;
|
---|
175 | boolean insertInSeg2 = false;
|
---|
176 | //find if the intersection point is at end point of one of the segments, if so use that point
|
---|
177 |
|
---|
178 | //segment 1
|
---|
179 | if (coordsComparator.compare(newNode, seg1Node1) == 0) {
|
---|
180 | intNode = seg1Node1;
|
---|
181 | } else if (coordsComparator.compare(newNode, seg1Node2) == 0) {
|
---|
182 | intNode = seg1Node2;
|
---|
183 | } else {
|
---|
184 | insertInSeg1 = true;
|
---|
185 | }
|
---|
186 |
|
---|
187 | //segment 2
|
---|
188 | if (coordsComparator.compare(newNode, seg2Node1) == 0) {
|
---|
189 | intNode = seg2Node1;
|
---|
190 | } else if (coordsComparator.compare(newNode, seg2Node2) == 0) {
|
---|
191 | intNode = seg2Node2;
|
---|
192 | } else {
|
---|
193 | insertInSeg2 = true;
|
---|
194 | }
|
---|
195 |
|
---|
196 | if (test) {
|
---|
197 | intersectionNodes.add(intNode);
|
---|
198 | return intersectionNodes;
|
---|
199 | }
|
---|
200 |
|
---|
201 | if (insertInSeg1) {
|
---|
202 | way1Nodes.add(seg1Pos +1, intNode);
|
---|
203 | changedWays[seg1Way] = true;
|
---|
204 |
|
---|
205 | //fix seg2 position, as indexes have changed, seg2Pos is always bigger than seg1Pos on the same segment.
|
---|
206 | if (seg2Way == seg1Way) {
|
---|
207 | seg2Pos++;
|
---|
208 | }
|
---|
209 | }
|
---|
210 |
|
---|
211 | if (insertInSeg2) {
|
---|
212 | way2Nodes.add(seg2Pos +1, intNode);
|
---|
213 | changedWays[seg2Way] = true;
|
---|
214 |
|
---|
215 | //Do not need to compare again to already split segment
|
---|
216 | seg2Pos++;
|
---|
217 | }
|
---|
218 |
|
---|
219 | intersectionNodes.add(intNode);
|
---|
220 |
|
---|
221 | if (intNode == newNode) {
|
---|
222 | cmds.add(new AddCommand(dataset, intNode));
|
---|
223 | }
|
---|
224 | }
|
---|
225 | } else if (test && !intersectionNodes.isEmpty())
|
---|
226 | return intersectionNodes;
|
---|
227 | }
|
---|
228 | }
|
---|
229 | }
|
---|
230 | }
|
---|
231 |
|
---|
232 |
|
---|
233 | for (int pos = 0; pos < ways.size(); pos++) {
|
---|
234 | if (!changedWays[pos]) {
|
---|
235 | continue;
|
---|
236 | }
|
---|
237 |
|
---|
238 | Way way = ways.get(pos);
|
---|
239 | Way newWay = new Way(way);
|
---|
240 | newWay.setNodes(newNodes[pos]);
|
---|
241 |
|
---|
242 | cmds.add(new ChangeCommand(dataset, way, newWay));
|
---|
243 | }
|
---|
244 |
|
---|
245 | return intersectionNodes;
|
---|
246 | }
|
---|
247 |
|
---|
248 | /**
|
---|
249 | * Tests if given point is to the right side of path consisting of 3 points.
|
---|
250 | *
|
---|
251 | * (Imagine the path is continued beyond the endpoints, so you get two rays
|
---|
252 | * starting from lineP2 and going through lineP1 and lineP3 respectively
|
---|
253 | * which divide the plane into two parts. The test returns true, if testPoint
|
---|
254 | * lies in the part that is to the right when traveling in the direction
|
---|
255 | * lineP1, lineP2, lineP3.)
|
---|
256 | *
|
---|
257 | * @param <N> type of node
|
---|
258 | * @param lineP1 first point in path
|
---|
259 | * @param lineP2 second point in path
|
---|
260 | * @param lineP3 third point in path
|
---|
261 | * @param testPoint point to test
|
---|
262 | * @return true if to the right side, false otherwise
|
---|
263 | */
|
---|
264 | public static <N extends INode> boolean isToTheRightSideOfLine(N lineP1, N lineP2, N lineP3, N testPoint) {
|
---|
265 | boolean pathBendToRight = angleIsClockwise(lineP1, lineP2, lineP3);
|
---|
266 | boolean rightOfSeg1 = angleIsClockwise(lineP1, lineP2, testPoint);
|
---|
267 | boolean rightOfSeg2 = angleIsClockwise(lineP2, lineP3, testPoint);
|
---|
268 |
|
---|
269 | if (pathBendToRight)
|
---|
270 | return rightOfSeg1 && rightOfSeg2;
|
---|
271 | else
|
---|
272 | return !(!rightOfSeg1 && !rightOfSeg2);
|
---|
273 | }
|
---|
274 |
|
---|
275 | /**
|
---|
276 | * This method tests if secondNode is clockwise to first node.
|
---|
277 | * @param <N> type of node
|
---|
278 | * @param commonNode starting point for both vectors
|
---|
279 | * @param firstNode first vector end node
|
---|
280 | * @param secondNode second vector end node
|
---|
281 | * @return true if first vector is clockwise before second vector.
|
---|
282 | */
|
---|
283 | public static <N extends INode> boolean angleIsClockwise(N commonNode, N firstNode, N secondNode) {
|
---|
284 | return angleIsClockwise(commonNode.getEastNorth(), firstNode.getEastNorth(), secondNode.getEastNorth());
|
---|
285 | }
|
---|
286 |
|
---|
287 | /**
|
---|
288 | * Finds the intersection of two line segments.
|
---|
289 | * @param p1 the coordinates of the start point of the first specified line segment
|
---|
290 | * @param p2 the coordinates of the end point of the first specified line segment
|
---|
291 | * @param p3 the coordinates of the start point of the second specified line segment
|
---|
292 | * @param p4 the coordinates of the end point of the second specified line segment
|
---|
293 | * @return EastNorth null if no intersection was found, the EastNorth coordinates of the intersection otherwise
|
---|
294 | */
|
---|
295 | public static EastNorth getSegmentSegmentIntersection(EastNorth p1, EastNorth p2, EastNorth p3, EastNorth p4) {
|
---|
296 |
|
---|
297 | CheckParameterUtil.ensure(p1, "p1", EastNorth::isValid);
|
---|
298 | CheckParameterUtil.ensure(p2, "p2", EastNorth::isValid);
|
---|
299 | CheckParameterUtil.ensure(p3, "p3", EastNorth::isValid);
|
---|
300 | CheckParameterUtil.ensure(p4, "p4", EastNorth::isValid);
|
---|
301 |
|
---|
302 | double x1 = p1.getX();
|
---|
303 | double y1 = p1.getY();
|
---|
304 | double x2 = p2.getX();
|
---|
305 | double y2 = p2.getY();
|
---|
306 | double x3 = p3.getX();
|
---|
307 | double y3 = p3.getY();
|
---|
308 | double x4 = p4.getX();
|
---|
309 | double y4 = p4.getY();
|
---|
310 |
|
---|
311 | //TODO: do this locally.
|
---|
312 | //TODO: remove this check after careful testing
|
---|
313 | if (!Line2D.linesIntersect(x1, y1, x2, y2, x3, y3, x4, y4)) return null;
|
---|
314 |
|
---|
315 | // solve line-line intersection in parametric form:
|
---|
316 | // (x1,y1) + (x2-x1,y2-y1)* u = (x3,y3) + (x4-x3,y4-y3)* v
|
---|
317 | // (x2-x1,y2-y1)*u - (x4-x3,y4-y3)*v = (x3-x1,y3-y1)
|
---|
318 | // if 0<= u,v <=1, intersection exists at ( x1+ (x2-x1)*u, y1 + (y2-y1)*u )
|
---|
319 |
|
---|
320 | double a1 = x2 - x1;
|
---|
321 | double b1 = x3 - x4;
|
---|
322 | double c1 = x3 - x1;
|
---|
323 |
|
---|
324 | double a2 = y2 - y1;
|
---|
325 | double b2 = y3 - y4;
|
---|
326 | double c2 = y3 - y1;
|
---|
327 |
|
---|
328 | // Solve the equations
|
---|
329 | double det = a1*b2 - a2*b1;
|
---|
330 |
|
---|
331 | double uu = b2*c1 - b1*c2;
|
---|
332 | double vv = a1*c2 - a2*c1;
|
---|
333 | double mag = Math.abs(uu)+Math.abs(vv);
|
---|
334 |
|
---|
335 | if (Math.abs(det) > 1e-12 * mag) {
|
---|
336 | double u = uu/det, v = vv/det;
|
---|
337 | if (u > -1e-8 && u < 1+1e-8 && v > -1e-8 && v < 1+1e-8) {
|
---|
338 | if (u < 0) u = 0;
|
---|
339 | if (u > 1) u = 1.0;
|
---|
340 | return new EastNorth(x1+a1*u, y1+a2*u);
|
---|
341 | } else {
|
---|
342 | return null;
|
---|
343 | }
|
---|
344 | } else {
|
---|
345 | // parallel lines
|
---|
346 | return null;
|
---|
347 | }
|
---|
348 | }
|
---|
349 |
|
---|
350 | /**
|
---|
351 | * Finds the intersection of two lines of infinite length.
|
---|
352 | *
|
---|
353 | * @param p1 first point on first line
|
---|
354 | * @param p2 second point on first line
|
---|
355 | * @param p3 first point on second line
|
---|
356 | * @param p4 second point on second line
|
---|
357 | * @return EastNorth null if no intersection was found, the coordinates of the intersection otherwise
|
---|
358 | * @throws IllegalArgumentException if a parameter is null or without valid coordinates
|
---|
359 | */
|
---|
360 | public static EastNorth getLineLineIntersection(EastNorth p1, EastNorth p2, EastNorth p3, EastNorth p4) {
|
---|
361 |
|
---|
362 | CheckParameterUtil.ensure(p1, "p1", EastNorth::isValid);
|
---|
363 | CheckParameterUtil.ensure(p2, "p2", EastNorth::isValid);
|
---|
364 | CheckParameterUtil.ensure(p3, "p3", EastNorth::isValid);
|
---|
365 | CheckParameterUtil.ensure(p4, "p4", EastNorth::isValid);
|
---|
366 |
|
---|
367 | // Basically, the formula from wikipedia is used:
|
---|
368 | // https://en.wikipedia.org/wiki/Line%E2%80%93line_intersection
|
---|
369 | // However, large numbers lead to rounding errors (see #10286).
|
---|
370 | // To avoid this, p1 is first subtracted from each of the points:
|
---|
371 | // p1' = 0
|
---|
372 | // p2' = p2 - p1
|
---|
373 | // p3' = p3 - p1
|
---|
374 | // p4' = p4 - p1
|
---|
375 | // In the end, p1 is added to the intersection point of segment p1'/p2'
|
---|
376 | // and segment p3'/p4'.
|
---|
377 |
|
---|
378 | // Convert line from (point, point) form to ax+by=c
|
---|
379 | double a1 = p2.getY() - p1.getY();
|
---|
380 | double b1 = p1.getX() - p2.getX();
|
---|
381 |
|
---|
382 | double a2 = p4.getY() - p3.getY();
|
---|
383 | double b2 = p3.getX() - p4.getX();
|
---|
384 |
|
---|
385 | // Solve the equations
|
---|
386 | double det = a1 * b2 - a2 * b1;
|
---|
387 | if (det == 0)
|
---|
388 | return null; // Lines are parallel
|
---|
389 |
|
---|
390 | double c2 = (p4.getX() - p1.getX()) * (p3.getY() - p1.getY()) - (p3.getX() - p1.getX()) * (p4.getY() - p1.getY());
|
---|
391 |
|
---|
392 | return new EastNorth(b1 * c2 / det + p1.getX(), -a1 * c2 / det + p1.getY());
|
---|
393 | }
|
---|
394 |
|
---|
395 | /**
|
---|
396 | * Check if the segment p1 - p2 is parallel to p3 - p4
|
---|
397 | * @param p1 First point for first segment
|
---|
398 | * @param p2 Second point for first segment
|
---|
399 | * @param p3 First point for second segment
|
---|
400 | * @param p4 Second point for second segment
|
---|
401 | * @return <code>true</code> if they are parallel or close to parallel
|
---|
402 | */
|
---|
403 | public static boolean segmentsParallel(EastNorth p1, EastNorth p2, EastNorth p3, EastNorth p4) {
|
---|
404 |
|
---|
405 | CheckParameterUtil.ensure(p1, "p1", EastNorth::isValid);
|
---|
406 | CheckParameterUtil.ensure(p2, "p2", EastNorth::isValid);
|
---|
407 | CheckParameterUtil.ensure(p3, "p3", EastNorth::isValid);
|
---|
408 | CheckParameterUtil.ensure(p4, "p4", EastNorth::isValid);
|
---|
409 |
|
---|
410 | // Convert line from (point, point) form to ax+by=c
|
---|
411 | double a1 = p2.getY() - p1.getY();
|
---|
412 | double b1 = p1.getX() - p2.getX();
|
---|
413 |
|
---|
414 | double a2 = p4.getY() - p3.getY();
|
---|
415 | double b2 = p3.getX() - p4.getX();
|
---|
416 |
|
---|
417 | // Solve the equations
|
---|
418 | double det = a1 * b2 - a2 * b1;
|
---|
419 | // remove influence of of scaling factor
|
---|
420 | det /= Math.sqrt(a1*a1 + b1*b1) * Math.sqrt(a2*a2 + b2*b2);
|
---|
421 | return Math.abs(det) < 1e-3;
|
---|
422 | }
|
---|
423 |
|
---|
424 | private static EastNorth closestPointTo(EastNorth p1, EastNorth p2, EastNorth point, boolean segmentOnly) {
|
---|
425 | CheckParameterUtil.ensureParameterNotNull(p1, "p1");
|
---|
426 | CheckParameterUtil.ensureParameterNotNull(p2, "p2");
|
---|
427 | CheckParameterUtil.ensureParameterNotNull(point, "point");
|
---|
428 |
|
---|
429 | double ldx = p2.getX() - p1.getX();
|
---|
430 | double ldy = p2.getY() - p1.getY();
|
---|
431 |
|
---|
432 | //segment zero length
|
---|
433 | if (ldx == 0 && ldy == 0)
|
---|
434 | return p1;
|
---|
435 |
|
---|
436 | double pdx = point.getX() - p1.getX();
|
---|
437 | double pdy = point.getY() - p1.getY();
|
---|
438 |
|
---|
439 | double offset = (pdx * ldx + pdy * ldy) / (ldx * ldx + ldy * ldy);
|
---|
440 |
|
---|
441 | if (segmentOnly && offset <= 0)
|
---|
442 | return p1;
|
---|
443 | else if (segmentOnly && offset >= 1)
|
---|
444 | return p2;
|
---|
445 | else
|
---|
446 | return p1.interpolate(p2, offset);
|
---|
447 | }
|
---|
448 |
|
---|
449 | /**
|
---|
450 | * Calculates closest point to a line segment.
|
---|
451 | * @param segmentP1 First point determining line segment
|
---|
452 | * @param segmentP2 Second point determining line segment
|
---|
453 | * @param point Point for which a closest point is searched on line segment [P1,P2]
|
---|
454 | * @return segmentP1 if it is the closest point, segmentP2 if it is the closest point,
|
---|
455 | * a new point if closest point is between segmentP1 and segmentP2.
|
---|
456 | * @see #closestPointToLine
|
---|
457 | * @since 3650
|
---|
458 | */
|
---|
459 | public static EastNorth closestPointToSegment(EastNorth segmentP1, EastNorth segmentP2, EastNorth point) {
|
---|
460 | return closestPointTo(segmentP1, segmentP2, point, true);
|
---|
461 | }
|
---|
462 |
|
---|
463 | /**
|
---|
464 | * Calculates closest point to a line.
|
---|
465 | * @param lineP1 First point determining line
|
---|
466 | * @param lineP2 Second point determining line
|
---|
467 | * @param point Point for which a closest point is searched on line (P1,P2)
|
---|
468 | * @return The closest point found on line. It may be outside the segment [P1,P2].
|
---|
469 | * @see #closestPointToSegment
|
---|
470 | * @since 4134
|
---|
471 | */
|
---|
472 | public static EastNorth closestPointToLine(EastNorth lineP1, EastNorth lineP2, EastNorth point) {
|
---|
473 | return closestPointTo(lineP1, lineP2, point, false);
|
---|
474 | }
|
---|
475 |
|
---|
476 | /**
|
---|
477 | * This method tests if secondNode is clockwise to first node.
|
---|
478 | *
|
---|
479 | * The line through the two points commonNode and firstNode divides the
|
---|
480 | * plane into two parts. The test returns true, if secondNode lies in
|
---|
481 | * the part that is to the right when traveling in the direction from
|
---|
482 | * commonNode to firstNode.
|
---|
483 | *
|
---|
484 | * @param commonNode starting point for both vectors
|
---|
485 | * @param firstNode first vector end node
|
---|
486 | * @param secondNode second vector end node
|
---|
487 | * @return true if first vector is clockwise before second vector.
|
---|
488 | */
|
---|
489 | public static boolean angleIsClockwise(EastNorth commonNode, EastNorth firstNode, EastNorth secondNode) {
|
---|
490 |
|
---|
491 | CheckParameterUtil.ensure(commonNode, "commonNode", EastNorth::isValid);
|
---|
492 | CheckParameterUtil.ensure(firstNode, "firstNode", EastNorth::isValid);
|
---|
493 | CheckParameterUtil.ensure(secondNode, "secondNode", EastNorth::isValid);
|
---|
494 |
|
---|
495 | double dy1 = firstNode.getY() - commonNode.getY();
|
---|
496 | double dy2 = secondNode.getY() - commonNode.getY();
|
---|
497 | double dx1 = firstNode.getX() - commonNode.getX();
|
---|
498 | double dx2 = secondNode.getX() - commonNode.getX();
|
---|
499 |
|
---|
500 | return dy1 * dx2 - dx1 * dy2 > 0;
|
---|
501 | }
|
---|
502 |
|
---|
503 | /**
|
---|
504 | * Returns the Area of a polygon, from its list of nodes.
|
---|
505 | * @param polygon List of nodes forming polygon
|
---|
506 | * @return Area for the given list of nodes (EastNorth coordinates)
|
---|
507 | * @since 6841
|
---|
508 | */
|
---|
509 | public static Area getArea(List<? extends INode> polygon) {
|
---|
510 | Path2D path = new Path2D.Double();
|
---|
511 |
|
---|
512 | boolean begin = true;
|
---|
513 | for (INode n : polygon) {
|
---|
514 | EastNorth en = n.getEastNorth();
|
---|
515 | if (en != null) {
|
---|
516 | if (begin) {
|
---|
517 | path.moveTo(en.getX(), en.getY());
|
---|
518 | begin = false;
|
---|
519 | } else {
|
---|
520 | path.lineTo(en.getX(), en.getY());
|
---|
521 | }
|
---|
522 | }
|
---|
523 | }
|
---|
524 | if (!begin) {
|
---|
525 | path.closePath();
|
---|
526 | }
|
---|
527 |
|
---|
528 | return new Area(path);
|
---|
529 | }
|
---|
530 |
|
---|
531 | /**
|
---|
532 | * Builds a path from a list of nodes
|
---|
533 | * @param polygon Nodes, forming a closed polygon
|
---|
534 | * @param path2d path to add to; can be null, then a new path is created
|
---|
535 | * @return the path (LatLon coordinates)
|
---|
536 | * @since 13638 (signature)
|
---|
537 | */
|
---|
538 | public static Path2D buildPath2DLatLon(List<? extends ILatLon> polygon, Path2D path2d) {
|
---|
539 | Path2D path = path2d != null ? path2d : new Path2D.Double();
|
---|
540 | boolean begin = true;
|
---|
541 | for (ILatLon n : polygon) {
|
---|
542 | if (begin) {
|
---|
543 | path.moveTo(n.lon(), n.lat());
|
---|
544 | begin = false;
|
---|
545 | } else {
|
---|
546 | path.lineTo(n.lon(), n.lat());
|
---|
547 | }
|
---|
548 | }
|
---|
549 | if (!begin) {
|
---|
550 | path.closePath();
|
---|
551 | }
|
---|
552 | return path;
|
---|
553 | }
|
---|
554 |
|
---|
555 | /**
|
---|
556 | * Calculate area in east/north space for given primitive. Uses {@link MultipolygonCache} for multipolygon relations.
|
---|
557 | * @param p the primitive
|
---|
558 | * @return the area in east/north space, might be empty if the primitive is incomplete or not closed or a node
|
---|
559 | * since 15938
|
---|
560 | */
|
---|
561 | public static Area getAreaEastNorth(IPrimitive p) {
|
---|
562 | if (p instanceof Way && ((Way) p).isClosed()) {
|
---|
563 | return Geometry.getArea(((Way) p).getNodes());
|
---|
564 | }
|
---|
565 | if (p.isMultipolygon() && !p.isIncomplete() && !((IRelation<?>) p).hasIncompleteMembers()) {
|
---|
566 | Multipolygon mp = MultipolygonCache.getInstance().get((Relation) p);
|
---|
567 | Path2D path = new Path2D.Double();
|
---|
568 | path.setWindingRule(Path2D.WIND_EVEN_ODD);
|
---|
569 | for (PolyData pd : mp.getCombinedPolygons()) {
|
---|
570 | path.append(pd.get(), false);
|
---|
571 | }
|
---|
572 | return new Area(path);
|
---|
573 | }
|
---|
574 | return new Area();
|
---|
575 | }
|
---|
576 |
|
---|
577 | /**
|
---|
578 | * Returns the Area of a polygon, from the multipolygon relation.
|
---|
579 | * @param multipolygon the multipolygon relation
|
---|
580 | * @return Area for the multipolygon (LatLon coordinates)
|
---|
581 | */
|
---|
582 | public static Area getAreaLatLon(Relation multipolygon) {
|
---|
583 | final Multipolygon mp = MultipolygonCache.getInstance().get(multipolygon);
|
---|
584 | Path2D path = new Path2D.Double();
|
---|
585 | path.setWindingRule(Path2D.WIND_EVEN_ODD);
|
---|
586 | for (Multipolygon.PolyData pd : mp.getCombinedPolygons()) {
|
---|
587 | buildPath2DLatLon(pd.getNodes(), path);
|
---|
588 | for (Multipolygon.PolyData pdInner : pd.getInners()) {
|
---|
589 | buildPath2DLatLon(pdInner.getNodes(), path);
|
---|
590 | }
|
---|
591 | }
|
---|
592 | return new Area(path);
|
---|
593 | }
|
---|
594 |
|
---|
595 | /**
|
---|
596 | * Tests if two polygons intersect.
|
---|
597 | * @param first List of nodes forming first polygon
|
---|
598 | * @param second List of nodes forming second polygon
|
---|
599 | * @return intersection kind
|
---|
600 | */
|
---|
601 | public static PolygonIntersection polygonIntersection(List<? extends INode> first, List<? extends INode> second) {
|
---|
602 | Area a1 = getArea(first);
|
---|
603 | Area a2 = getArea(second);
|
---|
604 | return polygonIntersection(a1, a2, INTERSECTION_EPS_EAST_NORTH);
|
---|
605 | }
|
---|
606 |
|
---|
607 | /**
|
---|
608 | * Tests if two polygons intersect. It is assumed that the area is given in East North points.
|
---|
609 | * @param a1 Area of first polygon
|
---|
610 | * @param a2 Area of second polygon
|
---|
611 | * @return intersection kind
|
---|
612 | * @since 6841
|
---|
613 | */
|
---|
614 | public static PolygonIntersection polygonIntersection(Area a1, Area a2) {
|
---|
615 | return polygonIntersection(a1, a2, INTERSECTION_EPS_EAST_NORTH);
|
---|
616 | }
|
---|
617 |
|
---|
618 | /**
|
---|
619 | * Tests if two polygons intersect.
|
---|
620 | * @param a1 Area of first polygon
|
---|
621 | * @param a2 Area of second polygon
|
---|
622 | * @param eps an area threshold, everything below is considered an empty intersection
|
---|
623 | * @return intersection kind
|
---|
624 | */
|
---|
625 | public static PolygonIntersection polygonIntersection(Area a1, Area a2, double eps) {
|
---|
626 | return polygonIntersectionResult(a1, a2, eps).a;
|
---|
627 | }
|
---|
628 |
|
---|
629 | /**
|
---|
630 | * Calculate intersection area and kind of intersection between two polygons.
|
---|
631 | * @param a1 Area of first polygon
|
---|
632 | * @param a2 Area of second polygon
|
---|
633 | * @param eps an area threshold, everything below is considered an empty intersection
|
---|
634 | * @return pair with intersection kind and intersection area (never null, but maybe empty)
|
---|
635 | * @since 15938
|
---|
636 | */
|
---|
637 | public static Pair<PolygonIntersection, Area> polygonIntersectionResult(Area a1, Area a2, double eps) {
|
---|
638 | Area inter = new Area(a1);
|
---|
639 | inter.intersect(a2);
|
---|
640 |
|
---|
641 | if (inter.isEmpty() || !checkIntersection(inter, eps)) {
|
---|
642 | return new Pair<>(PolygonIntersection.OUTSIDE, inter);
|
---|
643 | } else if (a2.getBounds2D().contains(a1.getBounds2D()) && inter.equals(a1)) {
|
---|
644 | return new Pair<>(PolygonIntersection.FIRST_INSIDE_SECOND, inter);
|
---|
645 | } else if (a1.getBounds2D().contains(a2.getBounds2D()) && inter.equals(a2)) {
|
---|
646 | return new Pair<>(PolygonIntersection.SECOND_INSIDE_FIRST, inter);
|
---|
647 | } else {
|
---|
648 | return new Pair<>(PolygonIntersection.CROSSING, inter);
|
---|
649 | }
|
---|
650 | }
|
---|
651 |
|
---|
652 | /**
|
---|
653 | * Check an intersection area which might describe multiple small polygons.
|
---|
654 | * Return true if any of the polygons is bigger than the given threshold.
|
---|
655 | * @param inter the intersection area
|
---|
656 | * @param eps an area threshold, everything below is considered an empty intersection
|
---|
657 | * @return true if any of the polygons is bigger than the given threshold
|
---|
658 | */
|
---|
659 | private static boolean checkIntersection(Area inter, double eps) {
|
---|
660 | PathIterator pit = inter.getPathIterator(null);
|
---|
661 | double[] res = new double[6];
|
---|
662 | Rectangle2D r = new Rectangle2D.Double();
|
---|
663 | while (!pit.isDone()) {
|
---|
664 | int type = pit.currentSegment(res);
|
---|
665 | switch (type) {
|
---|
666 | case PathIterator.SEG_MOVETO:
|
---|
667 | r = new Rectangle2D.Double(res[0], res[1], 0, 0);
|
---|
668 | break;
|
---|
669 | case PathIterator.SEG_LINETO:
|
---|
670 | r.add(res[0], res[1]);
|
---|
671 | break;
|
---|
672 | case PathIterator.SEG_CLOSE:
|
---|
673 | if (r.getWidth() > eps || r.getHeight() > eps)
|
---|
674 | return true;
|
---|
675 | break;
|
---|
676 | default:
|
---|
677 | break;
|
---|
678 | }
|
---|
679 | pit.next();
|
---|
680 | }
|
---|
681 | return false;
|
---|
682 | }
|
---|
683 |
|
---|
684 | /**
|
---|
685 | * Tests if point is inside a polygon. The polygon can be self-intersecting. In such case the contains function works in xor-like manner.
|
---|
686 | * @param polygonNodes list of nodes from polygon path.
|
---|
687 | * @param point the point to test
|
---|
688 | * @return true if the point is inside polygon.
|
---|
689 | */
|
---|
690 | public static boolean nodeInsidePolygon(INode point, List<? extends INode> polygonNodes) {
|
---|
691 | if (polygonNodes.size() < 2)
|
---|
692 | return false;
|
---|
693 |
|
---|
694 | //iterate each side of the polygon, start with the last segment
|
---|
695 | INode oldPoint = polygonNodes.get(polygonNodes.size() - 1);
|
---|
696 |
|
---|
697 | if (!oldPoint.isLatLonKnown()) {
|
---|
698 | return false;
|
---|
699 | }
|
---|
700 |
|
---|
701 | boolean inside = false;
|
---|
702 | INode p1, p2;
|
---|
703 |
|
---|
704 | for (INode newPoint : polygonNodes) {
|
---|
705 | //skip duplicate points
|
---|
706 | if (newPoint.equals(oldPoint)) {
|
---|
707 | continue;
|
---|
708 | }
|
---|
709 |
|
---|
710 | if (!newPoint.isLatLonKnown()) {
|
---|
711 | return false;
|
---|
712 | }
|
---|
713 |
|
---|
714 | //order points so p1.lat <= p2.lat
|
---|
715 | if (newPoint.getEastNorth().getY() > oldPoint.getEastNorth().getY()) {
|
---|
716 | p1 = oldPoint;
|
---|
717 | p2 = newPoint;
|
---|
718 | } else {
|
---|
719 | p1 = newPoint;
|
---|
720 | p2 = oldPoint;
|
---|
721 | }
|
---|
722 |
|
---|
723 | EastNorth pEN = point.getEastNorth();
|
---|
724 | EastNorth opEN = oldPoint.getEastNorth();
|
---|
725 | EastNorth npEN = newPoint.getEastNorth();
|
---|
726 | EastNorth p1EN = p1.getEastNorth();
|
---|
727 | EastNorth p2EN = p2.getEastNorth();
|
---|
728 |
|
---|
729 | if (pEN != null && opEN != null && npEN != null && p1EN != null && p2EN != null) {
|
---|
730 | //test if the line is crossed and if so invert the inside flag.
|
---|
731 | if ((npEN.getY() < pEN.getY()) == (pEN.getY() <= opEN.getY())
|
---|
732 | && (pEN.getX() - p1EN.getX()) * (p2EN.getY() - p1EN.getY())
|
---|
733 | < (p2EN.getX() - p1EN.getX()) * (pEN.getY() - p1EN.getY())) {
|
---|
734 | inside = !inside;
|
---|
735 | }
|
---|
736 | }
|
---|
737 |
|
---|
738 | oldPoint = newPoint;
|
---|
739 | }
|
---|
740 |
|
---|
741 | return inside;
|
---|
742 | }
|
---|
743 |
|
---|
744 | /**
|
---|
745 | * Returns area of a closed way in square meters.
|
---|
746 | *
|
---|
747 | * @param way Way to measure, should be closed (first node is the same as last node)
|
---|
748 | * @return area of the closed way.
|
---|
749 | */
|
---|
750 | public static double closedWayArea(Way way) {
|
---|
751 | return getAreaAndPerimeter(way.getNodes(), Projections.getProjectionByCode("EPSG:54008")).getArea();
|
---|
752 | }
|
---|
753 |
|
---|
754 | /**
|
---|
755 | * Returns area of a multipolygon in square meters.
|
---|
756 | *
|
---|
757 | * @param multipolygon the multipolygon to measure
|
---|
758 | * @return area of the multipolygon.
|
---|
759 | */
|
---|
760 | public static double multipolygonArea(Relation multipolygon) {
|
---|
761 | double area = 0.0;
|
---|
762 | final Multipolygon mp = MultipolygonCache.getInstance().get(multipolygon);
|
---|
763 | for (Multipolygon.PolyData pd : mp.getCombinedPolygons()) {
|
---|
764 | area += pd.getAreaAndPerimeter(Projections.getProjectionByCode("EPSG:54008")).getArea();
|
---|
765 | }
|
---|
766 | return area;
|
---|
767 | }
|
---|
768 |
|
---|
769 | /**
|
---|
770 | * Computes the area of a closed way and multipolygon in square meters, or {@code null} for other primitives
|
---|
771 | *
|
---|
772 | * @param osm the primitive to measure
|
---|
773 | * @return area of the primitive, or {@code null}
|
---|
774 | * @since 13638 (signature)
|
---|
775 | */
|
---|
776 | public static Double computeArea(IPrimitive osm) {
|
---|
777 | if (osm instanceof Way && ((Way) osm).isClosed()) {
|
---|
778 | return closedWayArea((Way) osm);
|
---|
779 | } else if (osm instanceof Relation && ((Relation) osm).isMultipolygon() && !((Relation) osm).hasIncompleteMembers()) {
|
---|
780 | return multipolygonArea((Relation) osm);
|
---|
781 | } else {
|
---|
782 | return null;
|
---|
783 | }
|
---|
784 | }
|
---|
785 |
|
---|
786 | /**
|
---|
787 | * Determines whether a way is oriented clockwise.
|
---|
788 | *
|
---|
789 | * Internals: Assuming a closed non-looping way, compute twice the area
|
---|
790 | * of the polygon using the formula {@code 2 * area = sum (X[n] * Y[n+1] - X[n+1] * Y[n])}.
|
---|
791 | * If the area is negative the way is ordered in a clockwise direction.
|
---|
792 | *
|
---|
793 | * See http://paulbourke.net/geometry/polyarea/
|
---|
794 | *
|
---|
795 | * @param w the way to be checked.
|
---|
796 | * @return true if and only if way is oriented clockwise.
|
---|
797 | * @throws IllegalArgumentException if way is not closed (see {@link Way#isClosed}).
|
---|
798 | */
|
---|
799 | public static boolean isClockwise(Way w) {
|
---|
800 | return isClockwise(w.getNodes());
|
---|
801 | }
|
---|
802 |
|
---|
803 | /**
|
---|
804 | * Determines whether path from nodes list is oriented clockwise.
|
---|
805 | * @param nodes Nodes list to be checked.
|
---|
806 | * @return true if and only if way is oriented clockwise.
|
---|
807 | * @throws IllegalArgumentException if way is not closed (see {@link Way#isClosed}).
|
---|
808 | * @see #isClockwise(Way)
|
---|
809 | */
|
---|
810 | public static boolean isClockwise(List<? extends INode> nodes) {
|
---|
811 | int nodesCount = nodes.size();
|
---|
812 | if (nodesCount < 3 || nodes.get(0) != nodes.get(nodesCount - 1)) {
|
---|
813 | throw new IllegalArgumentException("Way must be closed to check orientation.");
|
---|
814 | }
|
---|
815 | double area2 = 0.;
|
---|
816 |
|
---|
817 | for (int node = 1; node <= /*sic! consider last-first as well*/ nodesCount; node++) {
|
---|
818 | INode coorPrev = nodes.get(node - 1);
|
---|
819 | INode coorCurr = nodes.get(node % nodesCount);
|
---|
820 | area2 += coorPrev.lon() * coorCurr.lat();
|
---|
821 | area2 -= coorCurr.lon() * coorPrev.lat();
|
---|
822 | }
|
---|
823 | return area2 < 0;
|
---|
824 | }
|
---|
825 |
|
---|
826 | /**
|
---|
827 | * Returns angle of a segment defined with 2 point coordinates.
|
---|
828 | *
|
---|
829 | * @param p1 first point
|
---|
830 | * @param p2 second point
|
---|
831 | * @return Angle in radians (-pi, pi]
|
---|
832 | */
|
---|
833 | public static double getSegmentAngle(EastNorth p1, EastNorth p2) {
|
---|
834 |
|
---|
835 | CheckParameterUtil.ensure(p1, "p1", EastNorth::isValid);
|
---|
836 | CheckParameterUtil.ensure(p2, "p2", EastNorth::isValid);
|
---|
837 |
|
---|
838 | return Math.atan2(p2.north() - p1.north(), p2.east() - p1.east());
|
---|
839 | }
|
---|
840 |
|
---|
841 | /**
|
---|
842 | * Returns angle of a corner defined with 3 point coordinates.
|
---|
843 | *
|
---|
844 | * @param p1 first point
|
---|
845 | * @param common Common end point
|
---|
846 | * @param p3 third point
|
---|
847 | * @return Angle in radians (-pi, pi]
|
---|
848 | */
|
---|
849 | public static double getCornerAngle(EastNorth p1, EastNorth common, EastNorth p3) {
|
---|
850 |
|
---|
851 | CheckParameterUtil.ensure(p1, "p1", EastNorth::isValid);
|
---|
852 | CheckParameterUtil.ensure(common, "p2", EastNorth::isValid);
|
---|
853 | CheckParameterUtil.ensure(p3, "p3", EastNorth::isValid);
|
---|
854 |
|
---|
855 | double result = getSegmentAngle(common, p1) - getSegmentAngle(common, p3);
|
---|
856 | if (result <= -Math.PI) {
|
---|
857 | result += 2 * Math.PI;
|
---|
858 | }
|
---|
859 |
|
---|
860 | if (result > Math.PI) {
|
---|
861 | result -= 2 * Math.PI;
|
---|
862 | }
|
---|
863 |
|
---|
864 | return result;
|
---|
865 | }
|
---|
866 |
|
---|
867 | /**
|
---|
868 | * Get angles in radians and return it's value in range [0, 180].
|
---|
869 | *
|
---|
870 | * @param angle the angle in radians
|
---|
871 | * @return normalized angle in degrees
|
---|
872 | * @since 13670
|
---|
873 | */
|
---|
874 | public static double getNormalizedAngleInDegrees(double angle) {
|
---|
875 | return Math.abs(180 * angle / Math.PI);
|
---|
876 | }
|
---|
877 |
|
---|
878 | /**
|
---|
879 | * Compute the centroid/barycenter of nodes
|
---|
880 | * @param nodes Nodes for which the centroid is wanted
|
---|
881 | * @return the centroid of nodes
|
---|
882 | * @see Geometry#getCenter
|
---|
883 | */
|
---|
884 | public static EastNorth getCentroid(List<? extends INode> nodes) {
|
---|
885 | return getCentroidEN(nodes.stream().map(INode::getEastNorth).collect(Collectors.toList()));
|
---|
886 | }
|
---|
887 |
|
---|
888 | /**
|
---|
889 | * Compute the centroid/barycenter of nodes
|
---|
890 | * @param nodes Coordinates for which the centroid is wanted
|
---|
891 | * @return the centroid of nodes
|
---|
892 | * @since 13712
|
---|
893 | */
|
---|
894 | public static EastNorth getCentroidEN(List<EastNorth> nodes) {
|
---|
895 |
|
---|
896 | final int size = nodes.size();
|
---|
897 | if (size == 1) {
|
---|
898 | return nodes.get(0);
|
---|
899 | } else if (size == 2) {
|
---|
900 | return nodes.get(0).getCenter(nodes.get(1));
|
---|
901 | }
|
---|
902 |
|
---|
903 | BigDecimal area = BigDecimal.ZERO;
|
---|
904 | BigDecimal north = BigDecimal.ZERO;
|
---|
905 | BigDecimal east = BigDecimal.ZERO;
|
---|
906 |
|
---|
907 | // See https://en.wikipedia.org/wiki/Centroid#Of_a_polygon for the equation used here
|
---|
908 | for (int i = 0; i < size; i++) {
|
---|
909 | EastNorth n0 = nodes.get(i);
|
---|
910 | EastNorth n1 = nodes.get((i+1) % size);
|
---|
911 |
|
---|
912 | if (n0 != null && n1 != null && n0.isValid() && n1.isValid()) {
|
---|
913 | BigDecimal x0 = BigDecimal.valueOf(n0.east());
|
---|
914 | BigDecimal y0 = BigDecimal.valueOf(n0.north());
|
---|
915 | BigDecimal x1 = BigDecimal.valueOf(n1.east());
|
---|
916 | BigDecimal y1 = BigDecimal.valueOf(n1.north());
|
---|
917 |
|
---|
918 | BigDecimal k = x0.multiply(y1, MathContext.DECIMAL128).subtract(y0.multiply(x1, MathContext.DECIMAL128));
|
---|
919 |
|
---|
920 | area = area.add(k, MathContext.DECIMAL128);
|
---|
921 | east = east.add(k.multiply(x0.add(x1, MathContext.DECIMAL128), MathContext.DECIMAL128));
|
---|
922 | north = north.add(k.multiply(y0.add(y1, MathContext.DECIMAL128), MathContext.DECIMAL128));
|
---|
923 | }
|
---|
924 | }
|
---|
925 |
|
---|
926 | BigDecimal d = new BigDecimal(3, MathContext.DECIMAL128); // 1/2 * 6 = 3
|
---|
927 | area = area.multiply(d, MathContext.DECIMAL128);
|
---|
928 | if (area.compareTo(BigDecimal.ZERO) != 0) {
|
---|
929 | north = north.divide(area, MathContext.DECIMAL128);
|
---|
930 | east = east.divide(area, MathContext.DECIMAL128);
|
---|
931 | }
|
---|
932 |
|
---|
933 | return new EastNorth(east.doubleValue(), north.doubleValue());
|
---|
934 | }
|
---|
935 |
|
---|
936 | /**
|
---|
937 | * Compute center of the circle closest to different nodes.
|
---|
938 | *
|
---|
939 | * Ensure exact center computation in case nodes are already aligned in circle.
|
---|
940 | * This is done by least square method.
|
---|
941 | * Let be a_i x + b_i y + c_i = 0 equations of bisectors of each edges.
|
---|
942 | * Center must be intersection of all bisectors.
|
---|
943 | * <pre>
|
---|
944 | * [ a1 b1 ] [ -c1 ]
|
---|
945 | * With A = [ ... ... ] and Y = [ ... ]
|
---|
946 | * [ an bn ] [ -cn ]
|
---|
947 | * </pre>
|
---|
948 | * An approximation of center of circle is (At.A)^-1.At.Y
|
---|
949 | * @param nodes Nodes parts of the circle (at least 3)
|
---|
950 | * @return An approximation of the center, of null if there is no solution.
|
---|
951 | * @see Geometry#getCentroid
|
---|
952 | * @since 6934
|
---|
953 | */
|
---|
954 | public static EastNorth getCenter(List<? extends INode> nodes) {
|
---|
955 | int nc = nodes.size();
|
---|
956 | if (nc < 3) return null;
|
---|
957 | /**
|
---|
958 | * Equation of each bisector ax + by + c = 0
|
---|
959 | */
|
---|
960 | double[] a = new double[nc];
|
---|
961 | double[] b = new double[nc];
|
---|
962 | double[] c = new double[nc];
|
---|
963 | // Compute equation of bisector
|
---|
964 | for (int i = 0; i < nc; i++) {
|
---|
965 | EastNorth pt1 = nodes.get(i).getEastNorth();
|
---|
966 | EastNorth pt2 = nodes.get((i+1) % nc).getEastNorth();
|
---|
967 | a[i] = pt1.east() - pt2.east();
|
---|
968 | b[i] = pt1.north() - pt2.north();
|
---|
969 | double d = Math.sqrt(a[i]*a[i] + b[i]*b[i]);
|
---|
970 | if (d == 0) return null;
|
---|
971 | a[i] /= d;
|
---|
972 | b[i] /= d;
|
---|
973 | double xC = (pt1.east() + pt2.east()) / 2;
|
---|
974 | double yC = (pt1.north() + pt2.north()) / 2;
|
---|
975 | c[i] = -(a[i]*xC + b[i]*yC);
|
---|
976 | }
|
---|
977 | // At.A = [aij]
|
---|
978 | double a11 = 0, a12 = 0, a22 = 0;
|
---|
979 | // At.Y = [bi]
|
---|
980 | double b1 = 0, b2 = 0;
|
---|
981 | for (int i = 0; i < nc; i++) {
|
---|
982 | a11 += a[i]*a[i];
|
---|
983 | a12 += a[i]*b[i];
|
---|
984 | a22 += b[i]*b[i];
|
---|
985 | b1 -= a[i]*c[i];
|
---|
986 | b2 -= b[i]*c[i];
|
---|
987 | }
|
---|
988 | // (At.A)^-1 = [invij]
|
---|
989 | double det = a11*a22 - a12*a12;
|
---|
990 | if (Math.abs(det) < 1e-5) return null;
|
---|
991 | double inv11 = a22/det;
|
---|
992 | double inv12 = -a12/det;
|
---|
993 | double inv22 = a11/det;
|
---|
994 | // center (xC, yC) = (At.A)^-1.At.y
|
---|
995 | double xC = inv11*b1 + inv12*b2;
|
---|
996 | double yC = inv12*b1 + inv22*b2;
|
---|
997 | return new EastNorth(xC, yC);
|
---|
998 | }
|
---|
999 |
|
---|
1000 | /**
|
---|
1001 | * Tests if the {@code node} is inside the multipolygon {@code multiPolygon}. The nullable argument
|
---|
1002 | * {@code isOuterWayAMatch} allows to decide if the immediate {@code outer} way of the multipolygon is a match.
|
---|
1003 | * For repeated tests against {@code multiPolygon} better use {@link Geometry#filterInsideMultipolygon}.
|
---|
1004 | * @param node node
|
---|
1005 | * @param multiPolygon multipolygon
|
---|
1006 | * @param isOuterWayAMatch allows to decide if the immediate {@code outer} way of the multipolygon is a match
|
---|
1007 | * @return {@code true} if the node is inside the multipolygon
|
---|
1008 | */
|
---|
1009 | public static boolean isNodeInsideMultiPolygon(INode node, Relation multiPolygon, Predicate<Way> isOuterWayAMatch) {
|
---|
1010 | return isPolygonInsideMultiPolygon(Collections.singletonList(node), multiPolygon, isOuterWayAMatch);
|
---|
1011 | }
|
---|
1012 |
|
---|
1013 | /**
|
---|
1014 | * Tests if the polygon formed by {@code nodes} is inside the multipolygon {@code multiPolygon}. The nullable argument
|
---|
1015 | * {@code isOuterWayAMatch} allows to decide if the immediate {@code outer} way of the multipolygon is a match.
|
---|
1016 | * For repeated tests against {@code multiPolygon} better use {@link Geometry#filterInsideMultipolygon}.
|
---|
1017 | * <p>
|
---|
1018 | * If {@code nodes} contains exactly one element, then it is checked whether that one node is inside the multipolygon.
|
---|
1019 | * @param nodes nodes forming the polygon
|
---|
1020 | * @param multiPolygon multipolygon
|
---|
1021 | * @param isOuterWayAMatch allows to decide if the immediate {@code outer} way of the multipolygon is a match
|
---|
1022 | * @return {@code true} if the multipolygon is valid and the polygon formed by nodes is inside the multipolygon
|
---|
1023 | */
|
---|
1024 | public static boolean isPolygonInsideMultiPolygon(List<? extends INode> nodes, Relation multiPolygon, Predicate<Way> isOuterWayAMatch) {
|
---|
1025 | try {
|
---|
1026 | return isPolygonInsideMultiPolygon(nodes, MultipolygonBuilder.joinWays(multiPolygon), isOuterWayAMatch);
|
---|
1027 | } catch (MultipolygonBuilder.JoinedPolygonCreationException ex) {
|
---|
1028 | Logging.trace(ex);
|
---|
1029 | Logging.debug("Invalid multipolygon " + multiPolygon);
|
---|
1030 | return false;
|
---|
1031 | }
|
---|
1032 | }
|
---|
1033 |
|
---|
1034 | /**
|
---|
1035 | * Tests if the polygon formed by {@code nodes} is inside the multipolygon {@code multiPolygon}. The nullable argument
|
---|
1036 | * {@code isOuterWayAMatch} allows to decide if the immediate {@code outer} way of the multipolygon is a match.
|
---|
1037 | * For repeated tests against {@code multiPolygon} better use {@link Geometry#filterInsideMultipolygon}.
|
---|
1038 | * <p>
|
---|
1039 | * If {@code nodes} contains exactly one element, then it is checked whether that one node is inside the multipolygon.
|
---|
1040 | * @param nodes nodes forming the polygon
|
---|
1041 | * @param outerInner result of {@link MultipolygonBuilder#joinWays(Relation)}
|
---|
1042 | * @param isOuterWayAMatch allows to decide if the immediate {@code outer} way of the multipolygon is a match
|
---|
1043 | * @return {@code true} if the multipolygon is valid and the polygon formed by nodes is inside the multipolygon
|
---|
1044 | * @since 15069
|
---|
1045 | */
|
---|
1046 | public static boolean isPolygonInsideMultiPolygon(List<? extends INode> nodes, Pair<List<JoinedPolygon>,
|
---|
1047 | List<JoinedPolygon>> outerInner, Predicate<Way> isOuterWayAMatch) {
|
---|
1048 | Area a1 = nodes.size() == 1 ? null : getArea(nodes);
|
---|
1049 | // Test if object is inside an outer member
|
---|
1050 | for (JoinedPolygon out : outerInner.a) {
|
---|
1051 | if (a1 == null
|
---|
1052 | ? nodeInsidePolygon(nodes.get(0), out.nodes)
|
---|
1053 | : PolygonIntersection.FIRST_INSIDE_SECOND == polygonIntersection(a1, out.area)) {
|
---|
1054 | boolean insideInner = false;
|
---|
1055 | // If inside an outer, check it is not inside an inner
|
---|
1056 | for (JoinedPolygon in : outerInner.b) {
|
---|
1057 | if (a1 == null ? nodeInsidePolygon(nodes.get(0), in.nodes)
|
---|
1058 | : in.area.getBounds2D().contains(a1.getBounds2D())
|
---|
1059 | && polygonIntersection(a1, in.area) == PolygonIntersection.FIRST_INSIDE_SECOND
|
---|
1060 | && polygonIntersection(in.area, out.area) == PolygonIntersection.FIRST_INSIDE_SECOND) {
|
---|
1061 | insideInner = true;
|
---|
1062 | break;
|
---|
1063 | }
|
---|
1064 | }
|
---|
1065 | if (!insideInner) {
|
---|
1066 | // Final check using predicate
|
---|
1067 | if (isOuterWayAMatch == null || isOuterWayAMatch.test(out.ways.get(0)
|
---|
1068 | /* TODO give a better representation of the outer ring to the predicate */)) {
|
---|
1069 | return true;
|
---|
1070 | }
|
---|
1071 | }
|
---|
1072 | }
|
---|
1073 | }
|
---|
1074 | return false;
|
---|
1075 | }
|
---|
1076 |
|
---|
1077 | /**
|
---|
1078 | * Find all primitives in the given collection which are inside the given polygon.
|
---|
1079 | *
|
---|
1080 | * @param primitives the primitives
|
---|
1081 | * @param polygon the closed way or multipolygon relation
|
---|
1082 | * @return a new list containing the found primitives, empty if polygon is invalid or nothing was found.
|
---|
1083 | * @see Geometry#filterInsidePolygon
|
---|
1084 | * @see Geometry#filterInsideMultipolygon
|
---|
1085 | * @since 15730
|
---|
1086 | */
|
---|
1087 | public static List<IPrimitive> filterInsideAnyPolygon(Collection<IPrimitive> primitives, IPrimitive polygon) {
|
---|
1088 | if (polygon instanceof IWay<?>) {
|
---|
1089 | return filterInsidePolygon(primitives, (IWay<?>) polygon);
|
---|
1090 | } else if (polygon instanceof Relation && polygon.isMultipolygon()) {
|
---|
1091 | return filterInsideMultipolygon(primitives, (Relation) polygon);
|
---|
1092 | }
|
---|
1093 | return Collections.emptyList();
|
---|
1094 | }
|
---|
1095 |
|
---|
1096 | /**
|
---|
1097 | * Find all primitives in the given collection which are inside the given polygon.
|
---|
1098 | * Unclosed ways and multipolygon relations with unclosed outer rings are ignored.
|
---|
1099 | *
|
---|
1100 | * @param primitives the primitives
|
---|
1101 | * @param polygon the polygon
|
---|
1102 | * @return a new list containing the found primitives, empty if polygon is invalid or nothing was found.
|
---|
1103 | * @since 15069 (for {@link List} of {@code primitives}, 15730 for a {@link Collection} of {@code primitives})
|
---|
1104 | */
|
---|
1105 | public static List<IPrimitive> filterInsidePolygon(Collection<IPrimitive> primitives, IWay<?> polygon) {
|
---|
1106 | List<IPrimitive> res = new ArrayList<>();
|
---|
1107 | if (!polygon.isClosed() || polygon.getNodesCount() <= 3)
|
---|
1108 | return res;
|
---|
1109 | /** polygon area in east north space, calculated only when really needed */
|
---|
1110 | Area polygonArea = null;
|
---|
1111 | for (IPrimitive p : primitives) {
|
---|
1112 | if (p instanceof INode) {
|
---|
1113 | if (nodeInsidePolygon((INode) p, polygon.getNodes())) {
|
---|
1114 | res.add(p);
|
---|
1115 | }
|
---|
1116 | } else if (p instanceof IWay) {
|
---|
1117 | if (((IWay<?>) p).isClosed()) {
|
---|
1118 | if (polygonArea == null) {
|
---|
1119 | polygonArea = getArea(polygon.getNodes());
|
---|
1120 | }
|
---|
1121 | if (PolygonIntersection.FIRST_INSIDE_SECOND == polygonIntersection(getArea(((IWay<?>) p).getNodes()),
|
---|
1122 | polygonArea)) {
|
---|
1123 | res.add(p);
|
---|
1124 | }
|
---|
1125 | }
|
---|
1126 | } else if (p.isMultipolygon()) {
|
---|
1127 | if (polygonArea == null) {
|
---|
1128 | polygonArea = getArea(polygon.getNodes());
|
---|
1129 | }
|
---|
1130 | Multipolygon mp = new Multipolygon((Relation) p);
|
---|
1131 | boolean inside = true;
|
---|
1132 | // a (valid) multipolygon is inside the polygon if all outer rings are inside
|
---|
1133 | for (PolyData outer : mp.getOuterPolygons()) {
|
---|
1134 | if (!outer.isClosed()
|
---|
1135 | || PolygonIntersection.FIRST_INSIDE_SECOND != polygonIntersection(getArea(outer.getNodes()),
|
---|
1136 | polygonArea)) {
|
---|
1137 | inside = false;
|
---|
1138 | break;
|
---|
1139 | }
|
---|
1140 | }
|
---|
1141 | if (inside) {
|
---|
1142 | res.add(p);
|
---|
1143 | }
|
---|
1144 | }
|
---|
1145 | }
|
---|
1146 | return res;
|
---|
1147 | }
|
---|
1148 |
|
---|
1149 | /**
|
---|
1150 | * Find all primitives in the given collection which are inside the given multipolygon. Members of the multipolygon are
|
---|
1151 | * ignored. Unclosed ways and multipolygon relations with unclosed outer rings are ignored.
|
---|
1152 | * @param primitives the primitives
|
---|
1153 | * @param multiPolygon the multipolygon relation
|
---|
1154 | * @return a new list containing the found primitives, empty if multipolygon is invalid or nothing was found.
|
---|
1155 | * @since 15069
|
---|
1156 | */
|
---|
1157 | public static List<IPrimitive> filterInsideMultipolygon(Collection<IPrimitive> primitives, Relation multiPolygon) {
|
---|
1158 | List<IPrimitive> res = new ArrayList<>();
|
---|
1159 | if (primitives.isEmpty())
|
---|
1160 | return res;
|
---|
1161 |
|
---|
1162 | final Pair<List<JoinedPolygon>, List<JoinedPolygon>> outerInner;
|
---|
1163 | try {
|
---|
1164 | outerInner = MultipolygonBuilder.joinWays(multiPolygon);
|
---|
1165 | } catch (MultipolygonBuilder.JoinedPolygonCreationException ex) {
|
---|
1166 | Logging.trace(ex);
|
---|
1167 | Logging.debug("Invalid multipolygon " + multiPolygon);
|
---|
1168 | return res;
|
---|
1169 | }
|
---|
1170 |
|
---|
1171 | Set<OsmPrimitive> members = multiPolygon.getMemberPrimitives();
|
---|
1172 | for (IPrimitive p : primitives) {
|
---|
1173 | if (members.contains(p))
|
---|
1174 | continue;
|
---|
1175 | if (p instanceof Node) {
|
---|
1176 | if (isPolygonInsideMultiPolygon(Collections.singletonList((Node) p), outerInner, null)) {
|
---|
1177 | res.add(p);
|
---|
1178 | }
|
---|
1179 | } else if (p instanceof Way) {
|
---|
1180 | if (((IWay<?>) p).isClosed() && isPolygonInsideMultiPolygon(((Way) p).getNodes(), outerInner, null)) {
|
---|
1181 | res.add(p);
|
---|
1182 | }
|
---|
1183 | } else if (p.isMultipolygon()) {
|
---|
1184 | Multipolygon mp = new Multipolygon((Relation) p);
|
---|
1185 | boolean inside = true;
|
---|
1186 | // a (valid) multipolygon is inside multiPolygon if all outer rings are inside
|
---|
1187 | for (PolyData outer : mp.getOuterPolygons()) {
|
---|
1188 | if (!outer.isClosed() || !isPolygonInsideMultiPolygon(outer.getNodes(), outerInner, null)) {
|
---|
1189 | inside = false;
|
---|
1190 | break;
|
---|
1191 | }
|
---|
1192 | }
|
---|
1193 | if (inside) {
|
---|
1194 | res.add(p);
|
---|
1195 | }
|
---|
1196 | }
|
---|
1197 | }
|
---|
1198 | return res;
|
---|
1199 | }
|
---|
1200 |
|
---|
1201 | /**
|
---|
1202 | * Data class to hold two double values (area and perimeter of a polygon).
|
---|
1203 | */
|
---|
1204 | public static class AreaAndPerimeter {
|
---|
1205 | private final double area;
|
---|
1206 | private final double perimeter;
|
---|
1207 |
|
---|
1208 | /**
|
---|
1209 | * Create a new {@link AreaAndPerimeter}
|
---|
1210 | * @param area The area
|
---|
1211 | * @param perimeter The perimeter
|
---|
1212 | */
|
---|
1213 | public AreaAndPerimeter(double area, double perimeter) {
|
---|
1214 | this.area = area;
|
---|
1215 | this.perimeter = perimeter;
|
---|
1216 | }
|
---|
1217 |
|
---|
1218 | /**
|
---|
1219 | * Gets the area
|
---|
1220 | * @return The area size
|
---|
1221 | */
|
---|
1222 | public double getArea() {
|
---|
1223 | return area;
|
---|
1224 | }
|
---|
1225 |
|
---|
1226 | /**
|
---|
1227 | * Gets the perimeter
|
---|
1228 | * @return The perimeter length
|
---|
1229 | */
|
---|
1230 | public double getPerimeter() {
|
---|
1231 | return perimeter;
|
---|
1232 | }
|
---|
1233 | }
|
---|
1234 |
|
---|
1235 | /**
|
---|
1236 | * Calculate area and perimeter length of a polygon.
|
---|
1237 | *
|
---|
1238 | * Uses current projection; units are that of the projected coordinates.
|
---|
1239 | *
|
---|
1240 | * @param nodes the list of nodes representing the polygon
|
---|
1241 | * @return area and perimeter
|
---|
1242 | */
|
---|
1243 | public static AreaAndPerimeter getAreaAndPerimeter(List<? extends ILatLon> nodes) {
|
---|
1244 | return getAreaAndPerimeter(nodes, null);
|
---|
1245 | }
|
---|
1246 |
|
---|
1247 | /**
|
---|
1248 | * Calculate area and perimeter length of a polygon in the given projection.
|
---|
1249 | *
|
---|
1250 | * @param nodes the list of nodes representing the polygon
|
---|
1251 | * @param projection the projection to use for the calculation, {@code null} defaults to {@link ProjectionRegistry#getProjection()}
|
---|
1252 | * @return area and perimeter
|
---|
1253 | * @since 13638 (signature)
|
---|
1254 | */
|
---|
1255 | public static AreaAndPerimeter getAreaAndPerimeter(List<? extends ILatLon> nodes, Projection projection) {
|
---|
1256 | CheckParameterUtil.ensureParameterNotNull(nodes, "nodes");
|
---|
1257 | double area = 0;
|
---|
1258 | double perimeter = 0;
|
---|
1259 | Projection useProjection = projection == null ? ProjectionRegistry.getProjection() : projection;
|
---|
1260 |
|
---|
1261 | if (!nodes.isEmpty()) {
|
---|
1262 | boolean closed = nodes.get(0) == nodes.get(nodes.size() - 1);
|
---|
1263 | int numSegments = closed ? nodes.size() - 1 : nodes.size();
|
---|
1264 | EastNorth p1 = nodes.get(0).getEastNorth(useProjection);
|
---|
1265 | for (int i = 1; i <= numSegments; i++) {
|
---|
1266 | final ILatLon node = nodes.get(i == numSegments ? 0 : i);
|
---|
1267 | final EastNorth p2 = node.getEastNorth(useProjection);
|
---|
1268 | if (p1 != null && p2 != null) {
|
---|
1269 | area += p1.east() * p2.north() - p2.east() * p1.north();
|
---|
1270 | perimeter += p1.distance(p2);
|
---|
1271 | }
|
---|
1272 | p1 = p2;
|
---|
1273 | }
|
---|
1274 | }
|
---|
1275 | return new AreaAndPerimeter(Math.abs(area) / 2, perimeter);
|
---|
1276 | }
|
---|
1277 |
|
---|
1278 | /**
|
---|
1279 | * Get the closest primitive to {@code osm} from the collection of
|
---|
1280 | * OsmPrimitive {@code primitives}
|
---|
1281 | *
|
---|
1282 | * The {@code primitives} should be fully downloaded to ensure accuracy.
|
---|
1283 | *
|
---|
1284 | * Note: The complexity of this method is O(n*m), where n is the number of
|
---|
1285 | * children {@code osm} has plus 1, m is the number of children the
|
---|
1286 | * collection of primitives have plus the number of primitives in the
|
---|
1287 | * collection.
|
---|
1288 | *
|
---|
1289 | * @param <T> The return type of the primitive
|
---|
1290 | * @param osm The primitive to get the distances from
|
---|
1291 | * @param primitives The collection of primitives to get the distance to
|
---|
1292 | * @return The closest {@link OsmPrimitive}. This is not determinative.
|
---|
1293 | * To get all primitives that share the same distance, use
|
---|
1294 | * {@link Geometry#getClosestPrimitives}.
|
---|
1295 | * @since 15035
|
---|
1296 | */
|
---|
1297 | public static <T extends OsmPrimitive> T getClosestPrimitive(OsmPrimitive osm, Collection<T> primitives) {
|
---|
1298 | Collection<T> collection = getClosestPrimitives(osm, primitives);
|
---|
1299 | return collection.iterator().next();
|
---|
1300 | }
|
---|
1301 |
|
---|
1302 | /**
|
---|
1303 | * Get the closest primitives to {@code osm} from the collection of
|
---|
1304 | * OsmPrimitive {@code primitives}
|
---|
1305 | *
|
---|
1306 | * The {@code primitives} should be fully downloaded to ensure accuracy.
|
---|
1307 | *
|
---|
1308 | * Note: The complexity of this method is O(n*m), where n is the number of
|
---|
1309 | * children {@code osm} has plus 1, m is the number of children the
|
---|
1310 | * collection of primitives have plus the number of primitives in the
|
---|
1311 | * collection.
|
---|
1312 | *
|
---|
1313 | * @param <T> The return type of the primitive
|
---|
1314 | * @param osm The primitive to get the distances from
|
---|
1315 | * @param primitives The collection of primitives to get the distance to
|
---|
1316 | * @return The closest {@link OsmPrimitive}s. May be empty.
|
---|
1317 | * @since 15035
|
---|
1318 | */
|
---|
1319 | public static <T extends OsmPrimitive> Collection<T> getClosestPrimitives(OsmPrimitive osm, Collection<T> primitives) {
|
---|
1320 | double lowestDistance = Double.MAX_VALUE;
|
---|
1321 | TreeSet<T> closest = new TreeSet<>();
|
---|
1322 | for (T primitive : primitives) {
|
---|
1323 | double distance = getDistance(osm, primitive);
|
---|
1324 | if (Double.isNaN(distance)) continue;
|
---|
1325 | int comp = Double.compare(distance, lowestDistance);
|
---|
1326 | if (comp < 0) {
|
---|
1327 | closest.clear();
|
---|
1328 | lowestDistance = distance;
|
---|
1329 | closest.add(primitive);
|
---|
1330 | } else if (comp == 0) {
|
---|
1331 | closest.add(primitive);
|
---|
1332 | }
|
---|
1333 | }
|
---|
1334 | return closest;
|
---|
1335 | }
|
---|
1336 |
|
---|
1337 | /**
|
---|
1338 | * Get the furthest primitive to {@code osm} from the collection of
|
---|
1339 | * OsmPrimitive {@code primitives}
|
---|
1340 | *
|
---|
1341 | * The {@code primitives} should be fully downloaded to ensure accuracy.
|
---|
1342 | *
|
---|
1343 | * It does NOT give the furthest primitive based off of the furthest
|
---|
1344 | * part of that primitive
|
---|
1345 | *
|
---|
1346 | * Note: The complexity of this method is O(n*m), where n is the number of
|
---|
1347 | * children {@code osm} has plus 1, m is the number of children the
|
---|
1348 | * collection of primitives have plus the number of primitives in the
|
---|
1349 | * collection.
|
---|
1350 | *
|
---|
1351 | * @param <T> The return type of the primitive
|
---|
1352 | * @param osm The primitive to get the distances from
|
---|
1353 | * @param primitives The collection of primitives to get the distance to
|
---|
1354 | * @return The furthest {@link OsmPrimitive}. This is not determinative.
|
---|
1355 | * To get all primitives that share the same distance, use
|
---|
1356 | * {@link Geometry#getFurthestPrimitives}
|
---|
1357 | * @since 15035
|
---|
1358 | */
|
---|
1359 | public static <T extends OsmPrimitive> T getFurthestPrimitive(OsmPrimitive osm, Collection<T> primitives) {
|
---|
1360 | return getFurthestPrimitives(osm, primitives).iterator().next();
|
---|
1361 | }
|
---|
1362 |
|
---|
1363 | /**
|
---|
1364 | * Get the furthest primitives to {@code osm} from the collection of
|
---|
1365 | * OsmPrimitive {@code primitives}
|
---|
1366 | *
|
---|
1367 | * The {@code primitives} should be fully downloaded to ensure accuracy.
|
---|
1368 | *
|
---|
1369 | * It does NOT give the furthest primitive based off of the furthest
|
---|
1370 | * part of that primitive
|
---|
1371 | *
|
---|
1372 | * Note: The complexity of this method is O(n*m), where n is the number of
|
---|
1373 | * children {@code osm} has plus 1, m is the number of children the
|
---|
1374 | * collection of primitives have plus the number of primitives in the
|
---|
1375 | * collection.
|
---|
1376 | *
|
---|
1377 | * @param <T> The return type of the primitive
|
---|
1378 | * @param osm The primitive to get the distances from
|
---|
1379 | * @param primitives The collection of primitives to get the distance to
|
---|
1380 | * @return The furthest {@link OsmPrimitive}s. It may return an empty collection.
|
---|
1381 | * @since 15035
|
---|
1382 | */
|
---|
1383 | public static <T extends OsmPrimitive> Collection<T> getFurthestPrimitives(OsmPrimitive osm, Collection<T> primitives) {
|
---|
1384 | double furthestDistance = Double.NEGATIVE_INFINITY;
|
---|
1385 | TreeSet<T> furthest = new TreeSet<>();
|
---|
1386 | for (T primitive : primitives) {
|
---|
1387 | double distance = getDistance(osm, primitive);
|
---|
1388 | if (Double.isNaN(distance)) continue;
|
---|
1389 | int comp = Double.compare(distance, furthestDistance);
|
---|
1390 | if (comp > 0) {
|
---|
1391 | furthest.clear();
|
---|
1392 | furthestDistance = distance;
|
---|
1393 | furthest.add(primitive);
|
---|
1394 | } else if (comp == 0) {
|
---|
1395 | furthest.add(primitive);
|
---|
1396 | }
|
---|
1397 | }
|
---|
1398 | return furthest;
|
---|
1399 | }
|
---|
1400 |
|
---|
1401 | /**
|
---|
1402 | * Get the distance between different {@link OsmPrimitive}s
|
---|
1403 | * @param one The primitive to get the distance from
|
---|
1404 | * @param two The primitive to get the distance to
|
---|
1405 | * @return The distance between the primitives in meters
|
---|
1406 | * (or the unit of the current projection, see {@link Projection}).
|
---|
1407 | * May return {@link Double#NaN} if one of the primitives is incomplete.
|
---|
1408 | *
|
---|
1409 | * Note: The complexity is O(n*m), where (n,m) are the number of child
|
---|
1410 | * objects the {@link OsmPrimitive}s have.
|
---|
1411 | * @since 15035
|
---|
1412 | */
|
---|
1413 | public static double getDistance(OsmPrimitive one, OsmPrimitive two) {
|
---|
1414 | double rValue = Double.MAX_VALUE;
|
---|
1415 | if (one == null || two == null || one.isIncomplete()
|
---|
1416 | || two.isIncomplete()) return Double.NaN;
|
---|
1417 | if (one instanceof Node && two instanceof Node) {
|
---|
1418 | rValue = ((Node) one).getCoor().greatCircleDistance(((Node) two).getCoor());
|
---|
1419 | } else if (one instanceof Node && two instanceof Way) {
|
---|
1420 | rValue = getDistanceWayNode((Way) two, (Node) one);
|
---|
1421 | } else if (one instanceof Way && two instanceof Node) {
|
---|
1422 | rValue = getDistanceWayNode((Way) one, (Node) two);
|
---|
1423 | } else if (one instanceof Way && two instanceof Way) {
|
---|
1424 | rValue = getDistanceWayWay((Way) one, (Way) two);
|
---|
1425 | } else if (one instanceof Relation && !(two instanceof Relation)) {
|
---|
1426 | for (OsmPrimitive osmPrimitive: ((Relation) one).getMemberPrimitives()) {
|
---|
1427 | double currentDistance = getDistance(osmPrimitive, two);
|
---|
1428 | if (currentDistance < rValue) rValue = currentDistance;
|
---|
1429 | }
|
---|
1430 | } else if (!(one instanceof Relation) && two instanceof Relation) {
|
---|
1431 | rValue = getDistance(two, one);
|
---|
1432 | } else if (one instanceof Relation && two instanceof Relation) {
|
---|
1433 | for (OsmPrimitive osmPrimitive1 : ((Relation) one).getMemberPrimitives()) {
|
---|
1434 | for (OsmPrimitive osmPrimitive2 : ((Relation) two).getMemberPrimitives()) {
|
---|
1435 | double currentDistance = getDistance(osmPrimitive1, osmPrimitive2);
|
---|
1436 | if (currentDistance < rValue) rValue = currentDistance;
|
---|
1437 | }
|
---|
1438 | }
|
---|
1439 | }
|
---|
1440 | return rValue != Double.MAX_VALUE ? rValue : Double.NaN;
|
---|
1441 | }
|
---|
1442 |
|
---|
1443 | /**
|
---|
1444 | * Get the distance between a way and a node
|
---|
1445 | * @param way The way to get the distance from
|
---|
1446 | * @param node The node to get the distance to
|
---|
1447 | * @return The distance between the {@code way} and the {@code node} in
|
---|
1448 | * meters (or the unit of the current projection, see {@link Projection}).
|
---|
1449 | * May return {@link Double#NaN} if the primitives are incomplete.
|
---|
1450 | * @since 15035
|
---|
1451 | */
|
---|
1452 | public static double getDistanceWayNode(Way way, Node node) {
|
---|
1453 | if (way == null || node == null || way.getNodesCount() < 2 || !node.isLatLonKnown())
|
---|
1454 | return Double.NaN;
|
---|
1455 |
|
---|
1456 | double smallest = Double.MAX_VALUE;
|
---|
1457 | EastNorth en0 = node.getEastNorth();
|
---|
1458 | // go through the nodes as if they were paired
|
---|
1459 | Iterator<Node> iter = way.getNodes().iterator();
|
---|
1460 | EastNorth en1 = iter.next().getEastNorth();
|
---|
1461 | while (iter.hasNext()) {
|
---|
1462 | EastNorth en2 = iter.next().getEastNorth();
|
---|
1463 | double distance = getSegmentNodeDistSq(en1, en2, en0);
|
---|
1464 | if (distance < smallest)
|
---|
1465 | smallest = distance;
|
---|
1466 | en1 = en2;
|
---|
1467 | }
|
---|
1468 | return smallest != Double.MAX_VALUE ? Math.sqrt(smallest) : Double.NaN;
|
---|
1469 | }
|
---|
1470 |
|
---|
1471 | /**
|
---|
1472 | * Get the closest {@link WaySegment} from a way to a primitive.
|
---|
1473 | * @param way The {@link Way} to get the distance from and the {@link WaySegment}
|
---|
1474 | * @param primitive The {@link OsmPrimitive} to get the distance to
|
---|
1475 | * @return The {@link WaySegment} that is closest to {@code primitive} from {@code way}.
|
---|
1476 | * If there are multiple {@link WaySegment}s with the same distance, the last
|
---|
1477 | * {@link WaySegment} with the same distance will be returned.
|
---|
1478 | * May return {@code null} if the way has fewer than two nodes or one
|
---|
1479 | * of the primitives is incomplete.
|
---|
1480 | * @since 15035
|
---|
1481 | */
|
---|
1482 | public static WaySegment getClosestWaySegment(Way way, OsmPrimitive primitive) {
|
---|
1483 | if (way == null || primitive == null || way.isIncomplete()
|
---|
1484 | || primitive.isIncomplete()) return null;
|
---|
1485 | double lowestDistance = Double.MAX_VALUE;
|
---|
1486 | Pair<Node, Node> closestNodes = null;
|
---|
1487 | for (Pair<Node, Node> nodes : way.getNodePairs(false)) {
|
---|
1488 | Way tWay = new Way();
|
---|
1489 | tWay.addNode(nodes.a);
|
---|
1490 | tWay.addNode(nodes.b);
|
---|
1491 | double distance = getDistance(tWay, primitive);
|
---|
1492 | if (distance < lowestDistance) {
|
---|
1493 | lowestDistance = distance;
|
---|
1494 | closestNodes = nodes;
|
---|
1495 | }
|
---|
1496 | }
|
---|
1497 | if (closestNodes == null) return null;
|
---|
1498 | return lowestDistance != Double.MAX_VALUE ? WaySegment.forNodePair(way, closestNodes.a, closestNodes.b) : null;
|
---|
1499 | }
|
---|
1500 |
|
---|
1501 | /**
|
---|
1502 | * Get the distance between different ways. Iterates over the nodes of the ways, complexity is O(n*m)
|
---|
1503 | * (n,m giving the number of nodes)
|
---|
1504 | * @param w1 The first {@link Way}
|
---|
1505 | * @param w2 The second {@link Way}
|
---|
1506 | * @return The shortest distance between the ways in meters
|
---|
1507 | * (or the unit of the current projection, see {@link Projection}).
|
---|
1508 | * May return {@link Double#NaN}.
|
---|
1509 | * @since 15035
|
---|
1510 | */
|
---|
1511 | public static double getDistanceWayWay(Way w1, Way w2) {
|
---|
1512 | if (w1 == null || w2 == null || w1.getNodesCount() < 2 || w2.getNodesCount() < 2)
|
---|
1513 | return Double.NaN;
|
---|
1514 | double rValue = Double.MAX_VALUE;
|
---|
1515 | Iterator<Node> iter1 = w1.getNodes().iterator();
|
---|
1516 | Node w1N1 = iter1.next();
|
---|
1517 | while (iter1.hasNext()) {
|
---|
1518 | Node w1N2 = iter1.next();
|
---|
1519 | Iterator<Node> iter2 = w2.getNodes().iterator();
|
---|
1520 | Node w2N1 = iter2.next();
|
---|
1521 | while (iter2.hasNext()) {
|
---|
1522 | Node w2N2 = iter2.next();
|
---|
1523 | double distance = getDistanceSegmentSegment(w1N1, w1N2, w2N1, w2N2);
|
---|
1524 | if (distance < rValue)
|
---|
1525 | rValue = distance;
|
---|
1526 | w2N1 = w2N2;
|
---|
1527 | }
|
---|
1528 | w1N1 = w1N2;
|
---|
1529 | }
|
---|
1530 | return rValue != Double.MAX_VALUE ? rValue : Double.NaN;
|
---|
1531 | }
|
---|
1532 |
|
---|
1533 | /**
|
---|
1534 | * Get the distance between different {@link WaySegment}s
|
---|
1535 | * @param ws1 A {@link WaySegment}
|
---|
1536 | * @param ws2 A {@link WaySegment}
|
---|
1537 | * @return The distance between the two {@link WaySegment}s in meters
|
---|
1538 | * (or the unit of the current projection, see {@link Projection}).
|
---|
1539 | * May return {@link Double#NaN}.
|
---|
1540 | * @since 15035
|
---|
1541 | */
|
---|
1542 | public static double getDistanceSegmentSegment(WaySegment ws1, WaySegment ws2) {
|
---|
1543 | return getDistanceSegmentSegment(ws1.getFirstNode(), ws1.getSecondNode(), ws2.getFirstNode(), ws2.getSecondNode());
|
---|
1544 | }
|
---|
1545 |
|
---|
1546 | /**
|
---|
1547 | * Get the distance between different {@link WaySegment}s
|
---|
1548 | * @param ws1Node1 The first node of the first WaySegment
|
---|
1549 | * @param ws1Node2 The second node of the second WaySegment
|
---|
1550 | * @param ws2Node1 The first node of the second WaySegment
|
---|
1551 | * @param ws2Node2 The second node of the second WaySegment
|
---|
1552 | * @return The distance between the two {@link WaySegment}s in meters
|
---|
1553 | * (or the unit of the current projection, see {@link Projection}).
|
---|
1554 | * May return {@link Double#NaN}.
|
---|
1555 | * @since 15035
|
---|
1556 | */
|
---|
1557 | public static double getDistanceSegmentSegment(Node ws1Node1, Node ws1Node2, Node ws2Node1, Node ws2Node2) {
|
---|
1558 | EastNorth enWs1Node1 = ws1Node1.getEastNorth();
|
---|
1559 | EastNorth enWs1Node2 = ws1Node2.getEastNorth();
|
---|
1560 | EastNorth enWs2Node1 = ws2Node1.getEastNorth();
|
---|
1561 | EastNorth enWs2Node2 = ws2Node2.getEastNorth();
|
---|
1562 | if (enWs1Node1 == null || enWs1Node2 == null || enWs2Node1 == null || enWs2Node2 == null)
|
---|
1563 | return Double.NaN;
|
---|
1564 | if (getSegmentSegmentIntersection(enWs1Node1, enWs1Node2, enWs2Node1, enWs2Node2) != null)
|
---|
1565 | return 0;
|
---|
1566 |
|
---|
1567 | double dist1sq = getSegmentNodeDistSq(enWs1Node1, enWs1Node2, enWs2Node1);
|
---|
1568 | double dist2sq = getSegmentNodeDistSq(enWs1Node1, enWs1Node2, enWs2Node2);
|
---|
1569 | double dist3sq = getSegmentNodeDistSq(enWs2Node1, enWs2Node2, enWs1Node1);
|
---|
1570 | double dist4sq = getSegmentNodeDistSq(enWs2Node1, enWs2Node2, enWs1Node2);
|
---|
1571 | double smallest = Math.min(Math.min(dist1sq, dist2sq), Math.min(dist3sq, dist4sq));
|
---|
1572 | return smallest != Double.MAX_VALUE ? Math.sqrt(smallest) : Double.NaN;
|
---|
1573 | }
|
---|
1574 |
|
---|
1575 | /**
|
---|
1576 | * Calculate closest distance between a line segment s1-s2 and a point p
|
---|
1577 | * @param s1 start of segment
|
---|
1578 | * @param s2 end of segment
|
---|
1579 | * @param p the point
|
---|
1580 | * @return the square of the euclidean distance from p to the closest point on the segment
|
---|
1581 | */
|
---|
1582 | private static double getSegmentNodeDistSq(EastNorth s1, EastNorth s2, EastNorth p) {
|
---|
1583 | EastNorth c1 = closestPointTo(s1, s2, p, true);
|
---|
1584 | return c1.distanceSq(p);
|
---|
1585 | }
|
---|
1586 | }
|
---|