def plotOrs(*ptss):

    mlab.figure(34); mlab.clf()

    r = 1
    phi, theta = mgrid[0:pi:101j, 0:2*pi:101j]
    
    x = r*sin(phi)*cos(theta)
    y = r*sin(phi)*sin(theta)
    z = r*cos(phi)

    
    mlab.mesh(x,y,z, colormap='gray',opacity=.2)
    
    
    
    colors = [(1,0,0),(0,1,0),(0,0,1)]
    print len(colors)
    print len(ptss)
    for (pts,col) in izip(ptss,colors):
        print col        
        ors = normr(pts)
        # Create a sphere
        
        x,y,z = ors.T
        
        mlab.points3d(x,y,z,color=col)
        mlab.plot3d(x,y,z,color=col,tube_radius=.025)

def _plot_subnetwork_graph(Theta, coords, network_ix):
    (x, y, z) = coords
    p = Theta.shape[0]
    network_ix = network_ix * np.ones((p,))
    # 3D glass image of brain
    fig = mlab.figure(bgcolor=(1, 1, 1), size=(900, 769))
    mlab.clf()
    fig.scene.disable_render = True
    vmin = np.min(Theta[np.abs(Theta) != 0])
    vmax = np.max(np.abs(Theta))
    tubes, nodes = plot_graph(-Theta, x, y, z,
                              node_size=.6,
                              edge_vmin=vmin,
                              edge_vmax=vmax,
                              node_colormap='spectral',
                              node_color=(0.2, 0.2, 0.2),
                              node_scalar=network_ix,
                              tube_radius=.15)
    lut = tubes.module_manager.scalar_lut_manager.lut.table.to_array()
    lut = 255 * plt.cm.hot_r(np.linspace(0, 1, 256))
    tubes.module_manager.scalar_lut_manager.lut.table = lut
    tubes.update_pipeline()
    #nodes.module_manager.scalar_lut_manager.lut.table = color
    nodes.update_pipeline()

    viz3d.plot_anat_3d(outline_color=(0, 0, 0), gyri_opacity=0.15)
    fig.scene.disable_render = False
    return fig

    def update_plot(self,s):   
	mlab.clf()
        self.s = s
	self.s1 = self.vessel()
        self.src = mlab.pipeline.scalar_field(self.s, 
                                         scaling=(1, 1, 1), 
                                         origin=(1,1,1))
        
        self.src.spacing = [ 0.75,1, 1]
        
        #iso surface
	if ((self.n_min == 1) and (self.n_max== 1)):
	    mincon = 0
	    self.m=mlab.pipeline.iso_surface(self.src, vmin=0, vmax=1, contours=200, opacity=0.3) 
	else:
	    mincon = self.n_min/100
	    self.m=mlab.pipeline.iso_surface(self.src, vmin=self.n_min/100, vmax=self.n_max/100, contours=200, opacity=0.3)
	
	vol=s
	vol[nonzero(vol>=mincon)]=1
	vol[nonzero(vol<mincon)]=0
	volvox=sum(vol)*0.0156
	volvox = unicode(volvox)
	print sum(vol)
	print "volume objek = ", volvox, "cm3"
	
        self.m.module_manager.scalar_lut_manager.show_scalar_bar = True
	self.src1 = mlab.pipeline.scalar_field(self.s1,scaling=(1, 1, 1),origin=(1,1,1))
	self.src1.spacing = [0.75,1 , 1] 
	self.m1=mlab.pipeline.iso_surface(self.src1,contours=[self.s1.min()+0.01*self.s1.ptp(), ],opacity=0.1)         
        self.mlab.colorbar(orientation='vertical', title="capasitance") 
	cancer_vol = "Ca. Volume: " + volvox + " cm3"
	#self.mlab.text(0.1,0.1,unicode(cancer_vol))
        self.mlab.show()     

def run_mlab_examples():
    from mayavi import mlab
    from mayavi.tools.animator import Animator

    ############################################################
    # run all the "test_foobar" functions in the mlab module.
    for name, func in getmembers(mlab):
        if not callable(func) or not name[:4] in ('test', 'Test'):
            continue

        if sys.platform == 'win32' and name == 'test_mesh_mask_custom_colors':
            # fixme: This test does not seem to work on win32, disabling for now.
            continue

        mlab.clf()
        GUI.process_events()
        obj = func()

        if isinstance(obj, Animator):
            obj.delay = 10
            # Close the animation window.
            obj.close()
            while is_timer_running(obj.timer):
                GUI.process_events()
                sleep(0.05)

        # Mayavi has become too fast: the operator cannot see if the
        # Test function was succesful.
        GUI.process_events()
        sleep(0.1)

def plotvfonsph3D(theta_rad, phi_rad, E_th, E_ph, freq=0.0,
                     vcoord='sph', projection='equirectangular'):
    PLOT3DTYPE = "quiver"
    (x, y, z) = sph2crtISO(theta_rad, phi_rad)
    from mayavi import mlab
    
    mlab.figure(1, bgcolor=(1, 1, 1), fgcolor=(0, 0, 0), size=(400, 300))
    mlab.clf()
    if PLOT3DTYPE == "MESH_RADIAL" :
        r_Et = numpy.abs(E_th)
        r_Etmx = numpy.amax(r_Et)
        mlab.mesh(r_Et*(x)-1*r_Etmx, r_Et*y, r_Et*z, scalars=r_Et)
        r_Ep = numpy.abs(E_ph)
        r_Epmx = numpy.amax(r_Ep)
        mlab.mesh(r_Ep*(x)+1*r_Epmx , r_Ep*y, r_Ep*z, scalars=r_Ep)
    elif PLOT3DTYPE == "quiver":
        ##Implement quiver plot
        s2cmat = getSph2CartTransfMatT(numpy.array([x,y,z]))
        E_r = numpy.zeros(E_th.shape)
        E_fldsph = numpy.rollaxis(numpy.array([E_r, E_ph, E_th]), 0, 3)[...,numpy.newaxis]
        E_fldcrt = numpy.rollaxis(numpy.matmul(s2cmat, E_fldsph).squeeze(), 2, 0)
        #print E_fldcrt.shape
        mlab.quiver3d(x+1.5, y, z,
                      numpy.real(E_fldcrt[0]),
                      numpy.real(E_fldcrt[1]),
                      numpy.real(E_fldcrt[2]))
        mlab.quiver3d(x-1.5, y, z,
                      numpy.imag(E_fldcrt[0]),
                      numpy.imag(E_fldcrt[1]),
                      numpy.imag(E_fldcrt[2]))              
    mlab.show()

def plot_sphere_func(f, grid='Clenshaw-Curtis', theta=None, phi=None, colormap='jet', fignum=0):

    # Note: all grids except Clenshaw-Curtis have holes at the poles

    import matplotlib
    matplotlib.use('WxAgg')
    matplotlib.interactive(True)
    from mayavi import mlab

    if grid == 'Driscoll-Healy':
        b = f.shape[0] / 2
    elif grid == 'Clenshaw-Curtis':
        b = (f.shape[0] - 2) / 2
    elif grid == 'SOFT':
        b = f.shape[0] / 2
    elif grid == 'Gauss-Legendre':
        b = (f.shape[0] - 2) / 2

    if theta is None or phi is None:
        theta, phi = meshgrid(b=b, convention=grid)

    phi = np.r_[phi, phi[0, :][None, :]]
    theta = np.r_[theta, theta[0, :][None, :]]
    f = np.r_[f, f[0, :][None, :]]

    x = np.sin(theta) * np.cos(phi)
    y = np.sin(theta) * np.sin(phi)
    z = np.cos(theta)

    mlab.figure(fignum, bgcolor=(1, 1, 1), fgcolor=(0, 0, 0), size=(600, 400))
    mlab.clf()
    mlab.mesh(x, y, z, scalars=f, colormap=colormap)

    # mlab.view(90, 70, 6.2, (-1.3, -2.9, 0.25))
    mlab.show()

def plot3D(name, X, Y, Z, zlabel):
    """
    Plots a 3d surface plot of Z using the mayavi mlab.mesh function.

    Parameters
    ----------
    name: string
        The name of the figure.
    X: 2d ndarray
        The x-axis data.
    Y: 2d ndarray
        The y-axis data.
    Z: 2d nd array
        The z-axis data.
    zlabel: The title that appears on the z-axis.
    """
    mlab.figure(name)
    mlab.clf()
    plotData = mlab.mesh(X/(np.max(X) - np.min(X)),
                         Y/(np.max(Y) - np.min(Y)),
                         Z/(np.max(Z) - np.min(Z)))
    mlab.outline(plotData)
    mlab.axes(plotData, ranges=[np.min(X), np.max(X),
                                np.min(Y), np.max(Y),
                                np.min(Z), np.max(Z)])
    mlab.xlabel('Space ($x$)')
    mlab.ylabel('Time ($t$)')
    mlab.zlabel(zlabel)

    def show_volume(self, bbox):
        """
        show the volume with the given bounding box
        """
        # load the data
        with h5py.File(self.h5fn, 'r') as f:
            seg1 = f["segmentation"]["labels"][bbox[2]:bbox[5],bbox[1]:bbox[4],bbox[0]:bbox[3]]
            raw = f["raw"]["volume"][bbox[2]:bbox[5],bbox[1]:bbox[4],bbox[0]:bbox[3]]
        
        print "Drawing volume"
        t0 = time.time()
        
        #draw everything
        fig1 = mlab.figure(1, size=(500,450))
        mlab.clf(fig1)
        visCell.drawImagePlane(fig1, raw, 'gist_ncar')
        visCell.drawVolumeWithoutReferenceCell(fig1, seg1, np.array((-1,)), (0,0,1),0.5)
        with h5py.File(self.h5fn, 'r') as f:
            visCell.drawLabels(fig1, f, seg1, bbox)

        fig2 = mlab.figure(2, size=(500,450))
        mlab.clf(fig2)

        visCell.draw2DView(fig2, raw[20:-20,:,:], seg1[20:-20,:,:], -1)
        
        t = time.time() - t0
        print "Time for drawing:",t

def process_launch():
    '''Procédure reliant une fenetre graphique et le coeur du programme'''
    global nb_etapesIV
    nb_etapes=nb_etapesIV.get()#On récupère le nombre d'étapes
    fig=mlab.figure(1)
    mlab.clf()#La fenêtre de dessin est initialisée
    mlab.draw(terrain([(0,1,2),(2,3,4),(4,5,6)],[(Point(0,0,0),Point(1,0,0)),(Point(1,0,0),Point(1,1,0)),(Point(0,0,0),Point(1,1,0)),(Point(1,1,0),Point(0,1,0)),(Point(0,0,0),Point(0,1,0)),(Point(0,0,0),Point(-1,1,0)),(Point(-1,1,0),Point(0,1,0))],nb_etapes))#On affiche le dessin

	def _plotbutton1_fired(self):
		mlab.clf()
		self.loaddata()
		self.sregion[np.where(self.sregion<self.datamin)] = self.datamin
		self.sregion[np.where(self.sregion>self.datamax)] = self.datamax

		# The following codes from: http://docs.enthought.com/mayavi/mayavi/auto/example_atomic_orbital.html#example-atomic-orbital
		field = mlab.pipeline.scalar_field(self.sregion)     # Generate a scalar field
		colored = self.sregion
		vol=self.sregion.shape
		for v in range(0,vol[2]-1):
			colored[:,:,v] = self.extent[4] + v*(-1)*abs(self.hdr['cdelt3'])
		new = field.image_data.point_data.add_array(colored.T.ravel())
		field.image_data.point_data.get_array(new).name = 'color'
		field.image_data.point_data.update()

		field2 = mlab.pipeline.set_active_attribute(field, point_scalars='scalar')
		contour = mlab.pipeline.contour(field2)
		contour2 = mlab.pipeline.set_active_attribute(contour, point_scalars='color')

		mlab.pipeline.surface(contour2, colormap='jet', opacity=self.opacity)
		
		## Insert a continuum plot
		##im = pyfits.open('g28_SMA1.cont.image.fits')
		##dat = im[0].data
		##dat0 = dat[0]
		##channel = dat0[0]
		##region = np.swapaxes(channel[self.xstart:self.xend,self.ystart:self.yend]*1000.,0,1)
		##field = mlab.contour3d(region, colormap='gist_ncar')
		##field.contour.minimum_contour = 5

		self.field = field
		self.labels()
		mlab.view(azimuth=0, elevation=0, distance='auto')
		mlab.show()

def bigtest():
    from jds_image_proc.clouds import voxel_downsample
    from jds_image_proc.pcd_io import load_xyz
    import mayavi.mlab as mlab
    
    pts = load_xyz("/home/joschu/Data/scp/three_objs_ds.pcd")
    #pts = voxel_downsample(xyz, .03, False)


    mlab.clf()
    mlab.points3d(pts[:,0], pts[:,1], pts[:,2], color = (1,1,1), scale_factor=.01)
    
    clus = []
    
    labels = decompose(pts, .025)
    for i in xrange(labels.max()+1):
        clu = np.flatnonzero(labels == i)
        clus.append(clu)

    for clu in sorted(clus, key=len, reverse=True):
        if len(clu) < 10: break
        dirs = ss.get_sphere_points(1)
        sup_pd = np.dot(pts[clu,:], dirs.T)
        best_d = sup_pd.max(axis=0)
        print "max deficit",(sup_pd - best_d[None,:]).max(axis=1).min()
        mlab.points3d(pts[clu,0], pts[clu,1], pts[clu,2], color = (rand(),rand(),rand()), scale_factor=.01)
        raw_input()

def mark_gt_box3d( lidar_dir, gt_boxes3d_dir, mark_dir):

    os.makedirs(mark_dir, exist_ok=True)
    fig   = mlab.figure(figure=None, bgcolor=(0,0,0), fgcolor=None, engine=None, size=(500, 500))
    dummy = np.zeros((10,10,3),dtype=np.uint8)

    for file in sorted(glob.glob(lidar_dir + '/*.npy')):
        name = os.path.basename(file).replace('.npy','')

        lidar_file   = lidar_dir     +'/'+name+'.npy'
        boxes3d_file = gt_boxes3d_dir+'/'+name+'.npy'
        lidar   = np.load(lidar_file)
        boxes3d = np.load(boxes3d_file)

        mlab.clf(fig)
        draw_didi_lidar(fig, lidar, is_grid=1, is_axis=1)
        if len(boxes3d)!=0:
            draw_didi_boxes3d(fig, boxes3d)

        azimuth,elevation,distance,focalpoint = MM_PER_VIEW1
        mlab.view(azimuth,elevation,distance,focalpoint)
        mlab.show(1)
        imshow('dummy',dummy)
        cv2.waitKey(1)

        mlab.savefig(mark_dir+'/'+name+'.png',figure=fig)

def init_mlab_scene(size):
	fig = mlab.figure('Viz', size=size, bgcolor=(0,0,0))
	fig.scene.set_size(size)
	fig.scene.anti_aliasing_frames = 0

	mlab.clf()
	return fig

def showDsweep():
    k_s_sweep = [10**(x-5) for x in range(10)]
    sl_div_diam_sweep = [5.0*(x+1)/1000 for x in range(20)]
    vol_ratio_sweep = [5.0*(x+1)/100 for x in range(10)]    
    
    Dmsd = numpy.load(data_dir + "/D_numpy.npy")
    kDa = 1.660538921e-30;
    mass = 40.0*kDa; 
    viscosity = 8.9e-4; 
    diameter = 5e-9; 
    T = 300.0; 
    Dbase = k_b*T/(3.0*numpy.pi*viscosity*diameter);
    Dmsd = Dmsd/Dbase
    
    mlab.figure(1, size=(800, 800), fgcolor=(1, 1, 1),
                                    bgcolor=(0.5, 0.5, 0.5))
    mlab.clf()
    contours = numpy.arange(0.01,2,0.2).tolist()
    obj = mlab.contour3d(Dmsd,contours=contours,transparent=True,vmin=contours[0],vmax=contours[-1])
    outline = mlab.outline(color=(.7, .7, .7),extent=(0,10,0,20,0,10))
    axes = mlab.axes(outline, color=(.7, .7, .7),
            nb_labels = 5,
            ranges=(k_s_sweep[0], k_s_sweep[-1], sl_div_diam_sweep[0], sl_div_diam_sweep[-1], vol_ratio_sweep[0], vol_ratio_sweep[-1]), 
            xlabel='spring stiffness', 
            ylabel='step length',
            zlabel='volume ratio')
    mlab.colorbar(obj,title='D',nb_labels=5)

    mlab.show()

def show_contrasts(subject, contrasts, side, threshold):

    x, y, z, triangles = get_geometry(subject, side, "inflated")   ## inflated or white
    curv = get_curvature_sign(subject, side)

    f = mlab.figure()
    mlab.clf()

    # anatomical mesh
    mlab.triangular_mesh(x, y, z, triangles, transparent=False,
                         opacity=1., name=subject,
        scalars=curv, colormap="bone", vmin=-1, vmax=2)

    mlab.title(subject)

    cmaps = [colormaps[c.split("-")[0]]['colormap'] for c in contrasts]

    for contrast, colormap in zip(contrasts, cmaps):
        # functional mesh
        data = get_contrast(subject, contrast, side)
        func_mesh = mlab.pipeline.triangular_mesh_source(x, y, z, triangles,
                                                     scalars=data)
            # threshold
        thresh = mlab.pipeline.threshold(func_mesh, low=threshold)

        surf = mlab.pipeline.surface(thresh, colormap='hot', transparent=True,
                          opacity=.8) # diminuer pour avoir plus de transparence
        lut = (np.array([colormap(v) for v in np.linspace(.25, 1., 256)]) * 255
                       ).astype(int)

        surf.module_manager.scalar_lut_manager.lut.table = lut

    mlab.draw()

    return f

 def _render_model(self):
   figure = mlab.gcf()
   mlab.clf()
   s = mlab.triangular_mesh(self.coords[:,0], self.coords[:,1],
                            self.coords[:,2], self.tri_index,
                            color=(0.5,0.5,0.5))
   return s.scene

    def action(u, x, xv, y, yv, t, n):
        #print 'action, t=',t,'\nu=',u, '\Nx=',x, '\Ny=', y
        if plot == 1:
            mesh(xv, yv, u, title='t=%g' %t[n])
            time.sleep(0.2) # pause between frames

        elif plot == 2:
            # mayavi plotting
            mlab.clf()
            extent1 = (0, 20, 0, 20,-2, 2)
            s = mlab.surf(x , y, u, colormap='Blues', warp_scale=5,extent=extent1)
            mlab.axes(s, color=(.7, .7, .7), extent=extent1,
            ranges=(0, 10, 0, 10, -1, 1), xlabel='', ylabel='',
            zlabel='',
            x_axis_visibility=False, z_axis_visibility=False)
            mlab.outline(s, color=(0.7, .7, .7), extent=extent1)
            mlab.text(2, -2.5, '', z=-4, width=0.1)
            mlab.colorbar(object=None, title=None, orientation='horizontal', nb_labels=None, nb_colors=None, label_fmt=None)
            mlab.title('test 1D t=%g' % t[n])

            mlab.view(142, -72, 50)
            f = mlab.gcf()
            camera = f.scene.camera
            camera.yaw(0)
        
        
        if plot > 0:
            path = 'Figures_wave2D'
            time.sleep(0) # pause between frames
            if save_plot and plot != 2:
                filename = '%s/%08d.png' % (path, n)
                savefig(filename)  # time consuming!
            elif save_plot and plot == 2:
                filename = '%s/%08d.png' % (path,n)
                mlab.savefig(filename)  # time consuming!

def zoncaview(m):
    """
    m is a healpix sky map, such as provided by WMAP or Planck.
    """

    nside = hp.npix2nside(len(m))
    vmin = -1e3; vmax = 1e3

    # Set up some grids:
    xsize = ysize = 1000
    theta = np.linspace(np.pi, 0, ysize)
    phi   = np.linspace(-np.pi, np.pi, xsize)
    longitude = np.radians(np.linspace(-180, 180, xsize))
    latitude = np.radians(np.linspace(-90, 90, ysize))

    # Project the map to a rectangular matrix xsize x ysize:
    PHI, THETA = np.meshgrid(phi, theta)
    grid_pix = hp.ang2pix(nside, THETA, PHI)
    grid_map = m[grid_pix]

    # Create a sphere:
    r = 0.3
    x = r*np.sin(THETA)*np.cos(PHI)
    y = r*np.sin(THETA)*np.sin(PHI)
    z = r*np.cos(THETA)

    # The figure:
    mlab.figure(1, bgcolor=(1, 1, 1), fgcolor=(0, 0, 0), size=(400, 300))
    mlab.clf()

    mlab.mesh(x, y, z, scalars=grid_map, colormap="jet", vmin=vmin, vmax=vmax)

    mlab.draw()

    return

def draw_earth():

    mlab.figure(1, bgcolor=(0.48, 0.48, 0.48), fgcolor=(0, 0, 0),
                   size=(400, 400))
    mlab.clf()

    # Display continents outline, using the VTK Builtin surface 'Earth'
    continents_src = BuiltinSurface(source='earth', name='Continents')
    # The on_ratio of the Earth source controls the level of detail of the
    # continents outline.
    continents_src.data_source.on_ratio = 2 
    continents = mlab.pipeline.surface(continents_src, color=(0, 0, 0))

    # Display a semi-transparent sphere, for the surface of the Earth
    # We use a sphere Glyph, throught the points3d mlab function, rather than
    # building the mesh ourselves, because it gives a better transparent
    # rendering.
    sphere = mlab.points3d(0, 0, 0, scale_mode='none', scale_factor=2, color=(0.67, 0.77, 0.93),
                           resolution=50, opacity=0.7, name='Earth')

    # These parameters, as well as the color, where tweaked through the GUI,
    # with the record mode to produce lines of code usable in a script.
    sphere.actor.property.specular = 0.45
    sphere.actor.property.specular_power = 5
    # Backface culling is necessary for better transparent rendering.
    sphere.actor.property.backface_culling = True

 def test_test_backend(self):
     """Test if setting the backend to 'test' works."""
     mlab.options.backend = 'test'
     mlab.test_contour3d()
     mlab.clf()
     mlab.pipeline.open(get_example_data('cube.vti'))
     mlab.clf()