Estoy tratando de construir el diagrama voronoi a partir de la cuadrícula de puntos usando el código modificado de aquí . Esta es una consulta SQL después de mis modificaciones:
DROP TABLE IF EXISTS example.voronoi;
WITH
-- Sample set of points to work with
Sample AS (SELECT ST_SetSRID(ST_Union(geom), 0) geom FROM example."MeshPoints2d"),
-- Build edges and circumscribe points to generate a centroid
Edges AS (
SELECT id,
UNNEST(ARRAY['e1','e2','e3']) EdgeName,
UNNEST(ARRAY[
ST_MakeLine(p1,p2) ,
ST_MakeLine(p2,p3) ,
ST_MakeLine(p3,p1)]) Edge,
ST_Centroid(ST_ConvexHull(ST_Union(-- Done this way due to issues I had with LineToCurve
ST_CurveToLine(REPLACE(ST_AsText(ST_LineMerge(ST_Union(ST_MakeLine(p1,p2),ST_MakeLine(p2,p3)))),'LINE','CIRCULAR'),15),
ST_CurveToLine(REPLACE(ST_AsText(ST_LineMerge(ST_Union(ST_MakeLine(p2,p3),ST_MakeLine(p3,p1)))),'LINE','CIRCULAR'),15)
))) ct
FROM (
-- Decompose to points
SELECT id,
ST_PointN(g,1) p1,
ST_PointN(g,2) p2,
ST_PointN(g,3) p3
FROM (
SELECT (gd).Path id, ST_ExteriorRing((gd).geom) g -- ID andmake triangle a linestring
FROM (SELECT (ST_Dump(ST_DelaunayTriangles(geom))) gd FROM Sample) a -- Get Delaunay Triangles
)b
) c
)
SELECT ST_SetSRID((ST_Dump(ST_Polygonize(ST_Node(ST_LineMerge(ST_Union(v, ST_ExteriorRing(ST_ConvexHull(v)))))))).geom, 2180)
INTO example.voronoi
FROM (
SELECT -- Create voronoi edges and reduce to a multilinestring
ST_LineMerge(ST_Union(ST_MakeLine(
x.ct,
CASE
WHEN y.id IS NULL THEN
CASE WHEN ST_Within(
x.ct,
(SELECT ST_ConvexHull(geom) FROM sample)) THEN -- Don't draw lines back towards the original set
-- Project line out twice the distance from convex hull
ST_MakePoint(ST_X(x.ct) + ((ST_X(ST_Centroid(x.edge)) - ST_X(x.ct)) * 2),ST_Y(x.ct) + ((ST_Y(ST_Centroid(x.edge)) - ST_Y(x.ct)) * 2))
END
ELSE
y.ct
END
))) v
FROM Edges x
LEFT OUTER JOIN -- Self Join based on edges
Edges y ON x.id <> y.id AND ST_Equals(x.edge,y.edge)
) z
Abajo: resultado de mi consulta.
Como puede ver, obtengo un diagrama de voronoi "casi" correcto, excepto los puntos resaltados que no tienen un polígono único. A continuación se muestra lo que produce el algoritmo QGIS y lo que me gustaría obtener de la consulta. ¿Alguna sugerencia donde está el problema con mi código?
postgis
sql
voronoi-thiessen
DamnBack
fuente
fuente
Respuestas:
Si bien esto soluciona el problema inmediato con la consulta de los datos en cuestión, no estoy contento con esto como una solución para el uso general y volveré a visitar esta y la respuesta anterior cuando pueda.
El problema era que la consulta original estaba usando un casco convexo en los bordes voronoi para determinar el borde exterior del polígono voronoi. Esto significaba que algunos de los bordes del voronoi no estaban llegando al exterior cuando deberían haberlo estado. Miré el uso de la funcionalidad del casco cóncavo, pero realmente no funcionó como quería.
Como solución rápida, he cambiado el casco convexo para que se construya alrededor del conjunto cerrado de bordes voronoi más un búfer alrededor de los bordes originales. Los bordes de voronoi que no se cierran se extienden una gran distancia para tratar de asegurarse de que se cruzan con el exterior y se utilizan para construir polígonos. Es posible que desee jugar con los parámetros del búfer.
fuente
ST_Union(ST_Buffer(geom))), pero continuaré probando con otros conjuntos de puntos. Mientras tanto esperaré como dijiste: una solución más general. :)Siguiendo una sugerencia de @ LR1234567 para probar la nueva funcionalidad ST_Voronoi que ha desarrollado @dbaston, la sorprendente respuesta original de @MickyT (como se indica en la pregunta de OP) y el uso de los datos originales ahora se puede simplificar para:
lo que da como resultado esta salida, idéntica a la pregunta del OP.
Sin embargo, esto sufre el mismo problema, ahora corregido en la respuesta de MickyT, que los puntos en el casco cóncavo no obtienen un polígono cerrado (como se esperaría del algoritmo). Solucioné este problema con una consulta con la siguiente serie de pasos.
Diagrama 2 que muestra los puntos en el casco cóncavo (amarillo) y los puntos más cercanos al amortiguador en el casco (verde).
La consulta obviamente podría simplificarse / comprimirse, pero lo dejé es esta forma como una serie de CTE, ya que es más fácil seguir los pasos secuencialmente de esa manera. Esta consulta se ejecuta en el conjunto de datos original en milisegundos (promedio de 11 ms en un servidor de desarrollo), mientras que la respuesta de MickyT usando ST_Delauney se ejecuta en 4800 ms en el mismo servidor. DBaston afirma que se puede obtener otro aumento de velocidad de orden de magnitud construyendo contra el tronco actual de GEOS, 3.6dev, debido a mejoras en las rutinas de triangulación.
Diagrama 3 que muestra todos los puntos ahora encerrados en un polígono
Nota: Actualmente, ST_Voronoi implica construir Postgis desde la fuente (versión 2.3 o troncal) y vincularlo con GEOS 3.5 o superior.
Editar: Acabo de ver Postgis 2.3 ya que está instalado en Amazon Web Services, y parece que el nombre de la función ahora es ST_VoronoiPolygons.
Sin duda, esta consulta / algoritmo podría mejorarse. Sugerencias de bienvenida.
fuente
Si tiene acceso a PostGIS 2.3, pruebe la nueva función ST_Voronoi, recientemente comprometida:
http://postgis.net/docs/manual-dev/ST_Voronoi.html
Hay compilaciones precompiladas para Windows: http://postgis.net/windows_downloads/
fuente