How visualize 3D images without using paraview

[jgostick]

There are many external programs for viewing 3D images. We like to use Paraview since it can open a 3D tiff image with no issues, and offers ray-tracing which looks great. But it can be a pain to switch programs to see if an image is looking as expected. There is an exciting new project called Napari, which is a python package for viewing 3D images, and it can be called directly from the ipython command line. I have been playing with this a bit and cooked up the following recipe which may be helpful to people:

import napari
im = ps.generators.overlapping_spheres(shape=[300, 300, 300], r=10, porosity=0.75)
napari.view_image(~im, rendering='attenuated_mip', 
                  attenuation=1.0, 
                  ndisplay=3, 
                  gamma=2, 
                  rotate=[10, 35, 0]);

[jgostick]

I have just learned that matplotlib has a voxel function, which works almost exactly as we would like for 3D images. The only problem is that is unreasonably slow for image with any significant size. Nonetheless, it's a function worth keeping an eye on:

import porespy as ps
import matplotlib.pyplot as plt
im = ps.generators.blobs(shape=[50, 50, 50])
ax = plt.figure().add_subplot(projection='3d')
ax.voxels(im)
plt.show()

Produces something like this:

[jgostick]

This got me thinking....what about a 3D scatter plot? If the image is big enough, the points will all blend together. It's actually pretty fast, but doesn't look great:

import porespy as ps
import matplotlib.pyplot as plt
im = ps.generators.blobs(shape=[50, 50, 50])
fig = plt.figure()
ax = fig.add_subplot(projection='3d')
pts = np.where(1-im[30:60, 30:60, 30:60])
ax.scatter(*pts, marker='o');

[jgostick]

Just discovered Vedo, wow! There is a function for voxel display, with an example here.

[jgostick]

I have played with vedo a bit and it seems all our dreams have come true. Below is a basic code snippet...the vedo part is very fast!

import porespy as ps
im = ~ps.generators.rsa([100, 100, 100], r=10)
lt = ps.filters.local_thickness(im)

from vedo import Volume, show
vol = Volume(lt[:, :50, :])
lego = vol.legosurface(vmin=1, vmax=lt.max(), boundary=False)
lego.cmap('jet', on='cells', vmin=1, vmax=lt.max())
lego.addScalarBar()
show(lego).close()

[ma-sadeghi]

Wow, I'm impressed! It's quite fast, both rendering, and more importantly the interaction with the 3d image!

[jgostick]

Here is a recipe using open3d:

import open3d as o3d
pcd = o3d.geometry.PointCloud()
xyz = np.vstack(np.where(~im)).T
pcd.points = o3d.utility.Vector3dVector(xyz)
pcd.colors = o3d.utility.Vector3dVector(np.random.rand(*xyz.shape))
voxel_grid = o3d.geometry.VoxelGrid.create_from_point_cloud(pcd, voxel_size=1)
o3d.visualization.draw_geometries([voxel_grid])