Я пытаюсь построить диаграмму Вороного из сетки точек, используя модифицированный код отсюда . Это запрос SQL после моих изменений:
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
Ниже - результат моего запроса.
Как видите, я получаю «почти» правильную диаграмму вороной, за исключением выделенных точек, которые не имеют уникального многоугольника. Ниже описано, что производит алгоритм QGIS и что я хотел бы получить из запроса. Любые предложения, где проблема с моим кодом?
postgis
sql
voronoi-thiessen
DamnBack
источник
источник
Ответы:
Хотя это устраняет непосредственную проблему с запросом данных, я не доволен им как решением для общего пользования, и я вернусь к этому и предыдущему ответу, когда смогу.
Проблема заключалась в том, что исходный запрос использовал выпуклую оболочку по краям вороного, чтобы определить внешнее ребро для многоугольника вороного. Это означало, что некоторые из ребер вороного не достигли наружу, когда они должны были быть. Я смотрел на использование функциональности вогнутого корпуса, но на самом деле это не сработало, как я хотел.
В качестве быстрого исправления я изменил выпуклую оболочку, которая будет построена вокруг замкнутого набора ребер вороного, плюс буфер вокруг исходных ребер. Края voronoi, которые не закрываются, вытянуты на большое расстояние, чтобы попытаться убедиться, что они пересекают внешнюю сторону и используются для построения многоугольников. Вы можете поиграть с параметрами буфера.
источник
ST_Union(ST_Buffer(geom))
), но я продолжу тестирование с другими наборами точек. Пока что я буду ждать, как вы сказали, более общего решения. :)Следуя предложению @ LR1234567, чтобы опробовать новую функциональность ST_Voronoi , разработанную @dbaston, оригинальный удивительный ответ @MickyT (как указано в вопросе OP) и использование исходных данных теперь можно упростить до:
что приводит к такому выводу, идентично вопросу ОП.
Однако это страдает от той же проблемы, которая теперь исправлена в ответе MickyT, о том, что точки на вогнутом корпусе не получают вмещающий многоугольник (как и следовало ожидать от алгоритма). Я исправил эту проблему с помощью запроса с помощью следующей серии шагов.
Диаграмма 2, показывающая точки на вогнутом корпусе (желтый) и ближайшие точки для буферизации на корпусе (зеленый).
Очевидно, что запрос можно было бы упростить / сжать, но я оставил именно эту форму в виде серии CTE, поскольку таким образом легче последовательно выполнять шаги. Этот запрос выполняется для исходного набора данных в миллисекундах (в среднем 11 мс на сервере разработки), тогда как ответ MickyT с использованием ST_Delauney выполняется за 4800 мс на том же сервере. DBaston утверждает, что еще один порядок увеличения скорости может быть получен при построении по сравнению с текущей магистралью GEOS, 3.6dev, благодаря усовершенствованиям в процедурах триангуляции.
Диаграмма 3, показывающая все точки, теперь заключенные в многоугольник
Примечание. В настоящее время ST_Voronoi включает сборку Postgis из исходного кода (версия 2.3 или транк) и связывание с GEOS 3.5 или выше.
Изменить: я только что посмотрел на Postgis 2.3, как он установлен на Amazon Web Services, и кажется, что имя функции теперь ST_VoronoiPolygons.
Без сомнения, этот запрос / алгоритм может быть улучшен. Предложения приветствуются.
источник
Если у вас есть доступ к PostGIS 2.3, попробуйте новую функцию ST_Voronoi, недавно совершенную:
http://postgis.net/docs/manual-dev/ST_Voronoi.html
Для Windows есть скомпилированные сборки - http://postgis.net/windows_downloads/
источник