Le portail francophone de la géomatique

Bonjour,

Non, st_subdivide ne renvoie pas des polygones égaux:
elle découpe des polygones pour faire en sorte que le nombre de sommets du pg soit inférieure à la valeur passée avec subdivide (256 points par défaut)

Il n'y a pas de fonctions toutes faites pour cela dans PostGIS. Peut etre avec l'extension raster en rasterizant le polygone avec pixel = surface voulue. ?

Nicolas

Bonjour, oui j'ai regardé et installer une version 2.16 pour installer le plugin mais il ne fonctionne pas. Je ne sais pas comment sont geres les anciennes bibliotheque sous Qgis...c'est dommage. Merci pour l'idée

Ca peut facilement se mettre dans une fonction:

create or replace function st_dividePolygon(pg geometry, num_parts int default 10, num_points int default 4000)
returns geometry[]
language sql
as $$
    with pg_pts as (
        SELECT (ST_Dump(ST_GeneratePoints(pg, num_points))).geom AS geom
    ), pg_pts_clustered AS (
        SELECT geom, ST_ClusterKMeans(geom, num_parts) over () AS cluster
        FROM pg_pts
    ), pg_centers AS (
        SELECT cluster, ST_Centroid(ST_collect(geom)) AS geom
        FROM pg_pts_clustered
        GROUP BY cluster
    ), pg_voronoi AS (
        SELECT (ST_Dump(ST_VoronoiPolygons(ST_collect(geom)))).geom AS geom
        FROM pg_centers
    ) SELECT array_agg(ST_Intersection(pg, b.geom))
    FROM pg_voronoi b;
$$;

Pour calculer le nombre de morceaux, st_area(geom)/max_area
La méthode ne fabrique pas des morceaux qui ont tous la meme surface, mais ont des surfaces approchantes (ca peut qd meme varier du simple au double sur des pg a la forme complexe, quand meme )

Il existe des méthodes exactes (cf lien dans le blog de Paul Ramsey).

En image: les dept de france découpés pour avoir, pour chacun, des morceaux de surface ~ 55 km2:
(seuls les dept qui seraient découpés sont calculés)

create table deptcut_ea as
with tmp as (
    select gid, round(st_area(geom) / 550000000)::int as num_parts, geom
    from dept2154
) select gid, num_parts, unnest(st_dividePolygon(geom, num_parts)) as geom
from tmp
where num_parts > 1;

Nico