LargestEmptyCircle: is the boundary considered as obstacle as well?

[jorisvandenbossche]

I am running into a case where LargestEmptyCircle doesn't return the radius line of a circle that would be fully inside the envelope or convex hull of the input linestring.

Using geosop:

$ geosop -a 'LINESTRING (0 0, 0 10, 8 10, 8 0)' largestEmptyCircle 1
LINESTRING (4 4, 0 4)
$ geosop -a 'LINESTRING (0 0, 0 10, 8 10, 8 0, 9 0)' largestEmptyCircle 1
LINESTRING (3.9375 0.5625, 0 0.5625)

(the tolerance doesn't matter that much. Using a value of 1 gives a not super accurate radius in the second case, but what matters is the location of the center which is too much towards the bottom which causes the circle to be outside of the convex hull of the input, and this also happens for smaller tolerance values)

Visualizing the above:

vs

My expectation was that the second result would not be correct, although re-reading the documentation of LargestEmptyCircle now, I have to admit isn't explicit about this:

geos/capi/geos_c.h.in

Lines 3897 to 3902 in 386461d

	* The LEC is the largest circle which has its **center** inside the boundary,
	* and whose interior does not intersect with any obstacle. If no boundary is provided, the
	* convex hull of the obstacles is used as the boundary.
	* The LEC center is the point in the interior of the boundary which has the farthest distance from
	* the obstacles (up to tolerance). The LEC is determined by the center point and a point lying on an
	* obstacle indicating the circle radius.

Also in the second output above, the center is within the boundary (check). And the interior of the formed circle does not intersect with any of the obstacles (check). But I would have expected that the boundary is also treated as an additional "obstacle" (which, to be honest, is not said as such in the documentation). But so I suppose that is not the case?

I can manually add the envelope of the linestring to the input obstacles (creating a MultiLineString of the original input with the boundary of its envelope), and then I get the result I expected:

[dr-jts]

Correct, when computing the LEC the convex hull of the obstacles is not automatically included. However, the centre of the LEC always lies inside the convex hull of the input.

Here's a similar situation with point obstacles:

I'm going to call this a feature, since there may be situations where the convex hull is not desired as an obstacle? Here's an example where this might be required:

But agreed, the documentation needs to be more explicit about this.

[dr-jts]

A bigger issue is that there is currently no way to compute an LEC of a set of obstacles within a polygon, and constrain the LEC to lie within the polygon. The polygon boundary can be included as an additional obstacle, but it is treated as a linear constraint, not as an area.

This is most obvious for highly concave polygons, as in this example:

A use case for this is to compute a location within a given area which is farthest from a set of points within the area (for example, for siting a new facility farthest from existing ones).

This limitation may be relevant to your requirements as well. In the example you gave you could include the convex hull as an additional obstacle, and that works since there are no additional obstacles within the linework. But if the open linestring was not convex, the computed LEC might well lie outside it, and there is currently no way to ensure it lies "inside".

It's fairly easy to add an ability to constrain the LEC to lie within an area as well as respecting obstacles within the area. It just requires the desired area constraint to be provided as an additional argument, and to be considered during processing.

[jorisvandenbossche]

I'm going to call this a feature, since there may be situations where the convex hull is not desired as an obstacle? Here's an example where this might be required:

Yes, it will indeed by dependent on your actual use case whether that's something desired or not, and so you can quite easily add the envelope or convex hull to the obstacles manually, if desired, with the current interface.

A bigger issue is that there is currently no way to compute an LEC of a set of obstacles within a polygon, and constrain the LEC to lie within the polygon

Yes, I was exploring LEC today (to expose it in shapely, which also triggered me to open this issue), and encountered similar questions, see shapely/shapely#1307 (comment).

The first issue was not actually that the LEC might be "outside" of the polygon, but actually that it is within the interior of it, since LEC only considers the polygon's linework.
For example, copying from shapely/shapely#1307 (comment) (reusing the example of the PostGIS docs), this is the LEC, but is in the interior of the polygon:

I was trying to achieve the possible desired result ("inside" the polygon but not in its interior) by creating a custom boundary: subtracting the Polygon itself from its convex hull.
The problem is that this is a bit brittle: I did have to add back the the boundary of the input polygon to this custom boundary, to ensure it passes the constraint that the specified boundary fully covers the input geometries.

And in a similar way, it was somewhat possible to create the custom boundary in such a way that the resulting LEC is "inside" the polygon (in the case of the above example, just in the hole). But also here I needed to add back the boundary of the input polygon which caused some issues.

It would be nice to make this a bit easier to achieve as a user. Potential ideas (without knowing the details about how / constraints of the algorithm being used here):

• Would it be possible to flag that you actually want to include the interior of polygonal input as obstacles as well (and not just its linework)?
• Or, alternatively (to make it easier to achieve the previous point by specifying a custom boundary), would it be possible to relax the constraint that the boundary has to completely cover the input geometry? (i.e. the obstacles)

[dr-jts]

The first issue was not actually that the LEC might be "outside" of the polygon, but actually that it is within the interior of it, since LEC only considers the polygon's linework.

Right, that's an issue as well. Polygonal obstacles are not currently supported. It would be nice to have that capability. The algorithm is quite general, so it should be fairly easy to add this by testing distance to polygonal obstacles as well as linear and point ones. The computed LEC will be disjoint from all obstacle geometries, and the LEC centre will lie inside the boundary.

I suppose a further option would be to require the entire LEC to lie within the boundary. Have to think about this to see if that's feasible. (I suspect it is, since this is similar to the logic for the Maximum Inscribed Circle). EDIT: this can be done by including the linework of the boundary geometry as an obstacle.

I don't think there's a need to put any requirements on how the boundary relates to the obstacles. Although obviously if the obstacles cover the boundary no LEC exists. (I'm not sure why the constraint that the boundary cover the obstacles in in place).

[dr-jts]

I wonder if it's still feasible to have a single function for this (as in PostGIS)? If only the polygonal boundary is provided, the result is the MIC. If obstacles are provided then this determines the LEC within the boundary (if any).

[jorisvandenbossche]

Yes, but right now there is the requirement of the boundary covering the obstacles. So do I understand it correctly that you would be OK with removing this constraint?

I think so- but will have to test / analyze to confirm that this will work.

[dr-jts]

Yes, but right now there is the requirement of the boundary covering the obstacles. So do I understand it correctly that you would be OK with removing this constraint?

I think so- but will have to test / analyze to confirm that this will work.

In particular, there may be a failure mode if no LEC can be determined (e.g. when the obstacles cover the boundary completely).

[dr-jts]

I'm working on adding the ability to specify an LEC boundary polygon to JTS. It turns out that the GEOS restriction to have the boundary cover the obstacles is easily removed by a simple fix in the code. I'll port this to GEOS soon.

[dr-jts]

The final enhancement to fully handle polygonal obstacles has been added to JTS: locationtech/jts#988

This will be ported to GEOS soon.

[pramsey]

It's ported now.

[dr-jts]

I was trying to achieve the possible desired result ("inside" the polygon but not in its interior)

@jorisvandenbossche Here's the enhanced LargestEmptyCircle working on the circle-inside-a-hole case:

(Note that geosop now allows provides the largestEmptyCircleBdy operation to support a custom boundary - in this case, the shell of the polygon)

bin/geosop -a "POLYGON ((90 120, 120 310, 260 290, 344 346, 390 200, 260 160, 370 90, 200 10, 200 80, 90 120), (160 260, 220 250, 158 175, 160 260))" -b "POLYGON ((90 120, 120 310, 260 290, 344 346, 390 200, 260 160, 370 90, 200 10, 200 80, 90 120))" largestEmptyCircleBdy 0.001
====>
LINESTRING (180.47927856445312 235.24050903320312, 183.94083713840794 256.00986047693203)