def plot_inference_summary(inference_record):
     ll = []
     lp = []
     wlp_plus_ll=[]
     for step in inference_record.steps:
         ll += step.ll[1:] # start from 1 and not 0: to skip the initial guess
         try:
             lp += step.lp[1:] 
             
             wlp_plus_ll += list((step.wlp * np.asarray(step.lp[1:]) + 
                                         np.asarray(step.ll[1:])).tolist())
         except AttributeError:
             pass
                 
     plt.title('ll',fontsize=30)
     plt.plot(ll,lw=2)
     plt.plot(lp,lw=2)
     plt.plot(wlp_plus_ll,lw=2)
     
      
     counter = 0
     for i,step in enumerate(inference_record.steps):
         if i%2==1:
             facecolor = ".2"
         else:
             facecolor = ".5"           
         plt.axvspan(counter, counter+step.nAccepted, facecolor=facecolor, alpha=0.2)
         counter += step.nAccepted 

def plot_zipf(*freq):
	'''
	basic plotting using matplotlib and pylab
	'''
	ranks, frequencies = [], []
	langs, colors = [], []
	langs = ["English", "German", "Finnish"]
	colors = ['#FF0000', '#00FF00', '#0000FF']
	if bonus_part:
		colors.extend(['#00FFFF', '#FF00FF', '#FFFF00'])
		langs.extend(["English (Stemmed)", "German (Stemmed)", "Finnish (Stemmed)"])

	plt.subplot(111) # 1, 1, 1

	num = 6 if bonus_part else 3
	for i in xrange(num):
		ranks.append(range(1, len(freq[i]) + 1))
		frequencies.append([e[1] for e in freq[i]])

		# log x and y axi, both with base 10
		plt.loglog(ranks[i], frequencies[i], marker='', basex=10, color=colors[i], label=langs[i])

	plt.legend()
	plt.grid(True)
	plt.title("Zipf's law!")

	plt.xlabel('Rank')
	plt.ylabel('Frequency')

	plt.show()

def example_filterbank():
    from pylab import plt
    import numpy as np

    x = _create_impulse(2000)
    gfb = GammatoneFilterbank(density=1)

    analyse = gfb.analyze(x)
    imax, slopes = gfb.estimate_max_indices_and_slopes()
    fig, axs = plt.subplots(len(gfb.centerfrequencies), 1)
    for (band, state), imx, ax in zip(analyse, imax, axs):
        ax.plot(np.real(band))
        ax.plot(np.imag(band))
        ax.plot(np.abs(band))
        ax.plot(imx, 0, 'o')
        ax.set_yticklabels([])
        [ax.set_xticklabels([]) for ax in axs[:-1]]

    axs[0].set_title('Impulse responses of gammatone bands')

    fig, ax = plt.subplots()

    def plotfun(x, y):
        ax.semilogx(x, 20*np.log10(np.abs(y)**2))

    gfb.freqz(nfft=2*4096, plotfun=plotfun)
    plt.grid(True)
    plt.title('Absolute spectra of gammatone bands.')
    plt.xlabel('Normalized Frequency (log)')
    plt.ylabel('Attenuation /dB(FS)')
    plt.axis('Tight')
    plt.ylim([-90, 1])
    plt.show()

    return gfb

    def draw(cls, t_max, agents_proportions, eco_idx, parameters):

        color_set = ["green", "blue", "red"]

        for agent_type in range(3):
            plt.plot(np.arange(t_max), agents_proportions[:, agent_type],
                     color=color_set[agent_type], linewidth=2.0, label="Type-{} agents".format(agent_type))

            plt.ylim([-0.1, 1.1])

        plt.xlabel("$t$")
        plt.ylabel("Proportion of indirect exchanges")

        # plt.suptitle('Direct choices proportion per type of agents', fontsize=14, fontweight='bold')
        plt.legend(loc='upper right', fontsize=12)

        print(parameters)

        plt.title(
            "Workforce: {}, {}, {};   displacement area: {};   vision area: {};   alpha: {};   tau: {}\n"
            .format(
                parameters["x0"],
                parameters["x1"],
                parameters["x2"],
                parameters["movement_area"],
                parameters["vision_area"],
                parameters["alpha"],
                parameters["tau"]
                          ), fontsize=12)

        if not path.exists("../../figures"):
            mkdir("../../figures")

        plt.savefig("../../figures/figure_{}.pdf".format(eco_idx))
        plt.show()

def serve_css(name, length, keys, values):
    from pylab import plt, mpl
    mpl.rcParams['font.sans-serif'] = ['SimHei']
    mpl.rcParams['axes.unicode_minus'] = False
    from matplotlib.font_manager import FontProperties
    # font = FontProperties(fname="d:\Users\ll.tong\Desktop\msyh.ttf", size=12)
    font = FontProperties(fname="/usr/share/fonts/msyh.ttf", size=11)
    plt.xlabel(u'')
    plt.ylabel(u'出现次数',fontproperties=font)
    plt.title(u'词频统计',fontproperties=font)
    plt.grid()
    keys = keys.decode("utf-8").split(' ')
    values = values.split(' ')
    valuesInt = []
    for value in values:
        valuesInt.append(int(value))

    plt.xticks(range(int(length)), keys)
    plt.plot(range(int(length)), valuesInt)
    plt.xticks(rotation=defaultrotation, fontsize=9,fontproperties=font)
    plt.yticks(fontsize=10,fontproperties=font)
    name = name + str(datetime.now().date()).replace(':', '') + '.png'
    imgUrl = 'static/temp/' + name
    fig = matplotlib.pyplot.gcf()
    fig.set_size_inches(12.2, 2)
    plt.savefig(imgUrl, bbox_inches='tight', figsize=(20,4), dpi=100)
    plt.close()
    tempfile = static_file(name, root='./static/temp/')
    #os.remove(imgUrl)
    return tempfile

def plot_pre(fn):
    t = read_t(fn)
    y = read_y(fn)
    yres = read_yres(fn)

    plt.plot_date(t+_adj_dates, y, 'x', color='lightblue')
    plt.plot_date(t+_adj_dates, yres, 'x', color='lightgreen')

    linsec = read_linsec(fn)

    ch = cut_ts(t, linsec)
    
    plt.plot_date(t[ch]+_adj_dates, y[ch], 'x', color='blue', label='original')
    plt.plot_date(t[ch]+_adj_dates, yres[ch], 'x', color='green', label='residual')

    outliers = read_outlier(fn)
    idx = outlier_index(t, outliers)

    plt.plot_date(t[idx]+_adj_dates, y[idx], 'o', mec='red', mew=1, mfc='blue')
    plt.plot_date(t[idx]+_adj_dates, yres[idx], 'o', mec='red', mew=1, mfc='green')

    for jump in read_jumps(fn):
        plt.axvline(jump + _adj_dates, color='red', ls='--')

    plt.grid('on')
    site = basename(fn).split('.')[0]
    cmpt = basename(fn).split('.')[1]
    plt.title('%s - %s'%(site, cmpt))

def plot_response(data, plate_name, save_folder = 'Figures/'):
    """
    """
    if not os.path.isdir(save_folder):
        os.makedirs(save_folder)

    for block in data:
        #
        group = group_similar(data[block].keys())
        names = data[block].keys()
        names.sort()
        #
        plt.figure(figsize=(16, 4 + len(names)/8), dpi=300)
        #
        for i, name in enumerate(names):
            a, b, c = get_index(group, name)
            color, pattern = color_shade_pattern(a, b, c, group)
            mean = data[block][name]['mean'][0]
            std = data[block][name]['std'][0]

            plt.barh([i], [mean], height=1.0, color=color, hatch=pattern)
            plt.errorbar([mean], [i+0.5], xerr=[std], ecolor = [0,0,0], linestyle = '')

        plt.yticks([i+0.5 for i in xrange(len(names))], names, size = 8)
        plt.title(plate_name)
        plt.ylim(0, len(names))
        plt.xlabel('change')
        plt.tight_layout()

        plt.savefig(save_folder + 'response_' + str(block + 1))
    #
    return None

    def generate_start_time_figures(self):
        recording_time_grouped_by_patient = self.pain_data[["PatientID", "NRSTimeFromEndSurgery_mins"]].groupby("PatientID")
        recording_start_minutes = recording_time_grouped_by_patient.min()

        fig1 = "fig1.pdf"
        fig2 = "fig2.pdf"

        plt.figure(figsize=[8,4])
        plt.title("Pain score recording start times", fontsize=14).set_y(1.05) 
        plt.ylabel("Occurrences", fontsize=14)
        plt.xlabel("Recording Start Time (minutes)", fontsize=14)
        plt.hist(recording_start_minutes.values, bins=20, color="0.5")
        plt.savefig(os.path.join(self.tmp_directory, fig1), bbox_inches="tight")

        plt.figure(figsize=[8,4])
        plt.title("Pain score recording start times, log scale", fontsize=14).set_y(1.05) 
        plt.ylabel("Occurrences", fontsize=14)
        plt.xlabel("Recording Start Time (minutes)", fontsize=14)
        plt.hist(recording_start_minutes.values, bins=20, log=True, color="0.5")
        plt.savefig(os.path.join(self.tmp_directory, fig2), bbox_inches="tight")

        #save the figures in panel format
        f = open(os.path.join(self.tmp_directory, "tmp.tex"), 'w')
        f.write(r"""
            \documentclass[%
            ,float=false % this is the new default and can be left away.
            ,preview=true
            ,class=scrartcl
            ,fontsize=20pt
            ]{standalone}
            \usepackage[active,tightpage]{preview}
            \usepackage{varwidth}
            \usepackage{graphicx}
            \usepackage[justification=centering]{caption}
            \usepackage{subcaption}
            \usepackage[caption=false,font=footnotesize]{subfig}
            \renewcommand{\thesubfigure}{\Alph{subfigure}}
            \begin{document}
            \begin{preview}
            \begin{figure}[h]
                \begin{subfigure}{0.5\textwidth}
                        \includegraphics[width=\textwidth]{""" + fig1 + r"""}
                        \caption{Normal scale}
                \end{subfigure}\begin{subfigure}{0.5\textwidth}
                        \includegraphics[width=\textwidth]{""" + fig2 + r"""}
                        \caption{Log scale}
                \end{subfigure}
            \end{figure}
            \end{preview}
            \end{document}
        """)
        f.close()
        subprocess.call(["pdflatex", 
                            "-halt-on-error", 
                            "-output-directory", 
                            self.tmp_directory, 
                            os.path.join(self.tmp_directory, "tmp.tex")])
        shutil.move(os.path.join(self.tmp_directory, "tmp.pdf"), 
                    os.path.join(self.output_directory, "pain_score_start_times.pdf"))

    def plot(self, new_plot=False, xlim=None, ylim=None, title=None, figsize=None,
             xlabel=None, ylabel=None, fontsize=None, show_legend=True, grid=True):
        """
        Plot data using matplotlib library. Use show() method for matplotlib to see result or ::

            %pylab inline

        in IPython to see plot as cell output.

        :param bool new_plot: create or not new figure
        :param xlim: x-axis range
        :param ylim: y-axis range
        :type xlim: None or tuple(x_min, x_max)
        :type ylim: None or tuple(y_min, y_max)
        :param title: title
        :type title: None or str
        :param figsize: figure size
        :type figsize: None or tuple(weight, height)
        :param xlabel: x-axis name
        :type xlabel: None or str
        :param ylabel: y-axis name
        :type ylabel: None or str
        :param fontsize: font size
        :type fontsize: None or int
        :param bool show_legend: show or not labels for plots
        :param bool grid: show grid or not

        """
        xlabel = self.xlabel if xlabel is None else xlabel
        ylabel = self.ylabel if ylabel is None else ylabel
        figsize = self.figsize if figsize is None else figsize
        fontsize = self.fontsize if fontsize is None else fontsize
        self.fontsize_ = fontsize
        self.show_legend_ = show_legend
        title = self.title if title is None else title
        xlim = self.xlim if xlim is None else xlim
        ylim = self.ylim if ylim is None else ylim
        new_plot = self.new_plot or new_plot

        if new_plot:
            plt.figure(figsize=figsize)

        plt.xlabel(xlabel, fontsize=fontsize)
        plt.ylabel(ylabel, fontsize=fontsize)
        plt.title(title, fontsize=fontsize)
        plt.tick_params(axis='both', labelsize=fontsize)
        plt.grid(grid)

        if xlim is not None:
            plt.xlim(xlim)

        if ylim is not None:
            plt.ylim(ylim)

        self._plot()

        if show_legend:
            plt.legend(loc='best', scatterpoints=1)

 def plot_smoothed_alpha_comparison(self,rmsval,suffix=''):
     plt.plot(self.f,self.alpha,'ko',label='data set')
     plt.plot(self.f,self.salpha,'c-',lw=2,label='smoothed angle $\phi$')
     plt.xlabel('frequency in Hz')
     plt.ylabel('angle $\phi$ in coordinates of circle')
     plt.legend()
     ylims=plt.axes().get_ylim()
     plt.yticks((arange(9)-4)*0.5*pi, ['$-2\pi$','$-3\pi/2$','$-\pi$','$-\pi/2$','$0$','$\pi/2$','$\pi$','$3\pi/2$','$2\pi$'])
     plt.ylim(ylims)
     plt.title('RMS offset from smooth curve: {:.4f}'.format(rmsval))
     if self.show: plt.show()
     else: plt.savefig(join(self.sdc.plotpath,'salpha','c{}_salpha_on_{}_circle'.format(self.sdc.case,self.ZorY)+self.sdc.suffix+self.sdc.outsuffix+suffix+'.png'), dpi=240)
     plt.close()

def plot_cf(cf, color):
    t = cf.data.t
    y0 = cf.data.y0
    plt.plot_date(t+_adj_dates,y0,'x',label=cf.SITE+cf.CMPT,
              color='light'+color)
    plt.plot_date(cf.data._t+_adj_dates,
              cf.data._y0,'x',label=cf.SITE+cf.CMPT,
              color=color)
    t1 = min(t)
    t2 = max(t)
    ls=200
    plot_func(cf.func,linspace(t1,t2,ls))
    plt.title(cf.SITE+'-'+cf.CMPT)
    plt.gcf().autofmt_xdate()

def main(args=sys.argv[1:]):
    # there are some cases when this script is run on systems without DISPLAY variable being set
    # in such case matplotlib backend has to be explicitly specified
    # we do it here and not in the top of the file, as inteleaving imports with code lines is discouraged
    import matplotlib
    matplotlib.use('Agg')
    from pylab import plt, ylabel, grid, xlabel, array

    parser = argparse.ArgumentParser()
    parser.add_argument("rst_file", help="location of rst file in TRiP98 format", type=str)
    parser.add_argument("output_file", help="location of PNG file to save", type=str)
    parser.add_argument("-s", "--submachine", help="Select submachine to plot.", type=int, default=1)
    parser.add_argument("-f", "--factor", help="Factor for scaling the blobs. Default is 1000.", type=int, default=1000)
    parser.add_argument("-v", "--verbosity", action='count', help="increase output verbosity", default=0)
    parser.add_argument('-V', '--version', action='version', version=pt.__version__)
    args = parser.parse_args(args)

    file = args.rst_file
    sm = args.submachine
    fac = args.factor

    a = pt.Rst()
    a.read(file)

    # convert data in submachine to a nice array
    b = a.machines[sm]

    x = []
    y = []
    z = []
    for _x, _y, _z in b.raster_points:
        x.append(_x)
        y.append(_y)
        z.append(_z)

    title = "Submachine: {:d} / {:d} - Energy: {:.3f} MeV/u".format(sm, len(a.machines), b.energy)
    print(title)
    cc = array(z)

    cc = cc / cc.max() * fac

    fig = plt.figure()
    ax = fig.add_subplot(111)
    ax.scatter(x, y, c=cc, s=cc, alpha=0.75)
    ylabel("mm")
    xlabel("mm")
    grid(True)
    plt.title(title)
    plt.savefig(args.output_file)
    plt.close()

 def plot_overview(self,suffix=''):
     x=self.x; y=self.y; r=self.radius; cx,cy=self.center.real,self.center.imag
     ax=plt.axes()
     plt.scatter(x,y, marker='o', c='b', s=40)
     plt.axhline(y=0,color='grey', zorder=-1)
     plt.axvline(x=0,color='grey', zorder=-2)
     t=linspace(0,2*pi,201)
     circx=r*cos(t) + cx
     circy=r*sin(t) + cy
     plt.plot(circx,circy,'g-')
     plt.plot([cx],[cy],'gx',ms=12)
     if self.ZorY == 'Z':
         philist,flist=[self.phi_a,self.phi_p,self.phi_n],[self.fa,self.fp,self.fn]
     elif self.ZorY == 'Y':
         philist,flist=[self.phi_m,self.phi_s,self.phi_r],[self.fm,self.fs,self.fr]
     for p,f in zip(philist,flist):
         if f is not None:
             xpos=cx+r*cos(p); ypos=cy+r*sin(p); xos=0.2*(xpos-cx); yos=0.2*(ypos-cy)
             plt.plot([0,xpos],[0,ypos],'co-')
             ax.annotate('{:.3f} Hz'.format(f), xy=(xpos,ypos),  xycoords='data',
                         xytext=(xpos+xos,ypos+yos), textcoords='data', #textcoords='offset points',
                         arrowprops=dict(arrowstyle="->", shrinkA=0, shrinkB=10)
                         )
     #plt.xlim(0,0.16)
     #plt.ylim(-0.1,0.1)
     plt.axis('equal')
     if self.ZorY == 'Z':
         plt.xlabel(r'resistance $R$ in Ohm'); plt.ylabel(r'reactance $X$ in Ohm')
     if self.ZorY == 'Y':
         plt.xlabel(r'conductance $G$ in Siemens'); plt.ylabel(r'susceptance $B$ in Siemens')
     plt.title("fitting the admittance circle with Powell's method")
     tx1='best fit (fmin_powell):\n'
     tx1+='center at G+iB = {:.5f} + i*{:.8f}\n'.format(cx,cy)
     tx1+='radius = {:.5f};  '.format(r)
     tx1+='residue: {:.2e}'.format(self.resid)
     txt1=plt.text(-r,cy-1.1*r,tx1,fontsize=8,ha='left',va='top')
     txt1.set_bbox(dict(facecolor='gray', alpha=0.25))
     idxlist=self.to_be_annotated('triple')
     ofs=self.annotation_offsets(idxlist,factor=0.1,xshift=0.15)
     for i,j in enumerate(idxlist):
         xpos,ypos = x[j],y[j]; xos,yos = ofs[i].real,ofs[i].imag
         ax.annotate('{:.1f} Hz'.format(self.f[j]), xy=(xpos,ypos),  xycoords='data',
                     xytext=(xpos+xos,ypos+yos), textcoords='data', #textcoords='offset points',
                     arrowprops=dict(arrowstyle="->", shrinkA=0, shrinkB=10)
                     )
     if self.show: plt.show()
     else: plt.savefig(join(self.sdc.plotpath,'c{}_fitted_{}_circle'.format(self.sdc.case,self.ZorY)+suffix+'.png'), dpi=240)
     plt.close()

    def test_dep(self):
        xf = arange(0, 425)
        deps = self.fm.get_dep(xf)
        plt.plot(xf,deps)

        plt.gca().set_yticks(self.fm.DEP)
        plt.gca().set_xticks(self.fm.Y_PC)
        
        plt.grid('on')
        plt.title('Ground x versus depth')
        plt.xlabel('Ground X (km)')
        plt.ylabel('depth (km)')
        plt.axis('equal')
        plt.gca().invert_yaxis()
        plt.savefig(join(self.outs_dir, '~Y_PC_vs_deps.png'))
        plt.close()

def plot(site):
    tp = np.loadtxt('../post_offsets/%s.post'%site)

    t = dc.asmjd([ii[0] for ii in tp]) + dc.adjust_mjd_for_plot_date
    e = [ii[1] for ii in tp]
    n = [ii[2] for ii in tp]
    u = [ii[3] for ii in tp]

    plt.plot_date(t,e,'x-', label = 'eastings')
    plt.plot(t,n,'x-', label = 'northings')
    plt.plot(t,u,'x-', label = 'upings')
    plt.gcf().autofmt_xdate()
    plt.legend(loc=0)
    plt.title(site)
    plt.savefig('%s.png'%site)
    #plt.show()
    plt.close()

def plot_baseline(data, plate_name, save_folder = r'Figures/'):
    """
    """
    colors = ((0.2, 0.2, 0.2),
              (0.5, 0.5, 0.5),
              (0.7, 0.7, 0.7),
              (0.3, 0.3, 0.3))

    names = data.keys()
    names.sort()
    fig, axs = plt.subplots(figsize=(8,3))
    for index, name in enumerate(names):
        for value in data[name]['original_data']:
            plot_color = colors[index % len(colors)]
            
            if abs(value - data[name]['mean'][0]) > data[name]['std'][0] * 2.0:
                axs.plot([value], [index], 'ko', markerfacecolor = [1,1,1])
            else:
                axs.plot([value], [index], 'ko', color = plot_color)

        axs.plot([data[name]['mean'][0] for _ in xrange(2)],
                 [index-0.25, index+0.25],
                  'k-')
        axs.plot([data[name]['mean'][0] - data[name]['std'][0] for _ in xrange(2)],
                 [index-0.25, index+0.25],
                 'k--')
        axs.plot([data[name]['mean'][0] + data[name]['std'][0] for _ in xrange(2)],
                 [index-0.25, index+0.25],
                 'k--')

    plt.yticks([i for i in xrange(len(names))], names, size = 10)
    plt.title(plate_name)
    plt.ylim(-0.5,len(names)-0.5)
    plt.xlabel('Fluorescent intensity')
    plt.tight_layout()

    save_filename = save_folder + 'baseline_average'

    pdf = PdfPages(save_filename.split('.')[0] + '.pdf')
    pdf.savefig(fig)
    pdf.close()

    plt.savefig(save_filename)
    #
    return None

def convolve(arrays, melBank, genere, filter_idx):
  x = []
  melBank_time = np.fft.ifft(melBank) #need to transform melBank to time domain
  for eachClip in arrays:
    result = np.convolve(eachClip, melBank_time)
    x.append(result)
    plotBeforeAfterFilter(eachClip, melBank, melBank_time, result, genere, filter_idx)

  m = np.asmatrix(np.array(x))
  fig, ax = plt.subplots()
  ax.matshow(m.real) #each element has imaginary part. So just plot real part
  plt.axis('equal')
  plt.axis('tight')
  plt.title(genere)
  plt.tight_layout()
  # filename = "./figures/convolution/Convolution_"+"Filter"+str(filter_idx)+genere+".png"
  # plt.savefig(filename)
  plt.show()

def dynamic_img_show(img,title_str='',fig_size=[14,8],hide_axes=True):
    '''Show image <img>. If called repeatedly within a cycle will dynamically redraw image.
    #DEMO
    import time

    for i in range(10):
        img = np.zeros([50,50])
        img[:i*5]=1
        dynamic_img_show(img,'iter=%s'%i)
        time.sleep(0.1)
    '''
    plt.clf()
    plt.title(title_str)
    plt.imshow(img)
    plt.xticks([]); plt.yticks([]);
    plt.gcf().set_size_inches(fig_size)
    display.display(plt.gcf())
    display.clear_output(wait=True)   

	def plot(self):
		Zdelta = []
		for aa in self.Z:
			Zdelta += [aa-1]
		figure()
		imshow(Zdelta, interpolation='bilinear', origin='lower',
	                	cmap=cm.bone, extent=(self.Vzs, self.Vze, self.V1s, self.V1e))
		CS = contour(self.X, self.Y, Zdelta, [0], linewidths=4, colors='white')
		CS2 = contour(self.X, self.Y, Zdelta, 16, linewidths=1, colors='k')
		plt.clabel(CS2, fontsize=6, inline=1)
		plt.clabel(CS, fontsize=9, inline=1)
		if self.mode == 'stab':
			plt.title('Radial stability map %s V\n %s\n alpha=%s\n b=%s' % (str(self.setup.ener), str(self.setup.V), str(self.setup.alpha), str(self.setup.b)))
		else:
			plt.title('Synchronization stability map %s V\n %s\n alpha=%s\n b=%s' % (str(self.setup.ener), str(self.setup.V), str(self.setup.alpha), str(self.setup.b)))
		xlabel('Vz (V)')
		ylabel('V1 (V)')
		show()

def plotter(mode,Bc,Tc,Q):
    col = ['#000080','#0000FF','#4169E1','#6495ED','#00BFFF','#B0E0E6']
    plt.figure()
    ax = plt.subplot(111)
    for p in range(Bc.shape[1]):
        plt.plot(Tc[:,p],Bc[:,p],'-',color=str(col[p]))
    plt.xlabel('Tc [TW]')
    plt.ylabel('Bc normalised to total EU load')
    plt.title(str(mode)+' flow')
    
    # Shrink current axis by 25% to make room for legend
    box = ax.get_position()
    ax.set_position([box.x0, box.y0, box.width * 0.75, box.height])

    plt.legend(\
        ([str(Q[i]*100) for i in range(len(Q))]),\
        loc='center left', bbox_to_anchor=(1, 0.5),title='Quantiles')
    
    plt.savefig('figures/bctc_'+str(mode)+'.eps')