5

I have a polygon and a line that goes through it. That line might exit my polygon and enter it again at a later point.

So now I'd like to have a relative position inside my line, where those segments are.

As an example. I have the following polygon and a line. My result should be TWO segments, where segment one is located at about 0.1 - 0.4% of my line and segment two is located at about 0.6 - 0.9% of my line.

How can I do that using postgis?

What I've tried so far: ST_Intersections -> I do get the segments, but not the indexes. ST_LineLocatePoint -> I would get the index of a point (e.g. the start and end-points of my segments)

Problem with that is that if my line crosses itself, I get wrong results, as those start/end points would exist at least twice.

enter image description here

So what I can imagine would be nice is a function (I haven't found it and don't know if it exists) that would give me indexes of relative positions of matches as well. Am I just blind? Does such a function exist?

EDIT: So this here is a sample of a polygon I want to cut with a line. It gives me all parts of my line. Inside and outside, which would help me in my further analysis. But the order of ST_Split in this case is just not correct. All sub-segments I get from it are in correct order, expect the last two in that example. With other examples the order is reversed. But I NEED to know the proper order to tell where in my original line that intersection happened.

        "type": "LineString",
        "coordinates": [
          [
            14.965792894363402,
            47.526699241115935
          ],
          [
            14.970245361328125,
            47.52540968472776
          ],
          [
            14.967606067657472,
            47.52504019924539
          ],
          [
            14.968990087509153,
            47.52488081254758
          ],
          [
            14.966779947280882,
            47.52460550711063
          ],
          [
            14.972037076950071,
            47.52419254624608
          ],
          [
            14.97059941291809,
            47.523286918933046
          ],
          [
            14.973410367965698,
            47.52304783071382
          ]
        ]
      }'), st_geomfromgeojson('{
        "type": "Polygon",
        "coordinates": [
          [
            [
              14.96858239173889,
              47.52632976471708
            ],
            [
              14.965996742248535,
              47.52501121988183
            ],
            [
              14.96784210205078,
              47.52302609536715
            ],
            [
              14.971328973770142,
              47.52361294656684
            ],
            [
              14.96858239173889,
              47.52632976471708
            ]
          ]
        ]
      }')));
3
  • ST_Intersection will return a collection (a MultiLineString in the natural case of multi-overlay between areal and lineal features), and the underlying engine should keep order when creating the intersecting segments: ST_DumpPoints on that returned geometry returns a geometry_dump, having a path field of type NUMERIC[], where path[1] is the geometry sequence number within the collection, and path[2] the vertex sequence number of the given member. The same goes for ST_Dump, having only path[1] - with these, you should be able to identify the resulting segments.
    – geozelot
    Commented Feb 15, 2022 at 19:09
  • Thanks a lot for your help! I don't really understand it yet though. So I see that ST_DumpPoints gives me all the points each and that ST_Dump will give me the whole sequences as linestrings. My main problem: So if my input is a number of points or a linestring (the results of the dump, so basically a sub-line of my original line) how can I get the relative start/end inside the originating line? So how can I say, given that input line, and this sub-line-string -> where inside that input line is that sub-line-string? Which are the relative start/end positions?
    – Georg
    Commented Feb 16, 2022 at 6:43
  • @geozelot I updated my question with a concrete sample of ST_Split (but i've seen the same for ST_Intersection) where the order of the result is just not correct. But without that proper order, I can't find the exact places where inside my original line, my intersection lies. Would be very very glad for any further advice!
    – Georg
    Commented Mar 4, 2022 at 12:21

1 Answer 1

4

Run

SELECT ln.id,
       its.path[1] AS seq,
       ROUND(start_frac::NUMERIC, 3) AS start_frac,
       ROUND(end_frac::NUMERIC, 3) AS end_frac,
       its.geom
FROM   <polygons> AS ply
JOIN   <lines> AS ln
  ON   ST_Intersects(ply.geom, ln.geom)
CROSS JOIN LATERAL
       ST_Dump(ST_Intersection(ply.geom, ln.geom)) AS its
CROSS JOIN LATERAL
       ST_LineLocatePoint(ln.geom, ST_StartPoint(its.geom)) AS start_frac
CROSS JOIN LATERAL
       ST_LineLocatePoint(ln.geom, ST_EndPoint(its.geom)) AS end_frac
;

Here

  • we first calculate the ST_Intersection, and ST_Dump its result
  • we then find the fraction of line length at which both ST_StartPoint and ST_EndPoint of each dumped segment project on its parent line

I use [CROSS JOIN] LATERAL here since it let's us elegantly expand sets of rows returned by their statements.

The seq holds a sequence (INT starting at 1) referencing the relative segment position to their parent line - you can uniquely identify segments by ln.id & seq.


enter image description here Tick-marks correspond to 1%, 5% and 10% steps of line length.

3
  • Thank you very very much for your detailed help! This seems very promising!!! May I ask one more question. If my input line crosses itself (my input is a gps track and that quite often intersects itself), the results don't seem to be correct any more. I get a number of segments where e.g. end-fraction is way past start-fraction. Some where there are segments where start-section and end-section cover most of the line (which it does not) and then mans sections in-between. Would be super happy if you have an idea for that as well! Thanks a lot!!!!
    – Georg
    Commented Mar 7, 2022 at 8:24
  • For some reason I got a lot of null seq values with this query, but it seems to work out fine if just using order by start_frac, end_frac rather than relying on seq. Excellent answer.
    – Paul
    Commented Jul 2 at 14:00
  • @Paul yeah, that path member of the geometry_dump is an array of variable length, depending on the composition of the original geometry. Might need to switch the index for it to work in your case - or just do what you do ,)
    – geozelot
    Commented Jul 2 at 15:40

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.