def pie(codes,pounds,filename):

	plt.xlabel("Date")
	plt.ylabel("Pound")
	fig = plt.figure(facecolor='white')

	ax = fig.add_axes([0.05, 0.0, 0.80, 0.94])

	# Set graph label size
	matplotlib.rcParams['font.size'] = 7

	# Generate pie chart
	patches, texts, autotexts = ax.pie(pounds,labels=codes, autopct='%1.1f%%')
	
	# Don't show percentages on pie chart
	proptease = fm.FontProperties()
	proptease.set_size('0')
	plt.setp(autotexts, fontproperties=proptease)

	# Set graph size
	F = pylab.gcf()
	DefaultSize = F.get_size_inches()
	F.set_size_inches( (DefaultSize[0]*0.7, DefaultSize[1]*0.7) )
	
	# Save graph
	plt.savefig(filename) 
	plt.close()

    def plot_std_meshlines(self, step=0.1):
        '''
        plot mesh circles for stdv
        '''

        color = self.std_color

        nstdmax = self.stdmax
        if self.negative:
            axmin = -np.pi / 2.
        else:
            axmin = 0.

        th = np.arange(axmin, np.pi / 2, 0.01)

        for ra in np.arange(0, nstdmax + 0.1 * step, step):
            self.ax.plot(ra * np.sin(th), ra * np.cos(th), ':', color=color)

        if self.normalize:
            self.ax.set_ylabel('$\sigma / \sigma_{obs}$', color=color)
            self.ax.set_xlabel('$\sigma / \sigma_{obs}$', color=color)
        else:
            self.ax.set_ylabel('Standard Deviation', color=color)
            self.ax.set_xlabel('Standard Deviation', color=color)

        xticklabels = plt.getp(plt.gca(), 'xticklabels')
        plt.setp(xticklabels, color=color)
        yticklabels = plt.getp(plt.gca(), 'yticklabels')
        plt.setp(yticklabels, color=color)

    def __init__(self, ax, collection, mmc, img):
        self.colornormalizer = Normalize(vmin=0, vmax=1, clip=False)
        self.scat = plt.scatter(img[:, 0], img[:, 1], c=mmc.classvec)
        plt.gray()
        plt.setp(ax.get_yticklabels(), visible=False)
        ax.yaxis.set_tick_params(size=0)
        plt.setp(ax.get_xticklabels(), visible=False)
        ax.xaxis.set_tick_params(size=0)
        self.img = img
        self.canvas = ax.figure.canvas
        self.collection = collection
        #self.alpha_other = alpha_other
        self.mmc = mmc
        self.prevnewclazz = None

        self.xys = collection
        self.Npts = len(self.xys)
        
        self.lockedset = set([])

        self.lasso = LassoSelector(ax, onselect=self.onselect)#, lineprops = {:'prism'})
        self.lasso.disconnect_events()
        self.lasso.connect_event('button_press_event', self.lasso.onpress)
        self.lasso.connect_event('button_release_event', self.onrelease)
        self.lasso.connect_event('motion_notify_event', self.lasso.onmove)
        self.lasso.connect_event('draw_event', self.lasso.update_background)
        self.lasso.connect_event('key_press_event', self.onkeypressed)
        #self.lasso.connect_event('button_release_event', self.onrelease)
        self.ind = []
        self.slider_axis = plt.axes(slider_coords, visible = False)
        self.slider_axis2 = plt.axes(obj_fun_display_coords, visible = False)
        self.in_selection_slider = None
        newws = list(set(range(len(self.collection))) - self.lockedset)
        self.mmc.new_working_set(newws)
        self.lasso.line.set_visible(False)

def test2():
    mu = 10.0
    sigma = 2.0

    x = Variable("x", float)
    loggauss = -0.5 * math.log( 2.0 * math.pi * sigma * sigma ) - 0.5 * ((x - mu) ** 2) / (sigma * sigma)

    f = Function(
            name="foo",
            params=(x,),
            rettype=float,
            expr=loggauss)
    engine = create_execution_engine()
    module = f.compile(engine)
    func_ptr = engine.get_pointer_to_function(module.get_function("foo"))

    samples = metropolis_hastings(func_ptr, sigma, 0.0, 1000, 2)
    #plt.plot(np.arange(len(samples)), samples)

    n, bins, patches = plt.hist(samples[500:], 25, normed=1, histtype='stepfilled')
    plt.setp(patches, 'facecolor', 'g', 'alpha', 0.75)

    # add a line showing the expected distribution
    y = plt.normpdf(bins, mu, sigma)
    l = plt.plot(bins, y, 'k--', linewidth=1.5)
    plt.show()

def startAnim(x, m, th, Tsim, inter=1, Tstart=0, h=0.002):
#    fig, axM = subplo1ts()
#    axX = axM.twinx()
    fig = figure()
    axM = subplot(211)
    axX = subplot(212, sharex=axM)
    anim = MyAnim(x, m, th, Tsim, inter, Tstart/h, h)

    anim.line1, = axM.plot([], [], 'b')
    anim.line2, = axX.plot([], [], 'g')
    setp(axM.get_xticklabels(), visible=False)
    axM.set_xlim([-pi, pi])
    axX.set_xlim([-pi, pi])
#    axM.set_ylim([0., 3])
#    axX.set_ylim([0., 1.])
    axM.set_ylim([amin(m), amax(m)])
    axX.set_ylim([amin(x), amax(x)])

    axM.set_ylabel(r"$m$")
    axX.set_ylabel(r"$x$")
    axX.set_xlabel(r"$\theta$")
    hold(False)
    for tl in axM.get_yticklabels():
        tl.set_color('b')
    for tl in axX.get_yticklabels():
        tl.set_color('g')
    anim.axM = axM
    anim.axX = axX
    fig.canvas.mpl_connect('button_press_event', anim.onClick)
    a = FuncAnimation(fig, anim, frames=anim.dataGen, init_func=anim.init,
                 interval=0, blit=True, repeat=True)
    show()
    return a

def snIaspec_with_medbands(z = 1.2):
    """ plot a SNIa spectrum at the given z, overlaid with medium bands
    matching the SN Camille obs set.
    :param z:
    :return:
    """
    import stardust
    from demofig import w763,f763,w845,f845,w139,f139
    w1a, f1a = stardust.snsed.getsed( sedfile='/usr/local/SNDATA_ROOT/snsed/Hsiao07.dat', day=0 )

    w1az = w1a * (1+z)
    f1az = f1a / f1a.max() / 2.
    ax18 = pl.gca() # subplot(3,2,1)
    ax18.plot(w1az, f1az, ls='-', lw=0.7, color='0.5', label='_nolegend_')
    ax18.plot(w763, f763, ls='-', color='DarkOrchid',label='F763M')
    ax18.plot(w845, f845, ls='-',color='Teal', label='F845M')
    ax18.plot(w139, f139, ls='-',color='Maroon', label='F139M')
    #ax18.fill_between( w1az, np.where(f763>0.01,f1az,0), color='DarkOrchid', alpha=0.3 )
    #ax18.fill_between( w1az, f1az, where=((w1az>13500) & (w1az<14150)), color='teal', alpha=0.3 )
    #ax18.fill_between( w1az, f1az, where=((w1az>15000) & (w1az<15700)), color='Maroon', alpha=0.3 )
    ax18.text(0.95,0.4, 'SNIa\[email protected] z=%.1f'%(z), color='k',ha='right',va='bottom',fontweight='bold', transform=ax18.transAxes, fontsize='large')
    ax18.set_xlim( 6500, 16000 )
    pl.setp(ax18.get_xticklabels(), visible=False)
    pl.setp(ax18.get_yticklabels(), visible=False)
    ax18.text( 7630, 0.65, 'F763M', ha='right', va='center', color='DarkOrchid', fontweight='bold')
    ax18.text( 8450, 0.65, 'F845M', ha='center', va='center', color='Teal', fontweight='bold')
    ax18.text( 13900, 0.65, 'F139M', ha='left', va='center', color='Maroon', fontweight='bold')

    def OnButtonTidyButton(self, event):
        
        # for easy coding
        T = self.TreeCtrlMain
        s = T.GetSelection()
        f = self.GetTreeItemData(s, "figure") 
        w = self.GetTreeItemData(s, "window")
        
        # set the current figure
        pylab.figure(f.number)
        
        # first set the size of the window
        w.SetSize([500,500])
        
        # now loop over all the data and get the range
        lines = f.axes[0].get_lines()
        
        # we want thick lines
        f.axes[0].get_frame().set_linewidth(3.0)

        # get the tick lines in one big list
        xticklines = f.axes[0].get_xticklines()
        yticklines = f.axes[0].get_yticklines()
        
        # set their marker edge width
        pylab.setp(xticklines+yticklines, mew=2.0)
        
        # set what kind of tickline they are (outside axes)
        for l in xticklines: l.set_marker(matplotlib.lines.TICKDOWN)
        for l in yticklines: l.set_marker(matplotlib.lines.TICKLEFT)
        
        # get rid of the top and right ticks
        f.axes[0].xaxis.tick_bottom()
        f.axes[0].yaxis.tick_left()
        
        # we want bold fonts
        pylab.xticks(fontsize=20, fontweight='bold', fontname='Arial')
        pylab.yticks(fontsize=20, fontweight='bold', fontname='Arial')

        # we want to give the labels some breathing room (1% of the data range)
        for label in pylab.xticks()[1]:
            label.set_y(-0.02)
        for label in pylab.yticks()[1]:
            label.set_x(-0.01)
            
        # set the position/size of the axis in the window
        f.axes[0].set_position([0.1,0.1,0.8,0.8])
        
        # set the axis labels
        f.axes[0].set_title('')
        f.axes[0].set_xlabel('')
        f.axes[0].set_ylabel('')

        # set the position of the legend far away
        f.axes[0].legend(loc=[1.2,0])
        
        f.canvas.Refresh()
        
        # autoscale
        self.OnButtonAutoscaleButton(None)

    def plot(self,file):
        cds = CaseDataset(file, 'bson')
        data = cds.data.driver('driver').by_variable().fetch()
        cds2 = CaseDataset('../output/therm_mc_20141110173851.bson', 'bson')
        data2 = cds2.data.driver('driver').by_variable().fetch()
        
        #temp
        temp_boundary_k = data['hyperloop.temp_boundary']
        temp_boundary_k.extend(data2['hyperloop.temp_boundary'])
        temp_boundary = [((x-273.15)*1.8 + 32) for x in temp_boundary_k]
        #histogram
        n, bins, patches = plt.hist(temp_boundary, 100, normed=1, histtype='stepfilled')
        plt.setp(patches, 'facecolor', 'b', 'alpha', 0.75)

        #stats
        mean = np.average(temp_boundary)
        std = np.std(temp_boundary)
        percentile = np.percentile(temp_boundary,99.5)
        print "mean: ", mean, " std: ", std, " 99.5percentile: ", percentile
        x = np.linspace(50,170,150)
        plt.plot(x,mlab.normpdf(x,mean,std), color='black', lw=2)
        plt.xlim([60,160])
        plt.ylabel('Probability', fontsize=18)
        plt.xlabel(u'Equilibrium Temperature, \N{DEGREE SIGN}F', fontsize=18)
        #plt.show()
        plt.tight_layout()
        plt.savefig('../output/histo.pdf', dpi=300)

def interev_mag(times, mags):
    r"""Function to plot interevent times against magnitude for given times
    and magnitudes.

    :type times: list of datetime
    :param times: list of the detection times, must be sorted the same as mags
    :type mags: list of float
    :param mags: list of magnitudes
    """
    l = [(times[i], mags[i]) for i in xrange(len(times))]
    l.sort(key=lambda tup: tup[0])
    times = [x[0] for x in l]
    mags = [x[1] for x in l]
    # Make two subplots next to each other of time before and time after
    fig, axes = plt.subplots(1, 2, sharey=True)
    axes = axes.ravel()
    pre_times = []
    post_times = []
    for i in range(len(times)):
        if i > 0:
            pre_times.append((times[i] - times[i - 1]) / 60)
        if i < len(times) - 1:
            post_times.append((times[i + 1] - times[i]) / 60)
    axes[0].scatter(pre_times, mags[1:])
    axes[0].set_title('Pre-event times')
    axes[0].set_ylabel('Magnitude')
    axes[0].set_xlabel('Time (Minutes)')
    plt.setp(axes[0].xaxis.get_majorticklabels(), rotation=30)
    axes[1].scatter(pre_times, mags[:-1])
    axes[1].set_title('Post-event times')
    axes[1].set_xlabel('Time (Minutes)')
    plt.setp(axes[1].xaxis.get_majorticklabels(), rotation=30)
    plt.show()

    def draw(self):
        data_len = len(self.data)
        X = range(data_len)
        data = [float(i) for i in self.data]
        color = self.cp.get("Colors", item_name("Line", self.num))
        line = self.ax.plot(X, data, marker="o", markeredgecolor=color,
            markerfacecolor="white", markersize=13, linewidth=5,
            markeredgewidth=5, color=color, label=self.cp.get("Labels",
            item_name("Legend", self.num)))
        max_ax = max(data)*1.1 if data else 0
        if max_ax <= 10:
            self.ax.set_ylim(-0.1, max_ax)
        else:
            self.ax.set_ylim(-0.5, max_ax)
        self.ax.set_xlim(-1, data_len+1)

        if self.mode == "hourly":
            self.ax.xaxis.set_ticks((0, 8, 16, 24))
        elif self.mode == "daily":
            self.ax.xaxis.set_ticks((0, 15, 30))
        elif self.mode == "monthly":
            self.ax.xaxis.set_ticks((0, 4, 8, 12))
        else:
            raise Exception("Unknown time interval mode.")

        if self.legend:
            legend = self.ax.legend(loc=9, mode="expand",
                bbox_to_anchor=(0.25, 1.02, 1., .102))
            setp(legend.get_frame(), visible=False)
            setp(legend.get_texts(), size=int(self.cp.get("Sizes",
                "LegendSize")))

    def draw_plot(self):
        """ Redraws the plot
        """
        # when xmin is on auto, it "follows" xmax to produce a 
        # sliding window effect. therefore, xmin is assigned after
        # xmax.
        #
        xwin_size = 360
        if self.xmax_control.is_auto():
            xmax = len(self.data) if len(self.data) > xwin_size else xwin_size
        else:
            xmax = int(self.xmax_control.manual_value())
            
        if self.xmin_control.is_auto():            
            xmin = xmax - xwin_size
        else:
            xmin = int(self.xmin_control.manual_value())

        # for ymin and ymax, find the minimal and maximal values
        # in the data set and add a mininal margin.
        # 
        # note that it's easy to change this scheme to the 
        # minimal/maximal value in the current display, and not
        # the whole data set.
        # 
        if self.ymin_control.is_auto():
            ymin = round(min(self.data), 0) - 1
        else:
            ymin = int(self.ymin_control.manual_value())
        
        if self.ymax_control.is_auto():
            ymax = round(max(self.data), 0) + 1
        else:
            ymax = int(self.ymax_control.manual_value())

        self.axes.set_xbound(lower=xmin, upper=xmax)
        self.axes.set_ybound(lower=ymin, upper=ymax)
        
        # anecdote: axes.grid assumes b=True if any other flag is
        # given even if b is set to False.
        # so just passing the flag into the first statement won't
        # work.
        #
        if self.cb_grid.IsChecked():
            self.axes.grid(True, color='gray')
        else:
            self.axes.grid(False)

        # Using setp here is convenient, because get_xticklabels
        # returns a list over which one needs to explicitly 
        # iterate, and setp already handles this.
        #  
        pylab.setp(self.axes.get_xticklabels(),
            visible=self.cb_xlab.IsChecked())
        
        self.plot_data.set_xdata(np.arange(len(self.data)))
        self.plot_data.set_ydata(np.array(self.data))
        
        self.canvas.draw()

def dateticks(fmt='%Y-%m', **kwargs):
    '''setup the date ticks'''
    dateticker = ticker.FuncFormatter(lambda numdate, _: num2date(numdate).strftime(fmt))
    pylab.gca().xaxis.set_major_formatter(dateticker)
    # pylab.gcf().autofmt_xdate()
    tmp = dict(rotation=30, ha='right')
    tmp.update(kwargs)
    pylab.setp(pylab.xticks()[1], **tmp)

def heatmap(data, row_labels, col_labels, ax=None,
            cbar_kw={}, cbarlabel="", title = "",  **kwargs):
    """
    Create a heatmap from a numpy array and two lists of labels.
    Arguments:
        data       : A 2D numpy array of shape (N,M)
        row_labels : A list or array of length N with the labels
                     for the rows
        col_labels : A list or array of length M with the labels
                     for the columns
    Optional arguments:
        ax         : A matplotlib.axes.Axes instance to which the heatmap
                     is plotted. If not provided, use current axes or
                     create a new one.
        cbar_kw    : A dictionary with arguments to
                     :meth:`matplotlib.Figure.colorbar`.
        cbarlabel  : The label for the colorbar
    All other arguments are directly passed on to the imshow call.
    """

    if not ax:
        ax = plt.gca()

    # Plot the heatmap
    im = ax.imshow(data, **kwargs)
    ax.set_title(title, pad =50.0)
    # create an axes on the right side of ax. The width of cax will be 5%
    # of ax and the padding between cax and ax will be fixed at 0.05 inch.
    divider = make_axes_locatable(ax)
    cax = divider.append_axes("right", size="5%", pad=0.05)
    cbar = ax.figure.colorbar(im, ax=ax, cax=cax, **cbar_kw)
    cbar.ax.set_ylabel(cbarlabel, rotation=-90, va="bottom")

    # We want to show all ticks...
    ax.set_xticks(np.arange(data.shape[1]))
    ax.set_yticks(np.arange(data.shape[0]))
    # ... and label them with the respective list entries.
    ax.set_xticklabels(col_labels)
    ax.set_yticklabels(row_labels)

    # Let the horizontal axes labeling appear on top.
    ax.tick_params(top=True, bottom=False,
                   labeltop=True, labelbottom=False)

    # Rotate the tick labels and set their alignment.
    plt.setp(ax.get_xticklabels(), rotation=-30, ha="right",
             rotation_mode="anchor")

    # Turn spines off and create white grid.
    # for edge, spine in ax.spines.items():
    #    spine.set_visible(False)

    ax.set_xticks(np.arange(0, data.shape[1] + 1) - 0.5, minor=True)
    ax.set_yticks(np.arange(0, data.shape[0] + 1) - 0.5, minor=True)
    # ax.grid(which="minor", color="w", linestyle='-', linewidth=3)
    # ax.tick_params(which="minor", bottom=False, left=False)

    return im, cbar

    def create(self, images, x, y, mag, plot=False):
        """Using some input images and a catalog entry, define the pixel aperture for a star."""

        # keep track of the basics of this aperture
        self.x, self.y, self.mag = x, y, mag

        # figure out which label in the labeled image is relevant to this star
        label = images['labeled'][np.round(y), np.round(x)]
        if label == 0:
            self.row, self.col = np.array([np.round(y).astype(np.int)]), np.array([np.round(x).astype(np.int)])
        else:

            # identify and sort the pixels that might contribute
            ok = images['labeled'] == label
            okr, okc = ok.nonzero()
            sorted = np.argsort(images['stars'][ok]/images['noise'][ok])[::-1]

            signal = np.cumsum(images['stars'][ok][sorted])
            noise = np.sqrt(np.cumsum(images['noise'][ok][sorted]**2))
            snr = signal/noise
            mask = ok*0
            toinclude = sorted[:np.argmax(snr)+1]
            mask[okr[toinclude], okc[toinclude]] = 1

            self.row, self.col = okr[toinclude], okc[toinclude]

            if plot:
                fi = plt.figure('selecting an aperture', figsize=(10,3), dpi=100)
                fi.clf()
                gs = gridspec.GridSpec(3,2,width_ratios=[1,.1], wspace=0.1, hspace=0.01, top=0.9, bottom=0.2)
                ax_line = plt.subplot(gs[:,0])
                ax_image = plt.subplot(gs[0,1])
                ax_ok = plt.subplot(gs[1,1])
                ax_mask = plt.subplot(gs[2,1])

                shape = ok.shape
                row, col = ok.nonzero()
                left, right = np.maximum(np.min(col)-1, 0), np.minimum(np.max(col)+2, shape[1])
                bottom, top = np.maximum(np.min(row)-1, 0),  np.minimum(np.max(row)+2, shape[1])

                kwargs = dict( extent=[left, right, bottom, top],  interpolation='nearest')
                ax_image.imshow(np.log(images['median'][bottom:top, left:right]), cmap='gray_r', **kwargs)
                ax_ok.imshow(ok[bottom:top, left:right], cmap='Blues', **kwargs)
                ax_mask.imshow(mask[bottom:top, left:right], cmap='YlOrRd', **kwargs)

                for a in [ax_ok, ax_mask, ax_image]:
                    plt.setp(a.get_xticklabels(), visible=False)
                    plt.setp(a.get_yticklabels(), visible=False)

                ax_line.plot(signal.flatten(), signal.flatten()/noise.flatten(), marker='o', linewidth=1, color='black',alpha=0.5, markersize=10)
                ax_line.set_xlabel('Total Star Flux in Aperture')
                ax_line.set_ylabel('Total Signal-to-Noise Ratio')
                ax_line.set_title('{0:.1f} magnitude star'.format(mag))
                plt.draw()
                self.input('hmmm?')
        self.n = len(self.row)
        logger.info('created {0}'.format(self))

def plot_data(hist):
    X = np.arange(len(hist))
    plt.bar(X, hist.values(), align='center', width=0.5)
    plt.xticks(X, hist.keys())
    locs, labels = plt.xticks()
    plt.setp(labels, rotation=90)
    ymax = max(hist.values()) + 0.1
    plt.ylim(0, ymax)
    plt.show()

def print_histogram(y, num_bins):   
    # Prints a histogram of input array with equally spaced bins
    
    # Input
    # y: array
    # num_bins: number of bins in histogram
    
    _, _, patches = py.hist(y, num_bins, histtype='stepfilled')
    py.setp(patches, 'facecolor', 'g', 'alpha', 0.75)
    py.show()

 def make_chart(self,data=np.zeros(1),bins=np.arange(10)+1,polar=False):
     self.polar_chart = polar
     self.hist_bins = []
     self.hist_patches = []
     self.x = np.arange(17)
     self.means = np.zeros(self.x.size)
     self.stds = np.zeros(self.x.size)
     
     self.fitting_x = np.linspace(self.x[0], self.x[-1], 100, endpoint=True)
     self.fitting_y = np.zeros(self.fitting_x.size)
     self.fig.clear()
     
     grid = 18
     height = grid // 9
     gs = gridspec.GridSpec(grid, grid)
     # make tuning curve plot
     axes = self.fig.add_subplot(gs[:-height*2,height//2:-height//2],polar=polar)
     if polar:
         self.curve_data = axes.plot(self.x, self.means, 'ko-')[0]
     else:
         adjust_spines(axes,spines=['left','bottom','right'],spine_outward=['left','right','bottom'],xoutward=10,youtward=30,\
                       xticks='bottom',yticks='both',tick_label=['x','y'],xaxis_loc=5,xminor_auto_loc=2,yminor_auto_loc=2)
         axes.set_ylabel('Response(spikes/sec)',fontsize=12)
         self.curve_data = axes.plot(self.x, self.means, self.data_point_styles[0])[0]
     self.errbars = axes.errorbar(self.x, self.means, yerr=self.stds, fmt='k.') if self.showing_errbar else None
     self.curve_axes = axes
     
     self.fitting_data = axes.plot(self.fitting_x, self.fitting_y, self.fitting_curve_styles[0])[0]
     
     axes.set_ylim(0,100)
     axes.relim()
     axes.autoscale_view(scalex=True, scaley=False)
     axes.grid(b=True, which='major',axis='both',linestyle='-.')
     # make histgrams plot
     rows,cols = (grid-height,grid)
     for row in range(rows,cols)[::height]:
         for col in range(cols)[::height]:
             axes = self.fig.add_subplot(gs[row:row+height,col:col+height])
             axes.set_axis_bgcolor('white')
             #axes.set_title('PSTH', size=8)
             axes.set_ylim(0,100)
             if col == 0:
                 adjust_spines(axes,spines=['left','bottom'],xticks='bottom',yticks='left',tick_label=['y'],xaxis_loc=4,yaxis_loc=3)
                 axes.set_ylabel('Spikes',fontsize=11)
             elif col == cols-height:
                 adjust_spines(axes,spines=['right','bottom'],xticks='bottom',yticks='right',tick_label=['y'],xaxis_loc=4,yaxis_loc=3)
             else:
                 adjust_spines(axes,spines=['bottom'],xticks='bottom',yticks='none',tick_label=[],xaxis_loc=4,yaxis_loc=3)
             pylab.setp(axes.get_xticklabels(), fontsize=8)
             pylab.setp(axes.get_yticklabels(), fontsize=8)
             _n, bins, patches = axes.hist(data, bins, facecolor='black', alpha=1.0)
             self.hist_bins.append(bins)
             self.hist_patches.append(patches)
             
     self.fig.canvas.draw()

  def plotCorrEnh(self, parsList=None, **plotArgs):
    """
      Produces enhanced correlation plots.
      
      Parameters
      ----------
      parsList : list of string,  optional
          If not given, all available traces are used.
          Otherwise a list of at least two parameters
          has to be specified.
      plotArgs : dict, optional
          Keyword arguments handed to plot procedures of
          pylab. The following keywords are available:
          contour,bins,cmap,origin,interpolation,colors
    """
    if not ic.check["matplotlib"]:
      PE.warn(PE.PyARequiredImport("To use 'plotCorr', matplotlib has to be installed.", \
                                   solution="Install matplotlib."))
      return
    tracesDic = {}
    if parsList is not None:
      for parm in parsList:
        self._parmCheck(parm)
        tracesDic[parm] = self[parm]
      if len(tracesDic) < 2:
        raise(PE.PyAValError("For plotting correlations, at least two valid parameters are needed.", \
                             where="TraceAnalysis::plotCorr"))
    else:
      # Use all available traces
      for parm in self.availableParameters():
        tracesDic[parm] = self[parm]

    pars = tracesDic.keys()
    traces = tracesDic.values()

    fontmap = {1:10, 2:9, 3:8, 4:8, 5:8}
    if not len(tracesDic)-1 in fontmap:
      fontmap[len(tracesDic)-1] = 8

    k = 1
    for j in range(len(tracesDic)):
      for i in range(len(tracesDic)):
        if i>j:
          plt.subplot(len(tracesDic)-1,len(tracesDic)-1,k)
#          plt.title("Pearson's R: %1.5f" % self.pearsonr(pars[j],pars[i])[0], fontsize='x-small')
          plt.xlabel(pars[j], fontsize=fontmap[len(tracesDic)-1])
          plt.ylabel(pars[i], fontsize=fontmap[len(tracesDic)-1])
          tlabels = plt.gca().get_xticklabels()
          plt.setp(tlabels, 'fontsize', fontmap[len(tracesDic)-1])
          tlabels = plt.gca().get_yticklabels()
          plt.setp(tlabels, 'fontsize', fontmap[len(tracesDic)-1])
#          plt.plot(traces[j],traces[i],'.',**plotArgs)
          self.__hist2d(traces[j],traces[i],**plotArgs)
        if i!=j:
          k+=1

 def modifValZone(self, nameVar, ind, val):
     """modifier dans la zone nameVar la valeur
     ou liste valeur, attention le texte de la zone contient nom et valeur"""
     obj = self.listeZoneText[nameVar][ind].getObject()
     if type(obj) == type([2, 3]):
         for i in range(len(obj)):
             pl.setp(obj[i], text=str(val[i]))
     else:
         nom = pl.getp(obj, "text").split("\n")[0]
         pl.setp(obj, text=nom + "\n" + str(val)[:16])
     self.redraw()

def test_geometric_mechanism(alpha, n_bins=10):  
    X = geometric_mechanism(alpha, 10000) 
    # plot geometric distribution
    n, bins, patches = P.hist(X, 2*n_bins, range=(-n_bins,n_bins), normed=1, histtype='stepfilled')   # bins: ndarray of 51 elements
    P.setp(patches, 'facecolor', 'g', 'alpha', 0.75)
    print bins
    # add a line showing the expected distribution
    x = np.array(range(-n_bins, n_bins+1))          # x: -n_bins -> n_bins
    y = np.power(alpha, abs(x))*(1-alpha)/(1+alpha)    # y = (1-alpha)/(1+alpha)*alpha^|x|
    l = P.plot(x, y, 'k--', linewidth=1.5)

    P.show()