从使用 matplotlib 生成的 delaunay 三角剖分中获取外心

发布于 2024-10-31 08:33:06 字数 93 浏览 2 评论 0原文

如果我使用 matplotlib 为一组点生成 delaunay 三角剖分,那么获取已生成的三角形的外心的最合适方法是什么?我尚未在三角测量库中找到明显的方法来执行此操作。

If I use matplotlib to generate a delaunay triangulation for a group of points, what is the most appropraite way of getting the circumcentres of the triangles that have been geenrated? I haven't yet managed to find an obvious method in the Triangulation library to do this.

如果你对这篇内容有疑问,欢迎到本站社区发帖提问 参与讨论,获取更多帮助,或者扫码二维码加入 Web 技术交流群。

扫码二维码加入Web技术交流群

发布评论

需要 登录 才能够评论, 你可以免费 注册 一个本站的账号。

评论(3

溺孤伤于心 2024-11-07 08:33:16

我认为这个解决方案太过分了..您可以直接表示每个三角形的顶点,例如:

mid_points = tri.points[tri.vertices].mean(axis=1)

在此处输入图像描述

I think that solution is too overkill.. you can directly mean out the vertices of each triangle, like:

mid_points = tri.points[tri.vertices].mean(axis=1)

enter image description here

攀登最高峰 2024-11-07 08:33:13

这是计算它们的函数。它也可以用于其他三角测量结构,例如 scipy 的 Delaunay 三角剖分(见下文)。

def compute_triangle_circumcenters(xy_pts, tri_arr):
    """
    Compute the centers of the circumscribing circle of each triangle in a triangulation.
    :param np.array xy_pts : points array of shape (n, 2)
    :param np.array tri_arr : triangles array of shape (m, 3), each row is a triple of indices in the xy_pts array

    :return: circumcenter points array of shape (m, 2)
    """
    tri_pts = xy_pts[tri_arr]  # (m, 3, 2) - triangles as points (not indices)

    # finding the circumcenter (x, y) of a triangle defined by three points:
    # (x-x0)**2 + (y-y0)**2 = (x-x1)**2 + (y-y1)**2
    # (x-x0)**2 + (y-y0)**2 = (x-x2)**2 + (y-y2)**2
    #
    # becomes two linear equations (squares are canceled):
    # 2(x1-x0)*x + 2(y1-y0)*y = (x1**2 + y1**2) - (x0**2 + y0**2)
    # 2(x2-x0)*x + 2(y2-y0)*y = (x2**2 + y2**2) - (x0**2 + y0**2)
    a = 2 * (tri_pts[:, 1, 0] - tri_pts[:, 0, 0])
    b = 2 * (tri_pts[:, 1, 1] - tri_pts[:, 0, 1])
    c = 2 * (tri_pts[:, 2, 0] - tri_pts[:, 0, 0])
    d = 2 * (tri_pts[:, 2, 1] - tri_pts[:, 0, 1])

    v1 = (tri_pts[:, 1, 0] ** 2 + tri_pts[:, 1, 1] ** 2) - (tri_pts[:, 0, 0] ** 2 + tri_pts[:, 0, 1] ** 2)
    v2 = (tri_pts[:, 2, 0] ** 2 + tri_pts[:, 2, 1] ** 2) - (tri_pts[:, 0, 0] ** 2 + tri_pts[:, 0, 1] ** 2)

    # solve 2x2 system (see https://en.wikipedia.org/wiki/Invertible_matrix#Inversion_of_2_%C3%97_2_matrices)
    det = (a * d - b * c)
    detx = (v1 * d - v2 * b)
    dety = (a * v2 - c * v1)

    x = detx / det
    y = dety / det

    return (np.vstack((x, y))).T

来自 @JoshAdel 的数据上面的答案,添加以下代码:

cc = compute_triangle_circumcenters(np.vstack([tt.x, tt.y]).T, tt.triangle_nodes)
plt.plot(cc[:, 0], cc[:, 1], ".k")

我得到下图:

在此处输入图像描述

它也可以在 scipy.spatial.Delaunay 上使用,如下所示:

from scipy.spatial import Delaunay
xy_pts = np.vstack([x, y]).T
dt = Delaunay(xy_pts)
cc = compute_triangle_circumcenters(dt.points, dt.simplices)

Here is a function that computes them. It can also be used on other triangulation structures, e.g. scipy's Delaunay triangulation (see below).

def compute_triangle_circumcenters(xy_pts, tri_arr):
    """
    Compute the centers of the circumscribing circle of each triangle in a triangulation.
    :param np.array xy_pts : points array of shape (n, 2)
    :param np.array tri_arr : triangles array of shape (m, 3), each row is a triple of indices in the xy_pts array

    :return: circumcenter points array of shape (m, 2)
    """
    tri_pts = xy_pts[tri_arr]  # (m, 3, 2) - triangles as points (not indices)

    # finding the circumcenter (x, y) of a triangle defined by three points:
    # (x-x0)**2 + (y-y0)**2 = (x-x1)**2 + (y-y1)**2
    # (x-x0)**2 + (y-y0)**2 = (x-x2)**2 + (y-y2)**2
    #
    # becomes two linear equations (squares are canceled):
    # 2(x1-x0)*x + 2(y1-y0)*y = (x1**2 + y1**2) - (x0**2 + y0**2)
    # 2(x2-x0)*x + 2(y2-y0)*y = (x2**2 + y2**2) - (x0**2 + y0**2)
    a = 2 * (tri_pts[:, 1, 0] - tri_pts[:, 0, 0])
    b = 2 * (tri_pts[:, 1, 1] - tri_pts[:, 0, 1])
    c = 2 * (tri_pts[:, 2, 0] - tri_pts[:, 0, 0])
    d = 2 * (tri_pts[:, 2, 1] - tri_pts[:, 0, 1])

    v1 = (tri_pts[:, 1, 0] ** 2 + tri_pts[:, 1, 1] ** 2) - (tri_pts[:, 0, 0] ** 2 + tri_pts[:, 0, 1] ** 2)
    v2 = (tri_pts[:, 2, 0] ** 2 + tri_pts[:, 2, 1] ** 2) - (tri_pts[:, 0, 0] ** 2 + tri_pts[:, 0, 1] ** 2)

    # solve 2x2 system (see https://en.wikipedia.org/wiki/Invertible_matrix#Inversion_of_2_%C3%97_2_matrices)
    det = (a * d - b * c)
    detx = (v1 * d - v2 * b)
    dety = (a * v2 - c * v1)

    x = detx / det
    y = dety / det

    return (np.vstack((x, y))).T

On the data from @JoshAdel's answer above, adding the following code:

cc = compute_triangle_circumcenters(np.vstack([tt.x, tt.y]).T, tt.triangle_nodes)
plt.plot(cc[:, 0], cc[:, 1], ".k")

I get the following figure:

enter image description here

It can also be used on scipy.spatial.Delaunay like this:

from scipy.spatial import Delaunay
xy_pts = np.vstack([x, y]).T
dt = Delaunay(xy_pts)
cc = compute_triangle_circumcenters(dt.points, dt.simplices)
两仪 2024-11-07 08:33:11

您应该能够使用 matplotlib.delaunay.triangulate.Triangulation 来计算它:

三角测量(x, y)
x, y -- 作为一维浮点数数组的点的坐标



属性:(所有属性均应视为
只读以保持一致性)
x, y -- 作为一维浮点数数组的点的坐标。

circumcenters -- (ntriangles, 2) 浮点数数组,给出 (x,y)
    每个三角形的外心坐标(由 triangle_id 索引)。

改编自 matplotlib 示例之一(可能有一种更简洁的方法来执行此操作,但它应该有效):

import matplotlib.pyplot as plt
import matplotlib.delaunay
import matplotlib.tri as tri
import numpy as np
import math

# Creating a Triangulation without specifying the triangles results in the
# Delaunay triangulation of the points.

# First create the x and y coordinates of the points.
n_angles = 36
n_radii = 8
min_radius = 0.25
radii = np.linspace(min_radius, 0.95, n_radii)

angles = np.linspace(0, 2*math.pi, n_angles, endpoint=False)
angles = np.repeat(angles[...,np.newaxis], n_radii, axis=1)
angles[:,1::2] += math.pi/n_angles

x = (radii*np.cos(angles)).flatten()
y = (radii*np.sin(angles)).flatten()

tt = matplotlib.delaunay.triangulate.Triangulation(x,y)
triang = tri.Triangulation(x, y)

# Plot the triangulation.
plt.figure()
plt.gca().set_aspect('equal')
plt.triplot(triang, 'bo-')

plt.plot(tt.circumcenters[:,0],tt.circumcenters[:,1],'r.')
plt.show()

You should be able to calculate it using matplotlib.delaunay.triangulate.Triangulation:

Triangulation(x, y)
x, y -- the coordinates of the points as 1-D arrays of floats

.
.
.

Attributes: (all should be treated as
read-only to maintain consistency)
x, y -- the coordinates of the points as 1-D arrays of floats.

  circumcenters -- (ntriangles, 2) array of floats giving the (x,y)
    coordinates of the circumcenters of each triangle (indexed by a triangle_id).

Adapted from one of the matplotlib examples (there is probably a cleaner way to do this, but it should work):

import matplotlib.pyplot as plt
import matplotlib.delaunay
import matplotlib.tri as tri
import numpy as np
import math

# Creating a Triangulation without specifying the triangles results in the
# Delaunay triangulation of the points.

# First create the x and y coordinates of the points.
n_angles = 36
n_radii = 8
min_radius = 0.25
radii = np.linspace(min_radius, 0.95, n_radii)

angles = np.linspace(0, 2*math.pi, n_angles, endpoint=False)
angles = np.repeat(angles[...,np.newaxis], n_radii, axis=1)
angles[:,1::2] += math.pi/n_angles

x = (radii*np.cos(angles)).flatten()
y = (radii*np.sin(angles)).flatten()

tt = matplotlib.delaunay.triangulate.Triangulation(x,y)
triang = tri.Triangulation(x, y)

# Plot the triangulation.
plt.figure()
plt.gca().set_aspect('equal')
plt.triplot(triang, 'bo-')

plt.plot(tt.circumcenters[:,0],tt.circumcenters[:,1],'r.')
plt.show()
~没有更多了~
我们使用 Cookies 和其他技术来定制您的体验包括您的登录状态等。通过阅读我们的 隐私政策 了解更多相关信息。 单击 接受 或继续使用网站,即表示您同意使用 Cookies 和您的相关数据。
原文