def plotP2():
  '''A function to plot the two-angle Normalized Scattering
  Phase Function. For comparison see Myneni 1988c fig.1. 
  Output: plot of P2 function
  '''
  N = 32
  g_wt = np.array(gauss_wt[str(N)])
  mu_l = np.array(gauss_mu[str(N)])
  zen = np.arccos(mu_l)
  a = 180.*np.pi/180. # view azimuth
  view = []
  for z in zen:
    view.append((z,a))
  sun = (170./180.*np.pi,0./180.*np.pi) # sun zenith, azimuth
  arch = 'u'
  refl = 0.5
  trans = 0.5
  y = []
  for v in view:
    y.append(P2(v, sun, arch, refl, trans))
  fig, ax = plt.subplots(1)
  plt.plot(mu_l,y, 'r--')
  s = '''sun zenith:%.2f, sun azimuth:%.2f, arch:%s
  view azimuth:%.2f, refl:%.2f, trans:%.2f''' % \
      (sun[0]*180/np.pi,sun[1]*180/np.pi, arch,\
      a*180/np.pi, refl, trans)
  props = dict(boxstyle='round',facecolor='wheat',alpha=0.5)
  plt.text(.5,.5, s, bbox = props, transform=ax.transAxes,\
      horizontalalignment='center', verticalalignment='center')
  plt.xlabel(r"$\mu$ (cosine of exit zenith)")
  plt.ylabel(r"P($\Omega$, $\Omega$0)")
  plt.title("P (Normalized Scattering Phase Function)") 
  plt.show()

def _fig_density(sweight, surweight, pval, nlm):
    """
    Plot the histogram of sweight across the image
    and the thresholds implied by the surrogate model (surweight)
    """
    import matplotlib.pylab as mp
    # compute some thresholds
    nlm = nlm.astype('d')
    srweight = np.sum(surweight,1)
    srw = np.sort(srweight)
    nitem = np.size(srweight)
    thf = srw[int((1-min(pval,1))*nitem)]
    mnlm = max(1,nlm.mean())
    imin = min(nitem-1,int((1.-pval/mnlm)*nitem))
    
    thcf = srw[imin]
    h,c = np.histogram(sweight,100)
    I = h.sum()*(c[1]-c[0])
    h = h/I
    h0,c0 = np.histogram(srweight,100)
    I0 = h0.sum()*(c0[1]-c0[0])
    h0 = h0/I0
    mp.figure(1)
    mp.plot(c,h)
    mp.plot(c0,h0)
    mp.legend(('true histogram','surrogate histogram'))
    mp.plot([thf,thf],[0,0.8*h0.max()])
    mp.text(thf,0.8*h0.max(),'p<0.2, uncorrected')
    mp.plot([thcf,thcf],[0,0.5*h0.max()])
    mp.text(thcf,0.5*h0.max(),'p<0.05, corrected')
    mp.savefig('/tmp/histo_density.eps')
    mp.show()

def plot_confusion_matrix(cm, title='', cmap=plt.cm.Blues):
    #print cm
    #display vehicle, idle, walking accuracy respectively
    #display overall accuracy
    print type(cm)
   # plt.figure(index
    plt.imshow(cm, interpolation='nearest', cmap=cmap)
    #plt.figure("")
    plt.title("Confusion Matrix")
    plt.colorbar()
    tick_marks = [0,1,2]
    target_name = ["driving","idling","walking"]


    plt.xticks(tick_marks,target_name,rotation=45)

    plt.yticks(tick_marks,target_name,rotation=45)
    print len(cm[0])

    for i in range(0,3):
        for j in range(0,3):
         plt.text(i,j,str(cm[i,j]))
    plt.tight_layout()
    plt.ylabel("Actual Value")
    plt.xlabel("Predicted Outcome")

 def double_label(self, labels, position="horizontal", **kwargs):
     """
     Plots 2 colors in a square:
         
                   vertical
                +-------------+
                |            /|
                |   1st    /  |
                |        /    |
     horizontal |      /      | horizontal
                |    /        |
                |  /   2nd    |
                |/            |
                +-------------+
                    vertical
     """
     assert len(labels) == 2
     if position == 'horizontal':
         x1 = self.x - self.dx*0.25
         x2 = self.x + self.dx*1.25
         y1 = y2 = self.y + self.dy/2
         ha1, ha2 = 'right', 'left'
         va1 = va2 = 'center'
     elif position == 'vertical':
         x1 = x2 = self.x + self.dx/2
         y1 = self.y + self.dy*1.25
         y2 = self.y - self.dy*0.25
         ha1 = ha2 = 'center'
         va1, va2 = 'bottom', 'top'
     else:
         raise ValueError("`position` must be one of:"
                 "'horizontal', 'vertical'")
     fontdict = kwargs.get('fontdict', {})
     plt.text(x1, y1, labels[0], ha=ha1, va=va1, fontdict=fontdict) 
     plt.text(x2, y2, labels[1], ha=ha2, va=va2, fontdict=fontdict) 

    def plot_histogram(self, main="", numrows=1, numcols=1, fignum=1):
        """Plot a histogram of choices and probability sums. Expects probabilities as (at least) a 2D array.
        """
        from matplotlib.pylab import bar, xticks, yticks, title, text, axis, figure, subplot

        probabilities = self.get_probabilities()
        if probabilities.ndim < 2:
            raise StandardError, "probabilities must have at least 2 dimensions."
        alts = probabilities.shape[1]
        width_par = (1 / alts + 1) / 2.0
        choice_counts = self.get_choice_histogram(0, alts)
        sum_probs = self.get_probabilities_sum()

        subplot(numrows, numcols, fignum)
        bar(arange(alts), choice_counts, width=width_par)
        bar(arange(alts) + width_par, sum_probs, width=width_par, color="g")
        xticks(arange(alts))
        title(main)
        Axis = axis()
        text(
            alts + 0.5,
            -0.1,
            "\nchoices histogram (blue),\nprobabilities sum (green)",
            horizontalalignment="right",
            verticalalignment="top",
        )

    def plot_biggest(self,n_plot):
        # plots the n_plot biggest clusters
        s = self.cluster_membership.sum(0)
        order = np.argsort(s)
        
        for i in np.arange(s.size-1,s.size-n_plot-1,-1):
            cluster = order[0,i]
            peaks = np.nonzero(self.cluster_membership.getcol(cluster))[0]
            plt.figure(figsize=(8,8))
            plt.subplot(1,2,1)
            plt.plot(self.peak_data.mass[peaks],self.peak_data.rt[peaks],'ro')
            plt.plot(self.peak_data.transformed[peaks,cluster].toarray(),self.peak_data.rt[peaks],'ko')

            plt.subplot(1,2,2)
            for peak in peaks:
                plt.plot((self.peak_data.mass[peak], self.peak_data.mass[peak]),(0,self.peak_data.intensity[peak]))
                tr = self.peak_data.possible[peak,cluster]-1
                x = self.peak_data.mass[peak] + random.random()
                y = self.peak_data.intensity[peak] + random.random()
                plt.text(x,y,self.peak_data.transformations[tr].name)
                title_string = "Mean RT: " + str(self.cluster_model.cluster_rt_mean[cluster]) + "(" + \
                    str(1.0/self.cluster_model.cluster_rt_prec[cluster]) + ")"                    
                if hasattr(self.cluster_model, 'cluster_mass_mean'):
                    title_string += " Mean Mass: " + str(self.cluster_model.cluster_mass_mean[cluster]) + \
                        "(" + str(1.0/self.cluster_model.cluster_mass_prec[cluster]) + ")"
                plt.title(title_string) 
        plt.show()

def minj_mflux(*args,**kwargs):
    Qu = jana.quantities()
    Minj_List =[]
    Minj_MHD_List=[]
    for i in range(len(args[1])):
        Minj_List.append(Qu.Mflux(args[1][i],Mstar=args[0][i],scale=True)['Mfr']+Qu.Mflux(args[1][i],Mstar=args[0][i],scale=True)['Mfz'])
        
    MHD_30minj = Qu.Mflux(args[2],Mstar=30.0,scale=True)['Mfr']+Qu.Mflux(args[2],Mstar=30.0,scale=True)['Mfz']
    for i in args[0]:
        Minj_MHD_List.append(MHD_30minj*np.sqrt(i/30.0))
    
    f1 = plt.figure(num=1)
    ax1 = f1.add_subplot(211)
    plt.axis([10.0,70.0,2.0e-5,7.99e-5])
    plt.plot(args[0],Minj_MHD_List,'k*')
    plt.plot(args[0],Minj_List,'ko')
    plt.minorticks_on()
    locs,labels = plt.yticks()
    plt.yticks(locs, map(lambda x: "%.1f" % x, locs*1e5))
    plt.text(0.0, 1.03, r'$10^{-5}$', transform = plt.gca().transAxes)
    plt.xlabel(r'Stellar Mass : $M_{*} [M_{\odot}]$')
    plt.ylabel(r' $\dot{M}_{\rm vert} + \dot{M}_{\rm rad}\, [M_{\odot}\,\rm yr^{-1}]$')
    
    ax2 = f1.add_subplot(212)
    plt.axis([10.0,70.0,0.0,50.0])
    plt.plot(args[0],100*((np.array(Minj_List)-np.array(Minj_MHD_List))/np.array(Minj_MHD_List)),'k^')
    plt.minorticks_on()
    plt.xlabel(r'Stellar Mass : $M_{*} [M_{\odot}]$')
    plt.ylabel(r'$\%$ Change in Total Mass Outflow Rates')

 def segmentation(self, threshold):
     img = self.spectrogram["data"]
     mask = (img > threshold).astype(np.float)
     hist, bin_edges = np.histogram(img, bins=60)
     bin_centers = 0.5*(bin_edges[:-1] + bin_edges[1:])
     binary_img = mask > 0.5
     plt.figure(figsize=(11,8))
     plt.subplot(131)
     plt.imshow(img)
     plt.axis('off')
     plt.subplot(132)
     plt.plot(bin_centers, hist, lw=2)
     print(threshold)
     plt.axvline(threshold, color='r', ls='--', lw=2)
     plt.text(0.57, 0.8, 'histogram', fontsize=20, transform = plt.gca().transAxes)
     plt.text(0.45, 0.75, 'threshold = '+ str(threshold)[0:5], fontsize=15, transform = plt.gca().transAxes)
     plt.yticks([])
     plt.subplot(133)     
     plt.imshow(binary_img)
     plt.axis('off')
     plt.subplots_adjust(wspace=0.02, hspace=0.3, top=1, bottom=0.1, left=0, right=1)
     plt.show()
     print(img.max())
     print(binary_img.max())
     
     return mask

def plot_clustering(items, xs, ys, labels=None, texts=None, shapes=None):
    labels = [0] * len(xs) if labels is None else labels
    texts = [''] * len(xs) if texts is None else texts
    shapes = [0] * len(xs) if shapes is None else shapes
    for item, x, y, label, text, shape in zip(items, xs, ys, labels, texts, shapes):
        plt.plot(x, y, markers[shape], color=colors[label])
        plt.text(x, y, text)

def item_clustering(data, skill, cluster_number=3, plot=True):
    pk, level = data.get_skill_id(skill)
    items = data.get_items_df()
    items = items[items["skill_lvl_" + str(level)] == pk]
    items = items[items["visualization"] != "pairing"]

    corr = compute_corr(data)
    corr = pd.DataFrame(corr, index=items.index, columns=items.index)

    print("Corr ({}) contain total {} values and from that {} nans".format(corr.shape, corr.size, corr.isnull().sum().sum()))
    corr[corr.isnull()] = 0

    sc = SpectralClusterer(corr, kcut=corr.shape[0] / 2, mutual=True)
    # sc = SpectralClusterer(corr, kcut=30, mutual=True)
    labels = sc.run(cluster_number=cluster_number, KMiter=50, sc_type=2)

    if plot:
        colors = "rgbyk"
        visualizations = list(items["visualization"].unique())

        for i, p in enumerate(corr.columns):
            item = items.loc[p]
            plt.plot(sc.eig_vect[i,1], sc.eig_vect[i,2], "o", color=colors[visualizations.index(item["visualization"])])
            # plt.plot(sc.eig_vect[i, 1], sc.eig_vect[i, 2], "o", color=colors[labels[i]])
            plt.text(sc.eig_vect[i, 1], sc.eig_vect[i, 2], item["name"])

        for i, vis in enumerate(visualizations):
            plt.plot(0, 0, "o", color=colors[i], label=vis)
        plt.title(data)

        plt.legend(loc=3)

    return labels

def concept_clustering(data, skill, cluster_number=3, plot=True):
    pk, level = data.get_skill_id(skill)
    items = data.get_items_df()
    items = items[items["skill_lvl_" + str(level)] == pk]
    skills = data.get_skills_df()
    skill_ids = items[~items["skill_lvl_3"].isnull()]["skill_lvl_3"].unique()

    corr = compute_corr(data, merge_skills=True)
    corr = pd.DataFrame(corr, index=skill_ids, columns=skill_ids)
    print("Corr ({}) contain total {} values and from that {} nans".format(corr.shape, corr.size, corr.isnull().sum().sum()))
    corr[corr.isnull()] = 0

    try:
        sc = SpectralClusterer(corr, kcut=corr.shape[0] * 0.5, mutual=True)
        labels = sc.run(cluster_number=cluster_number, KMiter=50, sc_type=2)
    except np.linalg.linalg.LinAlgError:
        sc = SpectralClusterer(corr, kcut=corr.shape[0] * 0.5, mutual=False)
        labels = sc.run(cluster_number=cluster_number, KMiter=50, sc_type=2)

    if plot:
        colors = "rgbyk"
        for i, p in enumerate(corr.columns):
            skill = skills.loc[int(p)]
            plt.plot(sc.eig_vect[i, 1], sc.eig_vect[i, 2], "o", color=colors[labels[i]])
            plt.text(sc.eig_vect[i, 1], sc.eig_vect[i, 2], skill["name"])
        plt.title(data)

    return labels

    def test_fit_gaussian_image_l0(self):
        import matplotlib.pylab as plt
        from mpl_toolkits.axes_grid1.inset_locator import zoomed_inset_axes
        from mpl_toolkits.axes_grid1.inset_locator import mark_inset
        d1 = ap.catools.caget('LTB-BI:BD1{VF1}Img1:ArrayData')

        self.assertEqual(len(d1), 1220*1620)

        d2 = np.reshape(d1, (1220, 1620))
        params = ap.fitGaussianImage(d2)
        fit = self._gaussian(*params)

        plt.clf()
        extent=(500, 1620, 400, 1220)
        #plt.contour(fit(*np.indices(d2.shape)), cmap=plt.cm.copper)
        plt.imshow(d2, interpolation="nearest", cmap=plt.cm.bwr)
        ax = plt.gca()
        ax.set_xlim(750, 850)
        ax.set_ylim(550, 650)
        (height, y, x, width_y, width_x) = params
        #axins = zoomed_inset_axes(ax, 6, loc=1)
        #axins.contour(fit(*np.indices(d2.shape)), cmap=plt.cm.cool)
        #axins.set_xlim(759, 859)
        #axins.set_ylim(559, 660)
        #mark_inset(ax, axins, loc1=2, loc2=4, fc='none', ec='0.5')

        plt.text(0.95, 0.05, """
        x : %.1f
        y : %.1f
        width_x : %.1f
        width_y : %.1f""" %(x, y, width_x, width_y),
             fontsize=16, horizontalalignment='right',
             verticalalignment='bottom', transform=ax.transAxes)

        plt.savefig(figname("test2.png"))

    def plot_uv(self):
        """ Plot UV points. Just because its easy. """
        self.current_plot = 'uv'
        self.ax_zoomed = False
        
        uu = self.uv.d_uv_data['UU']*1e6
        vv = self.uv.d_uv_data['VV']*1e6
        xx = self.uv.d_array_geometry['STABXYZ'][:,0]
        yy = self.uv.d_array_geometry['STABXYZ'][:,1]
        pmax, pmin = np.max([uu, vv])*1.1, np.min([uu,vv])*1.1

        
        fig = self.sp_fig
        plt.subplot(121, aspect='equal')
        ax = plt.plot(xx, yy, 'bo')
        for i in range(len(xx)):
            plt.text(xx[i], yy[i], self.uv.d_array_geometry['ANNAME'][i].strip('Stand').strip('Tile'))

        plt.title("Antenna positions")
        plt.xlabel("X [m]")
        plt.ylabel("Y [m]")
        
        plt.subplot(122, aspect='equal')
        ax = plt.plot(uu, vv, 'b.')
        plt.title("UV data")
        plt.xlabel("UU [$\\mu s$]")
        plt.ylabel("VV [$\\mu s$]")
        plt.xlim(pmin, pmax)
        plt.ylim(pmin, pmax)
        return fig, ax

def lattice_plot(component_list, file_path):
    """
    Creates a lattice style plot of all graph components
    """
    graph_fig = plt.figure(figsize=(20, 10))  # Create figure

    # Set the number of rows in the plot based on an odd or
    # even number of components
    num_components = len(component_list)
    if num_components % 2 > 0:
        num_cols = (num_components / 2) + 1
    else:
        num_cols = num_components / 2

    # Plot subgraphs, with centrality annotation
    plot_count = 1
    for G in component_list:
        # Find actor in each component with highest degree
        in_cent = nx.degree(G)
        in_cent = [(b, a) for (a, b) in in_cent.items()]
        in_cent.sort()
        in_cent.reverse()
        high_in = in_cent[0][1]

        # Plot with annotation
        plt.subplot(2, num_cols, plot_count)
        nx.draw_spring(G, node_size=35, with_labels=False)
        plt.text(0, -0.1, "Highest degree: " + high_in, color="darkgreen")
        plot_count += 1

    plt.savefig(file_path)

def plot_prof(obj):
  '''Function that plots a vertical profile of scattering from
  TOC to BOC.
  Input: rt_layer instance object.
  Output: plot of scattering.
  '''
  I = obj.Inodes[:,1,:] \
      / -obj.mu_s
  y = np.linspace(obj.Lc, 0., obj.K+1)
  x = obj.views*180./np.pi
  xm = np.array([])
  for i in np.arange(0,len(x)-1):
    nx = x[i] + (x[i+1]-x[i])/2.
    xm = np.append(xm,nx)
  xm = np.insert(xm,0,0.)
  xm = np.append(xm,180.)
  xx, yy = np.meshgrid(xm, y)
  plt.pcolormesh(xx,yy,I)
  plt.colorbar()
  plt.title('Canopy Fluxes')
  plt.xlabel('Exit Zenith Angle')
  plt.ylabel('Cumulative LAI (0=soil)')
  plt.arrow(135.,3.5,0.,-3.,head_width=5.,head_length=.2,\
      fc='k',ec='k')
  plt.text(140.,2.5,'Downwelling Flux',rotation=90)
  plt.arrow(45.,.5,0.,3.,head_width=5.,head_length=.2,\
      fc='k',ec='k')
  plt.text(35.,2.5,'Upwelling Flux',rotation=270)
  plt.show()

def impact_turning_points(s, pivot_k=1, up_impact=0.01, down_impact=0.01, visible=False):

    import matplotlib.pylab as plt

    peaks, troughs = pivots(s, pivot_k)
    impacts_p_idx, impacts_p, impacts_t_idx, impacts_t = impact(s, peaks, troughs)

    new_peaks = []
    new_troughs = []

    for i, p in enumerate(impacts_p_idx):
        if impacts_p[i] >= down_impact :
            new_peaks.append(p)
    for i, t in enumerate(impacts_t_idx):
        if impacts_t[i] >= up_impact :
            new_troughs.append(t)
            
    peaks, troughs = alternate_pivots(s, new_peaks, new_troughs)
    
    fig = plt.figure(figsize=(8,6))
    s.plot()
    
    for p in peaks:
        plt.plot(s.index[p], s[p], 'ro', markersize=8)
        if visible == True:
            plt.text(s.index[p], s[p], str(p), fontsize=15)
    for t in troughs:
        plt.plot(s.index[t], s[t], 'go', markersize=8)
        if visible == True:
            plt.text(s.index[t], s[t], str(t), fontsize=15)

def difficulty_vs_time(data, the_skill, concepts=False):
    data.filter_data(0, 100)
    pk, level = data.get_skill_id(the_skill)
    data.trim_times()
    data.add_log_response_times()
    m = EloPriorCurrentModel(KC=2, KI=0.5)
    items = data.get_items_df()
    items = items[items["visualization"] != "pairing"]
    items = items.join(get_difficulty(data, m))
    items = items.join(pd.Series(data.get_dataframe_all().groupby(["item"])["log_response_time"].mean(), name="log_response_time_mean"))
    items = items[items["skill_lvl_"+str(level)] == pk]

    if concepts:
        skills = data.get_skills_df()
        skills = skills.join(items.groupby("skill_lvl_3")["difficulty"].mean())
        skills = skills.join(items.groupby("skill_lvl_3")["log_response_time_mean"].mean())
        skills = skills[skills.index.isin(items["skill_lvl_3"].unique())]
        for id, skill in skills.iterrows():
            plt.plot(skill["difficulty"], skill["log_response_time_mean"], "ok")
            plt.text(skill["difficulty"], skill["log_response_time_mean"], skill["name"])
    else:
        colors = "rgbyk"
        visualizations = list(items["visualization"].unique())
        for id, item in items.iterrows():
            plt.plot(item["difficulty"], item["log_response_time_mean"], "o", color=colors[visualizations.index(item["visualization"])])
            plt.text(item["difficulty"], item["log_response_time_mean"], item["name"])
        for i, vis in enumerate(visualizations):
            plt.plot(-1, 2, "o", color=colors[i], label=vis)
    plt.xlabel("difficulty according to " + str(m))
    plt.ylabel("mean of log time")
    plt.legend(loc=0)
    plt.title(the_skill)

    def test_fit_gaussian_2d(self):
        import matplotlib.pylab as plt
        # Create the gaussian data
        Xin, Yin = plt.mgrid[0:201, 0:201]
        data = self._gaussian(3, 100, 100, 20, 40)(Xin, Yin) + \
            np.random.random(Xin.shape)

        plt.clf()
        ax2 = plt.axes([0, 0.52, 1.0, 0.4])
        ax2.matshow(data, cmap=plt.cm.gist_earth_r)

        params = ap.fitGaussianImage(data)
        fit = self._gaussian(*params)

        ax2.contour(fit(*np.indices(data.shape)), cmap=plt.cm.cool)
        (height, x, y, width_x, width_y) = params

        plt.text(0.95, 0.05, """
        x : %.1f
        y : %.1f
        width_x : %.1f
        width_y : %.1f""" %(x, y, width_x, width_y),
             fontsize=16, horizontalalignment='right',
             verticalalignment='bottom', transform=ax2.transAxes)

        ax1 = plt.axes([0, 0.08, 1.0, 0.4])

def mass_flux_plot(*args,**kwargs):
    fltm = idls.read(args[0])
    injm = idls.read(args[1])
    f1 = plt.figure()

    ax1 = f1.add_subplot(211)
    plt.plot(injm.nt_sc,injm.nmf_rscale,'r')
    plt.plot(injm.nt_sc,injm.nmf_zscale,'b')
    plt.plot(injm.nt_sc,injm.nmf_z0scale,'k')
    plt.plot(injm.nt_sc,(injm.nmf_rscale+injm.nmf_zscale),'g')
    plt.axis([0.0,160.0,0.0,3.5e-5])
    plt.minorticks_on()
    locs,labels = plt.yticks()
    plt.yticks(locs, map(lambda x: "%.1f" % x, locs*1e5))
    plt.text(0.0, 1.03, r'$10^{-5}$', transform = plt.gca().transAxes)
    plt.xlabel(r'Time [yr]',labelpad=6)
    plt.ylabel(r'$\dot{\rm M}_{\rm out} [ \rm{M}_{\odot} \rm{yr}^{-1} ]$',labelpad=15)
    
    ax2 = f1.add_subplot(212)
    plt.plot(fltm.nt_sc,fltm.nmf_rscale,'r')
    plt.plot(fltm.nt_sc,fltm.nmf_zscale,'b')
    plt.plot(fltm.nt_sc,fltm.nmf_z0scale,'k')
    plt.plot(fltm.nt_sc,(fltm.nmf_rscale+fltm.nmf_zscale),'g')
    plt.axis([0.0,160.0,0.0,4.0e-5])
    plt.minorticks_on()
    locs,labels = plt.yticks()
    plt.yticks(locs, map(lambda x: "%.1f" % x, locs*1e5))
    plt.text(0.0, 1.03, r'$10^{-5}$', transform = plt.gca().transAxes)
    plt.xlabel(r'Time [yr]',labelpad=6)
    plt.ylabel(r'$\dot{\rm M}_{\rm out} [ \rm {M}_{\odot} \rm{yr}^{-1} ]$',labelpad=15)

    def draw(self, description=None, ofile="test.png"):

        plt.clf()

        for f in self.funcs:
            f()

        plt.axis("off")

        ax = plt.gca()
        ax.set_aspect("equal", "datalim")

        f = plt.gcf()
        f.set_size_inches(12.8, 7.2)

        if description is not None:
            plt.text(0.025, 0.05, description, transform=f.transFigure)

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

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

        plt.subplots_adjust(left=0, bottom=0, right=1, top=1, wspace=0, hspace=0)

        # dpi = 100 for 720p, 150 for 1080p
        plt.savefig(ofile, dpi=150)