Finding Station order number along way using PostGIS

Question

In QGIS and PostGIS I have two layers:

• way (id, name, geom (multistring))
• station (id, name, geom (polygon))

How do I find the station order along each way?

I tried using code of Finding all intersections of LineString and Polygon and order in which it intersects it using PostGIS

SELECT ROW_NUMBER() OVER() as id, station.nome as station, dp.path, dp.way_name, station.geom FROM s_602_ptrc.station, (
SELECT (ST_DumpPoints(geom)).path[2] as path,
       (ST_DumpPoints(geom)).geom as geom, 
       way.name as way_name from s_602_ptrc.way
) as dp 
WHERE ST_INTERSECTS(dp.geom, station.geom) 
ORDER BY dp.way_name, dp.path;

with this result.

Only the stations that are nodes of the polyline are considered - But how to consider them all?

UPDATE

the first part of the code of @Geozelot goes all null values ​​for the way_part:

Answer 1

Your query assumes that there is a vertex of your lines intersecting every polygon, but that seems not to be the case, and is a rather ineffective methodology at best; you'd need to ST_Segmentize your lines with a rather smallish max_segment_length, which effectively increases the vertex count of your LineStrings, to increase the chance that every polygon intersects a vertice.
However, you'd then also have to e.g. DISTINCT ON (station.id) to not count multiple intersections per polygon.

But what if a polygon happens to intersect vertices of two or more lines? Or none at all, as in your case here? A plain intersection predicate may very well be highly ambiguous (note that you will always have a certain ambiguity, but you can limit it to edge cases only).

---

First off, having circle Polygons for stations is likely a bad choice, at least in terms of network topology. I suggest to work with points, either the original points if available, or the Polygon centroids.

The same goes for MultiLineStrings; you will almost always have to dump them to work with them properly. And especially for large tables, the spatial indexes are much more effective (while dumping them on-the-fly will leave them having no index at all).

I would then find the nearest LineString to each station Point, and order the result by the fraction of line length of the shortest distance projection of each Point on that line. The two main tools are a LATERAL (K)NN subquery, and ST_LineLocatePoint.

I will use ST_Centroid on the Polygons, but a CTE to create the dumped LineStrings; note that these (within a CTE) will have no spatial index, and performance will suffer; I'd suggest to at least dump them into a table and index them properly:

WITH
  dumped_ways AS (
    SELECT "name",
           dmp.path[1] AS way_part,
           dmp.geom
    FROM   s_602_ptrc.way,
           LATERAL ST_Dump(geom) AS dmp
  )

SELECT ROW_NUMBER() OVER(PARTITION BY nway."name", nway.way_part ORDER BY nway._frac) AS station_id,
       station.nome AS station,
       nway."name",
       nway.way_part,
       station.geom
FROM   s_602_ptrc.station AS station
CROSS JOIN LATERAL (
  SELECT dumped_ways."name",
         dumped_ways.way_part,
         ST_LineLocatePoint(dumped_ways.geom, ST_Centroid(station.geom)) AS _frac
  FROM   dumped_ways
  ORDER BY
         dumped_ways.geom <-> ST_Centroid(station.geom)
  LIMIT  1
) AS nway
ORDER BY
       nway."name",
       nway.way_part,
       nway._frac
;

Update:

Since you followed my suggestion to dump your MultiLineStrings into a table of simple LineStrings, there's no need for either the CTE, or any reference to way_part.

Just run:

SELECT ROW_NUMBER() OVER(PARTITION BY nway."name" ORDER BY nway._frac) AS station_id,
       station.nome AS station,
       nway."name",
       station.geom
FROM   s_602_ptrc.station
CROSS JOIN LATERAL (
  SELECT way."name",
         ST_LineLocatePoint(way.geom, ST_Centroid(station.geom)) AS _frac
  FROM   s_602_ptrc.way AS way
  ORDER BY
         way.geom <-> ST_Centroid(station.geom)
  LIMIT  1
) AS nway
ORDER BY
       nway."name",
       nway._frac
;

This will find the nearest way.geom to every ST_Centroid(station.geom), calculate the fraction along the found line at which the given centroid can be projected on at shortest distance, and create an id in ORDER BY that fraction.

Note that

• this may yield wrong results in edge cases like intersections, where the closest edge to a centroid is semantically not the correct one
• the station_id increases in line direction, from ST_StartPoint to ST_EndPoint