Surface plots using a 4th parameter for colormap in Mayavi

Playing with Mayavi, I wanted to find an easy way to use the colormap in an mlab.surf() call to display of fourth dimension. It is a very nice way to explore datasets. The image here below shows an option valuation tool that overlays the Black & Scholes call delta value on top of the surface made using the strikes, time to expiry and Black & Scholes call value.

How can you do it ?

It is a bit tricky as there is no simple function call for that. The idea is pretty simple. You have to add a new point_scalars array to the datasource used to plot your surface and then make use of it. The following example presents a first way to do it by adding a second surface on top of the first one using the new point_scalars array for the colormap (some lines might be useless for the updates, I will correct that asap) :

src = surface.mlab_source
array_id = src.dataset.point_data.add_array(self.indicator.T.ravel())
src.dataset.point_data.get_array(array_id).name = 'indicator'    
surface.update_data()
# update the PolyDataNormals
vtkobject = surface.parent.parent
vtkobject.update_data()
updated_source = mlab.pipeline.set_active_attribute(vtkobject, point_scalars='indicator')	
updated_surface = mlab.pipeline.surface(updated_source, colormap='RdBu')

This works fine but has the drawback of creating two surfaces. Furthermore, it is not really convenient for updating lively the dataset.

The alternative to that is to build the entire pipeline from scratch adding a SetActiveAttribute filter at the right place (creating the datasource, then the pipeline : WarpScalar, PolyDataNormal, SetAttribute, Surface).

Here is the easy way to do it : you need to build the exact same pipeline as the mlab.surf() code but add a SetAttributeFilter for filtering onto a second array of data that will be used for the colormap (See Gael’s blog post for reference) :

# Create the data source
src = mlab.pipeline.array2d_source(x, y, z)

# Add the additional scalar information for the indicator
dataset = src.mlab_source.dataset
array_id = dataset.point_data.add_array(self.indicator.T.ravel())
dataset.point_data.get_array(array_id).name = ‘indicator’
dataset.point_data.update()

# Build the pipeline
warp = mlab.pipeline.warp_scalar(src)
normals = mlab.pipeline.poly_data_normals(warp)
active_attr = mlab.pipeline.set_active_attribute(normals,
                                            point_scalars=’indicator’)
surf = mlab.pipeline.surface(active_attr, vmin=0, vmax=1)

The pipeline now contains a SetActiveAttribute filter. Now you only need to update the datasource by adding/updating the array that you will use for the 4th dimension and then update the pipeline by selecting the array in the SetActiveAttribute filter. Something like the following snippet will work. It is the very same code used in the option valuation tool showed previously.

    def update_indicator(self):
        """
        Draw the selected indicator on top of the surface       

        FIXME : remove useless update_data() calls
        """
        # get the data to overlay
        self.indicator = self.indicator_fun(self.asset, self.strikes, self.d_yield, self.int_rate, self.sigma, self.expiry)		
		
        src = self.surface.mlab_source
	# remove previous array if existing
        src.dataset.point_data.remove_array('indicator')
	# add the dataset to the point_data of the surface
        array_id = src.dataset.point_data.add_array(self.indicator.T.ravel())
	# rename the new array 	
        src.dataset.point_data.get_array(array_id).name = 'indicator'    
        self.surface.update_data()
        # need to update the correct vtkobject, otherwise, we have an unstable behaviour		        
        filterobj = self.surface.parent.parent		
        filterobj.parent.update_data()
        filterobj.update_data()
        filterobj.update_pipeline() # must be called to have the list of available point_scalars updated
        filterobj.point_scalars_name = 'indicator'
        filterobj.update_data()

As soon as the pipeline is properly set, the updates are done smoothly and the user interaction runs fine.

The third wayd to do it is to use a mesh in place of a surface. The mlab.mesh command allows you to do exactly that but with the drawback of needing complete 3D array for the 4 inputs and thus has a cost in memory consumption.

About these ads

One Response to Surface plots using a 4th parameter for colormap in Mayavi

  1. [...] an additional information on a surface with Mayavi, using color not given eg by the elevation. A recent post on his blog by Didrik Pinte shows the problem quite [...]

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s

Follow

Get every new post delivered to your Inbox.

%d bloggers like this: