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

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

Пример один: Создание и Графическое изображение 2D триангуляции Делоне

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

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);
%
% Display the Vertex and Triangle labels on the plot
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

Пример два: Создание и Графическое изображение 3D триангуляции Делоне

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

X = rand(10,3)
X = 10×3

    0.6557    0.7060    0.4387
    0.0357    0.0318    0.3816
    0.8491    0.2769    0.7655
    0.9340    0.0462    0.7952
    0.6787    0.0971    0.1869
    0.7577    0.8235    0.4898
    0.7431    0.6948    0.4456
    0.3922    0.3171    0.6463
    0.6555    0.9502    0.7094
    0.1712    0.0344    0.7547

dt = delaunayTriangulation(X)
dt = 
  delaunayTriangulation with properties:

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

tetramesh(dt, 'FaceColor', 'cyan');

% To display large tetrahedral meshes use the convexHull method to
% compute the boundary triangulation and plot it using trisurf.
% For example;
% triboundary = convexHull(dt)
% trisurf(triboundary, X(:,1), X(:,2), X(:,3), 'FaceColor', 'cyan')

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

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

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

X = rand(10,2)
X = 10×2

    0.2760    0.7513
    0.6797    0.2551
    0.6551    0.5060
    0.1626    0.6991
    0.1190    0.8909
    0.4984    0.9593
    0.9597    0.5472
    0.3404    0.1386
    0.5853    0.1493
    0.2238    0.2575

dt = delaunayTriangulation(X)
dt = 
  delaunayTriangulation with properties:

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

% The triangulation datastructure is;
dt.ConnectivityList
ans = 12×3

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

% Indexing is a shorthand way to query the triangulation. The format is
% dt(i, j) where j is the j'th vertex of the i'th triangle, standard
% indexing rules apply.
% The triangulation datastructure is
dt(:,:)
ans = 12×3

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

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

dt(2,:)
ans = 1×3

     3     8     9

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

dt(2,3)
ans = 9

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

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

     4    10     1
     3     8     9
    10     4     5

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

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

% Construct a Delaunay triangulation from
% 10 random points within a unit square
x = rand(10,1);
y = rand(10,1);
dt = delaunayTriangulation(x,y)
dt = 
  delaunayTriangulation with properties:

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

% Insert 5 additional random points
dt.Points(end+(1:5),:) = rand(5,2)
dt = 
  delaunayTriangulation with properties:

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

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

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

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

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

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

              Points: [14x2 double]
    ConnectivityList: [19x3 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');
% Plot the constrained edges in red
hold on;
plot(X(C'),X(C'+size(X,1)),'-r', 'LineWidth', 2);
hold off;

% Now delete the constraints and plot the unconstrained Delaunay
dt.Constraints = [];
subplot(2,1,2);
triplot(dt);
axis([-1 17 -1 6]);
xlabel('Unconstrained Delaunay triangulation', 'fontweight','b');

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

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

clf
load usapolygon
% Define an edge constraint between two successive
% points that make up the polygonal boundary.
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 = dt.isInterior();
patch('faces',dt(io,:), 'vertices', dt.Points, 'FaceColor','r');
axis equal;
axis([-130 -60 20 55]);
xlabel('Constrained Delaunay Triangulation of usapolygon', 'fontweight','b');

Пример семь: Изгиб реконструкции от облака точек

Этот пример подсвечивает использование Триангуляции Делоне, чтобы восстановить многоугольный контур от облака точек. Реконструкция основана на изящном алгоритме Корки.

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

% Create a set of points representing the point cloud
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);
% Construct a Delaunay Triangulation of the point set.
dt = delaunayTriangulation(x,y);
tri = dt(:,:);
% Insert the location of the Voronoi vertices into the existing
% triangulation
V = dt.voronoiDiagram();
% Remove the infinite vertex
V(1,:) = [];
numv = size(V,1);
dt.Points(end+(1:numv),:) = V;
Warning: Duplicate data points have been detected and removed.
 The Triangulation indices are defined with respect to the unique set of points in delaunayTriangulation.
% The Delaunay edges that connect pairs of sample points represent the
% boundary.
delEdges = dt.edges();
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 a point cloud', 'fontweight','b');
hold off;

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

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

% Construct a constrained Delaunay triangulation of a sample of points
% on the domain boundary.
load trimesh2d
dt = delaunayTriangulation(x,y,Constraints);
inside = dt.isInterior();
% Construct a triangulation to represent the domain triangles.
tr = triangulation(dt(inside, :), dt.Points);

% Construct a set of edges that join the circumcenters of neighboring
% triangles; the additional logic constructs a unique set of such edges.
numt = size(tr,1);
T = (1:numt)';
neigh = tr.neighbors();
cc = tr.circumcenter();
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)]';
% Plot the domain triangles in green, the domain boundary in blue and the
% medial axis in red.
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 a Polygonal Domain', 'fontweight','b');
hold off;

Пример девять: Превращение 2D Mesh в измененный контур

В этом примере показано, как превратить сетку 2D области, чтобы вместить модификацию к доменному контуру.

Шаг 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');

Шаг 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');

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');

Шаг 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');