The illustration below shows the problem:

enter image description here

as in (a) I have a set of disjoint polygons, as geometries in PostGIS. I need something like (b), the "mosaic" of this set of polygons, building it by a "influence region" criteria... Is like a Voronoi construction (illustrated by (c)): in fact, if the polygons was points, the influence regions are Voronoi.

Summarizing: I need a SQL algorithm (or some specific for PostGIS) that generates the "mosaic" of a set of disjoint polygons. (perhaps a loop of little ST_Buffer and ST_Difference operations)

PS: I need, like Voronoi's, that space delimitation (a squared frame in (b)) is ignored.

This problem is similar to this other about lines.

EDIT (after @FelixIP comment)

I preffer to stay in vector universe, to not lost precision (ex. using ST_DelaunayTriangles and adding and subtracting interiors by the original polygons, them adapting a dual graph solution)... Some simple and automatic packages like pprepair (assisted like QGIS topological tools are not automatic). But raster is perhaps simpler and less CPU-consuming.

This "GRID process" illustration is also valid as solution, assuming that it can allowing same precision and "euclidean influence region growing".

enter image description here

In ARCGIS exists a spatial analysis tool known as Euclidean Allocation, so, perhaps there are a PostGIS similar solution, starting with the set of polygons (classifying, rasterizing and making-back the polygons).

  • Thanks @Nir, sorry the confusion with points, I am not using points, only polygons as items (a) and (b) of the illustration... Do you have a link of your clue about solution? – Peter Krauss Oct 11 '15 at 20:22
  • In arcgis there is a raster solution called eucledean allocation or proximity – FelixIP Oct 12 '15 at 5:10
  • Thanks @FelixIP, I edited to "well come raster solutions" ;-) – Peter Krauss Oct 13 '15 at 14:47
  • If you erase your polygons from their extent.aspolygon, you'll be left with polygon. Skeleton of it gis.stackexchange.com/questions/177/… is what you need I guess. Implementation is a biggie though – FelixIP Oct 22 '15 at 2:02
  • 1
    @Cyril, I can review your scripts... Next week. Posting something that seems a solution with modern PostGIS functions, will be a good first step. – Peter Krauss Apr 6 '19 at 16:28

So, I will prepare a cake for you - fruit platter, using PostGIS tools, as you requested, if I correctly understood the question, and as I mentioned, the responsibility for the operation of the PostGIS oven is borne by her creative team.

I will ask not to be offended by anyone in my humorous style and to understand it as a game!

The original file is sliced fruit and simple shapes (hereinafter referred to as fruit), see Figure 1 below.

enter image description here

Here is my recipe, and I will be helped in this by dear programmers, whom you will learn about later. Let's begin, and for this we will create a dough in which our fruits will be laid, for which run the script:

create table poly_extent as SELECT ST_SetSRID(ST_Buffer(ST_Envelope(ST_Extent(geom)),0.05),4326) as geom FROM poly;

See the result in Figure 2 below

enter image description here

Now, if there are few fruits, like in my picture, create the border of the external buffer on the fruit, or if there are many fruits, create the border of the negative buffer, for which run the script:

create table poly_buff_dump as SELECT ((ST_Dump(ST_Boundary(ST_Union(ST_Buffer((geom),0.01, 'join=mitre mitre_limit=5.0'))))).geom) geom FROM poly;

And slice the buffer lines around each fruit

UPDATE poly_buff_dump SET geom=ST_RemovePoint(geom, ST_NPoints(geom)-1) WHERE ST_IsClosed(geom)=true; See the result in Figure 3 below

enter image description here

(Actually, I thought that as a result I would get broken lines (such as in a circle), but if the figures are difficult, sometimes breaks are obtained, incorrect ones, for example, one side of the rectangle fell off, etc.)

Then you need to divide the obtained lines in a convenient way for you into equal segments and extract points from them

create table poly_buff_dump_pt as SELECT (ST_DumpPoints((geom))).geom geom FROM poly_buff_segm;

Result, see Figure 4 below

enter image description here

Now run the Voronoi tool, in this place I used the tool suggested by the link MickyT :https://gis.stackexchange.com/a/172246/120129 , as a result of which you will have created tables with the name “voronoi” for the fact that “my first assistant” is separate from the chef thanks from the chef! :-).

The second way in this step is to run the ST_VoronoiPolygons function.

Result, see Figure 5 below

enter image description here

Now, cut off the extra parts by running the script:

create table poly_voronoi_cut as SELECT ST_Intersection(a.geom, b.geom) geom FROM voronoi a INNER JOIN poly_extent b ON ST_Intersects(a.geom, b.geom); Result, see Figure 6 below.

enter image description here

Now run the script in order to align the geodata type in LineString:

create table poly_voronoi_dump as SELECT (ST_Dump(geom)).geom as geom FROM poly_voronoi_cut; And now I will ask "my second mate" to take up my duties and mix the cake wel (Jeff - https://gis.stackexchange.com/a/785/120129), leveling it in a single layer, and for that, thanks to me for it!

CREATE TABLE poly_overlay_cut AS SELECT geom FROM ST_Dump(( SELECT ST_Polygonize(geom) AS geom FROM ( SELECT ST_Union(geom) AS geom FROM ( SELECT ST_ExteriorRing(geom) AS geom FROM poly_voronoi_dump) AS lines ) AS noded_lines ) ); Now it’s time for me to get to work, for which I run the script:

create table poly_voronoi_union as SELECT b.id, (ST_ConvexHull(ST_Union(a.geom, b.geom))) geom FROM poly_overlay_cut a INNER JOIN poly_buff_dump b ON ST_Intersects(a.geom, b.geom) GROUP BY b.id, a.geom, b.geom; and another script:

create table poly_voronoi_union_area as SELECT ST_Union(ST_ConvexHull(ST_BuildArea(geom))) as geom FROM poly_voronoi_union GROUP BY id; see figure 7 below

enter image description here

As you can see in the picture, our cuts have small layers, which can be removed, as an option using ST_SnapToGrid (or in another way):

And finally, we will cut out our baked fruit from our pie, I even got a little tired standing by the oven, :-)

create table polygon_voronoi_result as SELECT (ST_Dump(ST_Difference(a.geom, b.geom))).geom as geom FROM poly_voronoi_union_area_snap as a JOIN poly b ON ST_Intersects(a.geom, b.geom); Result see figure 8 enter image description here

Everything from this day, now everyone will learn to bake delicious pies - fruit platter. Help yourself All, and choose the pieces you, likeenough for everyone.

(It’s a pity that I really can’t feed all people, not with electronic cakes, but with real cakes, perhaps hunger would end on Earth ...)

Edit: The cherry on the pie could look like this :-) :

tbla AS (SELECT (ST_DumpPoints(geom)).geom geom FROM poly),
tblb AS (SELECT ((ST_Dump(ST_VoronoiPolygons(ST_Collect(geom)))).geom) geom FROM tbla),
tblc AS (SELECT ST_Intersection(a.geom, b.geom) geom FROM tblb a JOIN poly_extent b ON ST_Intersects(a.geom,b.geom)),
tbld AS (SELECT id, ((ST_Dump(geom)).geom) geom FROM poly GROUP BY id, geom)
SELECT id, ST_Union(a.geom) as geom FROM tblc a JOIN tbld b ON ST_Intersects(a.geom, b.geom) GROUP BY id;


tbla AS (SELECT (ST_DumpPoints(geom)).geom geom FROM polygons),
tblb AS (SELECT ((ST_Dump(ST_VoronoiPolygons(ST_Collect(geom)))).geom) geom FROM tbla),
tblc AS (SELECT id, ((ST_Dump(geom)).geom) geom FROM polygons GROUP BY id, geom)
SELECT id, ST_Union(a.geom) geom FROM tblb a JOIN tblc b ON ST_Intersects(a.geom, b.geom) GROUP BY id;

Revise script 01.04.2020:

WITH tbla AS (
WITH atbl AS (SELECT id, (ST_ExteriorRing(((ST_Dump(geom)).geom))) geom FROM polygons),
intervals AS (SELECT generate_series (0, 501) as steps)
SELECT steps AS stp, ST_LineInterpolatePoint(geom, steps/(SELECT count(steps)::float-1 FROM intervals)) geom FROM atbl, intervals GROUP BY id, intervals.steps, geom),
tblb AS (SELECT ((ST_Dump(ST_VoronoiPolygons(ST_Collect(geom)))).geom) geom FROM tbla),
tblc AS (SELECT id, ((ST_Dump(geom)).geom) geom FROM polygons GROUP BY id, geom)
SELECT id, ST_Union(a.geom) geom FROM tblb a JOIN tblc b ON ST_Intersects(a.geom, b.geom) GROUP BY id;

With you was good and fair Mr.Baker, thank you all and good luck, :-)...

Original solutions.

This script is called: ST_VoronoiDiagramsFromPolygons.


SQL Function may look like this for now and only for one table and a small data set!!!

CREATE OR REPLACE FUNCTION ST_VoronoiDiagramsFromPolygons(
   geom GEOMETRY,
    n integer)
RETURNS TABLE (id bigint, geom GEOMETRY) AS  
tbla AS (SELECT row_number() over() AS id, ((ST_Dump(geom)).geom) geom FROM polygons),
tblb AS (WITH btbl AS (SELECT id, ST_ExteriorRing((ST_Dump(geom)).geom) geom FROM tbla GROUP BY id, geom),
    intervals AS (SELECT generate_series (0, n) as steps)
    SELECT steps AS stp, ST_LineInterpolatePoint(geom, steps/(SELECT count(steps)::float-1 FROM intervals)) geom FROM btbl, intervals GROUP BY id, stp, intervals.steps, geom),
tblc AS (SELECT (ST_Dump(ST_VoronoiPolygons(ST_Collect(geom)))).geom geom FROM tblb)
    SELECT b.id, ST_Union(a.geom) geom FROM tblc a JOIN tbla b ON ST_Intersects(a.geom, b.geom) GROUP BY id;


SELECT DISTINCT (ST_VoronoiDiagramsFromPolygons(geom, 300)).* geom FROM polygons ORDER BY id;

Maybe someday the PostgreSQL and PostGIS developers will bring its behavior up to date 🙂 ...

  • 1
    Hi, good illustrations, and seems good result (!). Suggestion to a brief discussion, see @geogeek solution, the steps in QGis seems the main steps of PostGIS here... Why e need ST_Buffer and ST_ConvexHull? There are an alternative algorithm? – Peter Krauss Apr 9 '19 at 10:47
  • 1) See, ST_Buffer and on polygons are created in order to reduce the separation lines between the original figures during the application of the Voronoi function, therefore, the larger the buffer size you can specify, the more smoothly the separation lines will be built as a Voronoi function ... 2) ST_ConvexHull from the multitude of figures on each polygon, the polygon we need ... 3) Today, this is my vision of solving your question, I don’t know what will happen to all of us and our decisions in the future ... – Cyril Mikhalchenko Apr 9 '19 at 11:45
  • Another important value of the buffer lines is that they will help us select fragments from the result of the function of the Voronoi to create combined polygons for each shape ... – Cyril Mikhalchenko Apr 9 '19 at 12:08
  • By the way, I welcome the improvement of the code and approach, so that they do not worsen the result ... – Cyril Mikhalchenko Apr 9 '19 at 16:17
  • 1
    Hi Cyril, I agree that ST_Buffer is a good option to normalize the contour lines, perfect (!). About ST_ConvexHull, is only a curiosity, what kind of results we can obtain after Figure 6, poly_voronoi_union without ST_ConvexHull. – Peter Krauss Apr 10 '19 at 20:08

Postgis do not have a dedicated function for voronoi, but Qgis contains vornoi function which could make voronoi polygons from points, so using qgis i've followed the following steps to have these result:

-make points from polygons using extract nodes functions.

-make vornoi polygons using voroi functions in Qgis.

-make a spatial join in Qgis.

-dissolve results.

enter image description here

  • 1
    Nice solution (!), and thanks the illustration. Can you enhance it? See the left rectangle below, it was 4 points at Voronoi, and the center of rectangle have no representantion, distorting its region (losting are to the triangle)... Perhaps enhance is possible doing a regular sampling of each polygon... and perhaps also some interior points. – Peter Krauss Oct 27 '15 at 17:50
  • I use densify tool in ArcGIS. This significantly improves proximity polygons. +1 – FelixIP Dec 1 '16 at 21:39

OK - Thought about this a little and found it was like something I have been looking at lately.

Take your starting polys:

enter image description here

Generate a new attribute with a number (100 in my case) Use the Vector-> Research tools -> Random points inside polygons tool this will generate (100) points inside each polygon: enter image description here

Then Vector-> Geometry tools -> Voronoi to generate polys based upon that point layer.

enter image description here

Now, you can use the Vector -> Spatial Query tool: Select the points that belong to one polygon (or one of the polygons) Use the spatial query tool to generate a selection of your voronoi polygons which apply to that polygon. Add an attribute to the voroni polygon which corresponds to the polygon of interest. (I just used 1,2,3,4)

Now you can Vector-> Geoprocessing Tools-> dissolve based upon your new attribute.

enter image description here

  • Thanks @AAmes (!), the illustrations are good, and the method is a good alternative... But the question is about "a SQL algorithm (or some specific for PostGIS)" and now, the bounty about the use of the function ST_VoronoiPolygons(), that perhaps solve all problem in one click ;-) – Peter Krauss Nov 30 '16 at 15:41
  • @AAmes A quick search reveals that this method can indeed be used in PostGIS as well. How to create random points in a polygon in postgis explains the point creation from polygons and from there it should be fairly straight forward. – Phil G Dec 1 '16 at 12:52
  • My bounty was to you, as incentive (welcome here!), but THIS ANSEWER IS WORNG for the question and bounty requirements. – Peter Krauss Dec 7 '16 at 13:42
  • Thanks Peter, I saw your notes but unfortunately I have not had the opportunity to come back and address it. I don't have any experience in PostGIS yet and it would have taken me more time than I had available to give any kind of good response. I hope that our exploration of this space is helpful to people in the future, maybe If I get a crack at PostGIS soon I will come back and take another whack at it for practice! – A.A Dec 7 '16 at 18:27

Random points are a good idea to generate a voronoi polygon from polygons, it works pretty well, but it's pretty bad for polygons close to each other: enter image description here enter image description here

ST_ApproximateMedialAxis is another good alternative if using PostGIS: Calculating Voronoi Diagrams for polygons


I posted a mini python package to make it - voronoi-diagram-for-polygons. It should be pointed out in advance that this package depends on v1.8.dev0 of shapely which is still in development. In other words, the dependency of voronoi-diagram-for-polygons cannot be installed by pip automatically. You have to install the latest Shapely firstly by the following:

pip install git+https://github.com/Toblerity/Shapely.



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.