Creating perpendicular line transects in PostGIS

Question

I have a river network (Line with arrow as shown in image below) and I need to create perpendicular lines to the river at every 15 meters as you walk along it (see image below).

Is there a function in PostGIS capable of doing this?

Answer 1

No clear existing function to do this, but working back from the great suggestions on this post: Making perpendicular lines through each polygon edge using PostGIS?

Assuming you are good with PostGIS syntax, this might help you:

1. Split the line into 15m segments, like this: https://www.northrivergeographic.com/qgis-split-a-line-at-an-interval
2. Loop every segment in the newly segmented line and for each segment:
3. Create parallel lines with ST_OffsetCurve with a desired distance D (can be as longer than the lake/sea as shown above, we will clip them later), on both sides of the line
4. For every 15m segment, make a new line from the centroid of that segment to the centroid of the corresponding parallel lines with ST_MakeLine
5. You should have perpendicular lines, spaced 15m extending beyond the river basin/lake
6. Clip the lines with the sea region

Answer 2

Unfortunately, I can't see your figure and as a result, I suggest to create the following pgSQL geo-function for implementation, as I perceived the content of the question.

Create the geo-function:

CREATE OR REPLACE FUNCTION ST_PerpendicularTransectsFromLine(
    geom GEOMETRY,
    sl integer)
    RETURNS TABLE (geom GEOMETRY) AS 
$BODY$
WITH
geodata AS (SELECT row_number() over() AS id, geom),
linecut AS (SELECT id, ST_LineSubstring(d.geom, substart, CASE WHEN subend > 1 THEN 1 ELSE subend END) geom
FROM (SELECT id, geom, ST_Length(((geom)::geometry)) len, sl sublen FROM geodata) AS d
CROSS JOIN LATERAL (SELECT i,  (sublen * i)/len AS substart, (sublen * (i+1))/len AS subend
FROM generate_series(0, floor( d.len / sublen)::integer) AS t(i) WHERE (sublen * i)/len <> 1.0) AS d2),
rotate AS (SELECT id, (ST_Rotate(ST_Collect(geom), -pi()/2, ST_Centroid(geom))) geom FROM linecut GROUP BY id, geom),
tbld AS (SELECT id, (ST_Dump(geom)).geom geom FROM rotate)
        SELECT ST_MakeLine(ST_StartPoint(geom), ST_EndPoint(geom)) geom FROM tbld;
$BODY$
LANGUAGE SQL

Run for your geodata (the function action is tested for geodata located in a flat rectangular coordinate system).

SELECT ST_PerpendicularTransectsFromLine(geom, 15) geom FROM <name_line_table>

As a result you should get a geo image of a river with almost (rivers are usually curves :-)) perpendicular transects of 15 meters long, as shown in the figure below.

Note: If necessary, crop them with the ST_Buffer function with the required size.

(-: FOGS :-)...

Answer 3

Here's another function to create transects, providing control over the transect length as well as the transect separation length

It operates as follows:

• split the input line into sections of length secLen
• for each section compute the transect angle and centre point of the section (the angle is Pi - azimuth(startPt, endPt))
• compute the transect lines from two points offset from the centre point

CREATE OR REPLACE FUNCTION ST_Transects(
    lineGeom geometry,
    secLen float,
    transectLen float)
    RETURNS geometry 
    LANGUAGE sql AS 
$BODY$
WITH
geodata AS (SELECT lineGeom, ST_Length(lineGeom) AS lineLen),
sections AS (SELECT ST_LineSubstring(lineGeom, secStart, CASE WHEN secEnd > 1 THEN 1 ELSE secEnd END) geom
    FROM geodata AS t
    CROSS JOIN LATERAL 
        (SELECT i * secLen/lineLen AS secStart, (i+1) * secLen/lineLen AS secEnd
            FROM generate_series(0, floor(lineLen / secLen)::integer) AS t(i) 
            WHERE (secLen * i)/lineLen <> 1.0) AS t2 ),
sectAnglePt AS (SELECT pi() - ST_Azimuth(ST_StartPoint(geom), ST_EndPoint(geom)) AS ang,
                  ST_LineInterpolatePoint(geom, 0.5) AS centre
                  FROM sections)
SELECT ST_Collect(ST_MakeLine(
                    ST_Point( ST_X(centre) - transectLen * cos(ang), 
                              ST_Y(centre) - transectLen * sin(ang)),
                    ST_Point( ST_X(centre) + transectLen * cos(ang), 
                              ST_Y(centre) + transectLen * sin(ang)) )) AS geom 
  FROM sectAnglePt;
$BODY$;

Call it like:

SELECT ST_Transects( line, 0.01, 0.005);

Answer 4

Another vector heavy solution to the issue, taking into account dr_jts (https://gis.stackexchange.com/users/14766/dr-jts?tab=profile) suggestions for curved river lines at this point in time (option 1 I posted above is better suited for straight sections).

Create an PostgreSQL/PostGIS geo-function:

CREATE OR REPLACE FUNCTION ST_PerpendicularTransectsFromRiverLine_v1(
    geom GEOMETRY,
    sl integer,
    radius_buffer float)
    RETURNS TABLE (geom GEOMETRY) AS 
$BODY$
WITH
geodata AS (SELECT row_number() over() AS id, geom),
linedivide AS (SELECT row_number() over() AS id, ST_LineSubstring(d.geom, substart, CASE WHEN subend > 1 THEN 1 ELSE subend END) geom
FROM (SELECT id, geom, ST_Length(((geom)::geometry)) len, sl sublen FROM geodata) AS d
CROSS JOIN LATERAL (SELECT i, (sublen * i)/len AS substart, (sublen * (i+1))/len AS subend
FROM generate_series(0, floor(d.len/sublen)::integer) AS t(i) WHERE (sublen * i)/len <> 1.0) AS d2),
tblb AS (SELECT id, ST_Buffer((geom), 0.5*(radius_buffer+0.000001),'endcap=round join=round') as geom FROM linedivide),
tblc_l AS (SELECT a.id, ST_LineMerge(ST_Difference(ST_ExteriorRing(a.geom), b.geom)) geom FROM tblb a LEFT JOIN tblb b ON (ST_Intersects(a.geom, b.geom) AND a.id < b.id)),
tblcpoly_l AS (SELECT id, ST_MakePolygon(ST_AddPoint(geom, ST_StartPoint(geom))) geom FROM tblc_l),
tbld_r AS (SELECT a.id, ST_LineMerge(ST_Difference(ST_ExteriorRing(a.geom), b.geom)) geom FROM tblb a LEFT JOIN tblb b ON (ST_Intersects(a.geom, b.geom) AND a.id > b.id)),
tbldpoly_r AS (SELECT id, ST_ExteriorRing(ST_MakePolygon(ST_AddPoint(geom, ST_StartPoint(geom)))) geom FROM tbld_r),
inter AS (SELECT ST_Intersection(a.geom, b.geom) as geom FROM tblcpoly_l a JOIN tbldpoly_r b ON ST_Intersects(a.geom, b.geom) AND a.id=b.id),
cut AS (SELECT ST_Intersection(ST_Buffer(a.geom, -0.000001), b.geom) geom FROM tblb a JOIN inter b ON ST_Intersects(a.geom, b.geom))
SELECT ST_Union(geom) geom FROM cut;
$BODY$
LANGUAGE SQL

Run:

SELECT ST_PerpendicularTransectsFromRiverLine_v1(geom,15,5) geom FROM <line_name_table>

Implementation sources:

1. https://postgis.net/docs/ST_LineSubstring.html;
2. https://gis.stackexchange.com/a/250496/120129.

Limitation: The segment length must be longer than the transect length!

Warning: Artifacts may appear in steep river bends!

Another vector solution:

Create an PostgreSQL/PostGIS geo-function:

CREATE OR REPLACE FUNCTION ST_PerpendicularTransectsFromRiverLine_v2(
    geom GEOMETRY,
    sl integer,
    radius_buffer float)
    RETURNS TABLE (geom GEOMETRY) AS 
$BODY$
WITH
geodata AS (SELECT row_number() over() AS id, geom),
linedivide AS (SELECT row_number() over() AS id, ST_LineSubstring(d.geom, substart, CASE WHEN subend > 1 THEN 1 ELSE subend END) geom
FROM (SELECT id, geom, ST_Length(((geom)::geometry)) len, sl sublen FROM geodata) AS d
CROSS JOIN LATERAL (SELECT i, (sublen * i)/len AS substart, (sublen * (i+1))/len AS subend
FROM generate_series(0, floor(d.len/sublen)::integer) AS t(i) WHERE (sublen * i)/len <> 1.0) AS d2),
bound AS (SELECT ST_Boundary(ST_Buffer(geom, radius_buffer, 'endcap=flat join=round')) geom FROM linedivide), 
segments AS (SELECT ST_MakeLine(pt1, pt2) geom FROM (SELECT ST_PointN(geom, generate_series(1, ST_NPoints(geom)-1)) pt1,
  ST_PointN(geom, generate_series(2, ST_NPoints(geom))) pt2 FROM bound) AS geom)
SELECT ST_Union(a.geom) geom FROM segments a JOIN geodata b ON ST_Intersects(a.geom,b.geom); 
$BODY$
LANGUAGE SQL

Run:

SELECT ST_PerpendicularTransectsFromRiverLine_v2(geom,15,5) geom FROM <line_name_table>

the choice is yours...

(-: FOGS :-)...