Extending a linestring to the edge of an enclosing polygon in PostGIS

Question

I have linestrings that are each entirely within a polygon. I want to extend the first and last segments to touch the first rings of the polygon in line with them.

Answer 1

I'm sure it could be prettier but I came up with getting the first and last segments, scaling them about the inner end of the segment so they are as long as the diagonal of the bounding box of the polygon then intersecting with the boundary of the polygon and selecting the points closest to the line endpoints then combining those new points with the existing linestring to give a new linestring.

select *, 
  (select ST_MakeLine(geom) from
  ((select array[0] as path, geom from ST_Dump(ST_Intersection(boundary, ST_Scale(first_segment, ST_MakePoint(dist/ST_length(first_segment), dist/ST_length(first_segment)), ST_EndPoint(first_segment)))) order by ST_Distance(geom, ST_StartPoint(first_segment)) limit 1)
  union all
  select * from ST_Dump(line)
  union all
  (select array [ST_NumPoints(line)+1] as path, geom from ST_Dump(ST_Intersection(boundary, ST_Scale(last_segment, ST_MakePoint(dist/ST_length(last_segment), dist/ST_length(last_segment)), ST_StartPoint(last_segment)))) order by ST_Distance(geom, ST_StartPoint(last_segment)) limit 1)) as points
)
from
  (select area, line, 
    ST_Length(ST_BoundingDiagonal(area)) as dist,
    (select st_makeline(geom) from ST_DumpPoints(line) where path in (array[1],array[2])) as first_segment,
    (select st_makeLine(geom) from ST_DumpPoints(line) where path in (array[ST_NPoints(line)],array[ST_NPoints(line)-1])) as last_segment,
    ST_Boundary(area) as boundary
    from (values ('POLYGON((-10 -25, -15 15, 20 20, 22 -10, -10 -25))'::geometry,'LINESTRING(-7 -10, -5 -2, 6 7, 15 8)'::geometry)) as t(area,line)) as parts

Answer 2

So let's sum it up,

1. My answer to your question is as follows:

Input data: table (layer) named line and table (layer) named polygon.

Run the script:

WITH
    tbla AS (SELECT (ST_Dump(geom)).geom FROM line), 
    tblb 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 tbla) geom),
    tblc AS (SELECT (a.geom) geom FROM tblb a JOIN tbla b ON ST_Intersects(a.geom, ST_StartPoint(b.geom))),
    tbld AS (SELECT (a.geom) geom FROM tblb a JOIN tbla b ON ST_Intersects(a.geom, ST_EndPoint(b.geom))),
    tble AS (SELECT (geom) geom FROM polygon),
    tblf AS (SELECT ST_Azimuth(ST_EndPoint(a.geom), ST_StartPoint(a.geom)) as azimuth, ST_Distance(ST_StartPoint(a.geom), ST_EndPoint(a.geom)) + ST_Length(ST_BoundingDiagonal(b.geom)) length FROM tblc a, tble b),
    tblg AS (SELECT ST_MakeLine(ST_StartPoint(a.geom), ST_Translate(ST_StartPoint(a.geom), sin(azimuth)*length, cos(azimuth)*length)) geom FROM tblc a, tblf b),
    tblh AS (SELECT ST_Azimuth(ST_StartPoint(a.geom), ST_EndPoint(a.geom)) as azimuth, ST_Distance(ST_StartPoint(a.geom), ST_EndPoint(a.geom)) + ST_Length(ST_BoundingDiagonal(b.geom)) length FROM tbld a, tble b),
    tbli AS (SELECT ST_MakeLine(ST_EndPoint(a.geom), ST_Translate(ST_EndPoint(a.geom), sin(azimuth)*length, cos(azimuth)*length)) geom FROM tbld a, tblh b),
    tblj AS (SELECT geom FROM tblg UNION SELECT geom FROM tbli),
    tblk AS (SELECT (ST_Dump(ST_Intersection(a.geom, b.geom))).geom geom FROM tblj a JOIN LATERAL (SELECT (geom) AS geom FROM tble) AS b ON true),
    tbll AS (SELECT (ST_Dump(a.geom)).geom geom FROM tblk a, tblk b, tblb c WHERE ST_Touches(a.geom, b.geom) OR ST_Intersects(a.geom, c.geom) UNION SELECT geom FROM tblb)
             SELECT ST_Union(geom) geom FROM tbll

Check the result.

The script is called: ST_ExpandEndsOfLineInDirectionsInToPolygon

2. As a result of my reflections, I created another custom function,

!Perhaps someday the function will work...🙂!

CREATE OR REPLACE FUNCTION ST_TransectsFromLine(
    geom GEOMETRY,
    length_left double precision,
    length_rigth double precision)
    RETURNS GEOMETRY AS
$BODY$
    WITH
        tbla AS (SELECT (ST_Dump(geom)).geom geom), 
        tblb 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 tbla) AS geom),
        tblc AS (SELECT (a.geom) geom FROM tblb a JOIN tbla b ON ST_Intersects(a.geom, ST_StartPoint(b.geom))),
        tbld AS (SELECT (a.geom) geom FROM tblb a JOIN tbla b ON ST_Intersects(a.geom, ST_EndPoint(b.geom))),
        tble AS (SELECT ST_Azimuth(ST_EndPoint(geom), ST_StartPoint(geom)) as azimuth, (((ST_Distance(ST_StartPoint(geom), ST_EndPoint(geom))/1855)/60)*pi()/180) + length_left length FROM tblc),
        tblf AS (SELECT ST_MakeLine(ST_StartPoint(a.geom), ST_Translate(ST_StartPoint(a.geom), sin(azimuth)*length, cos(azimuth)*length)) geom FROM tblc a, tble b LIMIT 1),
        tblg AS (SELECT ST_Azimuth(ST_StartPoint(geom), ST_EndPoint(geom)) as azimuth, (((ST_Distance(ST_StartPoint(geom), ST_EndPoint(geom))/1855)/60)*pi()/180) + length_rigth length FROM tbld),
        tblh AS (SELECT ST_MakeLine(ST_EndPoint(a.geom), ST_Translate(ST_EndPoint(a.geom), sin(azimuth)*length, cos(azimuth)*length)) geom FROM tbld a, tblg b LIMIT 1)
                SELECT geom FROM tbla UNION SELECT geom FROM tblf UNION SELECT geom FROM tblh
$BODY$
LANGUAGE SQL

Using:

SELECT ST_TransectsFromLine(geom, 1, 1) geom FROM <line_table>

The function is called: ST_TransectsFromLine

In the meantime, use a query with WITH

WITH
    tbla AS (SELECT (ST_Dump(geom)).geom geom FROM <name_table>), 
    tblb 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 tbla) AS geom),
    tblc AS (SELECT (a.geom) geom FROM tblb a JOIN tbla b ON ST_Intersects(a.geom, ST_StartPoint(b.geom))),
    tbld AS (SELECT (a.geom) geom FROM tblb a JOIN tbla b ON ST_Intersects(a.geom, ST_EndPoint(b.geom))),
    tble AS (SELECT ST_Azimuth(ST_EndPoint(geom), ST_StartPoint(geom)) as azimuth, (((ST_Distance(ST_StartPoint(geom), ST_EndPoint(geom))/1855)/60)*pi()/180) + 1 length FROM tblc),
    tblf AS (SELECT ST_MakeLine(ST_StartPoint(a.geom), ST_Translate(ST_StartPoint(a.geom), sin(azimuth)*length, cos(azimuth)*length)) geom FROM tblc a, tble b LIMIT 1),
    tblg AS (SELECT ST_Azimuth(ST_StartPoint(geom), ST_EndPoint(geom)) as azimuth, (((ST_Distance(ST_StartPoint(geom), ST_EndPoint(geom))/1855)/60)*pi()/180) + 1 length FROM tbld),
    tblh AS (SELECT ST_MakeLine(ST_EndPoint(a.geom), ST_Translate(ST_EndPoint(a.geom), sin(azimuth)*length, cos(azimuth)*length)) geom FROM tbld a, tblg b LIMIT 1)
            SELECT geom FROM tbla UNION SELECT geom FROM tblf UNION SELECT geom FROM tblh

3. And one more experimental custom function:

CREATE OR REPLACE FUNCTION ST_TransectsAndRotateFromLine(
   geom GEOMETRY,
   length_left double precision,
   length_rigth double precision,
   left_angle double precision,
   rigth_angle double precision)
RETURNS TABLE (id bigint, geom GEOMETRY, length double precision) AS 
$BODY$
     WITH
        tbla AS (SELECT (ST_Dump(geom)).geom), 
        tblb 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 tbla) AS geom),
        tblc AS (SELECT (a.geom) geom FROM tblb a JOIN tbla b ON ST_Intersects(a.geom, ST_StartPoint(b.geom))),
        tbld AS (SELECT (a.geom) geom FROM tblb a JOIN tbla b ON ST_Intersects(a.geom, ST_EndPoint(b.geom))),
        tble AS (SELECT ST_Azimuth(ST_EndPoint(geom), ST_StartPoint(geom)) as azimuth, (((ST_Distance(ST_StartPoint(geom), ST_EndPoint(geom))/1855)/60)*pi()/180) + length_left length FROM tblc),
        tblf AS (SELECT ST_MakeLine(ST_StartPoint(a.geom), ST_Translate(ST_StartPoint(a.geom), sin(azimuth)*length, cos(azimuth)*length)) geom FROM tblc a, tble b),
        tblg AS (SELECT ST_Azimuth(ST_StartPoint(geom), ST_EndPoint(geom)) as azimuth, (((ST_Distance(ST_StartPoint(geom), ST_EndPoint(geom))/1855)/60)*pi()/180) + length_rigth length FROM tbld),
        tblh AS (SELECT ST_MakeLine(ST_EndPoint(a.geom), ST_Translate(ST_EndPoint(a.geom), sin(azimuth)*length, cos(azimuth)*length)) geom FROM tbld a, tblg b),
        tbli AS (SELECT (ST_Rotate ((geom), azimuth-(azimuth+left_angle), ST_StartPoint(geom))) geom FROM tblf, tble LIMIT 1),
        tblj AS (SELECT (ST_Rotate ((geom), azimuth-(azimuth+rigth_angle), ST_StartPoint(geom))) geom FROM tblh, tblg LIMIT 1),
        tblk AS (SELECT ST_Union(geom) geom FROM (SELECT geom FROM tbla UNION SELECT geom FROM tbli UNION SELECT geom FROM tblj) foo),
        tbll AS (SELECT (ST_Dump(geom)).geom geom FROM tblk)
                (SELECT row_number() OVER (PARTITION BY (geom)) id, geom, ST_Length(geom::geography) length FROM tbll)              
$BODY$
LANGUAGE SQL

Using:

SELECT (ST_TransectsAndRotateFromLine(geom, 1, 1, -pi()/3, pi()/3)).* FROM <line_table>

OR

SELECT (ST_TransectsAndRotateFromLine(geom, 1, 1, degrees(ST_Azimuth(ST_EndPoint(geom),ST_StartPoint(geom))) * pi()/180, 
degrees(ST_Azimuth(ST_StartPoint(geom),ST_EndPoint(geom))) * pi()/180)).* FROM <line_table>)

The function is called: ST_TransectsAndRotateFromLine

I am responsible for my style, architecture and logic... 🙂.

Original spatial solutions ...