Triangulation Delaunay des points de la surface 2D en 3D avec python?
j'ai une collection de points 3D. Ces points sont échantillonnés à des niveaux constants (z=0,1,...,7). Une image doit préciser:
ces points sont dans un num PY ndarra de forme (N, 3)
X
. Au-dessus de la parcelle est créé à l'aide de:
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
X = load('points.npy')
fig = plt.figure()
ax = fig.gca(projection='3d')
ax.plot_wireframe(X[:,0], X[:,1], X[:,2])
ax.scatter(X[:,0], X[:,1], X[:,2])
plt.draw()
j'aimerais plutôt trianguler seulement la surface de cet objet, et tracer la surface. Je ne veux pas la coque convexe de cet objet, cependant, parce que cela perd la forme subtile des informations que je voudrais pouvoir inspecter.
j'ai essayé ax.plot_trisurf(X[:,0], X[:,1], X[:,2])
, mais cette résultats dans la suite de gâchis:
de l'aide?
exemple de données
voici un extrait pour générer des données 3D qui est représentatif du problème:
import numpy as np
X = []
for i in range(8):
t = np.linspace(0,2*np.pi,np.random.randint(30,50))
for j in range(t.shape[0]):
# random circular objects...
X.append([
(-0.05*(i-3.5)**2+1)*np.cos(t[j])+0.1*np.random.rand()-0.05,
(-0.05*(i-3.5)**2+1)*np.sin(t[j])+0.1*np.random.rand()-0.05,
i
])
X = np.array(X)
Exemple de données de l'image originale
Voici un pastebin à l'original données:
http://pastebin.com/YBZhJcsV
Voici les tranches le long de la constante z:
2 réponses
mise à jour 2
j'ai maintenant le faire comme suit:
- j'utilise le fait que les chemins d'accès dans chaque z-tranche sont fermés et simple d'utilisation
matplotlib.path
pour déterminer les points à l'intérieur et à l'extérieur du contour. En utilisant cette idée, je transforme les contours de chaque tranche en une image de valeur booléenne, qui est combinée en un volume de valeur booléenne. - Ensuite, j'utilise
skimage
marching_cubes
méthode pour obtenir une triangulation de la surface de visualisation.
Voici un exemple de la méthode. Je pense que les données sont légèrement différentes, mais vous pouvez certainement voir que les résultats sont beaucoup plus propres, et peut gérer les surfaces qui sont déconnectés ou ont des trous.
mise à jour de 1 (c'est toujours mauvais)
je devrais mettre à jour ceci pour les futures personnes qui le rencontrent. Alors que la méthode ci-dessus fonctionne la plupart du temps, il ne suppose (via sphérique coordonnée transformation) qu'il n'y a pas deux points le long du même rayon. Si vous remarquez l'artefact dans le milieu gauche de l'image ci-dessus, c'est pourquoi.
une meilleure approche consiste à faire de la" chirurgie " à la surface. En pensant à la surface comme une peau d'orange, on la Coupe d'un côté, puis on l'éclabousse, on l'étire. Vous avez alors un plan 2D que vous pouvez trianguler et interpoler. Vous devez simplement garder une trace de la façon de revenir à l'endroit correspondant en 3D. Pour mettre en œuvre idea prend un peu de travail, et la mise en œuvre nécessite également un soin particulier à la façon dont mes données sont représentées.
quoi qu'il en soit, ceci est juste pour donner une indication de la façon dont on pourrait aborder cela plus fermement.
réponse originale
Voici la solution que j'ai trouvée. Cela dépend fortement de mes données étant à peu près sphérique et échantillonné à uniformément dans z je pense. Certains des autres commentaires fournissent plus d'informations sur les solutions robustes. Depuis mes données grosso modo sphérique Je triangule les angles d'azimut et de Zénith à partir de la Transformée de coordonnées sphérique de mes points de données.
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
import matplotlib.tri as mtri
X = np.load('./mydatars.npy')
# My data points are strictly positive. This doesn't work if I don't center about the origin.
X -= X.mean(axis=0)
rad = np.linalg.norm(X, axis=1)
zen = np.arccos(X[:,-1] / rad)
azi = np.arctan2(X[:,1], X[:,0])
tris = mtri.Triangulation(zen, azi)
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.plot_trisurf(X[:,0], X[:,1], X[:,2], triangles=tris.triangles, cmap=plt.cm.bone)
plt.show()
en utilisant les données d'échantillon du pastebin ci-dessus, on obtient:
je me rends compte que vous avez mentionné dans votre question que vous ne vouliez pas utiliser la coque convexe parce que vous pourriez perdre quelques informations de forme. J'ai une solution simple qui fonctionne assez bien pour vos données d'exemple 'sphériques jittered', bien qu'elle utilise scipy.spatial.ConvexHull
. J'ai pensé le partager ici de toute façon, juste au cas où il serait utile pour d'autres:
from matplotlib.tri import triangulation
from scipy.spatial import ConvexHull
# compute the convex hull of the points
cvx = ConvexHull(X)
x, y, z = X.T
# cvx.simplices contains an (nfacets, 3) array specifying the indices of
# the vertices for each simplical facet
tri = Triangulation(x, y, triangles=cvx.simplices)
fig = plt.figure()
ax = fig.gca(projection='3d')
ax.hold(True)
ax.plot_trisurf(tri, z)
ax.plot_wireframe(x, y, z, color='r')
ax.scatter(x, y, z, color='r')
plt.draw()
cela se passe plutôt bien dans ce cas, puisque vos données d'exemple finissent par être une surface plus ou moins convexe. Peut-être pourriez-vous faire des exemples de données plus difficiles? Une surface toroïdale serait un bon cas d'essai que la méthode de la coque convexe échouerait évidemment.
cartographier une surface 3D arbitraire à partir d'un nuage de points est un vraiment problème difficile. Voici un question connexe contenant quelques liens qui pourraient vous être utiles.