我有一个 3 维 numpy 数组。我想（在 matplotlib 中）显示该数组的等值面的漂亮 3D 图（或者更严格地说，显示通过在采样点之间插值定义的 3D 标量场的等值面）。
matplotlib 的 mplot3D 部分提供了很好的 3D 绘图支持，但是（据我所知）它的 API 没有任何东西可以简单地采用 3D 标量值数组并显示等值面。但是，它确实支持显示多边形的集合，所以大概我可以实现行进立方体算法来生成这样的多边形。
似乎很可能已经在某个地方实现了一个对 scipy 友好的行进立方体，但我还没有找到它，或者我错过了一些简单的方法。或者，我欢迎任何指向其他工具的指针，这些工具可以从 Python/numpy/scipy 世界轻松使用可视化 3D 数组数据。
只是为了详细说明我上面的评论，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()
补充@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()
正如@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()
添加到我之前的答案中，我应该提一下，从那时起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()
结果如下
正如@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()
如果你想把你的图保存在 matplotlib 中（在我看来比 mayavi 更容易产生出版质量的图像），那么你可以使用在 skimage 中实现的 marching_cubes 函数，然后在 matplotlib 中使用
mpl_toolkits.mplot3d.art3d.Poly3DCollection
如上面的链接所示。Matplotlib 在渲染等值面方面做得很好。这是我用一些真实的断层扫描数据制作的一个例子：
