Создание и редактирование триангуляций Делоне

В этом примере показано, как создать, изменить и запросить триангуляции Делоне с помощью delaunayTriangulation класс. Триангуляция Делоне является наиболее широко используемой триангуляцией в научных вычислениях. Свойства, связанные с триангуляцией, обеспечивают базис для решения множества геометрических задач. Также показано конструкция ограниченных триангуляций Делоне вместе с приложениями, охватывающими расчет медиальной оси и морфирование сетки.

Пример первый: создайте и постройте график 2-D триангуляции Делоне

В этом примере показано, как вычислить 2-D триангуляцию Делоне, а затем построить график триангуляции вместе с метками вершин и треугольников.

rng default
x = rand(10,1);
y = rand(10,1);
dt = delaunayTriangulation(x,y)
dt = 
  delaunayTriangulation with properties:

              Points: [10x2 double]
    ConnectivityList: [11x3 double]
         Constraints: []

triplot(dt)

Отображение меток вершин и треугольников на графике.

hold on
vxlabels = arrayfun(@(n) {sprintf('P%d', n)}, (1:10)');
Hpl = text(x,y,vxlabels,'FontWeight','bold','HorizontalAlignment',...
   'center','BackgroundColor','none');
ic = incenter(dt);
numtri = size(dt,1);
trilabels = arrayfun(@(x) {sprintf('T%d',x)}, (1:numtri)');
Htl = text(ic(:,1),ic(:,2),trilabels,'FontWeight','bold', ...
   'HorizontalAlignment','center','Color','blue');
hold off

Figure contains an axes. The axes contains 22 objects of type line, text.

Пример второй: создайте и постройте график 3-D триангуляции Делоне

В этом примере показано, как вычислить и построить 3-D триангуляцию Делоне.

rng default
X = rand(10,3);
dt = delaunayTriangulation(X)
dt = 
  delaunayTriangulation with properties:

              Points: [10x3 double]
    ConnectivityList: [18x4 double]
         Constraints: []

tetramesh(dt)
view([10 20])

Figure contains an axes. The axes contains 18 objects of type patch.

Чтобы отобразить большие четырехгранные сетки, используйте convexHull метод для вычисления граничной триангуляции и построения графика с помощью trisurf. Для примера:

triboundary = convexHull(dt);
trisurf(triboundary, X(:,1), X(:,2), X(:,3),'FaceColor','cyan')

Пример третий: Доступ к структуре данных триангуляции

Существует два способа доступа к структуре данных триангуляции. Один из способов - через Triangulation свойство, другой способ использует индексацию.

Создайте 2-D триангуляцию Делоне из 10 случайных точек.

rng default
X = rand(10,2);
dt = delaunayTriangulation(X)
dt = 
  delaunayTriangulation with properties:

              Points: [10x2 double]
    ConnectivityList: [11x3 double]
         Constraints: []

Один из способов доступа к структуре данных триангуляции - это ConnectivityList свойство.

dt.ConnectivityList
ans = 11×3

     8     2     3
     6     7     3
     5     2     8
     7     8     3
     7     5     8
     7     6     1
     4     7     1
     9     5     4
     4     5     7
     9     2     5
      ⋮

Индексирование является сокращенным способом запроса триангуляции. Синтаксис следующий dt(i,j), где j является jвторая вершина i1й треугольник. Применяются стандартные правила индексации.

Запрос структуры данных триангуляции с индексацией.

dt(:,:)
ans = 11×3

     8     2     3
     6     7     3
     5     2     8
     7     8     3
     7     5     8
     7     6     1
     4     7     1
     9     5     4
     4     5     7
     9     2     5
      ⋮

Второй треугольник:

dt(2,:)
ans = 1×3

     6     7     3

Третья вершина второго треугольника:

dt(2,3)
ans = 3

Первые три треугольника:

dt(1:3,:)
ans = 3×3

     8     2     3
     6     7     3
     5     2     8

Пример четвертый: Редактирование триангуляции Делоне для вставки или удаления точек

В этом примере показано, как использовать основанную на индексе индексирование для вставки или удаления точек. Более эффективно редактировать delaunayTriangulation внести незначительные изменения в отличие от воссоздания нового delaunayTriangulation с нуля это особенно верно, если набор данных велик.

Создайте триангуляцию Делоне из 10 случайных точек в модуль квадрате.

rng default
x = rand(10,1);
y = rand(10,1);
dt = delaunayTriangulation(x,y)
dt = 
  delaunayTriangulation with properties:

              Points: [10x2 double]
    ConnectivityList: [11x3 double]
         Constraints: []

Вставьте 5 дополнительных случайных точек.

dt.Points(end+(1:5),:) = rand(5,2)
dt = 
  delaunayTriangulation with properties:

              Points: [15x2 double]
    ConnectivityList: [20x3 double]
         Constraints: []

Замените пятую точку.

dt.Points(5,:) = [0 0]
dt = 
  delaunayTriangulation with properties:

              Points: [15x2 double]
    ConnectivityList: [20x3 double]
         Constraints: []

Удалите четвертую точку.

dt.Points(4,:) = []
dt = 
  delaunayTriangulation with properties:

              Points: [14x2 double]
    ConnectivityList: [18x3 double]
         Constraints: []

Пример пятый: создайте ограниченную триангуляцию Делоне

Этот пример показывает, как создать ограниченную триангуляцию Делоне и иллюстрирует эффект ограничений.

Создайте и постройте график ограниченной триангуляции Делоне.

X = [0 0; 16 0; 16 2; 2 2; 2 3; 8 3; 8 5; 0 5];
C = [1 2; 2 3; 3 4; 4 5; 5 6; 6 7; 7 8; 8 1];
dt = delaunayTriangulation(X,C);
subplot(2,1,1)
triplot(dt)
axis([-1 17 -1 6])
xlabel('Constrained Delaunay triangulation','FontWeight','b')

Постройте график ребер с ограничениями красным цветом.

hold on
plot(X(C'),X(C'+size(X,1)),'-r','LineWidth',2)
hold off

Теперь удалите ограничения и постройте график триангуляции Делоне без ограничений.

dt.Constraints = [];
subplot(2,1,2)
triplot(dt)
axis([-1 17 -1 6])
xlabel('Unconstrained Delaunay triangulation','FontWeight','b')

Figure contains 2 axes. Axes 1 contains 9 objects of type line. Axes 2 contains an object of type line.

Пример шестой: создайте ограниченную триангуляцию Делоне географической карты

Загрузите карту периметра контерминозного США. Создайте ограниченную триангуляцию Делоне, представляющую многоугольник. Эта триангуляция охватывает область, которая ограничена выпуклой оболочкой множества точек. Отфильтровать треугольники, которые находятся в области многоугольника, и построить их. Примечание: Набор данных содержит дублирующие точки данных; то есть две или более точек данных имеют одно и то же местоположение. Повторяющиеся точки отклоняются, и delaunayTriangulation переформатирует ограничения соответственно.

clf
load usapolygon

Задайте ограничение ребра между двумя последовательными точками, которые образуют полигональные контуры, и создайте триангуляцию Делоне.

nump = numel(uslon);
C = [(1:(nump-1))' (2:nump)'; nump 1];
dt = delaunayTriangulation(uslon,uslat,C);
Warning: Duplicate data points have been detected and removed.
 The Triangulation indices and constraints are defined with respect to the unique set of points in delaunayTriangulation.
Warning: Intersecting edge constraints have been split, this may have added new points into the triangulation.
io = isInterior(dt);
patch('Faces',dt(io,:),'Vertices',dt.Points,'FaceColor','r')
axis equal
axis([-130 -60 20 55])
xlabel('Constrained Delaunay Triangulation of usapolygon','FontWeight','b')

Figure contains an axes. The axes contains an object of type patch.

Пример седьмой: реконструкция кривой из облака точек

В этом примере подсвечивается использование триангуляции Делоне для восстановления полигонального контура из облака точек. Реконструкция основана на элегантном алгоритме Crust.

Ссылка: Н. Амента, М. Берн и Д. Эппштейн. Корочка и бета-скелет: комбинаторная реконструкция кривой. Графические модели и обработка изображений, 60: 125-135, 1998.

Создайте набор точек, представляющих облако точек.

numpts = 192;
t = linspace( -pi, pi, numpts+1 )';
t(end) = [];
r = 0.1 + 5*sqrt( cos( 6*t ).^2 + (0.7).^2 );
x = r.*cos(t);
y = r.*sin(t);
ri = randperm(numpts);
x = x(ri);
y = y(ri);

Построение триангуляции Делоне набора точек.

dt = delaunayTriangulation(x,y);
tri = dt(:,:);

Вставьте положение вершин Вороного в существующую триангуляцию.

V = voronoiDiagram(dt);

Удалите бесконечную вершину и отфильтровайте повторяющиеся точки с помощью unique.

V(1,:) = [];
numv = size(V,1);
dt.Points(end+(1:numv),:) = unique(V,'rows');

Ребра Delaunay, которые соединяют пары точек выборки, представляют контур.

delEdges = edges(dt);
validx = delEdges(:,1) <= numpts;
validy = delEdges(:,2) <= numpts;
boundaryEdges = delEdges((validx & validy),:)';
xb = x(boundaryEdges);
yb = y(boundaryEdges);
clf
triplot(tri,x,y)
axis equal
hold on
plot(x,y,'*r')
plot(xb,yb,'-r')
xlabel('Curve reconstruction from point cloud','FontWeight','b')
hold off

Figure contains an axes. The axes contains 194 objects of type line.

Пример восемь: Вычисление аппроксимации медиальной оси полигональной области

В этом примере показано, как создать приблизительную Медиальную Ось полигональной области с помощью ограниченной триангуляции Делоне. Медиальная ось многоугольника определяется локусом центра максимального диска внутри многоугольника.

Создайте ограниченную триангуляцию Делоне для выборки точек на контуре области.

load trimesh2d
dt = delaunayTriangulation(x,y,Constraints);
inside = isInterior(dt);

Создайте триангуляцию, чтобы представлять треугольники области.

tr = triangulation(dt(inside,:),dt.Points);

Создайте набор ребер, которые соединяются с центрами описанной окружности соседних треугольников. Дополнительная логика создает уникальный набор таких ребер.

numt = size(tr,1);
T = (1:numt)';
neigh = neighbors(tr);
cc = circumcenter(tr);
xcc = cc(:,1);
ycc = cc(:,2);
idx1 = T < neigh(:,1);
idx2 = T < neigh(:,2);
idx3 = T < neigh(:,3);
neigh = [T(idx1) neigh(idx1,1); T(idx2) neigh(idx2,2); T(idx3) neigh(idx3,3)]';

Постройте график треугольников области зеленым цветом, контур области - синим цветом и медиальную ось - красным цветом.

clf
triplot(tr,'g')
hold on
plot(xcc(neigh), ycc(neigh), '-r','LineWidth',1.5)
axis([-10 310 -10 310])
axis equal
plot(x(Constraints'),y(Constraints'),'-b','LineWidth',1.5)
xlabel('Medial Axis of Polygonal Domain','FontWeight','b')
hold off

Figure contains an axes. The axes contains 364 objects of type line.

Пример девять: Морф 2-D Mesh в изменённый контур

В этом примере показано, как искривить mesh области 2-D для внесения изменений в контур области.

Шаг 1: Загрузите данные. Mesh, которая будет морфирована, определяется trife, xfe, и yfe, которая является триангуляцией в формате грани-вершины.

load trimesh2d
clf
triplot(trife,xfe,yfe)
axis equal
axis([-10 310 -10 310])
axis equal
xlabel('Initial Mesh','FontWeight','b')

Figure contains an axes. The axes contains an object of type line.

Шаг 2: Создайте триангуляцию фона - триангуляцию Ограниченного Делоне набора точек, представляющих контур mesh. Для каждой вершины mesh вычислите дескриптор, который задает ее местоположение относительно триангуляции фона. Дескриптор является окружающим треугольником вместе с барицентрическими координатами относительно этого треугольника.

dt = delaunayTriangulation(x,y,Constraints);
clf
triplot(dt)
axis equal
axis([-10 310 -10 310])
axis equal
xlabel('Background Triangulation','FontWeight','b')

Figure contains an axes. The axes contains an object of type line.

descriptors.tri = pointLocation(dt,xfe,yfe);
descriptors.baryCoords = cartesianToBarycentric(dt,descriptors.tri,[xfe yfe]);

Шаг 3: Отредактируйте триангуляцию фона, чтобы включить необходимую модификацию в контур области.

cc1 = [210 90];
circ1 = (143:180)';
x(circ1) = (x(circ1)-cc1(1))*0.6 + cc1(1);
y(circ1) = (y(circ1)-cc1(2))*0.6 + cc1(2);
tr = triangulation(dt(:,:),x,y);
clf
triplot(tr)
axis([-10 310 -10 310])
axis equal
xlabel('Edited Background Triangulation - Hole Size Reduced','FontWeight','b')

Figure contains an axes. The axes contains an object of type line.

Шаг 4: Преобразуйте дескрипторы обратно в Декартовы координаты, используя деформированную триангуляцию фона в качестве базиса для оценки.

Xnew = barycentricToCartesian(tr,descriptors.tri,descriptors.baryCoords);
tr = triangulation(trife,Xnew);
clf
triplot(tr)
axis([-10 310 -10 310])
axis equal
xlabel('Morphed Mesh','FontWeight','b')

Figure contains an axes. The axes contains an object of type line.