3

I am doing a assignment in which I created spatial tables for 3 warehouse location. These 3 warehouses has polygon geometry. I want add query to show which gemetry is not rectangle. Below is table I created: this table shows the name of the warehouse and its POLYGON coordinates

4
  • 2
    If you check basic definition of "a rectangle" you may see "A rectangle is a quadrilateral with four right angles." So, probably you will need to [1] calculate interior angles if they all are 90° and [2] check if opposite sides are equal in their length.
    – Taras
    Oct 14 at 13:21
  • 2
    What have you tried? Do you know how to tell if four lines (or four points) form a rectangle, and you just need to know how to implement that in SQL? eg check side lengths and diagonal lengths? onlinemath4all.com/…
    – Spacedman
    Oct 14 at 13:22
  • 1
    What is a rectangle when dealing with unprojected data? Let's not forget that a degree of longitude does not have the same ground length at various latitudes, and that angles are a bit more complex over a sphere than over a flat surface..
    – JGH
    Oct 14 at 14:31
  • Because the geodetic distances are very small, it's probably reasonable to compute as if the geometries were planar.
    – dr_jts
    Oct 14 at 18:02
3

If you are using projected data, you can check if the minimum oriented bounding box of the geometry is the same as the original geometry. Since a bounding box is a rectangle, you can deduce that the input is also a rectangle.

with src (id, geom) as ( values 
   (1,'polygon((0 0, 0 1, 1 1,1 0, 0 0))'::geometry),
   (2,'polygon((0 0, 0 1, 1 2,1 0, 0 0))'::geometry))
select id, ST_Equals(ST_OrientedEnvelope(geom),geom) is_rectangle
from src;

 id | is_rectangle
----+--------------
  1 | t
  2 | f
(2 rows)
9
  • 1
    This doesn't allow for slight inaccuracies in the input data. To handle this you could compute the Hausdorff distance between the rings of the input polygon and the oriented envelope and test if it is less than a desired amount (possibly a percentage of the perimeter or area?)
    – dr_jts
    Oct 14 at 17:57
  • 1
    Or just SELECT ..., ST_Area(geom) / ST_Area(ST_OrientedEnvelope(geom)) > 0.99 - no need to care for vertex count/order, simple and efficient.
    – geozelot
    Oct 14 at 18:30
  • @geozelot yes, that's simpler, but won't detect narrow "gores".
    – dr_jts
    Oct 14 at 22:11
  • @dr_jts I mean, neither does comparing diagonals - obviously being unreliable for anything but ST_NPoints = 5 (but those rectangular shapes do exist...); no chance to detect 2 diagonals (one, yes, e.g. ST_LongestLine) and it's vulnerable to anomalies, e.g. with 'POLYGON((0 0, 0.5 0, 1 0, 1.5 0, 0.75 1, 0 0))', which, I believe, is a triangle with a perfect ratio. Hausdorff is equally vulnerable to outliers, e.g. 'POLYGON((0 0, 0.5 0.4, 1 0, 1 1, 0 1, 0 0))'. The equality approach here leaves us with no tolerance at all, which is, as we know, not fitting for floating point precision.
    – geozelot
    Oct 15 at 11:25
  • So, for all those cases which cannot hold that very strict set of requirements, i.e. "has 4 corners / 5 vertices" - and this one requirement is needed for both approaches, comparing diagonals and making the equality check more likely to succeed - one of the most reliable, general purpose 'rectangularity' checks is comparing Polygon area with the minimum oriented envelope. And it's equally precise for the 5 vertices case...
    – geozelot
    Oct 15 at 11:32
2

If rectangles are required to have exactly 4 corners, then only polygons where ST_NPoints(geom) > 5 need to be considered.

To test quadrilaterals for rectangularity, compare the lengths of their diagonals. In a rectangle the diagonal lengths are (almost) equal. Since finite numerical precision means the lengths will rarely be exactly equal, a tolerance factor is needed To use a dimension-free tolerance the lengths can be normalized by the total length.

Here's a query computing the "rectangularity ratio" from your data and a perfect rectangle:

WITH data(id, geom) AS (VALUES
  (1, 'POLYGON((144.78116 -37.824855, 144.780843 -37.826916, 144.782018 -37.827019, 144.78232 -37.82496, 144.78116 -37.824855))')
 ,(3, 'POLYGON((153.193238 -27.682795, 153.19302 -27.68375, 153.193568 -27.683843, 153.193795 -27.682894,153.193238 -27.682795))')
 ,(4, 'POLYGON ((153.1931 -27.6828, 153.1937 -27.6828, 153.19370000000004 -27.6838, 153.1931 -27.6838, 153.1931 -27.6828))')
)
SELECT id,
     (ST_Distance( ST_PointN(ST_ExteriorRing(geom), 1), ST_PointN(ST_ExteriorRing(geom), 3))
   -  ST_Distance( ST_PointN(ST_ExteriorRing(geom), 2), ST_PointN(ST_ExteriorRing(geom), 4)))
   / (ST_Distance( ST_PointN(ST_ExteriorRing(geom), 1), ST_PointN(ST_ExteriorRing(geom), 3))
   +  ST_Distance( ST_PointN(ST_ExteriorRing(geom), 2), ST_PointN(ST_ExteriorRing(geom), 4))) AS rect_ratio
FROM data;

which gives the result:

 id |      rect_ratio      
----+-----------------------
  1 | -0.025764420345746754
  3 |  -0.02484385683541937
  4 | 6.269495874867416e-12
2
  • 1
    This is really helpful. I'm so grateful. I used ST_NPoints(geom) > and it worked Oct 14 at 18:58
  • As pointed out above, just checking the "diagonal" lengths does have a failure case (triangle with multiple vertices on one edge). So there needs to be an area check or an angle check as well to eliminate this.
    – dr_jts
    Oct 15 at 16:50

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.