postgis - Find start of river

Edit: Follow-up question here (= Find headwater polygons).

How can I determine the start of a river in PostGIS?

I have a river network (Multiline) and want to find the startpoints of the rivers.

I can select the startpoints (rectangles in the graphic) using

SELECT
  ST_StartPoint(ST_LineMerge(a.geom)) as stp, a.gid AS gid
FROM 

  spatial.stream AS a

and similarly the endpoints (crosses). But how can I find the start of the river - not of the line segments?

I tried something like (find the startpoints, that do not intersect with endpoints):

SELECT
  ST_StartPoint(ST_LineMerge(a.geom)) as stp, a.gid AS gid
FROM 
  spatial.stream AS a, spatial.stream AS b
WHERE   
  ST_Disjoint(ST_StartPoint(ST_LineMerge(a.geom)), ST_EndPoint(ST_LineMerge(b.geom)))

But this takes forever:

"Nested Loop  (cost=0.00..1019343418.86 rows=525839789 width=323)"
"  Join Filter: st_disjoint(st_startpoint(st_linemerge(a.geom)), st_endpoint(st_linemerge(b.geom)))"
"  ->  Seq Scan on stream_typ b  (cost=0.00..4498.18 rows=39718 width=323)"
"  ->  Materialize  (cost=0.00..6364.77 rows=39718 width=319)"
"        ->  Seq Scan on stream_typ a  (cost=0.00..4498.18 rows=39718 width=319)"

Edit: I did not run it till it finished, so I even don't know if this query returns the desired result. But looking at the returned rows, this does not look correct (too many)).

I the end I want to find the polygons (grey in the graphic) where a/any river starts (the green polygons), but not those where the rivers just passes (the red polygon, no starts inside). But having just the startpoints would be good start!

Any idea how this could be done in PostGis? (also other open source solutions like GRASS, R etc. are welcome).

Update: My idea just to extract the startpoints was not concise enough :( E.g. consinder the following situation:

The green polygons both have (true) startspoints inside and are those that I want. The red polygon has also startpoints inside, but the river flows through (so no headwater polygon). With John's solution below I get both. I only want the headwater streams. I am not sure if my idea with startpoints will lead to the solution (i make up a new question if desired). I though of two Where clauses:

• polygon contains a start point


• no flow through (= only 1 intersection of polygon with stream).

I tried this:

SELECT 
    polyg.*

FROM 
    polyg, start_points, stream
WHERE 
    st_contains(polyg.geom, start_points.geom)
    AND ST_Npoints(ST_Intersection(poly.geom, stream.geom)) < 2

But it never finishes:( Perhabs there is a better way to query the headwater polygons ?

Answer

ST_StartPoint is the correct function to find single nodes at the start of a Linestring, however it does not work with MultiLinestring, so you will need to use ST_Dump to get the constituent Linestrings. If I have understood you question correctly, you then want all start points which are not also end points for more than one line, ie, points where two rivers join.

As an example, the following MultiLinestring has 3 line segments with start points that meet at the point (0, 0), and one that flows out of it.

WITH start_nodes as 
  (SELECT ST_StartPoint((ST_Dump(ST_Geomfromtext('MULTILINESTRING((10 10, 0 0), (5 5, 0 0), 
    (1 1, 0 0), (0 0, 20 20))'))).geom) as geom),
  end_nodes as 
   (SELECT ST_EndPoint((ST_Dump(st_geomfromtext('MULTILINESTRING((10 10, 0 0), (5 5, 0 0), 
     (1 1, 0 0), (0 0, 20 20))'))).geom) as geom)
  SELECT st_astext(sn.geom) 
  FROM start_nodes sn 
  WHERE sn.geom NOT IN (SELECT geom FROM end_nodes);

which returns

 POINT(10 10)
 POINT(5 5)
 POINT(1 1)

excluding Point(0,0), as expected. The idea is to only include those start points that are not also end points.

In your case, the query could be written as,

WITH linestrings as 
  (SELECT (ST_DUMP(ST_LineMerge(geom))).geom as geom FROM spatial.stream),
start_points as

   (SELECT ST_StartPoint(geom) as geom, from linestrings),
end_points as 
   (SELECT ST_EndPoint(geom) as geom from linestrings)
SELECT sp.geom from start_points sp 
WHERE sp.geom not in (SELECT geom from end_points);

The idea it to first merge and split (st_dump) the lines making up the river, then grab the start and end points, and then finally select only those start points that are not also end points -- as this is where two, or more, rivers join.

Disclaimer: I have not tested this second query and I may have misunderstood the question or made a logical error, but I believe this is a good starting point.