我有一个三维的numpy数组。我想要显示(在matplotlib中)这个数组的等值面的一个漂亮的3D图(或者更严格地说,显示一个通过在采样点之间插值定义的3D标量场的等值面)。
matplotlib的mplot3D部件提供了很好的3D绘图支持,但(据我所知)它的API没有任何东西可以简单地获取标量值的3D数组并显示等值面。但是,它确实支持显示多边形的集合,所以我可能会实现marching cubes算法来生成这样的多边形。
似乎很可能已经在某个地方实现了对scipy友好的行进立方体,而我还没有找到它,或者我错过了一些简单的实现方法。或者,我欢迎任何指向可视化3D数组数据的其他工具的指针,这些工具很容易从Python/numpy/scipy世界中使用。
发布于 2011-05-17 21:55:47
只是为了详细说明我上面的评论,matplotlib的3D绘图真的不是为等值面这样复杂的东西而设计的。它的目的是为真正简单的3D绘图生成漂亮的、出版物质量的矢量输出。它不能处理复杂的3D多边形,所以即使你自己实现行进立方体来创建等值面,它也无法正确地渲染它。
但是,您可以使用mayavi (它的mlab API比直接使用mayavi更方便),它使用VTK来处理和可视化多维数据。
作为一个快速示例(从mayavi画廊示例之一修改而来):
import numpy as np
from enthought.mayavi import mlab
x, y, z = np.ogrid[-10:10:20j, -10:10:20j, -10:10:20j]
s = np.sin(x*y*z)/(x*y*z)
src = mlab.pipeline.scalar_field(s)
mlab.pipeline.iso_surface(src, contours=[s.min()+0.1*s.ptp(), ], opacity=0.3)
mlab.pipeline.iso_surface(src, contours=[s.max()-0.1*s.ptp(), ],)
mlab.show()
发布于 2016-02-18 11:14:44
补充@DanHickstein的答案,您还可以使用trisurf
来可视化在行进立方体阶段中获得的多边形。
import numpy as np
from numpy import sin, cos, pi
from skimage import measure
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
def fun(x, y, z):
return cos(x) + cos(y) + cos(z)
x, y, z = pi*np.mgrid[-1:1:31j, -1:1:31j, -1:1:31j]
vol = fun(x, y, z)
iso_val=0.0
verts, faces = measure.marching_cubes(vol, iso_val, spacing=(0.1, 0.1, 0.1))
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.plot_trisurf(verts[:, 0], verts[:,1], faces, verts[:, 2],
cmap='Spectral', lw=1)
plt.show()
更新: 2018年5月11日
正如@DrBwts所提到的,现在marching_cubes返回4个值。下面的代码可以工作。
import numpy as np
from numpy import sin, cos, pi
from skimage import measure
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
def fun(x, y, z):
return cos(x) + cos(y) + cos(z)
x, y, z = pi*np.mgrid[-1:1:31j, -1:1:31j, -1:1:31j]
vol = fun(x, y, z)
iso_val=0.0
verts, faces, _, _ = measure.marching_cubes(vol, iso_val, spacing=(0.1, 0.1, 0.1))
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.plot_trisurf(verts[:, 0], verts[:,1], faces, verts[:, 2],
cmap='Spectral', lw=1)
plt.show()
更新时间: 2020年2月2日
除了我之前的回答之外,我还应该提到,从那时起,PyVista已经发布了,它使得这类任务变得比较容易。
遵循与前面相同的示例。
from numpy import cos, pi, mgrid
import pyvista as pv
#%% Data
x, y, z = pi*mgrid[-1:1:31j, -1:1:31j, -1:1:31j]
vol = cos(x) + cos(y) + cos(z)
grid = pv.StructuredGrid(x, y, z)
grid["vol"] = vol.flatten()
contours = grid.contour([0])
#%% Visualization
pv.set_plot_theme('document')
p = pv.Plotter()
p.add_mesh(contours, scalars=contours.points[:, 2], show_scalar_bar=False)
p.show()
结果如下所示
更新时间: 2020年2月24日
正如@HenriMenke所提到的,marching_cubes
已被重命名为marching_cubes_lewiner
。“新”代码片段如下所示。
import numpy as np
from numpy import cos, pi
from skimage.measure import marching_cubes_lewiner
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
x, y, z = pi*np.mgrid[-1:1:31j, -1:1:31j, -1:1:31j]
vol = cos(x) + cos(y) + cos(z)
iso_val=0.0
verts, faces, _, _ = marching_cubes_lewiner(vol, iso_val, spacing=(0.1, 0.1, 0.1))
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.plot_trisurf(verts[:, 0], verts[:,1], faces, verts[:, 2], cmap='Spectral',
lw=1)
plt.show()
发布于 2014-10-01 22:06:31
如果您希望将绘图保存在matplotlib中(在我看来,生成出版质量的图像比mayavi容易得多),那么可以使用marching_cubes function implemented in skimage,然后在matplotlib中使用
mpl_toolkits.mplot3d.art3d.Poly3DCollection
如上面的链接所示。Matplotlib在渲染等值面方面做得很好。这是我用一些真实的断层扫描数据做的一个例子:
https://stackoverflow.com/questions/6030098
复制相似问题