GEOSOffsetCurve_r behaviour has changed on self-intersecting lines

[garci66]

Hello! We are using GEOS through QGIS and we noticed that when QGIS changed to GEOS 3.11, the "line offset" capability stopped working for self intersecting features.

Looking at the code, it seems to invoke GEOSOffsetCurve_r and theproblem can be seen here:

The green line is the "regular" self intersecting polyline. The Red line had a negative offset and thre black has a possitive offset.

I would argue that both results are not "right" . One becase its truncated earlier and the line no longer has this offset and the second one since it doesnt follow the "loop".

I saw some discussions on #477 and I think also related to #552. In particular, seeing that QGIS using the 3.10 and earlier releases was working fine, it would seem related to #530 and the new algortith breaking "something"

Sorry for the vague information, I can try to collect any needed information. There is an open issue on QGIS where I added som details: https://github.com/qgis/QGIS/issues/51583

[garci66]

I just wanted to leave a note to some further tests done with geos 3.10 vs 3.11. This is from an issue in QGIS where the changes in GEOS are causing a bit of breakage when it comes to offsets https://github.com/qgis/QGIS/issues/51583#issuecomment-1414330658

[dr-jts]

Agreed that the offset curve behaviour was changed in 3.11. This was semi-deliberate, since the algorithm is now very different, and it was considered harder to reproduce the previous behaviour, and also was not clear (at the time of design) that the old behaviour was better.

Can you elaborate on why the old behaviour is considered better? Are there any documented requirements for offset curves in QGIS (or elsewhere)?

I can see why it might be better to provide "more" linework that can be considered to be part of the offset. It might be possible to improve the 3.11 implementation to provide results that more closely match the old behaviour.

[dr-jts]

@garci66 The original report in https://github.com/qgis/QGIS/issues/51583 seems to be complaining that the behaviour for offsets of closed lines is wrong. Is this the case?

If so, that should be filed in a separate GEOS issue, with examples of the results that are considered incorrect, and the desired result.

[garci66]

Hi!

Sorry if we were mixing things up on the QGIS issue. At first I pegged it on a "offset is working weirdly in QGIS" issue and since in https://github.com/qgis/QGIS/issues/51583 one of the reportes realised it seemed to be related to GEOS, thats why its all in the same QGIS issue.

Relating specifically to GEOS, my issue here is that while indeed the algorithm has changed for self interesecting lines (avoiding the intersection in some parts (and, in my particular case, I don't think the new algorithm is "correct"), the biggest thing is that in some scenarios, the line created is plain wrong.

So, taking the following example curve: LINESTRING(0 0,100 0, 70 -30 ,70 100)

which looks like this

We ran the offset in GEOS 3.10 (offset 10) geosop -a "LINESTRING(0 0,100 0, 70 -30 ,70 100)" -f wkt offsetCurve 10

the output is: MULTILINESTRING ((80 10, 100 10, 101.95090322016128 9.807852804032304, 103.8268343236509 9.238795325112868, 105.55570233019603 8.314696123025453, 107.07106781186548 7.071067811865475, 108.31469612302546 5.555702330196022, 109.23879532511287 3.826834323650898, 109.80785280403231 1.9509032201612824, 110 0, 109.80785280403231 -1.9509032201612824, 109.23879532511287 -3.826834323650898, 108.31469612302546 -5.555702330196021, 107.07106781186548 -7.071067811865475, 77.07106781186548 -37.071067811865476, 75.55570233019601 -38.314696123025456, 73.8268343236509 -39.23879532511287, 71.95090322016128 -39.80785280403231, 70 -40, 68.04909677983872 -39.8078528040323, 66.1731656763491 -39.23879532511286, 64.44429766980397 -38.314696123025456, 62.928932188134524 -37.071067811865476, 61.685303876974544 -35.55570233019602, 60.76120467488713 -33.8268343236509, 60.19214719596769 -31.950903220161283, 60 -30, 60 -10), (0 10, 60 10, 60 100))

which looks like this:

the same command on GEOS 3.11 is the following: LINESTRING (0 10, 60 10, 60 100)

which renders as follows:

One can argue whether the self intersection part missing is correct or not. But the offset starts and ends with the polygon, so its roughly correct

but with a negative offset, things get even wackier geosop -a "LINESTRING(0 0,100 0, 70 -30 ,70 100)" -f wkt offsetCurve N-10

this is the result for 3.10: MULTILINESTRING ((80 100, 80 10), (60 -10, 0 -10))

and this is the result for 3.11:

LINESTRING (0 -10, 60 -10)

which I dont really agree is correct under any criteria. The line just stops at the first intersection

Now, in my particular case, where I use the offset in QGIS to be able to show multiple lines running through the same path (fiber optic networks following the same route, the same set of telephone poles), the 3.10 and earlier behavior was ideal. I could apply an offset to each object and be able to display multiple objects, even if they were self-intersecting without much fear of missing any part of the network. With 3.11, this goes away as any self intersecting loop, no matter its relative size to the offset applied, is removed. And if the offset is negative, its particularly bad since the offset line just ends at the first intersection.

Is there any way both the old and new algorithms could be kept? maybe as an optional parameter to the GEOSOffsetCurve_r function so that, potentially, the older algorithm could be called from QGIS? (This would require changes in QGIS as well but i'd try my best to push those)

thanks!

[dr-jts]

Now, in my particular case, where I use the offset in QGIS to be able to show multiple lines running through the same path (fiber optic networks following the same route, the same set of telephone poles), the 3.10 and earlier behavior was ideal. I could apply an offset to each object and be able to display multiple objects, even if they were self-intersecting without much fear of missing any part of the network. With 3.11, this goes away as any self intersecting loop, no matter its relative size to the offset applied, is removed. And if the offset is negative, its particularly bad since the offset line just ends at the first intersection.

That's a good use case.

Is there any way both the old and new algorithms could be kept? maybe as an optional parameter to the GEOSOffsetCurve_r function so that, potentially, the older algorithm could be called from QGIS? (This would require changes in QGIS as well but i'd try my best to push those)

The old algorithm had some issues as well (e.g. #477). So I'd prefer to upgrade the new implementation to provide behaviour closer to (or identical with) the old code.

[garci66]

if a way to find a solution to the issues raised is found.. i'd be more than happy to keep the new algorithm

The "missing loop" I though was by-design, looking at the comments on #477 is doesnt seem easy. Personally, I think an output such as this one:

would be ideal

but I imagine that it could lead to examples such as

Or it can get even trickier if the offset is larger than the "radius" of a self intersecting feature:

For example, if the offset o is greater than the "diameter" d of the self intersection shown below:

what is the "correct" way to do it? something like this?

im not sure... but open to ideas.

[pramsey]

You hit on exactly the issues that informed the current behaviour, namely loops that are near the offset size (or, conversely, offsets that are near or greater than the loop size). Basically the logic that makes sense for one case does not make sense for the other, there is almost a thresholding, which I think is hard to detect in the midst of the assembly of the offset curve.

[pramsey]

I wonder if rather than trying to jimmy with the outputs (increasingly baroque attacks on the raw offset segments) jimmying with the inputs might result in more predictable/desireable results? Like, if step one was to collapse the input geometry using the offset tolerance, so small loops go away and narrow corridors, spikes and gores close up.

[garci66]

But I have issues even if the offset is much much smaller than the "loop size" . such as the following:

you can see that the line to the "right" of the line ends at the first intersection. Even if the loop is huge compared to the offset

[dr-jts]

The hardest aspect for computing offsets is determining the linework for concave corners which are shallower than the offset. (Taking this to the limit, what is the offset curve inside a circle at a distance greater than the radius of the circle?). This includes the case you mention, of a curve at a distance greater than the width of a self-intersecting loop.

The buffer algorithm has to deal with the same problem, but it can resolve it by discarding all linework which lies inside the buffer. An approach to offsetting which seemed tractable is to create an offset curve as a subset of the buffer outline. This sidesteps the issue, and should still provide "enough" linework to produce results which seem reasonable. (The old code was doing something similar, but operating on the raw buffer curve rather than the final "clean" one.)

The problem is that the new implementation was a bit too aggressive at discarding parts of the buffer curve which matched the offset curve. I think this can be improved.

[dr-jts]

I wonder if rather than trying to jimmy with the outputs (increasingly baroque attacks on the raw offset segments) jimmying with the inputs might result in more predictable/desireable results? Like, if step one was to collapse the input geometry using the offset tolerance, so small loops go away and narrow corridors, spikes and gores close up.

I can't say a priori that this approach wouldn't work. But I'm not seeing a path to an implementation this way. I'm not sure what a "collapsed" line would look like (in particular, for large loops and long corridors with loops at the end)? It would have to be something which allowed generating the offset linework correctly.

[dr-jts]

But I have issues even if the offset is much much smaller than the "loop size" .

Yes, this is a result of the current code "short-circuiting" as soon as it finds a single part of the offset. It should be possible to have it determine all parts of the offset lines.

[dr-jts]

Here's a question: if a line closes on itself so that an offset curve contains rings (or equivalently the buffer contains holes), is it always the case that the entire hole is either part of the offset curve or not in it? To put this another way, is the hole always formed by an offset on the same side of the input curve?

I can't think of a counter-example. @garci66 @pramsey any ideas?

If this is the case, it means that holes in the buffer just need to be checked to see which side of the input line they lie on, and then included in or excluded from the offset curve result in their entirety.

[dr-jts]

Here's a question: if a line closes on itself so that an offset curve contains rings (or equivalently the buffer contains holes), is it always the case that the entire hole is either part of the offset curve or not in it? To put this another way, is the hole always formed by an offset on the same side of the input curve?

Darn, here's a counter-example:

The middle hole and the one top-right are bounded by both sides of the input line.

So I guess the offset curve should be just the sections on the required side of the input line? Like this:

(Incidentally, that image come from the JTS TestBuilder. The code generating it just uses raw offset curves. It looks fine for this case, but it has the usual artifacts for narrow concave corners. No magic!)

[dr-jts]

What might be nice is if the offset curve was split into sections which were ordered to lie along the original line. This adheres more closely to the (often-stated) desire to have the offset curve be a single unbroken line along the input line. It would allow a simple labelling strategy of "place labels at intervals along the line".

As per example above, this would sometimes require splitting closed rings into separate sections (since they can be generated by input line sections very far apart).

It would be possible to take this one step further and link the ordered offset sections into a continuous line. This should provide the effect of creating a single, possibly-self-crossing, offset line. In the example above this would produce a very good quality result. Other more complex inputs might produce less quality results, however. Some experimenting will be required to see if this approach is generally acceptable.

[pramsey]

I'm going to refer you to this change 7561b50f943305b1419bd655e10e7dcad2cb9492, and hope it covers most cases. I think rational people can disagree on what the right behaviour of this function is.

[jorisvandenbossche]

We just got a report in Shapely (not with self-intersections, but just part of the lines coming closer together than the offset distance, interrupting the offset line) that regressed with GEOS 3.11 and now looks more correct again after the latest changes on GEOS main: https://github.com/shapely/shapely/issues/1789#issuecomment-1475252464

[dr-jts]

We just got a report in Shapely (not with self-intersections, but just part of the lines coming closer together than the offset distance, interrupting the offset line) that regressed with GEOS 3.11 and now looks more correct again after the latest changes on GEOS main:

Thanks, @jorisvandenbossche , good to hear that this is the desired behaviour.

[lexelby]

What might be nice is if the offset curve was split into sections which were ordered to lie along the original line. This adheres more closely to the (often-stated) desire to have the offset curve be a single unbroken line along the input line. It would allow a simple labelling strategy of "place labels at intervals along the line".

As per example above, this would sometimes require splitting closed rings into separate sections (since they can be generated by input line sections very far apart).

It would be possible to take this one step further and link the ordered offset sections into a continuous line. This should provide the effect of creating a single, possibly-self-crossing, offset line. In the example above this would produce a very good quality result. Other more complex inputs might produce less quality results, however. Some experimenting will be required to see if this approach is generally acceptable.

We did exactly this in Ink/Stitch here. Our goal was to convert an arbitrary (user-provided) linestring into two equal parallel offsets, as a precursor to doing machine embroidery stitching along the line. Imagine turning a possibly self-intersecting linestring into a railroad track, with the railroad ties being the stitches we want the machine to sew.

My algorithm attempts a parallel offset on the whole linestring and rejects if it gets back something that's not a simple linestring. When it rejects, it splits the path in "half" (by splitting the coordinate list in half) and recursively runs the same algorithm on each half. At the end it rejoins all the parts.

I know this is closed, and I'm probably shouting into the void here, but I wanted to chime in and say that this is definitely a use case that folks care about! It'd be great if geos (and/or Shapely) did all this for us. No other library I've found can.

[dr-jts]

My algorithm attempts a parallel offset on the whole linestring and rejects if it gets back something that's not a simple linestring. When it rejects, it splits the path in "half" (by splitting the coordinate list in half) and recursively runs the same algorithm on each half. At the end it rejoins all the parts.

@lexelby That's an interesting approach! Was there any problem with too much depth of recursion for complex inputs?

I know this is closed, and I'm probably shouting into the void here, but I wanted to chime in and say that this is definitely a use case that folks care about! It'd be great if geos (and/or Shapely) did all this for us. No other library I've found can.

If you are able to try Shapely with the latest GEOS it would be very interesting to see a comparison between the algorithms.

[lexelby]

Well, I set an arbitrary depth limit of 20, but I've personally never hit it, and we haven't had a user report hitting it. It's possible to hit it though, of course. We'll let you know if we try the new shapely/geos.

[dr-jts]

@lexelby are there any example inputs or unit test cases that I can try against the GEOS algorithm?

[lexelby]

Honestly, when I test this, I just draw a looping path in Inkscape. Anything with a self-crossing should do the job.