Opened 7 years ago

Closed 4 years ago

#3176 closed enhancement (fixed)

Add ST_MinimumRectangle

Reported by: Mike Taves Owned by: dbaston
Priority: medium Milestone: PostGIS 2.5.0
Component: postgis Version: 2.1.x
Keywords: Cc:

Description

Other GIS systems have a ST_MinimumRectangle function, for example H2GIS. A similar function for PostGIS would be desirable.

See https://trac.osgeo.org/geos/ticket/729 to port the getMinimumRectangle() JTS method to GEOS.

Change History (6)

comment:1 by nw, 7 years ago

I think the following works. I think but haven't proved that a side of a minimum rectangle must contain one of the sides of the convex hull of the point geometries. So this algorithm iterates through those sides and rotates the hull to be parallel to the y axis, calculates a bounding box, and counter rotates it. This gives a number of candidate rectangles, from which I choose the smallest, then the "first" based on the sides counting out from the first point in the case of ties.

This is just to make the algorithm stable. There isn't in general a unique minimum rectangle (think of a hexagon). There could be other choice mechanisms, since the rectangle chosen (potentially) depends on the ordering of the points. The rectangle could alternatively be chosen based on the angle of the major side of the rectangle to the X axis, which I think is then independent of the point ordering.

This also doesn't handle degenerate rectangles. If the approach seem sound, once a rectangle choice method is chosen, I will implement that and take care of the degenerate cases.

create or replace function ST_MinimumRectangle(g geometry) returns geometry
language 'plpgsql' as $$
declare
        hull geometry;
begin
        hull = ST_ExteriorRing(ST_ConvexHull(g));
        -- one side must lie along the rectangle.
        -- so, for each side, rotate, bbox, counter-rotate bbox
        with sides as (
                select ST_PointN(hull, n) as a, ST_PointN(hull, n+1) as b,
                        n as side
                        from generate_series(1,ST_NPoints(hull)-1) n
                ),
        angles as (
                select side, a, b, st_azimuth(a, b) as angle from sides
        ),
        boxes as (
                select ST_Rotate(ST_Envelope(ST_Rotate(hull, -angle)),angle) as rect, side, angle from angles
        )
        select rect into hull from boxes order by ST_Area(rect), side limit 1
        ;
        return hull;
end;
$$ immutable strict;

comment:2 by nw, 7 years ago

Owner: changed from pramsey to nw

comment:3 by dbaston, 6 years ago

Any updates on this? MinimumRectangle is now available in GEOS and could be called directly.

comment:4 by robe, 5 years ago

Milestone: PostGIS FuturePostGIS Fund Me

Milestone renamed

comment:5 by dbaston, 4 years ago

Milestone: PostGIS Fund MePostGIS 2.5.0
Owner: changed from nw to dbaston

comment:6 by komzpa, 4 years ago

Resolution: fixed
Status: newclosed

In trunk since r16369

Note: See TracTickets for help on using tickets.