Python matplotlib.pylab 模块，title() 实例源码

def plot_x_y_yhat(x, y, y_hat, xsz, ysz, binz=False):
    """Plot x, y and y_hat side by side."""
    plt.close("all")
    f = plt.figure(figsize=(15, 10.8), dpi=300)
    gs = gridspec.GridSpec(1, 3)
    if binz:
        y_hat = (y_hat > 0.5) * 1.
    ims = [x, y, y_hat]
    tils = [
        "x:" + str(xsz) + "x" + str(xsz),
        "y:" + str(ysz) + "x" + str(ysz),
        "yhat:" + str(ysz) + "x" + str(ysz)]
    for n, ti in zip([0, 1, 2], tils):
        f.add_subplot(gs[n])
        if n == 0:
            plt.imshow(ims[n], cmap=cm.Greys_r)
        else:
            plt.imshow(ims[n], cmap=cm.Greys_r)
        plt.title(ti)

    return f

def plot_x_x_yhat(x, x_hat):
    """Plot x, y and y_hat side by side."""
    plt.close("all")
    f = plt.figure()  # figsize=(15, 10.8), dpi=300
    gs = gridspec.GridSpec(1, 2)
    ims = [x, x_hat]
    tils = [
        "xin:" + str(x.shape[0]) + "x" + str(x.shape[1]),
        "xout:" + str(x.shape[1]) + "x" + str(x_hat.shape[1])]
    for n, ti in zip([0, 1], tils):
        f.add_subplot(gs[n])
        plt.imshow(ims[n], cmap=cm.Greys_r)
        plt.title(ti)
        ax = f.gca()
        ax.set_axis_off()

    return f

def _clear(self, title=True, xlabel=True, ylabel=True):
        '''
        Method for removing the title and/or xlabel and/or Ylabel.

        Parameters
        ----------
        Title : boolean, optional
            Boolean of if title will be cleared.  The default is True.
        xlabel : boolean, optional
            Boolean of if xlabel will be cleared.  The default is True.
        ylabel : boolean, optional
            Boolean of if ylabel will be cleared.  The default is True.

        '''
        if title:
            pyl.title('')
        if xlabel:
            pyl.xlabel('')
        if ylabel:
            pyl.ylabel('')

    # From mesa.py

def plot(l, x1, x2, y, e):
    # Plot
    time_range = numpy.arange(0, l)
    pl.figure(1)
    pl.subplot(221)
    pl.plot(time_range, x1)
    pl.title("Input signal")
    pl.subplot(222)
    pl.plot(time_range, x2, c="r")
    pl.plot(time_range, y, c="b")
    pl.title("Reference signal")
    pl.subplot(223)
    pl.plot(time_range, e, c="r")
    pl.title("Noise")
    pl.xlabel("time")
    pl.show()

def plot_prediction_MM(model, y_train, y_test, plot_title=''):
    T = y_test.shape[0]
    mx, vx, my, vy_noiseless, vy = model.predict_forward(T, prop_mode=PROP_MM)
    T_train = y_train.shape[0]
    fig = plt.figure()
    ax = fig.add_subplot(111)
    ax.plot(np.arange(T_train), y_train[:, 0], 'k+-')
    ttest = np.arange(T_train, T_train+T)
    # pdb.set_trace()
    ax.plot(ttest, my[:, 0], '-', color='b')
    ax.fill_between(
        ttest, 
        my[:, 0] + 2*np.sqrt(vy_noiseless[:, 0]),
        my[:, 0] - 2*np.sqrt(vy_noiseless[:, 0]),
        alpha=0.3, edgecolor='b', facecolor='b')
    ax.fill_between(
        ttest, 
        my[:, 0] + 2*np.sqrt(vy[:, 0]),
        my[:, 0] - 2*np.sqrt(vy[:, 0]),
        alpha=0.1, edgecolor='b', facecolor='b')
    ax.plot(ttest, y_test, 'ro')
    ax.set_xlim([T_train-5, T_train + T])
    plt.title(plot_title)
    plt.savefig('/tmp/kink_pred_MM_'+plot_title+'.pdf')
    # plt.savefig('/tmp/kink_pred_MM_'+plot_title+'.png')

def plot_confusion_matrix(cm, label_list, title='Confusion matrix', cmap=None):
    from matplotlib import pylab
    cm = np.asarray(cm, dtype=np.float32)
    for i, row in enumerate(cm):
        cm[i] = cm[i] / np.sum(cm[i])
    #import matplotlib.pyplot as plt
    #plt.ion()
    pylab.clf()
    pylab.matshow(cm, fignum=False, cmap='Blues', vmin=0, vmax=1.0)
    ax = pylab.axes()
    ax.set_xticks(range(len(label_list)))
    ax.set_xticklabels(label_list, rotation='vertical')
    ax.xaxis.set_ticks_position('bottom')
    ax.set_yticks(range(len(label_list)))
    ax.set_yticklabels(label_list)
    pylab.title(title)
    pylab.colorbar()
    pylab.grid(False)
    pylab.xlabel('Predicted class')
    pylab.ylabel('True class')
    pylab.grid(False)
    pylab.savefig('test.jpg')
    pylab.show()

def plotRes(pre, real, test_x,l):
    s = set(pre)
    col = ['r','b','g','y','m']
    fig = plt.figure()

    ax = fig.add_subplot(111)
    for i in range(0, len(s)):
        index1 = pre == i
        index2 = real == i
        x1 = test_x[index1, :]
        x2 = test_x[index2, :]
        ax.scatter(x1[:,0],x1[:,1],color=col[i],marker='v',linewidths=0.5)
        ax.scatter(x2[:,0],x2[:,1],color=col[i],marker='.',linewidths=12)
    plt.title('learning rating='+str(l))
    plt.legend(('c1:predict','c1:true',\
                'c2:predict','c2:true',
                'c3:predict','c3:true',
                'c4:predict','c4:true',
                'c5:predict','c5:true'), shadow = True, loc = (0.01, 0.4))
    plt.show()

def ramachandran(PDB_file, title="Ramachandran Plot", AA_list=None,
    pymol_selection=None, engine=None):
    """PROCHECK style Ramachandran Plot
    A wrapper around ramachandran_tkinter and ramachandran_matplotlib

    engine (graphic engine for plotting)
        None - (Default) Use ramachandran_matplotlib if matplotlib is present
               Use ramachandran_tkinter if matplotlib is not importable
        "matplotlib" - Use ramachandran_matplotlib
        "tkinter" - Use ramachandran_tkinter
    """
    if not engine:
        engine="tkinter"

    if engine.lower().startswith("matplotlib"):
        ramachandran_matplotlib(PDB_file=PDB_file,title=title,AA_list=AA_list)
    elif engine.lower().startswith("tk"):
        ramachandran_tkinter(PDB_file=PDB_file,title=title,AA_list=AA_list,
            pymol_selection=pymol_selection)

def plotLine(self, x_vals, y_vals, x_label, y_label, title, filename=None):
        plt.clf()

        plt.xlabel(x_label)
        plt.xlim(((min(x_vals) - 0.5), (max(x_vals) + 0.5)))
        plt.ylabel(y_label)
        plt.ylim(((min(y_vals) - 0.5), (max(y_vals) + 0.5)))

        plt.title(title)
        plt.plot(x_vals, y_vals, c='k', lw=2)
        #plt.plot(x_vals, len(x_vals) * y_vals[0], c='r', lw=2)

        if filename == None:
            plt.show()
        else:
            plt.savefig(self.outputPath + filename)

def plot_entropy():
    pylab.clf()
    pylab.figure(num=None, figsize=(5, 4))

    title = "Entropy $H(X)$"
    pylab.title(title)
    pylab.xlabel("$P(X=$coin will show heads up$)$")
    pylab.ylabel("$H(X)$")

    pylab.xlim(xmin=0, xmax=1.1)
    x = np.arange(0.001, 1, 0.001)
    y = -x * np.log2(x) - (1 - x) * np.log2(1 - x)
    pylab.plot(x, y)
    # pylab.xticks([w*7*24 for w in [0,1,2,3,4]], ['week %i'%(w+1) for w in
    # [0,1,2,3,4]])

    pylab.autoscale(tight=True)
    pylab.grid(True)

    filename = "entropy_demo.png"
    pylab.savefig(os.path.join(CHART_DIR, filename), bbox_inches="tight")

def plot_clustering(x, y, title, mx=None, ymax=None, xmin=None, km=None):
    pylab.figure(num=None, figsize=(8, 6))
    if km:
        pylab.scatter(x, y, s=50, c=km.predict(list(zip(x, y))))
    else:
        pylab.scatter(x, y, s=50)

    pylab.title(title)
    pylab.xlabel("Occurrence word 1")
    pylab.ylabel("Occurrence word 2")

    pylab.autoscale(tight=True)
    pylab.ylim(ymin=0, ymax=1)
    pylab.xlim(xmin=0, xmax=1)
    pylab.grid(True, linestyle='-', color='0.75')

    return pylab

def plot_feat_importance(feature_names, clf, name):
    pylab.figure(num=None, figsize=(6, 5))
    coef_ = clf.coef_
    important = np.argsort(np.absolute(coef_.ravel()))
    f_imp = feature_names[important]
    coef = coef_.ravel()[important]
    inds = np.argsort(coef)
    f_imp = f_imp[inds]
    coef = coef[inds]
    xpos = np.array(list(range(len(coef))))
    pylab.bar(xpos, coef, width=1)

    pylab.title('Feature importance for %s' % (name))
    ax = pylab.gca()
    ax.set_xticks(np.arange(len(coef)))
    labels = ax.set_xticklabels(f_imp)
    for label in labels:
        label.set_rotation(90)
    filename = name.replace(" ", "_")
    pylab.savefig(os.path.join(
        CHART_DIR, "feat_imp_%s.png" % filename), bbox_inches="tight")

def plot_confusion_matrix(cm, genre_list, name, title):
    pylab.clf()
    pylab.matshow(cm, fignum=False, cmap='Blues', vmin=0, vmax=1.0)
    ax = pylab.axes()
    ax.set_xticks(range(len(genre_list)))
    ax.set_xticklabels(genre_list)
    ax.xaxis.set_ticks_position("bottom")
    ax.set_yticks(range(len(genre_list)))
    ax.set_yticklabels(genre_list)
    pylab.title(title)
    pylab.colorbar()
    pylab.grid(False)
    pylab.show()
    pylab.xlabel('Predicted class')
    pylab.ylabel('True class')
    pylab.grid(False)
    pylab.savefig(
        os.path.join(CHART_DIR, "confusion_matrix_%s.png" % name), bbox_inches="tight")

def plot_roc(auc_score, name, tpr, fpr, label=None):
    pylab.clf()
    pylab.figure(num=None, figsize=(5, 4))
    pylab.grid(True)
    pylab.plot([0, 1], [0, 1], 'k--')
    pylab.plot(fpr, tpr)
    pylab.fill_between(fpr, tpr, alpha=0.5)
    pylab.xlim([0.0, 1.0])
    pylab.ylim([0.0, 1.0])
    pylab.xlabel('False Positive Rate')
    pylab.ylabel('True Positive Rate')
    pylab.title('ROC curve (AUC = %0.2f) / %s' %
                (auc_score, label), verticalalignment="bottom")
    pylab.legend(loc="lower right")
    filename = name.replace(" ", "_")
    pylab.savefig(
        os.path.join(CHART_DIR, "roc_" + filename + ".png"), bbox_inches="tight")

def plot_feat_importance(feature_names, clf, name):
    pylab.clf()
    coef_ = clf.coef_
    important = np.argsort(np.absolute(coef_.ravel()))
    f_imp = feature_names[important]
    coef = coef_.ravel()[important]
    inds = np.argsort(coef)
    f_imp = f_imp[inds]
    coef = coef[inds]
    xpos = np.array(range(len(coef)))
    pylab.bar(xpos, coef, width=1)

    pylab.title('Feature importance for %s' % (name))
    ax = pylab.gca()
    ax.set_xticks(np.arange(len(coef)))
    labels = ax.set_xticklabels(f_imp)
    for label in labels:
        label.set_rotation(90)
    filename = name.replace(" ", "_")
    pylab.savefig(os.path.join(
        CHART_DIR, "feat_imp_%s.png" % filename), bbox_inches="tight")

def plotKChart(self, misClassDict, saveFigPath):
        kList = []
        misRateList = []
        for k, misClassNum in misClassDict.iteritems():
            kList.append(k)
            misRateList.append(1.0 - 1.0/k*misClassNum)

        fig = plt.figure(saveFigPath)
        plt.plot(kList, misRateList, 'r--')
        plt.title(saveFigPath)
        plt.xlabel('k Num.')
        plt.ylabel('Misclassified Rate')
        plt.legend(saveFigPath)
        plt.grid(True)
        plt.savefig(saveFigPath)
        plt.show()

################################### PART3 TEST ########################################
# ??

def show_feature_importance(gbdt, feature_names=None):
    importance = gbdt.get_fscore(fmap='xgb.fmap')
    importance = sorted(importance.items(), key=operator.itemgetter(1))

    df = pd.DataFrame(importance, columns=['feature', 'fscore'])
    df['fscore'] = df['fscore'] / df['fscore'].sum()
    print "feature importance", df

    if feature_names is not None:
        used_features = df['feature']
        unused_features = [f for f in feature_names if f not in used_features]
        print "[IDF]Unused features:", str(unused_features)

    plt.figure()
    df.plot()
    df.plot(kind='barh', x='feature', y='fscore', legend=False, figsize=(6, 10))
    plt.title('XGBoost Feature Importance')
    plt.xlabel('relative importance')
    plt.gcf().savefig('feature_importance_xgb.png')

def plot_penalty_vl(debug, tag, fold_exp):
    plt.close("all")
    vl = np.array(debug["penalty"])
    fig = plt.figure(figsize=(15, 10.8), dpi=300)
    names = debug["names"]
    for i in range(vl.shape[1]):
        if vl.shape[1] > 1:
            plt.plot(vl[:, i], label="layer_"+str(names[i]))
        else:
            plt.plot(vl[:], label="layer_"+str(names[i]))
    plt.xlabel("mini-batchs")
    plt.ylabel("value of penlaty")
    plt.title(
        "Penalty value over layers:" + "_".join([str(k) for k in names]) +
        ". tag:" + tag)
    plt.legend(loc='upper right', fancybox=True, shadow=True, prop={'size': 8})
    plt.grid(True)
    fig.savefig(fold_exp+"/penalty.png", bbox_inches='tight')
    plt.close('all')
    del fig

def plot_x_y_yhat(x, y, y_hat, xsz, ysz, binz=False):
    """Plot x, y and y_hat side by side."""
    plt.close("all")
    f = plt.figure(figsize=(15, 10.8), dpi=300)
    gs = gridspec.GridSpec(1, 3)
    if binz:
        y_hat = (y_hat > 0.5) * 1.
    ims = [x, y, y_hat]
    tils = [
        "x:" + str(xsz) + "x" + str(xsz),
        "y:" + str(ysz) + "x" + str(ysz),
        "yhat:" + str(ysz) + "x" + str(ysz)]
    for n, ti in zip([0, 1, 2], tils):
        f.add_subplot(gs[n])
        if n == 0:
            plt.imshow(ims[n], cmap=cm.Greys_r)
        else:
            plt.imshow(ims[n], cmap=cm.Greys_r)
        plt.title(ti)

    return f

def plot_roc(y_test, y_pred, label=''):
    """Compute ROC curve and ROC area"""

    fpr, tpr, _ = roc_curve(y_test, y_pred)
    roc_auc = auc(fpr, tpr)

    # Plot of a ROC curve for a specific class
    plt.figure()
    plt.plot(fpr, tpr, label='ROC curve (area = %0.2f)' % roc_auc)
    plt.plot([0, 1], [0, 1], 'k--')
    plt.xlim([0.0, 1.0])
    plt.ylim([0.0, 1.05])
    plt.xlabel('False Positive Rate')
    plt.ylabel('True Positive Rate')
    plt.title('Receiver operating characteristic' + label)
    plt.legend(loc="lower right")
    plt.show()

def plot_confusion_matrix(cm, plot_title, filename, genres=None):
    if not genres:
        genres = GENRES

    pylab.clf()
    pylab.matshow(cm, fignum=False, cmap='Blues', vmin=0, vmax=100.0)

    axes = pylab.axes()
    axes.set_xticks(range(len(genres)))
    axes.set_xticklabels(genres, rotation=45)

    axes.set_yticks(range(len(genres)))
    axes.set_yticklabels(genres)
    axes.xaxis.set_ticks_position("bottom")

    pylab.title(plot_title, fontsize=14)
    pylab.colorbar()
    pylab.xlabel('Predicted class', fontsize=12)
    pylab.ylabel('Correct class', fontsize=12)
    pylab.grid(False)
    #pylab.show()
    pylab.savefig(os.path.join(PLOTS_DIR, "cm_%s.eps" % filename), bbox_inches="tight")

def plot_confusion_matrix(cm, label_list, title='Confusion matrix', cmap=None):
    from matplotlib import pylab
    cm = np.asarray(cm, dtype=np.float32)
    for i, row in enumerate(cm):
        cm[i] = cm[i] / np.sum(cm[i])
    #import matplotlib.pyplot as plt
    #plt.ion()
    pylab.clf()
    pylab.matshow(cm, fignum=False, cmap='Blues', vmin=0, vmax=1.0)
    ax = pylab.axes()
    ax.set_xticks(range(len(label_list)))
    ax.set_xticklabels(label_list, rotation='vertical')
    ax.xaxis.set_ticks_position('bottom')
    ax.set_yticks(range(len(label_list)))
    ax.set_yticklabels(label_list)
    pylab.title(title)
    pylab.colorbar()
    pylab.grid(False)
    pylab.xlabel('Predicted class')
    pylab.ylabel('True class')
    pylab.grid(False)
    pylab.savefig('test.jpg')
    pylab.show()

def plot_feat_importance(feature_names, clf, name):
    pylab.figure(num=None, figsize=(6, 5))
    coef_ = clf.coef_
    important = np.argsort(np.absolute(coef_.ravel()))
    f_imp = feature_names[important]
    coef = coef_.ravel()[important]
    inds = np.argsort(coef)
    f_imp = f_imp[inds]
    coef = coef[inds]
    xpos = np.array(list(range(len(coef))))
    pylab.bar(xpos, coef, width=1)

    pylab.title('Feature importance for %s' % (name))
    ax = pylab.gca()
    ax.set_xticks(np.arange(len(coef)))
    labels = ax.set_xticklabels(f_imp)
    for label in labels:
        label.set_rotation(90)
    filename = name.replace(" ", "_")
    pylab.savefig(os.path.join(
        CHART_DIR, "feat_imp_%s.png" % filename), bbox_inches="tight")

def plot_feat_importance(feature_names, clf, name):
    pylab.figure(num=None, figsize=(6, 5))
    coef_ = clf.coef_
    important = np.argsort(np.absolute(coef_.ravel()))
    f_imp = feature_names[important]
    coef = coef_.ravel()[important]
    inds = np.argsort(coef)
    f_imp = f_imp[inds]
    coef = coef[inds]
    xpos = np.array(list(range(len(coef))))
    pylab.bar(xpos, coef, width=1)

    pylab.title('Feature importance for %s' % (name))
    ax = pylab.gca()
    ax.set_xticks(np.arange(len(coef)))
    labels = ax.set_xticklabels(f_imp)
    for label in labels:
        label.set_rotation(90)
    filename = name.replace(" ", "_")
    pylab.savefig(os.path.join(
        CHART_DIR, "feat_imp_%s.png" % filename), bbox_inches="tight")

def nmf(fdoc, fvocab):
    T = 100

    nmf = NMF(fdoc, fvocab)
    nmf.train(T)
    nmf.get_words()
    # print(mf.R)

    plt.figure()
    plt.plot(range(1,T+1),nmf.objective)
    plt.xticks(np.linspace(1,T,10))
    plt.xlabel('Iterations')
    plt.ylabel('Objective')
    plt.title('Variation of objective with iterations')
    plt.savefig('hw5_2a.png')
    plt.show()

def gp_partd(Xtrain,ytrain,Xtest,ytest):
    gp = gaussian_process(Xtrain[:,3],ytrain,Xtrain[:,3],ytrain)


    gp.init_kernel_matrices(b=5,var=2)
    gp.predict_test()

    x = np.asarray(Xtrain[:,3]).flatten()
    xsortind = np.argsort(x)
    y1 = np.asarray(ytrain).flatten()
    y2 = np.asarray(gp.test_predictions).flatten()
    plt.figure()
    plt.scatter(x[xsortind],y1[xsortind])
    plt.plot(x[xsortind],y2[xsortind],'b-')
    plt.xlabel('Car Weight (Dimension 4)')
    plt.ylabel('Outcome')
    plt.title('Visualizing model through single dimension')
    plt.savefig('hw3_gaussian_dim4_viz')
    plt.show()

def plot_tree_data(data, indicies_x, indicies_y, model):
    plt.subplot(3, 1, 1)
    data, indicies_x, indicies_y, model = load_tree_data()
    data_line, = plt.plot(data, color="blue", label="data")
    data_indicies_line, = plt.plot(
        indicies_x,
        indicies_y,
        "o",
        color="green",
        label="fitness predictors"
    )
    model_line, = plt.plot(model, color="red", label="model")
    plt.title("Data and Model Output")
    plt.legend()

    return data_line, data_indicies_line, model_line

def plot_tree_data(data, indicies_x, indicies_y, model, plot_indicies=False):
    plt.subplot(3, 1, 1)
    plt.plot(data, "o", color="blue", label="data")
    plt.plot(model, color="red", label="model")
    plt.ylim([-10, 10])

    if plot_indicies:
        plt.plot(
            indicies_x,
            indicies_y,
            "o",
            color="green",
            label="fitness predictors"
        )

    plt.title("Data and Model Output")
    plt.legend()

def gen_aperture(imgsize,ypos,xpos,radius,pixval=1,showaperture=False,verbose=True):
    """
    Generating an aperture image

    --- INPUT ---
    imgsize       The dimensions of the array to return. Expects [y-size,x-size].
                  The aperture will be positioned in the center of a (+/-x-size/2., +/-y-size/2) sized array
    ypos          Pixel position in the y direction
    xpos          Pixel position in the x direction
    radius        Radius of aperture in pixels
    showaperture  Display image of generated aperture
    verbose       Toggle verbosity

    --- EXAMPLE OF USE ---
    import tdose_utilities as tu
    apertureimg  = tu.gen_aperture([20,40],10,5,10,showaperture=True)
    apertureimg  = tu.gen_aperture([2000,4000],900,1700,150,showaperture=True)

    """
    if verbose: print ' - Generating aperture in image (2D array)'
    y , x    = np.ogrid[-ypos:imgsize[0]-ypos, -xpos:imgsize[1]-xpos]
    mask     = x*x + y*y <= radius**2.
    aperture = np.zeros(imgsize)

    if verbose: print ' - Assigning pixel value '+str(pixval)+' to aperture'
    aperture[mask] = pixval

    if showaperture:
        if verbose: print ' - Displaying resulting image of aperture'
        plt.imshow(aperture,interpolation='none')
        plt.title('Generated aperture')
        plt.show()

    return aperture
# = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =

def gen_overview_plot_image(ax,imagefile,imgext=0,cubelayer=1,title='Img Title?',fontsize=6,lthick=2,alpha=0.5,
                            cmap='coolwarm'):
    """
    Plotting commands for image (cube layer) overview plotting

    --- INPUT ---

    cubelayer     If the content of the file is a cube, provide the cube layer to plot. If
                    cubelayer = 'fmax' the layer with most flux is plotted

    """

    ax.set_title(title,fontsize=fontsize)
    if os.path.isfile(imagefile):
        imgdata = pyfits.open(imagefile)[imgext].data

        if len(imgdata.shape) == 3: # it is a cube
            imgdata = imgdata[cubelayer,:,:]

        ax.imshow(imgdata, interpolation='None',cmap=cmap,aspect='equal', origin='lower')

        ax.set_xlabel('x-pixel')
        ax.set_ylabel('y-pixel ')
        ax.set_xticks([])
        ax.set_yticks([])

    else:
        textstr = 'No image\nfound'
        ax.text(1.0,22,textstr,horizontalalignment='center',verticalalignment='center',fontsize=fontsize)

        ax.set_ylim([28,16])
        ax.plot([0.0,2.0],[28,16],'r--',lw=lthick)
        ax.plot([2.0,0.0],[28,16],'r--',lw=lthick)

        ax.set_xlabel(' ')
        ax.set_ylabel(' ')
        ax.set_xticks([])
        ax.set_yticks([])

# = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =

def residual_multigauss(param, dataimage, nonfinite = 0.0, ravelresidual=True, showimages=False, verbose=False):
    """
    Calculating the residual bestween the multigaussian model with the paramters 'param' and the data.

    --- INPUT ---
    param         Parameters of multi-gaussian model to generate. See modelimage_multigauss() header for details
    dataimage     Data image to take residual
    nonfinite     Value to replace non-finite entries in residual with
    ravelresidual To np.ravel() the residual image set this to True. Needed by scipy.optimize.leastsq()
                  optimizer function
    showimages    To show model and residiual images set to True
    verbose       Toggle verbosity

    --- EXAMPLE OF USE ---
    import tdose_model_FoV as tmf
    param      = [18,31,1*0.3,2.1*0.3,1.2*0.3,30*0.3,    110,90,200*0.5,20.1*0.5,15.2*0.5,0*0.5]
    dataimg    = pyfits.open('/Users/kschmidt/work/TDOSE/mock_cube_sourcecat161213_tdose_mock_cube.fits')[0].data[0,:,:]
    residual   = tmf.residual_multigauss(param, dataimg, showimages=True)

    """
    if verbose: ' - Estimating residual (= model - data) between model and data image'
    imgsize      = dataimage.shape
    xgrid, ygrid = tu.gen_gridcomponents(imgsize)
    modelimg     = tmf.modelimage_multigauss((xgrid, ygrid),param,imgsize,showmodelimg=showimages, verbose=verbose)

    residualimg  = modelimg - dataimage

    if showimages:
        plt.imshow(residualimg,interpolation='none', vmin=1e-5, vmax=np.max(residualimg), norm=mpl.colors.LogNorm())
        plt.title('Resdiaul (= model - data) image')
        plt.show()

    if nonfinite is not None:
        residualimg[~np.isfinite(residualimg)] = 0.0

    if ravelresidual:
        residualimg = np.ravel(residualimg)

    return residualimg
# = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =

def _do_title_string(self,title_items,cycle):
        '''
        Create title string

        Private method that creates a title string for a cycle plot
        out of a list of title_items that are cycle attributes and can
        be obtained with self.get

        Parameters
        ----------
        title_items : list
            A list of cycle attributes.
        cycle : scalar
            The cycle for which the title string should be created.

        Returns
        -------
        title_string: string
            Title string that can be used to decorate plot.

        '''

        title_string=[]
        form_str='%4.1F'

        for item in title_items:
            num=self.get(item,fname=cycle)
            if num > 999:
                num=log10(num)
                prefix='log '
            else:
                prefix=''
            title_string.append(prefix+item+'='+form_str%num)
        tt=''
        for thing in title_string:
            tt = tt+thing+", "
        return tt.rstrip(', ')

def test_abu_evolution(self):
        from nugridpy import ppn, utils
        import matplotlib
        matplotlib.use('agg')
        import matplotlib.pylab as mpy
        import os

        # Perform tests within temporary directory
        with TemporaryDirectory() as tdir:
            # wget the data for a ppn run from the CADC VOspace
            os.system("wget -q --content-disposition --directory '" + tdir +  "' "\
                          + "'http://www.canfar.phys.uvic.ca/vospace/synctrans?TARGET="\
                          + "vos%3A%2F%2Fcadc.nrc.ca%21vospace%2Fnugrid%2Fdata%2Fprojects%2Fppn%2Fexamples%2F"\
                          + "ppn_Hburn_simple%2Fx-time.dat&DIRECTION=pullFromVoSpace&PROTOCOL"\
                          + "=ivo%3A%2F%2Fivoa.net%2Fvospace%2Fcore%23httpget'")

            #nugrid_dir= os.path.dirname(os.path.dirname(ppn.__file__))
            #NuPPN_dir= nugrid_dir + "/NuPPN"
            #test_data_dir= NuPPN_dir + "/examples/ppn_Hburn_simple/RUN_MASTER"

            symbs=utils.symbol_list('lines2')
            x=ppn.xtime(tdir)
            specs=['PROT','HE  4','C  12','N  14','O  16']
            i=0
            for spec in specs:
                x.plot('time',spec,logy=True,logx=True,shape=utils.linestyle(i)[0],show=False,title='')
                i += 1
            mpy.ylim(-5,0.2)
            mpy.legend(loc=0)
            mpy.xlabel('$\log t / \mathrm{min}$')
            mpy.ylabel('$\log X \mathrm{[mass fraction]}$')
            abu_evol_file = 'abu_evolution.png'
            mpy.savefig(abu_evol_file)
            self.assertTrue(os.path.exists(abu_evol_file))

def plot1D_mat(a, b, M, title=''):
    """ Plot matrix M  with the source and target 1D distribution

    Creates a subplot with the source distribution a on the left and
    target distribution b on the tot. The matrix M is shown in between.


    Parameters
    ----------
    a : np.array, shape (na,)
        Source distribution
    b : np.array, shape (nb,)
        Target distribution
    M : np.array, shape (na,nb)
        Matrix to plot
    """
    na, nb = M.shape

    gs = gridspec.GridSpec(3, 3)

    xa = np.arange(na)
    xb = np.arange(nb)

    ax1 = pl.subplot(gs[0, 1:])
    pl.plot(xb, b, 'r', label='Target distribution')
    pl.yticks(())
    pl.title(title)

    ax2 = pl.subplot(gs[1:, 0])
    pl.plot(a, xa, 'b', label='Source distribution')
    pl.gca().invert_xaxis()
    pl.gca().invert_yaxis()
    pl.xticks(())

    pl.subplot(gs[1:, 1:], sharex=ax1, sharey=ax2)
    pl.imshow(M, interpolation='nearest')
    pl.axis('off')

    pl.xlim((0, nb))
    pl.tight_layout()
    pl.subplots_adjust(wspace=0., hspace=0.2)

def plot_latent(model, y, plot_title=''):
    # make prediction on some test inputs
    N_test = 300
    C = model.get_hypers()['C_emission'][0, 0]
    x_test = np.linspace(-10, 8, N_test) / C
    x_test = np.reshape(x_test, [N_test, 1])
    if isinstance(model, aep.SGPSSM) or isinstance(model, vfe.SGPSSM):
        zu = model.dyn_layer.zu
    else:
        zu = model.sgp_layer.zu
    mu, vu = model.predict_f(zu)
    # mu, Su = model.dyn_layer.mu, model.dyn_layer.Su
    mf, vf = model.predict_f(x_test)
    my, vy = model.predict_y(x_test)
    # plot function
    fig = plt.figure()
    ax = fig.add_subplot(111)
    # ax.plot(x_test[:,0], kink_true(x_test[:,0]), '-', color='k')
    ax.plot(C*x_test[:,0], my[:,0], '-', color='r', label='y')
    ax.fill_between(
        C*x_test[:,0], 
        my[:,0] + 2*np.sqrt(vy[:, 0]), 
        my[:,0] - 2*np.sqrt(vy[:, 0]), 
        alpha=0.2, edgecolor='r', facecolor='r')
    ax.plot(
        y[0:model.N-1], 
        y[1:model.N], 
        'r+', alpha=0.5)
    mx, vx = model.get_posterior_x()
    ax.set_xlabel(r'$x_{t-1}$')
    ax.set_ylabel(r'$x_{t}$')
    plt.title(plot_title)
    plt.savefig('/tmp/lincos_'+plot_title+'.png')

# generate a dataset from the lincos function above

def plot_latent(model, y, plot_title=''):
    # make prediction on some test inputs
    N_test = 200
    C = model.get_hypers()['C_emission'][0, 0]
    x_test = np.linspace(-4, 6, N_test) / C
    x_test = np.reshape(x_test, [N_test, 1])
    zu = model.dyn_layer.zu
    mu, vu = model.predict_f(zu)
    # mu, Su = model.dyn_layer.mu, model.dyn_layer.Su
    mf, vf = model.predict_f(x_test)
    my, vy = model.predict_y(x_test)
    # plot function
    fig = plt.figure()
    ax = fig.add_subplot(111)
    # ax.plot(x_test[:,0], kink_true(x_test[:,0]), '-', color='k')
    ax.plot(C*x_test[:,0], my[:,0], '-', color='r', label='y')
    ax.fill_between(
        C*x_test[:,0], 
        my[:,0] + 2*np.sqrt(vy[:, 0]), 
        my[:,0] - 2*np.sqrt(vy[:, 0]), 
        alpha=0.2, edgecolor='r', facecolor='r')
    # ax.plot(zu, mu, 'ob')
    # ax.errorbar(zu, mu, yerr=3*np.sqrt(vu), fmt='ob')
    # ax.plot(x_test[:,0], mf[:,0], '-', color='b')
    # ax.fill_between(
    #     x_test[:,0], 
    #     mf[:,0] + 2*np.sqrt(vf[:,0]), 
    #     mf[:,0] - 2*np.sqrt(vf[:,0]), 
    #     alpha=0.2, edgecolor='b', facecolor='b')
    ax.plot(
        y[0:model.N-1], 
        y[1:model.N], 
        'r+', alpha=0.5)
    mx, vx = model.get_posterior_x()
    ax.set_xlabel(r'$x_{t-1}$')
    ax.set_ylabel(r'$x_{t}$')
    ax.set_xlim([-4, 6])
    # ax.set_ylim([-7, 7])
    plt.title(plot_title)
    # plt.savefig('/tmp/kink_'+plot_title+'.pdf')
    plt.savefig('/tmp/kink_'+plot_title+'.png')

def plot_prediction_MC(model, y_train, y_test, plot_title=''):
    T = y_test.shape[0]
    x_samples, my, vy = model.predict_forward(T, prop_mode=PROP_MC)
    T_train = y_train.shape[0]
    fig = plt.figure()
    ax = fig.add_subplot(111)
    ax.plot(np.arange(T_train), y_train[:, 0], 'k+-')
    ttest = np.arange(T_train, T_train+T)
    ttest = np.reshape(ttest, [T, 1])
    loglik, ranks = compute_log_lik(np.exp(2*model.sn), y_test, my[:, :, 0].T)
    red = 0.1
    green = 0. * red
    blue = 1. - red
    color = np.array([red, green, blue]).T
    for k in np.argsort(ranks):
        ax.plot(ttest, my[:, k, 0], '-', color=color*ranks[k], alpha=0.5)
    # ax.plot(np.tile(ttest, [1, my.shape[1]]), my[:, :, 0], '-x', color='r', alpha=0.3)
    # ax.plot(np.tile(ttest, [1, my.shape[1]]), x_samples[:, :, 0], 'x', color='m', alpha=0.3)
    ax.plot(ttest, y_test, 'ro')
    ax.set_xlim([T_train-5, T_train + T])
    plt.title(plot_title)
    plt.savefig('/tmp/kink_pred_MC_'+plot_title+'.pdf')
    # plt.savefig('/tmp/kink_pred_MC_'+plot_title+'.png')


# generate a dataset from the kink function above

def dispersion_plot(text, words, ignore_case=False, title="Lexical Dispersion Plot"):
    """
    Generate a lexical dispersion plot.

    :param text: The source text
    :type text: list(str) or enum(str)
    :param words: The target words
    :type words: list of str
    :param ignore_case: flag to set if case should be ignored when searching text
    :type ignore_case: bool
    """

    try:
        from matplotlib import pylab
    except ImportError:
        raise ValueError('The plot function requires matplotlib to be installed.'
                     'See http://matplotlib.org/')

    text = list(text)
    words.reverse()

    if ignore_case:
        words_to_comp = list(map(str.lower, words))
        text_to_comp = list(map(str.lower, text))
    else:
        words_to_comp = words
        text_to_comp = text

    points = [(x,y) for x in range(len(text_to_comp))
                    for y in range(len(words_to_comp))
                    if text_to_comp[x] == words_to_comp[y]]
    if points:
        x, y = list(zip(*points))
    else:
        x = y = ()
    pylab.plot(x, y, "b|", scalex=.1)
    pylab.yticks(list(range(len(words))), words, color="b")
    pylab.ylim(-1, len(words))
    pylab.title(title)
    pylab.xlabel("Word Offset")
    pylab.show()

def plot_word_freq_dist(text):
    fd = text.vocab()

    samples = [item for item, _ in fd.most_common(50)]
    values = [fd[sample] for sample in samples]
    values = [sum(values[:i+1]) * 100.0/fd.N() for i in range(len(values))]
    pylab.title(text.name)
    pylab.xlabel("Samples")
    pylab.ylabel("Cumulative Percentage")
    pylab.plot(values)
    pylab.xticks(range(len(samples)), [str(s) for s in samples], rotation=90)
    pylab.show()

def tabulate(self, *args, **kwargs):
        """
        Tabulate the given samples from the frequency distribution (cumulative),
        displaying the most frequent sample first.  If an integer
        parameter is supplied, stop after this many samples have been
        plotted.

        :param samples: The samples to plot (default is all samples)
        :type samples: list
        :param cumulative: A flag to specify whether the freqs are cumulative (default = False)
        :type title: bool
        """
        if len(args) == 0:
            args = [len(self)]
        samples = [item for item, _ in self.most_common(*args)]

        cumulative = _get_kwarg(kwargs, 'cumulative', False)
        if cumulative:
            freqs = list(self._cumulative_frequencies(samples))
        else:
            freqs = [self[sample] for sample in samples]
        # percents = [f * 100 for f in freqs]  only in ProbDist?

        width = max(len("%s" % s) for s in samples)
        width = max(width, max(len("%d" % f) for f in freqs))

        for i in range(len(samples)):
            print("%*s" % (width, samples[i]), end=' ')
        print()
        for i in range(len(samples)):
            print("%*d" % (width, freqs[i]), end=' ')
        print()

def plot_velocity(self, timestamps, vel_true, vel_est):
        N = vel_est.shape[1]
        t = timestamps[:N]
        vel_true = vel_true[:, :N]
        vel_est = vel_est[:, :N]

        # Figure
        plt.figure()
        plt.suptitle("Velocity")

        # X axis
        plt.subplot(311)
        plt.plot(t, vel_true[0, :], color="red", label="Ground_truth")
        plt.plot(t, vel_est[0, :], color="blue", label="Estimate")

        plt.title("x-axis")
        plt.xlabel("Date Time")
        plt.ylabel("ms^-1")
        plt.legend(loc=0)

        # Y axis
        plt.subplot(312)
        plt.plot(t, vel_true[1, :], color="red", label="Ground_truth")
        plt.plot(t, vel_est[1, :], color="blue", label="Estimate")

        plt.title("y-axis")
        plt.xlabel("Date Time")
        plt.ylabel("ms^-1")
        plt.legend(loc=0)

        # Z axis
        plt.subplot(313)
        plt.plot(t, vel_true[2, :], color="red", label="Ground_truth")
        plt.plot(t, vel_est[2, :], color="blue", label="Estimate")

        plt.title("z-axis")
        plt.xlabel("Date Time")
        plt.ylabel("ms^-1")
        plt.legend(loc=0)

def plot_attitude(self, timestamps, att_true, att_est):
        # Setup
        N = att_est.shape[1]
        t = timestamps[:N]
        att_true = att_true[:, :N]
        att_est = att_est[:, :N]

        # Figure
        plt.figure()
        plt.suptitle("Attitude")

        # X axis
        plt.subplot(311)
        plt.plot(t, att_true[0, :], color="red", label="Ground_truth")
        plt.plot(t, att_est[0, :], color="blue", label="Estimate")

        plt.title("x-axis")
        plt.legend(loc=0)
        plt.xlabel("Date Time")
        plt.ylabel("rad s^-1")

        # Y axis
        plt.subplot(312)
        plt.plot(t, att_true[1, :], color="red", label="Ground_truth")
        plt.plot(t, att_est[1, :], color="blue", label="Estimate")

        plt.title("y-axis")
        plt.legend(loc=0)
        plt.xlabel("Date Time")
        plt.ylabel("rad s^-1")

        # Z axis
        plt.subplot(313)
        plt.plot(t, att_true[2, :], color="red", label="Ground_truth")
        plt.plot(t, att_est[2, :], color="blue", label="Estimate")

        plt.title("z-axis")
        plt.legend(loc=0)
        plt.xlabel("Date Time")
        plt.ylabel("rad s^-1")

def plot_velocity(self, timestamps, vel_true, vel_est):
        N = vel_est.shape[1]
        t = timestamps[:N]
        vel_true = vel_true[:, :N]
        vel_est = vel_est[:, :N]

        # Figure
        plt.figure()
        plt.suptitle("Velocity")

        # X axis
        plt.subplot(311)
        plt.plot(t, vel_true[0, :], color="red", label="Ground_truth")
        plt.plot(t, vel_est[0, :], color="blue", label="Estimate")

        plt.title("x-axis")
        plt.xlabel("Date Time")
        plt.ylabel("ms^-1")
        plt.legend(loc=0)

        # Y axis
        plt.subplot(312)
        plt.plot(t, vel_true[1, :], color="red", label="Ground_truth")
        plt.plot(t, vel_est[1, :], color="blue", label="Estimate")

        plt.title("y-axis")
        plt.xlabel("Date Time")
        plt.ylabel("ms^-1")
        plt.legend(loc=0)

        # Z axis
        plt.subplot(313)
        plt.plot(t, vel_true[2, :], color="red", label="Ground_truth")
        plt.plot(t, vel_est[2, :], color="blue", label="Estimate")

        plt.title("z-axis")
        plt.xlabel("Date Time")
        plt.ylabel("ms^-1")
        plt.legend(loc=0)

def plot_attitude(self, timestamps, att_true, att_est):
        # Setup
        N = att_est.shape[1]
        t = timestamps[:N]
        att_true = att_true[:, :N]
        att_est = att_est[:, :N]

        # Figure
        plt.figure()
        plt.suptitle("Attitude")

        # X axis
        plt.subplot(311)
        plt.plot(t, att_true[0, :], color="red", label="Ground_truth")
        plt.plot(t, att_est[0, :], color="blue", label="Estimate")

        plt.title("x-axis")
        plt.legend(loc=0)
        plt.xlabel("Date Time")
        plt.ylabel("rad s^-1")

        # Y axis
        plt.subplot(312)
        plt.plot(t, att_true[1, :], color="red", label="Ground_truth")
        plt.plot(t, att_est[1, :], color="blue", label="Estimate")

        plt.title("y-axis")
        plt.legend(loc=0)
        plt.xlabel("Date Time")
        plt.ylabel("rad s^-1")

        # Z axis
        plt.subplot(313)
        plt.plot(t, att_true[2, :], color="red", label="Ground_truth")
        plt.plot(t, att_est[2, :], color="blue", label="Estimate")

        plt.title("z-axis")
        plt.legend(loc=0)
        plt.xlabel("Date Time")
        plt.ylabel("rad s^-1")

def plot_storage(self, storage):
        plt.figure()
        plt.plot(range(len(storage)), storage)
        plt.title("Num of tracks over time")
        plt.xlabel("Frame No.")
        plt.ylabel("Num of Tracks")