3D contour plot from data using Mayavi / Python

Question

I would like to make a three-dimensional contour graph using Mayavi, just like the third figure on this page (model of a cloud of hydrogen electrons):

http://www.sethanil.com/python-for-reseach/5

I have a set of data points that I created using my own model that I would like to use. Data points are stored in a multidimensional array of dummies, for example:

XYZV = [[1, 2, 3, 4], [6, 7, 8, 9], ... [4, 5, 6, 7]] 

Data points are not evenly distributed in XYZ space and are not stored in any particular order. I think the example uses a mesh grid to generate data points - I looked up, but I don't understand it at all. Any help would be much appreciated?

_{(source: sethanil.com )}

Answer 1

The trick is to interpolate over the grid before plotting - I would use scipy for this. Below, R is an array of (500.3) XYZ values, and V is the “magnitude” at each XYZ point.

 from scipy.interpolate import griddata import numpy as np # Create some test data, 3D gaussian, 200 points dx, pts = 2, 100j N = 500 R = np.random.random((N,3))*2*dx - dx V = np.exp(-( (R**2).sum(axis=1)) ) # Create the grid to interpolate on X,Y,Z = np.mgrid[-dx:dx:pts, -dx:dx:pts, -dx:dx:pts] # Interpolate the data F = griddata(R, V, (X,Y,Z)) 

From here it is easy to display our data:

 from mayavi.mlab import * contour3d(F,contours=8,opacity=.2 ) 

This gives a good (lumpy) Gauss.

See the docs for griddata , note that you can change the interpolation method. If you have more points (both on the interpolated grid and on the dataset), the interpolation is getting better and better represents the main function that you are trying to illustrate. Here is the above example with 10K points and a more accurate grid:

Answer 2

You can use the delaunay3d filter to create cells from dots. Then you can create iso_surface () to output UnstructuredGrid delaunay3d. If you want ImageData, you can use the image_data_probe filter.

 import numpy as np from tvtk.api import tvtk from mayavi import mlab points = np.random.normal(0, 1, (1000, 3)) ug = tvtk.UnstructuredGrid(points=points) ug.point_data.scalars = np.sqrt(np.sum(points**2, axis=1)) ug.point_data.scalars.name = "value" ds = mlab.pipeline.add_dataset(ug) delaunay = mlab.pipeline.delaunay3d(ds) iso = mlab.pipeline.iso_surface(delaunay) iso.actor.property.opacity = 0.1 iso.contour.number_of_contours = 10 mlab.show()