def show_or_save(plt, fig, use_x11, filename):
	if use_x11:
		plt.show()
	else:
		pp = PdfPages(filename)
		pp.savefig(fig)
		pp.close()

def main():
  data = scipy.io.loadmat('data.mat')
  x1 = data['x1'][0]
  x2 = data['x2'][0]
  n = len(x1)
  kl = [1, 7, 14, 28] # k = 14 may be an optimal
  x = np.arange(-6, 6.05, 0.05)

  fig = plt.figure()
  plt.rcParams['font.size'] = 10
  for i in range(len(kl)):
    k = kl[i]
    p1 = np.zeros(len(x))
    p2 = np.zeros(len(x))
    for j in range(len(x)):
      r1 = sorted(abs(x1 - x[j]))
      r2 = sorted(abs(x2 - x[j]))
      p1[j] = float(k) / (n * 2 * r1[k-1])
      p2[j] = float(k) / (n * 2 * r2[k-1])
    plt.subplot(2, 2, i+1)
    plt.plot(x, p1, label=r'$p(\mathbf{x} \mid c_1)$')
    plt.plot(x, p2, label=r'$p(\mathbf{x} \mid c_2)$')
    plt.legend(framealpha=0, fontsize=7)
    plt.title(r'$k = %d$' % k)
    plt.xlabel(r'$x$')
    plt.ylabel(r'$p(\mathbf{x} \mid c_i)$')
  plt.tight_layout()
  pp = PdfPages('knn.pdf')
  pp.savefig(fig)
  pp.close()
  plt.clf()

def make_lick_individual(targetSN, w1, w2):
    """ Make maps for the kinematics. """
    filename = "lick_corr_sn{0}.tsv".format(targetSN)
    binimg = pf.getdata("voronoi_sn{0}_w{1}_{2}.fits".format(targetSN, w1, w2))
    intens = "collapsed_w{0}_{1}.fits".format(w1, w2)
    extent = calc_extent(intens)
    bins = np.loadtxt(filename, usecols=(0,), dtype=str).tolist()
    bins = np.array([x.split("bin")[1] for x in bins]).astype(int)
    data = np.loadtxt(filename, usecols=np.arange(25)+1).T
    labels = [r'Hd$_A$', r'Hd$_F$', r'CN$_1$', r'CN$_2$', r'Ca4227', r'G4300',
             r'Hg$_A$', r'Hg$_F$', r'Fe4383', r'Ca4455', r'Fe4531', r'C4668',
             r'H$_\beta$', r'Fe5015', r'Mg$_1$', r'Mg$_2$', r'Mg$_b$', r'Fe5270',
             r'Fe5335', r'Fe5406', r'Fe5709', r'Fe5782', r'Na$_D$', r'TiO$_1$',
             r'TiO$_2$']
    mag = "[mag]"
    ang = "[\AA]"
    units = [ang, ang, mag, mag, ang, ang,
             ang, ang, ang, ang, ang, ang,
             ang, ang, mag, mag, ang, ang,
             ang, ang, ang, ang, ang, mag,
             mag]
    lims = [[None, None], [None, None], [None, None], [None, None],
            [None, None], [None, None], [None, None], [None, None],
            [None, None], [None, None], [None, None], [None, None],
            [None, None], [None, None], [None, None], [None, None],
            [None, None], [None, None], [None, None], [None, None],
            [None, None], [None, None], [None, None], [None, None],
            [None, None], [None, None], [None, None], [None, None]]
    pdf = PdfPages("figs/lick_sn{0}.pdf".format(targetSN))
    fig = plt.figure(1, figsize=(6.25,5))
    plt.subplots_adjust(bottom=0.12, right=0.97, left=0.09, top=0.96)
    plt.minorticks_on()
    ax = plt.subplot(111)
    ax.minorticks_on()
    plot_indices = np.arange(12,22)
    for i, vector in enumerate(data):
        if i not in plot_indices:
            continue
        print "Making plot for {0}...".format(labels[i])
        kmap = np.zeros_like(binimg)
        kmap[:] = np.nan
        for bin,v in zip(bins, vector):
            idx = np.where(binimg == bin)
            kmap[idx] = v
        vmin = lims[i][0] if lims[i][0] else np.median(vector) - 2 * vector.std()
        vmax = lims[i][1] if lims[i][1] else np.median(vector) + 2 * vector.std()
        m = plt.imshow(kmap, cmap="inferno", origin="bottom", vmin=vmin,
                   vmax=vmax, extent=extent, aspect="equal")
        make_contours()
        plt.minorticks_on()
        plt.xlabel("X [kpc]")
        plt.ylabel("Y [kpc]")
        plt.xlim(extent[0], extent[1])
        plt.ylim(extent[2], extent[3])
        cbar = plt.colorbar(m)
        cbar.set_label("{0} {1}".format(labels[i], units[i]))
        pdf.savefig()
        plt.clf()
    pdf.close()
    return

    def plot_data_comb_2D(self, results_path, file_n, data, fit, timepoints):

        pp = PdfPages(results_path+'/'+file_n)
        cc = 0
        for tp in timepoints:
            xmin, xmax = -3, 3
            ymin, ymax = -3, 3

            xx, yy = mgrid[xmin:xmax:100j, ymin:ymax:100j]
            positions = vstack([xx.ravel(), yy.ravel()])
            values = vstack([ log10(1+data[tp][:, 0]), log10(1+data[tp][:, 1])])
            kernel = st.gaussian_kde(values)
            f = reshape(kernel(positions).T, xx.shape)

            xxf, yyf = mgrid[xmin:xmax:100j, ymin:ymax:100j]
            positions_f = vstack([xxf.ravel(), yyf.ravel()])
            values_f = vstack([log10(1+fit[tp][:, 0]), log10(1+fit[tp][:, 1])])
            kernel_f = st.gaussian_kde(values_f)
            ff = reshape(kernel_f(positions_f).T, xxf.shape)

            ax = plt.subplot(4, 5, cc+1)
            ax.contourf(xx, yy, f, cmap='Blues')
            ax.contourf(xxf, yyf, ff, cmap='Reds')

            ax.set_xlim([-1, 3])
            ax.set_ylim([-1, 3])
            cc += 1
        pp.savefig()
        plt.close()
        pp.close()

def print_pdf_graph(file_f, regulon, conn):
  pdf = PdfPages(file_f)
  edgesLimits = [50, 100, 200, 300, 400, 500, 600, 700, 800, 900, 1000]
  #CRP = regulon_set['LexA']
  for lim in edgesLimits:
    print lim
    g = buildSimilarityGraph_top_10_v2(conn, lim)

    # Here the node is motif, eg 87878787_1, the first 8 digits represent gi
    node_color = [ 1 if node[0:8] in regulon else 0 for node in g ]

    pos = nx.graphviz_layout(g, prog="neato")
    plt.figure(figsize=(10.0, 10.0))
    plt.axis("off")
    nx.draw(g,
        pos,
        node_color = node_color,
        node_size = 20,
        alpha=0.8,
        with_labels=False,
        cmap=plt.cm.jet,
        vmax=1.0,
        vmin=0.0
        )
    pdf.savefig()
    plt.close()

  pdf.close()

    def pdfdiagnostics(self,what='specs',n_subplot = 5):
        print 'creating a diagnostic pdf of '+what
        from matplotlib.backends.backend_pdf import PdfPages
        exec('data = self.%s'%what)

        data.sort_index(axis=1,inplace=True)# arrange alphabetically
        pp = PdfPages('%s.pdf'%self.group)
        
        for i in xrange(0, len(data.columns), n_subplot+1):
            Axes = data[data.columns[i:i+n_subplot]].plot(subplots=True)  
            tick_params(labelsize=6)

            #y ticklabels
            [setp(item.yaxis.get_majorticklabels(), 'size', 7) for item in Axes.ravel()]
            #x ticklabels
            [setp(item.xaxis.get_majorticklabels(), 'size', 5) for item in Axes.ravel()]
            #y labels
            [setp(item.yaxis.get_label(), 'size', 10) for item in Axes.ravel()]
            #x labels
            [setp(item.xaxis.get_label(), 'size', 10) for item in Axes.ravel()]

            tight_layout() 
            ylabel('mix ratio')

            #plt.locator_params(axis='y',nbins=2)
            print '%.03f'%(float(i) / float(len(data.columns)) ) , '% done'  
            savefig(pp, format='pdf')
            close('all') 
                            
        pp.close()
        print 'PDF out'
        close('all')     

    def make_comp_plot_1D(self, results_path, file_n, data, sims, timepoints, ind=0):

        pp = PdfPages(results_path+'/'+file_n)
        cc = 0
        xmin, xmax = -1, 7
        x_grid = linspace(xmin, xmax, 1000)

        def kernel_est(d, ind, x_grid):
            dl = log10(1+d[:, ind])
            #dl[isneginf(dl)] = 0
            dl = dl[isfinite(dl)]
            kde = st.gaussian_kde(dl, bw_method=0.2)
            pdf = kde.evaluate(x_grid)
            return pdf

        for tp in timepoints:

            pdf_data = kernel_est(data[tp], ind, x_grid)
            pdf_sim = kernel_est(sims[tp], ind, x_grid)
            ax = plt.subplot(4, 5, cc + 1)
            ax.plot(x_grid, pdf_data, color='blue', alpha=0.5, lw=3)
            ax.plot(x_grid, pdf_sim, color='red', alpha=0.5, lw=3)
            cc += 1

        pp.savefig()
        plt.close()
        pp.close()

    def write_pressure_graph(self):
        P.xlabel("Temperature")
        P.ylabel("Pressure")
        P.title("Pressure per Temperature")
        P.axis([0.0,self.temp_points[-1]+10., 0.0,self.pressure_points[-1]+1])
        ax = P.gca()
        ax.set_autoscale_on(False)

        popt,pcov = curve_fit(fit,self.temp_points,self.pressure_points)
        y_fit = [popt[0]*x+popt[1] for x in self.temp_points]
        y_fit.insert(0, popt[1])
        fit_x_points = self.temp_points[:]
        fit_x_points.insert(0,0.)

        pp = PdfPages(str(SCREEN_SIZE)+"_"+str(N)+".pdf")
        P.plot(self.temp_points, self.pressure_points, "o", fit_x_points, y_fit, "--")
        #P.savefig()
        #P.plot(
        #pp.savefig()
        #print(self.pressure_points)
        #print(y_fit)
        #print(fit_x_points)
        #print(y_fit)
        pp.savefig()
        pp.close()

	def compare_board_estimations(esti_extrinsics, board, board_dim, \
								actual_boards, save_name=None):
		"""
		Plots true and estimated boards on the same figure
		Args:
			esti_extrinsics: dictionary, keyed by image number, values are Extrinsics
			board:
			board_dim: (board_width, board_height)
			actual_boards: list of dictionaries
			save_name: filename, string
		"""
		if save_name:
			pp = PdfPages(save_name)
		plt.clf()

		for i in xrange(len(actual_boards)):
			fig = plt.figure()
			ax = fig.add_subplot(111, projection='3d')

			act_board = actual_boards[i]
			aX, aY, aZ = util.board_dict2array(act_board, board_dim)
			ax.plot_wireframe(aX, aY, aZ, color='b')

			if i in esti_extrinsics:
				esti_loc = esti_extrinsics[i].trans_vec
				esti_board = util.move_board(board, esti_loc)
				eX, eY, eZ = util.board_dict2array(esti_board, board_dim)
				ax.plot_wireframe(eX, eY, eZ, color='r')

			if pp:
				pp.savefig()
			else:
				plt.show()
		if pp:
			pp.close()

def readCurvesFromFileCallback():

    # Ask for input file
    _filename = askForInputFile(fileFilter="*.cur")
    if _filename == "":
        return

    # Load data from file
    [_readCurves, _readVoltages] = loadDataFromFile(_filename)
    for i in range(0, len(_readVoltages)):
        _voltages = []
        _currents = []
        for j in range(0, len(_readCurves[i])):
            _voltages.extend([_readCurves[i][j][0]])
            _currents.extend([_readCurves[i][j][1]])

    # Plot read data
    onlyFilename = _filename.split("/")[-1]
    plotCurves(_readCurves, _readVoltages, onlyFilename, interpolate=True)

    # Ask if wants to save plot to PDF file
    d = YesNoDialog(rootWindowHandler, 'Save plot to PDF file?', 'Yes', 'No')
    rootWindowHandler.wait_window(d.top)
    if not yesNoReturnedValue:
        return
    _fileTypes = '*.*'
    _filename = askForOutputFilename(_fileTypes)
    if _filename == "":
        return
    _filename = fixExtensionOfFilename(_filename, 'pdf')
    pp = PdfPages(_filename)
    plot.savefig(pp, format='pdf')
    pp.close()
    print 'Saved curves to PDF file: "%s"' % _filename

 def complexAll(self, f1=0., f2=0., amax=.16, nrows=1, ncols=1, antList=allAnts ) :
   pyplot.ioff()
   pp = PdfPages( 'ComplexLeaks.pdf' )
   scale = 10./math.sqrt(ncols*nrows)
   if f1 == 0. :
     [f1, f2] = LkSet.xlimits( self )          # default is to find freq limits in the data
   print "frequency limits: %.3f - %.3f GHz" % (f1,f2)
   ymin = -1.*amax
   ymax = amax
   npanel = 0
   for ant in antList :
     npanel = npanel + 1
     if npanel > nrows * ncols :
       npanel = 1
       pyplot.clf()
     p = pyplot.subplot(nrows, ncols, npanel, aspect='equal')    # DL,DR in one panel
     p.tick_params( axis='both', which='major', labelsize=scale )
     p.axis( [ymin, ymax, ymin, ymax] )
     p.grid(True)
     for Leak in self.LeakList :
       if Leak.ant == ant :
         print "plotting DR and DL for antenna %d" % ant
         Leak.plotComplex( p, f1, f2 ) 
     #pyplot.title("C%d DR (circles, solid) and DL (diamonds, dashed)" % ant, fontdict={'fontsize': scale})
     if (npanel == nrows*ncols) or (ant == antList[-1] ) :
       pyplot.savefig( pp, format='pdf' )
   pp.close()

def plot_miri_comparison():

    inst = webbpsf.MIRI()
    filtlist_W = [f for f in inst.filter_list if f[-1] == 'W']
    filtlist_C = [f for f in inst.filter_list if f[-1] != 'W']

    from matplotlib.backends.backend_pdf import PdfPages
    pdf=PdfPages('weights_miri_comparison.pdf')


    for filts in [filtlist_W, filtlist_C]:

        try:
            os.unlink('/Users/mperrin/software/webbpsf/data/MIRI/filters')
        except: 
            pass
        os.symlink('/Users/mperrin/software/webbpsf/data/MIRI/real_filters', '/Users/mperrin/software/webbpsf/data/MIRI/filters')
        plotweights('miri', filtlist=filts)

        os.unlink('/Users/mperrin/software/webbpsf/data/MIRI/filters')
        os.symlink('/Users/mperrin/software/webbpsf/data/MIRI/fake_filters', '/Users/mperrin/software/webbpsf/data/MIRI/filters')
        plotweights('miri', filtlist=filts, overplot=True, ls='--')
        P.draw()
        pdf.savefig()

    pdf.close()

def multipage(filename, figs=None, dpi=200):
    pp = PdfPages(filename)
    if figs is None:
        figs = [plt.figure(n) for n in plt.get_fignums()]
    for fig in figs:
        fig.savefig(pp, format='pdf')
    pp.close()

def heat_map_single(data, file = "heat_map_plate.pdf", *args, **kwargs):
    """ Create a heat_map for a single readout

    Create a heat_map for a single readout

    ..todo:: Share code between heat_map_single and heat_map_multiple
    """

    np_data = data.data
    pp = PdfPages(os.path.join(PATH, file))

    fig, ax = plt.subplots()

    im = ax.pcolormesh(np_data, vmin=np_data.min(), vmax=np_data.max()) # cmap='RdBu'
    fig.colorbar(im)

    # put the major ticks at the middle of each cell
    ax.set_xticks(np.arange(np_data.shape[1]) + 0.5, minor=False)
    ax.set_yticks(np.arange(np_data.shape[0]) + 0.5, minor=False)

    # Invert the y-axis such that the data is displayed as it appears on the plate.
    ax.invert_yaxis()
    ax.xaxis.tick_top()

    ax.set_xticklabels(data.axes['x'], minor=False)
    ax.set_yticklabels(data.axes['y'], minor=False)

    pp.savefig(fig)
    pp.close()
    fig.clear()

    return ax

def anal2pdf():


    pp = PdfPages('../../datafiles/jul14/analplots.pdf')
   
    
    fnames = [
    '../../datafiles/jul14/vsweep_10_1.h5',
    '../../datafiles/jul14/vsweep_10_1b.h5',
    '../../datafiles/jul14/vsweep_10_2.h5',
    '../../datafiles/jul14/vsweep_10_3.h5',
    '../../datafiles/jul14/vsweep_10_4.h5',
    '../../datafiles/jul14/vsweep_10_5.h5',
    '../../datafiles/jul14/vsweep_10_6.h5']
   
    
    for fn in fnames:
       anal_vsweep(fn)
       figure(1)
       suptitle('Voltage sweep, 5096MHz, raw phase(Y), time_samples(X)')
       f=gcf()
       f.savefig(pp,format='pdf')
    
       figure(3)
       suptitle('Voltage sweep, 5096MHz, radians(Y) vs mV(X)')
       f=gcf()
       f.savefig(pp,format='pdf')
    
    pp.close()

def plotLoad(dataFolder, srv, TS, trace_type, tsInterval):
	srv_file = srv + "_" + TS + "_" + trace_type + ".json"
	srv_load = json.load(open(dataFolder + srv_file))

	fig, ax = plt.subplots()
	print "Ploting QoE evaluation for server :", srv
	ts = [int(x) for x in srv_load.keys() if int(x) > tsInterval[0] and int(x) < tsInterval[1]]
	sorted_ts = sorted(ts)
	tr_vals = [srv_load[str(cur_ts)] for cur_ts in sorted_ts]

	srvName = '$S_{' + str(int(srv.split('-')[1])) + '}$'
	plt.plot(sorted_ts, tr_vals, 'b-', label=srvName, linewidth=2.0, markersize=8)


	## Change the time stamp ticks
	num_intvs = int((tsInterval[1] - tsInterval[0])/900) + 1
	ts_labels = [tsInterval[0] + x*900 for x in range(num_intvs)]
	str_ts = [datetime.datetime.fromtimestamp(x*900 + tsInterval[0]).strftime('%H:%M') for x in range(num_intvs)]
	plt.xticks(ts_labels, str_ts, fontsize=15)

	box = ax.get_position()
	#ax.legend(loc='lower center', bbox_to_anchor=(0.5, 0),
	#	          fancybox=True, shadow=True, ncol=4, prop={'size':20})
	#params = {'legend.fontsize': 20,
	#          'legend.linewidth': 2}
	#plt.rcParams.update(params)
	# ax.set_title('Server QoE Score Observed on $S_{10}$', fontsize=20)
	# ax.set_xlabel("Time in a day", fontsize=20)
	ax.set_ylabel("Server Load", fontsize=20)
	# ax.set_ylim([0,5])
	plt.show()

	pdf = PdfPages(dataFolder + '/imgs/' + cache_agent + '_' + 'load.pdf')
	pdf.savefig(fig)
	pdf.close()

def plot_distributions(distributions, bucket_pct, axes, out_file):
	x_max, y_max = axes
	pp = PdfPages(out_file)
	variances = sorted(distributions)
	subsamples = sorted(distributions[variances[-1]])
	if x_max == 0.:
		x_max = distributions[0][-1]
	bucket_width = x_max * bucket_pct
	bucket_boundaries = np.arange(0, x_max + bucket_width / 2., bucket_width)
	x_axis_points = np.arange(bucket_width / 2., x_max, bucket_width)
	for i,v in enumerate(variances):
		plt.figure(i)
		plt.xlabel("regret distribution")
		if v == 0:
			plt.title("true game")
			cum_dist = np.array([bisect(distributions[0], b) for b in \
								bucket_boundaries])
			plt.plot(x_axis_points, (cum_dist[1:] - cum_dist[:-1]) / \
					float(cum_dist[-1]), label="true game")
		else:
			plt.title("$\sigma \\approx$" +str(v))
			for s in subsamples:
				cum_dist = np.array([bisect(distributions[v][s], b) for b in \
									bucket_boundaries])
				plt.plot(x_axis_points, (cum_dist[1:] - cum_dist[:-1]) / \
						float(cum_dist[-1]), label=str(s)+" samples")
		plt.legend(loc="upper right", prop={'size':6})
		if y_max != 0.:
			plt.axis([0, x_max, 0, y_max])
		pp.savefig()
	pp.close()

def plot_pnas_vs_tics(pnas_dir, tic_dir, pnas_names, directory, scale = 7.14, refcoords_file = None):
	pnas = np.concatenate(load_file(pnas_dir))
	pnas[:,0] *= scale
	print(np.shape(pnas))
	print(len(pnas_names))
	if("ktICA" in tic_dir):
		tics = load_dataset(tic_dir)
	else:
		tics = verboseload(tic_dir)
	print(np.shape(tics))
	tics = np.concatenate(tics)
	print(np.shape(tics))
	if len(pnas_names) != np.shape(pnas)[1]:
		print("Invalid pnas names")
		return

	for i in range(0,np.shape(pnas)[1]):
		for j in range(0,np.shape(tics)[1]):
			tic = tics[:,j]
			pnas_coord = pnas[:,i]
			plt.hexbin(tic, pnas_coord, bins = 'log', mincnt=1)
			coord_name = pnas_names[i]
			tic_name = "tIC.%d" %(j+1)
			plt.xlabel(tic_name)
			plt.ylabel(coord_name)
			pp = PdfPages("%s/%s_%s_hexbin.pdf" %(directory, tic_name, coord_name))
			pp.savefig()
			pp.close()
			plt.clf()

	return

    def plot_data_dict_1D(self, results_path, file_n, data, timepoints):
        print 'plotting now...'
        pp = PdfPages(results_path+'/'+file_n)
        #nBins = 50
        cc = 0
        xmin, xmax = -1, 7
        x_grid = linspace(xmin, xmax, 1000)

        for tp in timepoints:
            dat = log10(1+data[tp][:,0])
            dat[isneginf(dat)] = 0
            print dat
            kde = st.gaussian_kde(dat, bw_method=0.2)
            pdf = kde.evaluate(x_grid)

            ax = plt.subplot(4, 5, cc + 1)
            ax.plot(x_grid, pdf, color='blue', alpha=0.5, lw=3)
            ax.set_xlim([-1, 7])

            #plt.hist( data[tp], nBins )
            #plt.xscale('log')
            cc += 1
        pp.savefig()
        plt.close()
        pp.close()

def plot_tica(transformed_data_dir, lag_time):
	transformed_data = verboseload(transformed_data_dir)
	trajs = np.concatenate(transformed_data)
	plt.hexbin(trajs[:,0], trajs[:,1], bins='log', mincnt=1)
	pp = PdfPages("/scratch/users/enf/b2ar_analysis/tica_phi_psi_chi2_t%d.pdf" %lag_time)
	pp.savefig()
	pp.close()