11

I have a rather special requirement for a PostGIS geometry manipulation.

Consider the following case:

enter image description here

I have a polygon and a line feature. If I do a ST_ConvexHull() or ST_ConcaveHull() I would get something like on the left side. However, what I need is the area that encloses the polygon and the line where its perpendicular lines touch the polygon. Is there a function or another multiple step way to achieve this in PostGIS?

I think with ST_ShortestLine() from each vertex of the polygon to the line it may be possible. However, I still don't know, how I can get the area from the set of shortest lines.

4
  • Given your line is always horizontal as in your example, just take the min/max x and max y from the bbox of the polygon, plus the overall y from the line, and build a new polygon.
    – geozelot
    Aug 12 '21 at 19:44
  • 2
    And if the line is not horizontal, it should be possible to rotate it and the polygon to make it horizontal, compute the area, and rotate that back.
    – dr_jts
    Aug 15 '21 at 0:29
  • 1
    @geozelot Even with the min/max X It still seems quite fiddly to "build a new polygon". Can you elaborate?
    – dr_jts
    Aug 15 '21 at 0:31
  • @dr_jts done. Specific case handling only, but not as clumsy as you may think. And fast. I added my preferred generic solution as well.
    – geozelot
    Aug 16 '21 at 7:59
11

OK. I figured it out myself. By doing a ST_ShortestLine() on each vertex of the polygon and then a ST_Convexhull() on this set of lines and a final ST_Difference() between that area and the area of the polygon I will get what I need.

0
4

In the specific case (i.e. your example) where

  • the line follows the Y axis
  • the opposite (topmost) polygon border equals the (topmost) bounding box segment

you can simply extent the polygon bbox to get the desired area:

SELECT ST_MakeBox2D(
        ST_MakePoint(ST_XMin(ply.box), ST_YMin(ln.box)),
        ST_MakePoint(ST_XMax(ply.box), ST_YMax(ply.box))
       )::GEOMETRY(POLYGON) AS geom
FROM   (
    SELECT ply.geom::BOX2D AS box
    FROM   <polygon>
) AS ply, (
    SELECT ln.geom::BOX2D AS box
    FROM   <line>
) AS ln
;

This allows for input set relations (creating this area for closest pair of geomertries only etc.) and adjustments to the resulting area.


In all other cases it would likely be more efficient to simply dump all polygon vertices, get their projection on the line, and create the convex hull from the point set directly - either

  • by pre-collecting pairs:

    SELECT ST_ConvexHull(
             ST_Collect(geom)
           ) AS geom
    FROM   (
      SELECT ST_Collect(
               dmp.geom,
               ST_LineinterpolatePoint(
                 ln.geom,
                 ST_LineLocatePoint(
                   ln.geom,
                   dmp.geom
                 )
               )
             ) AS geom
      FROM   <polygon> AS ply,
             <line> AS ln,
             LATERAL ST_DumpPoints(ply.geom) AS dmp
    ) q
    ;
    
  • or a simple UNION of all individual points:

    WITH
      pts AS (
        SELECT dmp.geom
        FROM   ply,
               LATERAL ST_DumpPoints(geom) AS dmp
      )
    SELECT ST_ConvexHull(
             ST_Collect(geom)
           ) AS geom
    FROM   (
      SELECT geom
      FROM   pts
      UNION ALL
      SELECT ST_LineinterpolatePoint(
               ln.geom,
               ST_LineLocatePoint(
                 ln.geom,
                 pts.geom
               )
             ) AS geom
      FROM   ln, pts
    ) q
    ;
    

I'd prefer the former, except for when there is already a list of CTEs.

2

Another implementation.

So, if your shapes do have vertical lines perpendicular to the horizontal line, as shown in your illustration,

run the script:

WITH
tbla AS (SELECT ST_ShortestLine(a.geom, b.geom) geom FROM (SELECT (ST_DumpPoints(geom)).geom FROM polygon) a JOIN line b ON true),
tblb AS (SELECT geom, ST_Length(geom) length FROM tbla GROUP BY geom ORDER BY length DESC LIMIT 2)
SELECT ST_PointsInStarPolygon(geom) geom FROM (SELECT ST_Union(geom) geom FROM
        (SELECT (ST_DumpPoints(geom)).geom FROM tblb) a) b;

In addition to the basic functions, don't forget about the existence of custom functions, such as: ST_PointsInStarPolygon(https://gis.stackexchange.com/a/350004/120129) - developed by Martin Davis - (user:14766).

Original spatial solutions :-)...

Update

dr-jts, this is how new simple custom spatial functions written in GEO SQL are born :-),

I won't let good geo-tools written in SQL programming language go to waste :-)...

CREATE OR REPLACE FUNCTION ST_SnapPolygonToLine(
    geom1 GEOMETRY,
    geom2 GEOMETRY)
RETURNS GEOMETRY AS  
$BODY$
WITH
tbla AS (SELECT (ST_Dump($1)).geom),
tblb AS (SELECT (ST_Dump($2)).geom),
tblc AS (SELECT ST_ShortestLine(a.geom, b.geom) geom FROM (SELECT (ST_DumpPoints(geom)).geom FROM tbla) a JOIN tblb b ON true)
SELECT ST_PointsInStarPolygon(geom) geom FROM (SELECT ST_Union(geom) geom FROM
        (SELECT (ST_DumpPoints(geom)).geom FROM tblc) a) b;
$BODY$
LANGUAGE SQL

RUN

SELECT ST_SnapPolygonToLine(a.geom, b.geom) geom FROM polygon a JOIN line b ON true

1

Here's an outline of another solution:

  1. For each vertex in the polygon shell, use ST_LineLocatePoint to find its fractional index along the line.
  2. Construct the line segments from each vertex to the line, and pick the shortest ones with max and min index
  3. Construct the subline along the baseline between the max and min indices
  4. Extract all the line segments from the polygon
  5. Polygonize the extracted and constructed line segments. This should produce two polygons
  6. Union the polygonization results

This is going to be complex to implement in SQL. A function might be a better way to implement this.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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