def plot_3D( self ):
     x = self.compute3D_plot[0]
     y = self.compute3D_plot[1]
     z = self.compute3D_plot[2]
     # print x_axis, y_axis, z_axis
     if self.autowarp_bool:
         x = x / x[-1]
         y = y / y[-1]
         z = z / z[-1] * self.z_scale
     mlab.surf( x, y , z, representation = 'wireframe' )
     engine = Engine()
     engine.start()
     if len( engine.scenes ) == 75.5:
         engine.new_scene()
         surface = engine.scenes[0].children[0].children[0].children[0].children[0].children[0]
         surface.actor.mapper.scalar_range = np.array( [ 6.97602671, 8.8533387 ] )
         surface.actor.mapper.scalar_visibility = False
         scene = engine.scenes[0]
         scene.scene.background = ( 1.0, 1.0, 1.0 )
         surface.actor.property.specular_color = ( 0.0, 0.0, 0.0 )
         surface.actor.property.diffuse_color = ( 0.0, 0.0, 0.0 )
         surface.actor.property.ambient_color = ( 0.0, 0.0, 0.0 )
         surface.actor.property.color = ( 0.0, 0.0, 0.0 )
         surface.actor.property.line_width = 1.
         scene.scene.isometric_view()
         mlab.xlabel( self.x_name3D )
         mlab.ylabel( self.y_name3D )
     mlab.outline()
     mlab.show()

    def showPath(self):
        px, py, pz, dirx,diry,dirz = [],[],[],[],[],[]

        print self.path
        print len(self.z)

        for node in self.path:
            px.append(node[0])
            py.append(node[1])
            pz.append(self.z[node[0]][node[1]])
            dirx.append(cos(radians(node[2])))
            diry.append(sin(radians(node[2])))
            dirz.append(0)

        px = np.array(px)
        py = np.array(py)
        pz = np.array(pz)
        dirx = np.array(dirx)
        diry = np.array(diry)
        dirz = np.array(dirz)
        mlab.quiver3d(self.x,self.y,self.z,self.gx,self.gy,self.Gt,color=(1,0,0),scale_factor=1)
        mlab.quiver3d(px,py,pz,dirx,diry,dirz,color=(0,0,0),scale_factor=1)

        mlab.surf(self.x, self.y, self.z, representation='wireframe')
        mlab.show()

 def render(self, colour=(1, 0, 0), line_width=2, step=None,
            marker_style='2darrow', marker_resolution=8, marker_size=0.05,
            alpha=1.0):
     from mayavi import mlab
     warp_scale = kwargs.get('warp_scale', 'auto')
     mlab.surf(self.values, warp_scale=warp_scale)
     return self

 def plot_mayavi(self):
     """Use mayavi to plot a phenotype phase plane in 3D.
     The resulting figure will be quick to interact with in real time,
     but might be difficult to save as a vector figure.
     returns: mlab figure object"""
     from mayavi import mlab
     figure = mlab.figure(bgcolor=(1, 1, 1), fgcolor=(0, 0, 0))
     figure.name = "Phenotype Phase Plane"
     max = 10.0
     xmax = self.reaction1_fluxes.max()
     ymax = self.reaction2_fluxes.max()
     zmax = self.growth_rates.max()
     xgrid, ygrid = meshgrid(self.reaction1_fluxes, self.reaction2_fluxes)
     xgrid = xgrid.transpose()
     ygrid = ygrid.transpose()
     xscale = max / xmax
     yscale = max / ymax
     zscale = max / zmax
     mlab.surf(xgrid * xscale, ygrid * yscale, self.growth_rates * zscale,
               representation="wireframe", color=(0, 0, 0), figure=figure)
     mlab.mesh(xgrid * xscale, ygrid * yscale, self.growth_rates * zscale,
               scalars=self.shadow_prices1 + self.shadow_prices2,
               resolution=1, representation="surface", opacity=0.75,
               figure=figure)
     # draw axes
     mlab.outline(extent=(0, max, 0, max, 0, max))
     mlab.axes(opacity=0, ranges=[0, xmax, 0, ymax, 0, zmax])
     mlab.xlabel(self.reaction1_name)
     mlab.ylabel(self.reaction2_name)
     mlab.zlabel("Growth rates")
     return figure

def main():

    parser = argparse.ArgumentParser()
    parser.add_argument("point_x", help="x coordinate of point.", type=int)
    parser.add_argument("point_y", help="y coordinate of point.", type=int)
    parser.add_argument("file", help="The bathymetry file.")
    parser.add_argument("--halo", help="Size of halo.", type=int, default=50)
    args = parser.parse_args()

    f = nc.Dataset(args.file)
    topo = f.variables['depth'][:]

    # Calculate the extents of the topo array to show. We don't just add halo to (point_x, point_y)
    # because it may run off the edge of the map.
    north_ext = min(args.point_y + args.halo, topo.shape[0])
    south_ext = max(args.point_y - args.halo, 0)
    east_ext = min(args.point_x + args.halo, topo.shape[1])
    west_ext = max(args.point_x - args.halo, 0)

    width = east_ext - west_ext
    height = north_ext - south_ext

    # The origin of the figure in global coordinates.
    origin_x = west_ext + width / 2
    origin_y = south_ext + height / 2

    point_y = args.point_y - origin_y
    point_x = args.point_x - origin_x

    mlab.surf(topo[south_ext:north_ext, west_ext:east_ext], warp_scale=0.005)
    mlab.points3d([point_y], [point_x], [0], color=(1, 0, 0), scale_factor=1.0)
    mlab.show()

def k_COV_plots(k_arr, COV_arr):
    sig_max_arr = np.zeros((len(k_arr), len(COV_arr)))
    wmax_arr = np.zeros((len(k_arr), len(COV_arr)))
    mu_r, mu_tau = 0.01, 0.1
    for i, k in enumerate(k_arr):
        for j, cov in enumerate(COV_arr):
            Vf = Vf_k(k)
            loc_r = mu_r * (1 - cov * np.sqrt(3.0))
            scale_r = cov * 2 * np.sqrt(3.0) * mu_r
            loc_tau = mu_tau * (1 - cov * np.sqrt(3.0))
            scale_tau = cov * 2 * np.sqrt(3.0) * mu_tau
            reinf = ContinuousFibers(r=RV('uniform', loc=loc_r, scale=scale_r),
                                  tau=RV('uniform', loc=loc_tau, scale=scale_tau),
                                  xi=WeibullFibers(shape=7.0, sV0=0.003),
                                  V_f=Vf, E_f=200e3, n_int=100)
            ccb_view.model.reinforcement_lst = [reinf]
            sig_max, wmax = ccb_view.sigma_c_max
            sig_max_arr[i, j] = sig_max/Vf
            wmax_arr[i, j] = wmax
    ctrl_vars = orthogonalize([np.arange(len(k_arr)), np.arange(len(COV_arr))])
    print sig_max_arr
    print wmax_arr
    mlab.surf(ctrl_vars[0], ctrl_vars[1], sig_max_arr / np.max(sig_max_arr))
    mlab.surf(ctrl_vars[0], ctrl_vars[1], wmax_arr / np.max(wmax_arr))
    mlab.show()

 def draw_working_area(self, width, height):
     color = (0.9, 0.9, 0.7)
     x, y = np.mgrid[0:height+1:10, 0:width+1:10]
     z = np.zeros(x.shape)
     mlab.surf(x, y, z, color=color, line_width=1.0,
               representation='wireframe', warp_scale=self.scale)
     mlab.axes()

def ripple():
    """Plot z = sin(10(x^2 + y^2))/10 on [-1,1]x[-1,1] using Mayavi."""
    X, Y = np.mgrid[-1:1:.01, -1:1:0.01]
    Z = np.sin(10*(X**2+Y**2))/10
    mlab.surf(X,Y,Z, colormap='RdYlGn') 
    mlab.show()
    pass

 def plotSurf(self,det,axisX,axisY,start=None,stop=None):
     
     bin = self.detector[det]
     try:
         bin.data
     except AttributeError:
         self.arrayRead(det)
     sumaxis = [2,1,0]
     sumaxis.remove(axisX)
     sumaxis.remove(axisY)
     low,high = self.__ranges(det,start,stop)
     subdata = (slice(low[0],high[0]),slice(low[1],high[1]),slice(low[2],high[2]))
     plotdata = bin.data[subdata]
     plotdata = sum(plotdata,sumaxis[0])
     #Renormalize data
     if (bin.type in [1,7,11,17]): #Cylinder bins
         plotdata = plotdata/((pi*high[0]**2*high[2]*high[1]/(2*pi))-(pi*low[0]**2*low[2]*low[1]/(2*pi)))
     elif (bin.type in [0,3,4,5,10,13,14,15,16]): #Cartesian bins 
         plotdata = plotdata/((high[0]-low[0])*(high[1]-low[1])*(high[1]-low[1]))
         
     x = linspace(bin.start[axisX]+bin.delta[axisX]/2,bin.stop[axisX]+bin.delta[axisX]/2, num=bin.num[axisX], endpoint=False)[slice(low[axisX],high[axisX])]
     y = linspace(bin.start[axisY]+bin.delta[axisY]/2,bin.stop[axisY]+bin.delta[axisY]/2, num=bin.num[axisY], endpoint=False)[slice(low[axisY],high[axisY])]
     mlab.figure(1, fgcolor=(0, 0, 0), bgcolor=(1, 1, 1))
     z =plotdata
     # Visualize the points
     mlab.surf(x,y,z)
     mlab.show()

    def plot_3D_spectrum(self, xmin=None, xmax=None, xN=None, ymin=None,
                         ymax=None, yN=None, trajectory=False, tube_radius=1e-2,
                         part='imag'):
        """Plot the Riemann sheet structure around the EP.

            Parameters:
            -----------
                xmin, xmax: float
                    Dimensions in x-direction.
                ymin, ymax: float
                    Dimensions in y-direction.
                xN, yN: int
                    Number of sampling points in x and y direction.
                trajectory: bool
                    Whether to include a projected trajectory of the eigenbasis
                    coefficients.
                part: str
                    Which function to apply to the eigenvalues before plotting.
                tube_radius: float
                    Trajectory tube thickness.
        """
        from mayavi import mlab

        X, Y, Z = self.sample_H(xmin, xmax, xN, ymin, ymax, yN)
        Z0, Z1 = [Z[..., n] for n in (0, 1)]

        def get_min_and_max(*args):
            data = np.concatenate(*args)
            return data.min(), data.max()

        surf_kwargs = dict(colormap='Spectral', mask=np.diff(Z0.real) > 0.015)

        mlab.figure(0)
        Z_min, Z_max = get_min_and_max([Z0.real, Z1.real])
        mlab.surf(X.real, Y.real, Z0.real, vmin=Z_min, vmax=Z_max, **surf_kwargs)
        mlab.surf(X.real, Y.real, Z1.real, vmin=Z_min, vmax=Z_max, **surf_kwargs)
        mlab.axes(zlabel="Re(E)")

        mlab.figure(1)
        Z_min, Z_max = get_min_and_max([Z0.imag, Z1.imag])
        mlab.mesh(X.real, Y.real, Z0.imag, vmin=Z_min, vmax=Z_max, **surf_kwargs)
        mlab.mesh(X.real, Y.real, Z1.imag, vmin=Z_min, vmax=Z_max, **surf_kwargs)
        mlab.axes(zlabel="Im(E)")

        if trajectory:
            x, y = self.get_cycle_parameters(self.t)
            _, c1, c2 = self.solve_ODE()

            for i, part in enumerate([np.real, np.imag]):
                e1, e2 = [part(self.eVals[:, n]) for n in (0, 1)]
                z = map_trajectory(c1, c2, e1, e2)
                mlab.figure(i)
                mlab.plot3d(x, y, z, tube_radius=tube_radius)
                mlab.points3d(x[0], y[0], z[0],
                              # color=line_color,
                              scale_factor=1e-1,
                              mode='sphere')

        mlab.show()

    def mlab3Dview():

        e_arr = orthogonalize( [x, w] )
        n_e_arr = [ e / np.max( np.fabs( e ) ) for e in e_arr ]

        n_mu_q_arr = mu_q_arr / np.max( np.fabs( mu_q_arr ) )
        m.surf( n_e_arr[0], n_e_arr[1], n_mu_q_arr )
        m.show()

def plot_3d_mayavi(X, Y, Z):
    f = figure(bgcolor=(1, 1, 1))

    surf(X.T, Y.T, Z, figure=f, warp_scale='auto')
    axes(xlabel='N Samples', ylabel='Sample', zlabel='Gradient',
         z_axis_visibility=False, nb_labels=10, line_width=1)

    show()

def problem8():
    opener = urllib.URLopener()
    opener.retrieve('https://s3.amazonaws.com/storage.enthought.com/www/sample_data/N36W113.hgt.zip', 'N36W113.hgt.zip')
    data = np.fromstring(zipfile.ZipFile('N36W113.hgt.zip').read('N36W113.hgt'), '>i2').reshape((3601, 3601)).astype('float32')
    data = data[:1000, 900:1900]
    data[data == -32768] = data[data>0].min()
    mlab.figure(size=(400, 320), bgcolor = (.16, .28, .46))
    mlab.surf(data, colormap="gist_earth", warp_scale=.2, vmin=1200, vmax=1610)
    mlab.view(-5.9, 83, 570, [5.3, 20, 238])
    return mlab.gcf()

    def __init__(self, hex_list, np_file='/tmp/magnetic_ground_truth.np', robot_height=40, width=800, height=600,
                 start_point=(0, 0, 0), message='experiment default message...'):
        self.debug = False
        self.animator = None
        self.movement_mode = 0
        self.start_point = start_point
        self.robot_height = robot_height
        self.message = message

        self.start_time = int(time.time() * 1000)

        self.width = width
        self.height = height

        self.f = mlab.figure(size=(self.width, self.height))
        visual.set_viewer(self.f)

        v = mlab.view(270, 180)
        #print v

        engine = mlab.get_engine()
        self.s = engine.current_scene
        self.s.scene.interactor.add_observer('KeyPressEvent', self.keypress_callback)

        self.robots = []

        colors = list(PathBatterySimulator.color_codes)

        for key, local_hex_list in sorted(hex_list['internal_routes'].items()):
            color = colors.pop(0)

            ball = visual.sphere(color=color, radius=PathBatterySimulator.ball_radius)
            ball.x = self.start_point[0]
            ball.y = self.start_point[1]
            ball.z = self.start_point[2]

            r, g, b = color
            rt = r + (0.25 * (1 - r))
            gt = g + (0.25 * (1 - g))
            bt = b + (0.25 * (1 - b))
            curve_color = (rt, gt, bt)

            curve = visual.curve(color=curve_color, radius=PathBatterySimulator.curve_radius)

            r_ball = RobotBall(key, local_hex_list, hex_list['external_routes'][key], ball, curve)
            self.robots.append(r_ball)

        x = np.linspace(0, self.width, 1)
        y = np.linspace(0, self.height, 1)

        z = np.loadtxt(np_file)
        z *= 255.0/z.max()
        mlab.surf(x, y, z)

        self.master_cmd = MasterCommand(self.robots)

    def create_surf_map(self, zip_path, hgt_path, map_offset_x = 0.0,
                        map_offset_y = 0.0, map_offset_z = 0.0):

                self.parse_srtm_data(zip_path, hgt_path)

                mlab.surf(self._data, name=hgt_path, colormap='gist_earth', warp_scale=0.2,
                            vmin=self.vmin, vmax=self.vmax)

                self.adjust_offsets(map_offset_x, map_offset_y, map_offset_z)

                del self._data

    def __init__(self, **traits):
        # Has to be done since it's a Traits class
        super(MayaviTraitsApp, self).__init__(**traits)

        x, y = mgrid[-5:5:100j, -5:5:100j]
        r = sqrt(x**2 + y**2)
        z = 5*sin(r)/r

        # Build the visualization of the data: specifying which Mayavi scene to
        # use is necessary
        mlab.surf(x, y, z, figure=self.scene.mayavi_scene)

def plotData2D(a):
    """
    Reduce one dimension of the 3D array
    """
    n = np.empty((a.shape[0], a.shape[1]), dtype=a.dtype)
    for i in range(a.shape[0]): 
        for j in range(a.shape[1]): 
            s = np.sum(a[i, j, :])
            n[i, j] = np.round(s/20)
    
    mlab.surf(n)
    mlab.show()

    def load_3D_plot( self ):
        config = 0
        os.chdir( 'DATA' )
        # os.chdir( self.load_folder3D )
        os.chdir( self.dropdownls3D )
        x = np.load( 'x.npy' )
        y = np.load( 'y.npy' )
        z = np.load( 'z.npy' )
        print 'min_x', np.min( x )
        print 'max_x', np.max( x )
        print 'min_y', np.min( y )
        print 'max_y', np.max( y )
        print 'min_z', np.min( z )
        print 'max_z', np.max( z )
        # print x_axis, y_axis, z_axis
        engine = Engine()
        engine.start()
        if len( engine.scenes ) == 0:
            engine.new_scene  # print os.curdir
            # os.chdir( '.' )()
        if self.autowarp_bool:
            x = x / x[-1]
            y = y / y[0][-1]
            z = z / z[-1] * self.z_scale
        
        # print x, y, z
        mlab.surf( x, y , z, representation = 'wireframe', line_width = 10 )  # ,warp_scale = "auto"
        if 'config.py' in os.listdir( os.curdir ):
            execfile( 'config.py' )
            os.chdir( os.pardir )
            os.chdir( os.pardir )
        else:
            os.chdir( os.pardir )
            os.chdir( os.pardir )
            surface = engine.scenes[0].children[0].children[0].children[0].children[0].children[0]
            surface.actor.mapper.scalar_range = np.array( [ 6.97602671, 8.8533387 ] )
            surface.actor.mapper.scalar_visibility = False
            scene = engine.scenes[0]
            scene.scene.background = ( 1.0, 1.0, 1.0 )
            surface.actor.property.specular_color = ( 0.0, 0.0, 0.0 )
            surface.actor.property.diffuse_color = ( 0.0, 0.0, 0.0 )
            surface.actor.property.ambient_color = ( 0.0, 0.0, 0.0 )
            surface.actor.property.color = ( 0.0, 0.0, 0.0 )
            surface.actor.property.line_width = 1.
            scene.scene.isometric_view()
        mlab.outline()
        outline = engine.scenes[0].children[0].children[0].children[0].children[0].children[1]
        outline.actor.property.specular_color = ( 0.0, 0.0, 0.0 )
        outline.actor.property.diffuse_color = ( 0.0, 0.0, 0.0 )
        outline.actor.property.ambient_color = ( 0.0, 0.0, 0.0 )
        outline.actor.property.color = ( 0.0, 0.0, 0.0 )

        mlab.show()

    def plot_data(cls, geomdim, p, data, style="-", opacity=1.0):
        if geomdim == 1:
            # 2d plot
            plt.plot(p[0], data, style)
#            plt.ylim([-1, 1]) 
            plt.show()
        else:
            # 3d plot
            assert geomdim == 2
            if style == "-":
                mlab.surf(p[0], p[1], data, opacity=opacity)
            else:
                mlab.points3d(p[0], p[1], data)

def surf3d(df, **surfargs):
    """ Mayavi mlab surf plot from a dataframe.  For now, just posting here because it works, not sure
    how necessary it will be."""
    from mayavi import mlab

    warp_scale=surfargs.pop('warp_scale', 'auto')

    ### So this is how I need to overcome glitches to pass x,y in if they are needed.  If no axis needed,
    ### see below.
    #mlab.surf(np.asarray(list(df.columns)), np.asarray(list(df.index)) , np.asarray(df), warp_scale='auto')

    mlab.surf(np.asarray(df), warp_scale=warp_scale, **surfargs)
    mlab.show()