4

The problem I am trying to solve is to divide a given polygon into smaller polygons based on a given input number (proportion of the polygon).

Example polygon

In this example, the desired output would be a MultiPolygon composed of the coordinates (geometry) of the orange polygon and the green polygon or something similar.

Also, in this case, the input was 3x. The input is always a given side of the polygon (in this case, the height, x) multiplied or divided by a number (in this case, 3). The remaining area will vary depending on the polygon (here it is 2x).

I have thought about calculating the area of the polygon and then getting the area of the smaller polygons but I do not know how to get the geometry (polygon coordinates) of those smaller polygons. (Also, I do not know if it is possible to do this with PostGIS)

Is there any function or script capable of achieving this?

Taras
  • 32,823
  • 4
  • 66
  • 137
Belén F.
  • 49
  • 2
  • 1
    The polygons can be any shape? Does it have to be a straight split line? – BERA Sep 05 '23 at 05:49
  • Hi, thank you for your comment! Yes, actually the polygon will represent an agricultural field (so usually they are usually irregular polygons) and the line will be one side of the given polygon multiplied or divided by a number! – Belén F. Sep 05 '23 at 06:22
  • 1
    What is the allowable area error of the parts? – Comrade Che Sep 05 '23 at 07:54
  • 2
    Does this answer your question? Splitting polygon into equal area polygons using QGIS You could split it into 5 parts and then merge 3х and 2х together. – Comrade Che Sep 05 '23 at 07:56
  • 3
    For irregular polygonal areas (and with computational completeness in mind) the only way is a heuristic approach, i.e. Voronoi or Delaunay dissection. For near-regular areas a grid based space partitioning might be used, i.e. create the oriented envelope, rotate to axis colinearity, create a grid with cell edge length of x, rotate back and merge the intersecting parts of the Polygon with each cell according to your needs. – geozelot Sep 05 '23 at 08:21
  • You might want to use a subset of this answer https://gis.stackexchange.com/questions/449199/filling-polygon-with-multiple-unequally-different-colors-and-label-3-4-5-and/449241#449241 – Kasper Sep 06 '23 at 09:23

1 Answers1

1

Here's a query that might help you - it is based on the assumption that the Polygons in question are near-rectangular, using their ST_OrientedEnvelope to create a blade at the desired ratio to split them:

SELECT
    id,
    CASE WHEN ST_Area(poly_split.geom) / ST_Area(pd.geom) < 0.5
        THEN 'smaller'
        ELSE 'larger'
    END AS part,
    poly_split.geom
FROM
    <polygon_layer> AS pd,
    LATERAL ST_Boundary(ST_OrientedEnvelope(ST_Buffer(geom, 0.00000000001))) AS obb_bd,
    LATERAL ST_Distance(ST_PointN(obb_bd, 1), ST_PointN(obb_bd, 2)) AS len_min,
    LATERAL ST_Distance(ST_PointN(obb_bd, 2), ST_PointN(obb_bd, 3)) AS len_maj,
    LATERAL ST_MakeLine(ST_PointN(obb_bd, 2), ST_PointN(obb_bd, 3)) AS maj_1,
    LATERAL ST_MakeLine(ST_PointN(obb_bd, 5), ST_PointN(obb_bd, 4)) AS maj_2,
    LATERAL LEAST((len_min * <RATIO>) / len_maj, 1.0) AS rat_frac,
    LATERAL ST_LineInterpolatePoint(maj_1, rat_frac) AS blade_sp,
    LATERAL ST_LineInterpolatePoint(maj_2, rat_frac) AS blade_ep,
    LATERAL ST_MakeLine(blade_sp, blade_ep) AS blade_ln,
    LATERAL ST_Dump(ST_Split(geom, blade_ln)) AS poly_split
;

I kept it verbose, using LATERAL expressions to run a sequence of geometric computations needed to create the blade, where

  • we introduce a tiny ST_Buffer around the original Polygons to avoid precision issues during the split operation
  • we extract the ST_Boundary to access its vertices - note the reversed vertex order for the second ST_MakeLine in the next step
  • we create the major axii maj_1 & maj_2 - this is where you can decide along what sides of the Polygons the ratio is calculated; create the minor axii if that is what you want
  • we calculate the ratio_frac at which we ST_LineInterpolatePoint the blade start and end points of the blade - this is where you need to come up with a <RATIO>; note that we always fall back to 1.0 for cases where the ratio exceeds the axis lengths!
  • we ST_Split the Polygons with the created blade_ln and ST_Dump its results into individual Polygon parts
geozelot
  • 30,050
  • 4
  • 32
  • 56