¿Cómo mostrar un gráfico 3D de una matriz 3D isosurface en matplotlib mplot3D o similar?

Question

Tengo una matriz numpy tridimensional. Me gustaría mostrar (en matplotlib) una parcela 3D de un buen isosuperficie de esta matriz (o más estrictamente, mostrar un isosuperficie del campo escalar definido en 3D mediante interpolación entre los puntos de muestra).

La parte mplot3D de matplotlib proporciona un buen soporte para el trazado 3D, pero (hasta donde puedo ver) su API no tiene nada que simplemente tome una matriz tridimensional de valores escalares y muestre una isosuperficie. Sin embargo, sí es compatible con la visualización de una colección de polígonos, por lo que presumiblemente podría implementar el algoritmo de los cubos de marcha para generar dichos polígonos.

Parece bastante probable que ya se hayan implementado unos cubos de marcha fáciles de usar y que no los he encontrado, o que me falta una forma sencilla de hacerlo. De forma alternativa, agradecería cualquier apuntador a otras herramientas para visualizar datos de matriz 3D fácilmente utilizables del mundo de Python/numpy/scipy.

Answer 1

Sólo para explicar en detalle mi comentario anterior, 3D de matplotlib ploteo realmente no está destinado a algo tan complejo como isosuperficies. Está destinado a producir resultados vectoriales bonitos y de calidad de publicación para gráficos en 3D realmente simples. No puede manejar polígonos 3D complejos, por lo que incluso si se implementa marchas cubos para crear la isosuperficie, no se procesaría correctamente.

Sin embargo, lo que puede hacer en su lugar es utilizar mayavi (se mlab API es un poco más conveniente que usar directamente Mayavi), que utiliza VTK para procesar y visualizar datos multidimensionales.

Como un ejemplo rápido (modificado a partir de uno de los ejemplos galería 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() 

Answer 2

Si desea mantener sus parcelas en matplotlib (mucho más fácil de producir imágenes con calidad de publicación de Mayavi en mi opinión), entonces se puede utilizar la marching_cubes function implemented in skimage y luego graficar los resultados en matplotlib usando

mpl_toolkits.mplot3d.art3d.Poly3DCollection 

como se muestra en el enlace de arriba Matplotlib hace un buen trabajo al renderizar la isosuperficie. Aquí hay un ejemplo que hice de algunos datos de tomografía real:

Answer 3

Como complemento de la respuesta de @DanHickstein, también se puede utilizar para visualizar trisurf los polígonos obtenidos en la fase de cubos de marcha.

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) 
verts, faces = measure.marching_cubes(vol, 0, 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()