Python中的三角网格查询

debugcn 发表于 Dev

贾尼斯·埃德曼尼斯

我正在寻找一个支持网状查询的Python库。现在，我已经看过了openmesh，但是我有点担心这对于我的小型硕士论文项目来说是一个过大的杀伤力。我需要的功能是：

遍历给定顶点的顶点
遍历所有边缘，面，顶点
轻松将函数值与每个顶点，面，边关联（我认为这些几何实体已被索引）

如果我真的很成功，我可能还需要：

更改网格的拓扑，例如添加或删除顶点

有可能这样做，numpy这样我就可以减少我的出行清单吗？现在，我计划将使用distmesh（pydistmesh）生成初始网格。它是否包含对我的网格查询有用的零件？

贾尼斯·埃德曼尼斯

随着改进的基于人脸的数据结构被CGAL使用，这些查询变得非常容易和有效。在这里，我实现了围绕一个特定顶点的代码：

# The demonstration of improved face based data structure

from numpy import array

triangles = array([[ 5,  7, 10],
                  [ 7,  5,  6],
                  [ 4,  0,  3],
                  [ 0,  4,  6],
                  [ 4,  7,  6],
                  [ 4,  9, 10],
                  [ 7,  4, 10],
                  [ 0,  2,  1],
                  [ 2,  0,  6],
                  [ 2,  5,  1],
                  [ 5,  2,  6],
                  [ 8,  4,  3],
                  [ 4, 11,  9],
                  [ 8, 11,  4],
                  [ 9, 11,  3],
                  [11,  8,  3]], dtype=int)

points = array([[ 0.95448092,  0.45655774],
                [ 0.86370317,  0.02141752],
                [ 0.53821089,  0.16915935],
                [ 0.97218064,  0.72769053],
                [ 0.55030382,  0.70878147],
                [ 0.34692982,  0.08765148],
                [ 0.46289581,  0.29827649],
                [ 0.21159925,  0.39472549],
                [ 0.61679844,  0.79488884],
                [ 0.4272861 ,  0.93375762],
                [ 0.12451604,  0.54267654],
                [ 0.45974728,  0.91139648]])

import pylab as plt

fig = plt.figure()
pylab.triplot(points[:,0],points[:,1],triangles)

for i,tri in enumerate(triangles):
    v1,v2,v3 = points[tri]
    vavg = (v1 + v2 + v3)/3
    plt.text(vavg[0],vavg[1],i)

#plt.show()

## constructing improved face based data structure

def edge_search(v1,v2,skip):
    """
    Which triangle has edge with verticies i and j and aren't triangle <skip>?
    """
    neigh = -1
    for i,tri in enumerate(triangles):
        if (v1 in tri) and (v2 in tri):
            if i is skip:
                continue
            else:
                neigh = i
                break

    return(neigh)


def triangle_search(i):
    """
    For given vertex with index i return any triangle from neigberhood
    """
    for i,tri in enumerate(triangles):
        if i in tri:
            return(i)

neighberhood = []
for i,tri in enumerate(triangles):

    v1, v2, v3 = tri

    t3 = edge_search(v1,v2,i)
    t1 = edge_search(v2,v3,i)
    t2 = edge_search(v3,v1,i)

    neighberhood.append([t1,t2,t3])

neighberhood = array(neighberhood,dtype=int)

faces = []
for vi,_ in enumerate(points):
    faces.append(triangle_search(vi))

## Now walking over first ring can be implemented
def triangle_ring(vertex):

    tri_start = faces[vertex]
    tri = tri_start

    ## with asumption that vertex is not on the boundary
    for i in range(10):

        yield tri

        boolindx = triangles[tri]==vertex

        # permutating to next and previous vertex
        w = boolindx[[0,1,2]]
        cw = boolindx[[2,0,1]]
        ccw = boolindx[[1,2,0]]

        ct = neighberhood[tri][cw][0]

        if ct==tri_start:
            break
        else:
            tri=ct

for i in triangle_ring(6):
    print(i)

## Using it for drawing lines on plot

vertex = 6
ring_points = []

for i in triangle_ring(vertex):
    vi = triangles[i]
    cw = (vi==vertex)[[2,0,1]] 
    print("v={}".format(vi[cw][0]))
    ring_points.append(vi[cw][0])

data = array([points[i] for i in ring_points])
plt.plot(data[:,0],data[:,1],"ro")
#plt.savefig("topology.png")
plt.show()

input("Press Enter to continue...")
plt.close("all")