def draw_artifact(artifact):
    draw_artifact_ports(artifact)

    center = artifact.pos
    mlab.points3d(center[0],center[1],center[2],color=(1,1,1))
    wire_box(center=artifact.pos)

    func_points = np.array([func.pos for func in artifact.subfunctions])
    mlab.points3d(func_points[:,0],
                  func_points[:,1],
                  func_points[:,2],
                  color=(0,0,1),scale_factor=0.4)

    for port in artifact.in_ports:
        fcn = port['parent']
        src = fcn.pos
        for output in fcn.outputs:
            snk = output['flow'].sink().pos
            print src,snk
            line(src,snk,color=(1,0,0))
    for port in artifact.out_ports:
        fcn = port['parent']
        snk = fcn.pos
        for input in fcn.inputs:
            src = input['flow'].source().pos
            line(src,snk,color=(0,1,0))

def plot_pos_ori(pos, ori, color=(0., 0., 0.)):
    mlab.points3d(pos[:, 0], pos[:, 1], pos[:, 2], scale_factor=0.005,
                  color=color)
    mlab.quiver3d(pos[:, 0], pos[:, 1], pos[:, 2],
                  ori[:, 0], ori[:, 1], ori[:, 2],
                  scale_factor=0.03,
                  color=color)

def plot_predicted_labels(points, labels):
    print '[plot_points] Plotting points!'
    xs = np.array([int(point._x) for point in points])
    ys = np.array([int(point._y) for point in points])
    zs = np.array([int(point._z) for point in points])
    mlab.points3d(xs, ys, zs, labels, scale_factor = .4, mode='cube')
    mlab.show()

 def _sphere_viz_default(self):
     # different between wx and qt
     try:
         color_tuple = self.sphere_color.toTuple()
     except:
         color_tuple = self.sphere_color
         
     try:
         pts = mlab.points3d(
             self.sphere_coords[:,0],
             self.sphere_coords[:,1],
             self.sphere_coords[:,2],
             mode='cube',
             scale_factor=1,
             figure = self.scene3d.mayavi_scene,
             color  = (color_tuple[0]/255.,
                      color_tuple[1]/255.,
                      color_tuple[2]/255.)
             )
     except:
         pts = mlab.points3d(
             self.sphere_coords[:,0],
             self.sphere_coords[:,1],
             self.sphere_coords[:,2],
             mode='cube',
             scale_factor=1,
             figure = self.scene3d.mayavi_scene,
             color  = (1.,0.,0.)
             )
     return pts

    def xmas_balls(connectivity, node_data=None, edge_data=False):
        """
        Plots coloured balls at the region centres of connectivity, colour and
        size is determined by a vector of length number of regions (node_data).
        
        Optional: adds the connections between pair of nodes.
        
        """
        centres = connectivity.centres
        edges = numpy.array(numpy.nonzero(connectivity.weights))
        edges = numpy.array([(start, stop) for (start, stop) in edges.T if start != stop])

        if node_data is not None:
            data_scale = 13.0 / node_data.max()
            pts = mlab.points3d(centres[:, 0], centres[:, 1], centres[:, 2],
                                node_data, transparent=True,
                                scale_factor=data_scale,
                                colormap='Blues')

            mlab.colorbar(orientation="vertical")
        else:
            #NOTE: the magic numbers are used to align region centers and surface representation. 
            #Do not ask ... 
            pts = mlab.points3d(centres[:, 0] * 1.13, centres[:, 1] * 1.13 + 15, centres[:, 2] - 25)

        if edge_data:
            pts.mlab_source.dataset.lines = edges
            tube = mlab.pipeline.tube(pts, tube_radius=0.5)
            mlab.pipeline.surface(tube, colormap='binary', opacity=0.142)

def mayaView(
             atoms,
             figure=None,
             colormap='hot',
             resolution=32,
             vmax=None,
             depth=3,
             ):
    from mayavi import mlab
    numbers = atoms.get_atomic_numbers()
    u_numbers = np.unique(numbers)
    if vmax is None:
        vmax = atoms.positions[:,2].max()
    vmin = vmax - depth
    my_colors = {1:(1,0,0), 6:(0.5,0.5,1)}

    for number in u_numbers:
        take = numbers == number
        element_atoms = atoms[take]
        points = element_atoms.positions
        points = points.transpose()
        radius = my_radii[number]*2
        mlab.points3d(points[0], points[1], points[2], points[2], scale_mode='none', scale_factor=radius, colormap=colormap, vmin=vmin, vmax=vmax,
                    resolution=resolution,
                    reset_zoom=False,
                    figure=figure,
                    )
    return figure

    def render(self, scale_factor=1.0, text_scale=1.0, **kwargs):
        import mayavi.mlab as mlab
        # disabling the rendering greatly speeds up this for loop
        self.figure.scene.disable_render = True
        positions = []
        for label in self.lmark_group:
            p = self.lmark_group[label]
            for i, p in enumerate(p.points):
                positions.append(p)
                l = '%s_%d' % (label, i)
                # TODO: This is due to a bug in mayavi that won't allow
                # rendering text to an empty figure
                mlab.points3d(p[0], p[1], p[2], scale_factor=scale_factor)
                mlab.text3d(p[0], p[1], p[2], l, figure=self.figure,
                            scale=text_scale)
        positions = np.array(positions)
        os = np.zeros_like(positions)
        os[:, 2] = 1
        mlab.quiver3d(positions[:, 0], positions[:, 1], positions[:, 2],
                      os[:, 0], os[:, 1], os[:, 2], figure=self.figure)
        self.figure.scene.disable_render = False

        # Ensure everything fits inside the camera viewport
        mlab.get_engine().current_scene.scene.reset_zoom()

        return self

def plot_static(color):
    global max_p, min_p, ary_q, ary_p
    # --- Set labels and ticks. ---
#    ax.set_xlabel(r'$k$', fontsize=16)
#    ax.set_ylabel(r'$\theta$', fontsize=16)
#    ax.set_zlabel(r'$p_{\theta}$', fontsize=16)
#    plt.title('Standard map')
#    unit   = 0.5
#    q_tick = np.arange(0, 2.0+unit, unit)
#    q_label = [r"$0$",r"$\frac{1}{2}\pi$",r"$\pi$",r"$\frac{3}{2}\pi$",r"$2\pi$"]
#    ax.set_ylim(0, 2*pi)
#    ax.set_yticks(q_tick*pi)
#    ax.set_yticklabels(q_label, fontsize=18)
#    ax.set_xlim(0, kmax)
#    if max(ary_p[ntransient:]) > max_p:
#        max_p = max(ary_p[ntransient:])
#    if min(ary_p[ntransient:]) < min_p:
#        min_p = min(ary_p[ntransient:])
#    margin =  (max_p-min_p)*0.05
#    ax.set_zlim(min_p-margin, max_p+margin)

    # --- Plot the standard map. ---
    mlab.points3d([0.1, ] * (nmax - ntransient + 1),
                  ary_q[ntransient:], ary_p[ntransient:])

    plt.show()

      def generar(self, nivel, cursor, distancia, dim):

          if nivel == 0:
                  mlab.points3d([cursor.x, cursor.x+3, cursor.x, cursor.x+3],[cursor.y, cursor.y, cursor.y -3,cursor.y -3],
                                [cursor.z, cursor.z, cursor.z -3,cursor.z -3],colormap="copper", scale_factor=.25)
                                
                  cursor.x = cursor.x + distancia*np.cos(cursor.ang)
                  cursor.y = cursor.y + distancia*np.sin(cursor.ang)
                  cursor.z = cursor.z + distancia*np.sin(cursor.ang)
                  
                  #a = [(cursor.x, cursor.x, cursor.x)] #para trabajar con 1D
                  #b = (cursor.x + distancia*np.cos(cursor.ang),cursor.y + distancia*np.sin(cursor.ang))
                  #mlab.points3d(a, a, a, colormap="copper", scale_factor=.25)
                  
                  #cursor.x = b[0]
                  #cursor.y = b[1]
          else:
              print "nivel = ", nivel
              for i in self.gramatica[self.inicial]:
                  if i == 'A':
                      self.generar(nivel-1, cursor, distancia/3.0, dim)
                  elif i == 'B':
                      cursor.x += 3*distancia
                      
                      print "la x de cursor tiene: ",cursor.x

    def drawgranules(gids, gcolor, *arg):
      
      from params import granule_diam as diam
      x = []
      y = []
      z = []
      if self.colors:
        s = []
        
      if self.colors and not not_in_weights:
          gids = gids.intersection(self.colors.keys())
          
      for gid in gids:
        if self.projected:
          if scale_mode == 'scalar':
            u, p = gpo(gid)[-2:]

            if self.colors:
              try:
                s.append(self.colors[gid] / 20 * 20 * 0.9 + 10)
              except KeyError:
                continue

              dep = -(100 - s[-1]) / 20 * params.grid_dim # depth from colors
              p = [ dep * u[0] + p[0], dep * u[1] + p[1], dep * u[2] + p[2]]
              diam = 0.42
          elif scale_mode == 'none':
            p = gpo(gid)[2]
        else:
          p = gpo(gid)[0]
      
        if self.colors and not (self.projected and scale_mode == 'scalar'):
          
          try:
            color = self.colors[gid]
          except KeyError:
            color = 0
            
          if gthreshold:
            if color >= gthreshold:
              s.append(1)
            else:
              s.append(0)
          else:
            s.append(color)
          
        
        r = params.ranstream(gid, 0)
        r.uniform(-params.grid_dim * .5, params.grid_dim * .5)          
        x.append(p[0] + r.repick())
        y.append(p[1] + r.repick())
        z.append(p[2] + r.repick())

      #print 'drawn mitral:',len(x)

      if self.colors:
        if self.vmin != None and self.vmax != None:
          return points3d(x, y, z, s, scale_factor=diam, vmin=self.vmin, vmax=self.vmax, scale_mode=scale_mode, colormap=gpalette)
        return points3d(x, y, z, s, scale_factor=diam, scale_mode=scale_mode,colormap=gpalette)
      return points3d(x, y, z, scale_factor=diam, color=gcolor,colormap=gpalette)

 def plot_intersections(self,):
     '''Plot the intersection points on the plasma surface
     SRH: 12July2013
     '''
     from mayavi import mlab
     mlab.points3d(self.intersection1[self.valid_channels,0], self.intersection1[self.valid_channels,1], self.intersection1[self.valid_channels,2],scale_factor=0.02, color=(0,0,1))
     mlab.points3d(self.intersection2[self.valid_channels,0], self.intersection2[self.valid_channels,1], self.intersection2[self.valid_channels,2],scale_factor=0.02,color=(1,0,0))

def plot_matrix(connectmat_file, centers_file, threshold_pct=5, weight_edges=False,
                node_scale_factor=2, edge_radius=.5, resolution=8, name_scale_factor=1,
                names_file=None, node_indiv_colors=[], highlight_nodes=[], fliplr=False):
    """
    Given a connectivity matrix and a (x,y,z) centers file for each region, plot the 3D network
    """
    matrix = core.file_reader(connectmat_file)
    nodes = core.file_reader(centers_file)
    if names_file:
        names = core.file_reader(names_file,1)
    num_nodes = len(nodes)
    edge_thresh_pct = threshold_pct / 100.0
    matrix_flat = np.array(matrix).flatten()
    edge_thresh = np.sort(matrix_flat)[len(matrix_flat)-int(len(matrix_flat)*edge_thresh_pct)]
    
    matrix = core.file_reader(connectmat_file)
    ma = np.array(matrix)

    thresh = scipy.stats.scoreatpercentile(ma.ravel(),100-threshold_pct)
    ma_thresh = ma*(ma > thresh)
    
    if highlight_nodes:
        nr = ma.shape[0]
        subset_mat = np.zeros((nr, nr))
        for i in highlight_nodes:
            subset_mat[i,:] = 1
            subset_mat[:,i] = 1
        ma_thresh = ma_thresh * subset_mat
        
    if fliplr:
        new_nodes = []
        for node in nodes:
            new_nodes.append([45-node[0],node[1],node[2]]) # HACK
        nodes = new_nodes
    
    mlab.figure(bgcolor=(1, 1, 1), size=(400, 400))
    for count,(x,y,z) in enumerate(nodes):
        if node_indiv_colors:
            mlab.points3d(x,y,z, color=colors[node_indiv_colors[count]], scale_factor=node_scale_factor, resolution=resolution)
        else:
            mlab.points3d(x,y,z, color=(0,1,0), scale_factor=node_scale_factor, resolution=resolution)
        if names_file:
            width = .025*name_scale_factor*len(names[count])
            print width
            print names[count]
            mlab.text(x, y,names[count], z=z,width=.025*len(names[count]),color=(0,0,0))
    for i in range(num_nodes-1):    
        x0,y0,z0 = nodes[i]
        for j in range(i+1, num_nodes):
            #if matrix[i][j] > edge_thresh:
            if ma_thresh[i][j] > edge_thresh:
                x1,y1,z1 = nodes[j]
                if weight_edges:
                    mlab.plot3d([x0,x1], [y0,y1], [z0,z1],
                            tube_radius=matrix[i][j]/matrix_flat.max(),
                            color=(1,1,1))
                else:
                    mlab.plot3d([x0,x1], [y0,y1], [z0,z1],
                            tube_radius=edge_radius,
                            color=(1,1,1))

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 testVisually(self):
        '''blocks selected visually.'''
#        if self.comm.rank == 0:
        g2z,zvec = PETSc.Scatter().toZero(self.bnd.gindBlockWBand)
        g2z.scatter(self.bnd.gindBlockWBand,zvec, PETSc.InsertMode.INSERT)
        x = self.bnd.BlockSub2CenterCarWithoutBand(\
                                                   self.bnd.BlockInd2SubWithoutBand(zvec.getArray()) )
        lx = self.bnd.BlockSub2CenterCarWithoutBand(\
                                                    self.bnd.BlockInd2SubWithoutBand(self.bnd.gindBlockWBand.getArray()))
        try:
            try:
                from mayavi import mlab
            except ImportError:
                from enthought.mayavi import mlab

            if self.comm.rank == 0:
                mlab.figure()
                mlab.points3d(x[:,0],x[:,1],x[:,2])
            mlab.figure()
            mlab.points3d(lx[:,0],lx[:,1],lx[:,2])
            mlab.show()
            #fig.add(pts1)
            #fig.add(pts2)
        except ImportError:
            import pylab as pl
            from mpl_toolkits.mplot3d import Axes3D #@UnusedImport
            fig = pl.figure()
            ax = fig.add_subplot(111, projection='3d')
            ax.scatter3D(x[:,0],x[:,1],x[:,2],c='blue',marker='o')
            ax.scatter3D(lx[:,0],lx[:,1],lx[:,2],c='red',marker='D')
            pl.savefig('testVis{0}.png'.format(self.comm.rank))
            pl.show()

def plot(ply):
    '''
    Plot vertices and triangles from a PlyData instance. Assumptions:
        `ply' has a 'vertex' element with 'x', 'y', and 'z'
            properties;

        `ply' has a 'face' element with an integral list property
            'vertex_indices', all of whose elements have length 3.

    '''
    vertex = ply['vertex']

    (x, y, z) = (vertex[t] for t in ('x', 'y', 'z'))

    mlab.points3d(x, y, z, color=(1, 1, 1), mode='point')

    if 'face' in ply:
        tri_idx = ply['face']['vertex_indices']
        idx_dtype = tri_idx[0].dtype

        triangles = numpy.fromiter(tri_idx, [('data', idx_dtype, (3,))],
                                   count=len(tri_idx))['data']

        mlab.triangular_mesh(x, y, z, triangles,
                             color=(1, 0, 0.4), opacity=0.5)

 def plot_sr(self, mlab, lc=0, sr=0):
     lc = self.lc_arr[lc]
     gd = self.geo_data_dict
     mlab.points3d(gd['X'], gd['Y'], gd['Z'], lc[:, sr],
                   # colormap = "YlOrBr",
                   mode="cube",
                   scale_factor=0.1)

def make_figures(coor, fun, interp_results, error):
    # Figure of harmoinc function on sphere in fine cordinates
    # Points3d showing interpolation training points coloured to their value
    mlab.figure()
    vmax, vmin = np.max(fun.fine), np.min(fun.fine)
    mlab.mesh(coor.fine.x, coor.fine.y, coor.fine.z,
              scalars=fun.fine, vmax=vmax, vmin=vmin)
    mlab.points3d(coor.coarse.x, coor.coarse.y, coor.coarse.z, fun.coarse,
                  scale_factor=0.1, scale_mode='none', vmax=vmax, vmin=vmin)
    mlab.colorbar(title='Spherical Harmonic', orientation='vertical')
    # mlab.savefig('interppointssphere.png')

    # Figure showing results of rbf interpolation
    mlab.figure()
    mlab.mesh(coor.fine.x, coor.fine.y, coor.fine.z,
              scalars=interp_results, vmax=vmax, vmin=vmin)
    mlab.points3d(coor.coarse.x, coor.coarse.y, coor.coarse.z, fun.coarse,
                  scale_factor=0.1, scale_mode='none', vmax=vmax, vmin=vmin)
    mlab.colorbar(title='Interpolation', orientation='vertical')
    # mlab.savefig('interpsphere.png')

    mlab.figure()
    mlab.mesh(coor.fine.x, coor.fine.y, coor.fine.z,
              scalars=error.errors, vmax=error.max, vmin=-error.max)
    mlab.colorbar(title='Error', orientation='vertical')
    # mlab.points3d(coor.coarse.x, coor.coarse.y, coor.coarse.z, scalars, scale_factor=0.1, scale_mode='none',vmax=vmax, vmin=vmin)
    # mlab.savefig('interpsphere.png')

    mlab.show()

    def show3d(self, style={'A':1}):
        from mayavi import mlab
        scalar = []
        pos = []
        def rep(lst):
            for (x,y,z) in lst:
                yield x,y,z
                if x == 0:
                    yield 1,y,z
                    if y == 0:
                        yield 1,1,z
                if y == 0:
                    yield x,1,z
                    if z == 0:
                        yield x,1,1
                if z == 0:
                    yield x,y,1
                    if x == 0:
                        yield 1,y,1
                if x == y == z == 0:
                    yield 1,1,1

        for atom, site in self.sites:
            site = list(rep(site))
            scalar.extend([style[atom]]*len(site))
            pos.extend(self.v(site))
        scalar = np.asarray(scalar)
        pos = np.asarray(pos)


        mlab.points3d(pos[:,0],pos[:,1],pos[:,2],scalar,scale_factor=1, resolution=10)
        mlab.show()

def plot_element_values(n, nodes, values, resample_n=None,
        node_tuples=None, show_nodes=False):
    dims = len(nodes)

    orig_nodes = nodes
    orig_values = values

    if resample_n is not None:
        import modepy as mp
        basis = mp.simplex_onb(dims, n)
        fine_nodes = mp.equidistant_nodes(dims, resample_n)

        values = np.dot(mp.resampling_matrix(basis, fine_nodes, nodes), values)
        nodes = fine_nodes
        n = resample_n

    from pytools import generate_nonnegative_integer_tuples_summing_to_at_most \
            as gnitstam

    if dims == 1:
        import matplotlib.pyplot as pt
        pt.plot(nodes[0], values)
        if show_nodes:
            pt.plot(orig_nodes[0], orig_values, "x")
        pt.show()
    elif dims == 2:
        import mayavi.mlab as mlab
        mlab.triangular_mesh(
                nodes[0], nodes[1], values, submesh(list(gnitstam(n, 2))))
        if show_nodes:
            mlab.points3d(orig_nodes[0], orig_nodes[1], orig_values,
                    scale_factor=0.05)
        mlab.show()
    else:
        raise RuntimeError("unsupported dimensionality %d" % dims)

def make_figures(coor, fun, interp_results, error):
    '''Produce MayaVi figures for interpolation results'''
    mlab.figure()
    vmax, vmin = np.max(fun.fine), np.min(fun.fine)
    mlab.mesh(coor.fine.x, coor.fine.y, coor.fine.z,
              scalars=fun.fine, vmax=vmax, vmin=vmin)
    mlab.points3d(coor.coarse.x, coor.coarse.y, coor.coarse.z, fun.coarse,
                  scale_factor=0.1, scale_mode='none', vmax=vmax, vmin=vmin)
    mlab.colorbar(title='Spherical Harmonic', orientation='vertical')
    mlab.savefig('Figures/poorfunctionsphere.png')

    # Figure showing results of rbf interpolation
    mlab.figure()
    mlab.mesh(coor.fine.x, coor.fine.y, coor.fine.z,
              scalars=interp_results, vmax=vmax, vmin=vmin)
    mlab.points3d(coor.coarse.x, coor.coarse.y, coor.coarse.z, fun.coarse,
                  scale_factor=0.1, scale_mode='none', vmax=vmax, vmin=vmin)
    mlab.colorbar(title='Interpolation', orientation='vertical')
    mlab.savefig('Figures/poorinterpsphere.png')

    mlab.figure()
    mlab.mesh(coor.fine.x, coor.fine.y, coor.fine.z,
              scalars=error.errors, vmax=error.max, vmin=-error.max)
    mlab.colorbar(title='Error', orientation='vertical')
    # mlab.points3d(coor.coarse.x, coor.coarse.y, coor.coarse.z, scalars, scale_factor=0.1, scale_mode='none',vmax=vmax, vmin=vmin)
    mlab.savefig('Figures/poorerrorsphere.png')

    mlab.show()