Calculate median width of road polygons

Question

I'd like to calculate the median width of a polygon in PostGIS.

The problem is almost identical to the one described here, which has an excellent answer. My question is how to calculate the median, not the mean.

The polygons represent roads, and due to irregularities, the median is a more robust measure than the mean. I have the road centerlines if that would help, although due to positional inaccuracies they are not always precisely in the center.

--UPDATE--

Here's an example of where the mean would be different from the median - a cul-de-sac where the mean would be wider than the median due to the turning circle at the end. The yellow polygons are the street rights of way.

For polygons that represent roads, the median width of the road right-of-way. The polygons are roughly rectangular. — amball, Jun 08 '20 at 16:52
Yeah, but how do you define the median ? What's the difference between the mean and the median in this case ? — obchardon, Jun 08 '20 at 20:06
Run ST_ApproximateMedialAxis for a better centerline, get a bunch of points on it, for each point compute the distance to the street border ×2. — CL., Jun 09 '20 at 09:44

geozelot · Answer 1 · 2020-06-10T10:47:28.523

One dead simple method I used to use on non-curvy segments, which does have parts of the linked answers math baked in:

ST_Subdivide the polygons into minimum vertex parts, and get the largest by area (to erase outliers)
get minimum oriented envelope, its area and smaller side length (i.e. its width)
calculate fraction of minimum vertex polygon area from envelope area
get percentage length of envelope width with that fraction

Run

WITH
   env AS (
    SELECT id,
           a / ST_Area(geom) AS afrac,
           ST_ExteriorRing(geom) AS geom
    FROM   (
      SELECT id,
             ROW_NUMBER() OVER(PARTITION BY id ORDER BY ST_Area(geom) DESC) AS seq,
             ST_Area(geom) AS a,
             ST_OrientedEnvelope(geom) AS geom
      FROM   (
        SELECT  <id> AS id,
                ST_Subdivide(<geom>, 5) AS geom
        FROM    <polygons>
      ) q
    ) q
    WHERE  seq = 1
  )

SELECT id,
       LEAST(ST_Distance(ST_PointN(geom, 1), ST_PointN(geom, 2)), ST_Distance(ST_PointN(geom, 1), ST_PointN(geom, 4))) * afrac AS avg_width
FROM   env
;

to get the <id> and approximated avg_width.

Naturally, robustness of this approach is limited to similar shapes as seen below.

Compare these calculated widths to their approximated positions along example polygons:
[green: base polygons, subdivided; orange: largest part; red line: line of calculated length along polygon]

score 0 · Accepted Answer · answered Jun 10 '20 at 18:28

The basic approach is:

Create points on the street centerline (in the code below, I had already generated these)
Split the exterior ring of the street polygon into 5 segments, so that for the median point, the nearest two segments will be on opposite sides of the street)
Calculate the distance from each point to each segment, and for each point, calculate the width as the average of the smallest two distances
Take the median across points, and multiply by 2

WITH frcs AS
  (SELECT generate_series(0, 0.9, 0.2) AS frc1, 
          generate_series(0.2, 1, 0.2) AS frc2),
     rings AS -- exterior ring of polygon, split into 5 substrings
  (SELECT id, ST_LineSubstring(ST_ExteriorRing(geom), frc1, frc2) AS geom 
      FROM frcs, street_polygons)

SELECT id, percentile_cont(0.5) within group (order by avgdist) *2 AS width -- median distance * 2
FROM 
    (SELECT id, AVG(dist) avgdist -- average distance to the two closest segments of the polygon ring
     FROM
        (SELECT ROW_NUMBER() OVER (PARTITION BY id, ptid ORDER BY dist) AS r, *
         FROM 
            (SELECT id, ptid, ST_Distance(pt.geom, rings.geom) dist 
                FROM rings 
                JOIN street_points pt USING (id)) t0) t1
          WHERE r<3 -- limit to 2 closest segments
          GROUP BY id, ptid) t2 
    GROUP BY id;

Thanks to @CL for the suggestion to compute the distance from a bunch of points on the medial axis to the street border. Since I already had the centerlines, it was more efficient to calculate the distance from points on the centerlines to the street border, rather than calculate the medial axis.

However, since the centerlines are not exact centerlines, the approach here is to average the distance to the nearest two segments on the street border (which, for the median point, will be on opposite sides of the street)

Calculate median width of road polygons

2 Answers2