def testplot(cats, catsavg, xy, data, levels=[0, 1, 2, 3, 4, 5, 10, 15, 20, 25, 30, 40, 50, 100], title=""):
    """Quick test plot layout for this example file
    """
    colors = plt.cm.spectral(np.linspace(0, 1, len(levels)))
    mycmap, mynorm = from_levels_and_colors(levels, colors, extend="max")

    radolevels = [0, 1, 2, 3, 4, 5, 10, 15, 20, 25, 30, 40, 50, 100]
    radocolors = plt.cm.spectral(np.linspace(0, 1, len(radolevels)))
    radocmap, radonorm = from_levels_and_colors(radolevels, radocolors, extend="max")

    fig = plt.figure(figsize=(14, 8))

    # Average rainfall sum
    ax = fig.add_subplot(121, aspect="equal")
    coll = PatchCollection(cats, array=catsavg, cmap=mycmap, norm=mynorm, edgecolors='white', lw=0.5)
    ax.add_collection(coll)
    ax.autoscale()
    plt.colorbar(coll, ax=ax, shrink=0.5)
    plt.xlabel("GK2 Easting")
    plt.ylabel("GK2 Northing")
    plt.title(title)
    plt.draw()

    # Original radar data
    ax1 = fig.add_subplot(122, aspect="equal")
    pm = plt.pcolormesh(xy[:, :, 0], xy[:, :, 1], np.ma.masked_invalid(data), cmap=radocmap, norm=radonorm)
    coll = PatchCollection(cats, facecolor='None', edgecolor='white', lw=0.5)
    ax1.add_collection(coll)
    cb = plt.colorbar(pm, ax=ax1, shrink=0.5)
    cb.set_label("(mm/h)")
    plt.xlabel("GK2 Easting")
    plt.ylabel("GK2 Northing")
    plt.title("Original radar rain sums")
    plt.draw()
    plt.tight_layout()

    def geometry_plot(self,Npoints = 256):
        x = np.linspace(-self.a, self.a)
        y = np.linspace(-self.b, self.b)
        X,Y = np.meshgrid(x,y)
        n_plot = np.zeros(np.shape(X))
        k_plot = np.zeros(np.shape(X))
        for i,xx in enumerate(x):
            for j,yy in enumerate(y):
                n_plot[i,j] = ref([xx,yy])*10000000
                k_plot[i,j] = extinction([xx,yy])*10000000

                cmap1, norm1 = from_levels_and_colors([1,ref([r_core,0])*10000000,ref([r_clad,0])*10000000], ['blue', 'green','red'],extend='max')
                cmap2, norm2 = from_levels_and_colors([0,extinction([r_core,0])*10000000,extinction([r_clad,0])*10000000],  ['blue', 'green','red'],extend='max')

                fig = plt.figure(figsize=(20.0, 20.0))
                ax1 = fig.add_subplot(221)
                ax1.pcolormesh(X,Y,n_plot, cmap=cmap1, norm=norm1)
                ax1.set_title('real part profile')
                ax1.axis('equal')
                
                ax2 = fig.add_subplot(222)
                ax2.pcolormesh(X,Y,k_plot, cmap=cmap2, norm=norm2)
                ax2.set_title('Imaginary part  profile')
                ax2.axis('equal')
                return n_plot,k_plot

def test_cmap_and_norm_from_levels_and_colors2():
    levels = [-1, 2, 2.5, 3]
    colors = ['red', (0, 1, 0), 'blue', (0.5, 0.5, 0.5), (0.0, 0.0, 0.0, 1.0)]
    clr = mcolors.to_rgba_array(colors)
    bad = (0.1, 0.1, 0.1, 0.1)
    no_color = (0.0, 0.0, 0.0, 0.0)
    masked_value = 'masked_value'

    # Define the test values which are of interest.
    # Note: levels are lev[i] <= v < lev[i+1]
    tests = [('both', None, {-2: clr[0],
                             -1: clr[1],
                             2: clr[2],
                             2.25: clr[2],
                             3: clr[4],
                             3.5: clr[4],
                             masked_value: bad}),

             ('min', -1, {-2: clr[0],
                          -1: clr[1],
                          2: clr[2],
                          2.25: clr[2],
                          3: no_color,
                          3.5: no_color,
                          masked_value: bad}),

             ('max', -1, {-2: no_color,
                          -1: clr[0],
                          2: clr[1],
                          2.25: clr[1],
                          3: clr[3],
                          3.5: clr[3],
                          masked_value: bad}),

             ('neither', -2, {-2: no_color,
                              -1: clr[0],
                              2: clr[1],
                              2.25: clr[1],
                              3: no_color,
                              3.5: no_color,
                              masked_value: bad}),
             ]

    for extend, i1, cases in tests:
        cmap, norm = mcolors.from_levels_and_colors(levels, colors[0:i1],
                                                    extend=extend)
        cmap.set_bad(bad)
        for d_val, expected_color in cases.items():
            if d_val == masked_value:
                d_val = np.ma.array([1], mask=True)
            else:
                d_val = [d_val]
            assert_array_equal(expected_color, cmap(norm(d_val))[0],
                               'Wih extend={0!r} and data '
                               'value={1!r}'.format(extend, d_val))

    with pytest.raises(ValueError):
        mcolors.from_levels_and_colors(levels, colors)

def make_plot(a, b, z, title, maxterms):
    approx, terms = optimal_terms.asymptotic_series(a, b, z, maxterms)
    ref = np.float64(mpmath.hyp1f1(a, b, z))

    cd = correct_digits(approx, ref)
    termsize = np.abs(terms/terms[0])

    fig, ax1 = plt.subplots()
    ax1.plot(cd, '-', linewidth=2, color=GREEN)
    ax1.set_ylim(0, 17)
    ax1.set_xlabel('term number')
    # Make the y-axis label and tick labels match the line color.
    ax1.set_ylabel('correct digits', color=GREEN)
    for tl in ax1.get_yticklabels():
        tl.set_color(GREEN)

    ax2 = ax1.twinx()
    cmap, norm = colors.from_levels_and_colors([-np.inf, 0, np.inf],
                                               [BLUE, RED])
    ax2.scatter(np.arange(termsize.shape[0]), termsize,
                c=terms, cmap=cmap, norm=norm, edgecolors='')
    # ax2.semilogy(termsize, 'r--', linewidth=2)
    ax2.set_yscale('log')
    ax2.set_ylabel('relative term size', color=RED)
    # Set the limits, with a little margin.
    ax2.set_xlim(0, termsize.shape[0])
    ax2.set_ylim(np.min(termsize), np.max(termsize))
    for tl in ax2.get_yticklabels():
        tl.set_color(RED)

    ax1.set_title("a = {:.2e}, b = {:.2e}, z = {:.2e}".format(a, b, z))

    plt.savefig("{}.png".format(title))

    def animate(self, steps=10, start_fire_power=0):
        fig = plt.figure()

        colors = ['#990000', '#CC0000', '#FF0000', '#FF8000', '#FFFF00',
                  '#AFCB96',
                  '#66FF66', '#33FF33', '#00FF00', '#00CC00', '#009900']
        scale = [-1000, -80, -60, -40, -20, -0.1, 0.1, 50, 100, 200, 400, 1000]
        cmap, norm = mcolors.from_levels_and_colors(scale, colors)

        df = self.prepare_data_to_draw()
        images = [[plt.pcolor(df, cmap=cmap, norm=norm)]]

        self.start_fire(start_fire_power)
        df = self.prepare_data_to_draw()
        images.append([plt.pcolor(df, cmap=cmap, norm=norm)])
        drzewa = [sum([x.wood for x in sum(self.forest.forest, [])])]
        ogien = [sum([x.fire for x in sum(self.forest.forest, [])])]
        for i in range(steps):
            print i
            self.next_step()
            self.burn()
            drzewa.append(sum([x.wood for x in sum(self.forest.forest, [])]))
            ogien.append(sum([x.fire for x in sum(self.forest.forest, [])]))
            df = self.prepare_data_to_draw()
            images.append([plt.pcolor(df, cmap=cmap, norm=norm,)])

        ani = animation.ArtistAnimation(fig, images, interval=500, repeat_delay=0, repeat=True)
        plt.show()

        plt.plot(range(0, steps + 1), drzewa)
        plt.show()

        plt.plot(range(0, steps + 1), ogien)
        plt.show()

def smart_colormap(clevs, name='jet', extend='both', minval=0.0, maxval=1.0):
    '''
    Automatically grabs the colors to extend the colorbar from the colormap.
    '''
    
    # Define number of colors
    if extend == 'both':
        nrColors = len(clevs)+1
    elif (extend == 'min') | (extend == 'max'):
        nrColors = len(clevs)
    elif (extend == 'neither'):
        nrColors = len(clevs)-1
    else:
        nrColors = len(clevs)-1
        extend = 'neither'
    
    # Get colormap
    cmap = get_cmap(name, nrColors)
    
    # Truncate colormap if asked
    if (minval != 0.0) or (maxval != 1.0):
        cmap = truncate_colormap(cmap, minval=minval, maxval=maxval, n=nrColors/2)
    
    # Get the list of colors
    colors = []
    for i in range(0, nrColors):
        colors.append(cmap(i/(nrColors-1)))
    
    # Use utility function to get cmap and norm at the same time
    cmap, norm = from_levels_and_colors(clevs, colors, extend=extend)

    return(cmap, norm)

def create_wind_color_map(levels):
    # replicate XWS colours for wind speed
    #
    cmap = ["#ffffff", "#fff4e8", "#ffe1c1", "#ffaa4e", "#ff6d00",
    		"#d33100", "#890000"] #, "#650000", "#390000"]
    ccmap, norm = col.from_levels_and_colors(levels, cmap, 'neither')
    return ccmap, norm

    def generer_image_transp(self):
        """Fonction graphique.

        Elle est identique à la fonction précédente mais permet de générer
        une image transparente sauf pour les individus qui possèdent une
        couleur par tribu. Les deux images sont ensuite superposée pour
        donner l'illusion que la simulation se déroule sur la carte du
        dessous.

        Arguments:
            verbose (bool):
                Active ou désactive la sortie console.

        Retour:
            Modifie img_transp par référence interne.

        """
        nombre_tribus = len(self.liste_tribus)

        levels = list(range(0, 6)) + self.liste_tribus_int
        colors_void = [(0., 0., 0., 0.) for i in range(0, 6)]
        jet = plt.get_cmap('gist_rainbow')
        colors_tribu = [jet(x/(max(self.liste_tribus_int) - 10)) for x in range(0, nombre_tribus)]
        colors = colors_void + colors_tribu

        self.colormap_indiv, self.norm = from_levels_and_colors(levels, colors, extend='max')

        self.img_transp = self.matrice_tribus.astype(str).astype(int)

def heatmap(data_, num_levels_=20, mode=None, xlabel='', ylabel='', xlabels=None, ylabels=None, rotation=90, changeTicks=True):
    '''generates a nice heatmap'''
    if mode == 'special':
        vmin, vmax = data_.min(), data_.max()
        midpoint = data_.mean()
    else:
        vmin, vmax = -1.0001, 1.0001
        midpoint = 0.0
    levels = np.linspace(vmin, vmax, num_levels_)
    midp = np.mean(np.c_[levels[:-1], levels[1:]], axis=1)
    vals = np.interp(midp, [vmin, midpoint, vmax], [0, 0.5, 1])
    colors = plt.cm.seismic(vals)
    cmap, norm = from_levels_and_colors(levels, colors)
    fig, ax = plt.subplots()
    im = ax.imshow(data_, cmap=cmap, norm=norm, interpolation='none')
    fig.colorbar(im)

    ax.xaxis.tick_top()
    ax.set_xlabel(xlabel)
    ax.set_ylabel(ylabel)
    ax.xaxis.set_label_position('top')
    if xlabels is not None:
        ax.set_xticklabels(xlabels, minor=False)
    if ylabels is not None:
        ax.set_yticklabels(ylabels, minor=False)
    if changeTicks:
        plt.xticks(np.arange(len(data_)), rotation=rotation, fontsize=10)
        plt.yticks(np.arange(len(data_)), fontsize=10)

    plt.show()

    def view(self, original, segmented, reduced):
        r"""
        This method is implemented for test purposes, it takes as arguments
        an untreated slice, a segmented slice and a reduced and segmented
        slice showing its differences on screen using a matplotlib figure

            view(original, segmented, reduced)
        """

        f = pyplot.figure()
        levels = [0, 1, 2]
        colores = ['red', 'white', 'blue', 'red']
        cmap, norm = colors.from_levels_and_colors(levels, colores, extend='both')

       # org = self.sample_mask(original) #Mask irelevant data

        f.add_subplot(1, 3, 1)  # Original image
        pyplot.imshow(original, cmap='binary', interpolation='nearest')
        pyplot.colorbar()

        f.add_subplot(1, 3, 2)  # Segmented image  by K-means
        pyplot.imshow(segmented, interpolation='nearest', cmap=cmap, norm=norm)
        pyplot.colorbar()

        f.add_subplot(1, 3, 3)  # Reduced image
        pyplot.imshow(reduced, interpolation='nearest',
                      origin='lower', cmap = cmap, norm = norm)
        pyplot.colorbar()
        pyplot.show()

def plotheat(data_, num_levels_=20, mode=None, xlabel='', ylabel='',
             xlabels=None, ylabels=None, rotation=90,
             changeTicks=True, annotation=None, vmin=None, vmax=None,
             showplot=True, removebar=False):
    '''generates a nice heatmap'''
    if mode == 'special':
        if vmin is None or vmax is None:
            vmin, vmax = data_.min() - 0.0001, data_.max() + 0.0001
            midpoint = data_.mean()
        else:
            midpoint = data_.mean()
    else:
        vmin, vmax = -1.0001, 1.0001
        midpoint = 0.0
    levels = np.linspace(vmin, vmax, num_levels_)
    midp = np.mean(np.c_[levels[:-1], levels[1:]], axis=1)
    vals = np.interp(midp, [vmin, midpoint, vmax], [0, 0.5, 1])
    colors = plt.cm.seismic(vals)
    # == 6 option is photocopy friendly, == 7 is almost photocopy friendly
    if num_levels_ == 6:
        colors = [[215, 25, 28, 255], [253, 174, 97, 255], [255, 255, 191, 255],
                  [171, 221, 164, 255], [43, 131, 186, 255]]
        colors = [[x / 255.0 for x in c] for c in colors]
    if num_levels_ == 7:
        colors = [[213, 62, 79, 255], [252, 141, 89, 255], [254, 224, 139, 255],
                  [230, 245, 152, 255], [153, 213, 148, 255], [50, 136, 189, 255]]
        colors = [[x / 255.0 for x in c] for c in colors]
    # print(colors)
    cmap, norm = from_levels_and_colors(levels, colors)
    fig, ax = plt.subplots()
    im = ax.imshow(data_, cmap=cmap, norm=norm, interpolation='none')
    if not removebar:
        fig.colorbar(im)

    ax.xaxis.tick_top()
    ax.set_xlabel(xlabel)
    ax.set_ylabel(ylabel)
    ax.xaxis.set_label_position('top')
    if xlabels is not None:
        ax.set_xticklabels(xlabels, minor=False)
    if ylabels is not None:
        ax.set_yticklabels(ylabels, minor=False)
    if changeTicks:
        plt.xticks(np.arange(len(data_)), rotation=rotation, fontsize=14)
        plt.yticks(np.arange(len(data_)), fontsize=14)

    if annotation is not None:
        for t in annotation:
            if len(t[1]) == 1:
                x = -50
            else:
                x = -35
            ax.annotate(t[1], xy=(0, t[0]), xytext=(x, t[0]),
                        arrowprops=dict(facecolor='black', shrink=0.05),
                        annotation_clip=False
                        )

    if showplot:
        plt.show()

def create_color_map():
    # this is the color scale from Neu et al 2013 (BAMS)
    levels = numpy.array([0,0.0001,0.001,0.01,0.1,1.0,2.0,3.0,4.0,5.0,6.0], 'f')
    cmap = ["#ffffff", "#f8e73d", "#e9a833", "#009a4a",
            "#00aeca", "#0088b7", "#295393", "#2e3065", 
            "#973683", "#ff0000", "#999999"]
    ccmap, norm = col.from_levels_and_colors(levels, cmap, 'max')
    return ccmap, norm, levels

def testplot(cats, catsavg, xy, data, levels = np.arange(0,320,20), title=""):
    """Quick test plot layout for this example file
    """
    colors = plt.cm.spectral(np.linspace(0,1,len(levels)) )    
    mycmap, mynorm = from_levels_and_colors(levels, colors, extend="max")

    radolevels = levels.copy()#[0,1,2,3,4,5,10,15,20,25,30,40,50,100]
    radocolors = plt.cm.spectral(np.linspace(0,1,len(radolevels)) )    
    radocmap, radonorm = from_levels_and_colors(radolevels, radocolors, extend="max")

    fig = plt.figure(figsize=(12,12))
    # Average rainfall sum
    ax = fig.add_subplot(111, aspect="equal")
    wradlib.vis.add_lines(ax, cats, color='black', lw=0.5)
    catsavg[np.isnan(catsavg)]=-9999
    coll = PolyCollection(cats, array=catsavg, cmap=mycmap, norm=mynorm, edgecolors='none')
    ax.add_collection(coll)
    ax.autoscale()
    cb = plt.colorbar(coll, ax=ax, shrink=0.75)
    plt.xlabel("UTM 51N Easting")
    plt.ylabel("UTM 51N Northing")
#    plt.title("Subcatchment rainfall depth")
    plt.grid()
    plt.draw()
    plt.savefig(r"E:\docs\_projektantraege\SUSTAIN_EU_ASIA\workshop\maik\pampanga_cat.png")
#    plt.close()
    # Original RADOLAN data
    fig = plt.figure(figsize=(12,12))
    # Average rainfall sum
    ax1 = fig.add_subplot(111, aspect="equal")
    pm = plt.pcolormesh(xy[:, :, 0], xy[:, :, 1], np.ma.masked_invalid(data), cmap=radocmap, norm=radonorm)
    wradlib.vis.add_lines(ax1, cats, color='black', lw=0.5)
    plt.xlim(ax.get_xlim())
    plt.ylim(ax.get_ylim())
    cb = plt.colorbar(pm, ax=ax1, shrink=0.75)
    cb.set_label("(mm/h)")
    plt.xlabel("UTM 51N Easting")
    plt.ylabel("UTM 51N Northing")
#    plt.title("Composite rainfall depths")
    plt.grid()
    plt.draw()
    plt.savefig(r"E:\docs\_projektantraege\SUSTAIN_EU_ASIA\workshop\maik\pampanga_comp.png")
    plt.close()

def setup_colormap_with_zeroval(vmin, vmax, nlevs=5,
                                cmap=plt.get_cmap('Greens'),
                                extend='both'):
    """create a discrete colormap with reserved level for small values

    setup a colormap based on a existing colormap with a specified
    number N of levels reserving the lowest colormap level for
    extremely small values.

    Extremely small are currently defined as [0.0, 1e-8].

    Returns the N-level colormap and a normalizer to plot arbitrary
    data using the colormap.  The normalizing places the N levels at
    constant intervals between the vmin and vmax.

    ARGS:
        vmin (float): value to map to the lowest color level
        vmax (float): value to map to the highest color level
        nlevs (integer): number of levels for the colormap (default is 5)
        cmap (:class:`matplotlib.colors.Colormap` instance): colormap to use for the N-level colormap
        extend (string): ({'both'} | 'min' | 'max' | 'neither') should the top or bottom of the colormap use an arrow to indicate "and larger" or "and smaller"

    RETURNS:
        tuple (cmap, norm) containing a
        :class:`matplotlib.colors.Colormap` object and a
        :class:`matplotlib.colors.Normalize object`.

    EXAMPLE:
        >>> import matplotlib.pyplot as plt
        >>> import numpy as np
        >>> from timutils.colormap_nlevs import setup_colormap_with_zeroval
        >>> data = np.random.rand(100, 100)
        >>> mycmap, mynorm = setup_colormap_with_zeroval(vmin=data.min(), vmax=data.max(), nlevs=7, extend='neither')
        >>> fig, ax = plt.subplots()
        >>> cm = ax.pcolormesh(data, cmap=mycmap, norm=mynorm)
        >>> plt.colorbar(cm)
        >>> plt.title(('setup_colormap_with_zeroval example\\n data values 0 to 1; nlevs=5'))
        >>> plt.show()

    """

    # Pick some of the nicer colors from the palette...
    if extend is "neither":
        ncolors = nlevs
    elif (extend is "min") or (extend is "max"):
        ncolors = nlevs + 1
    elif extend is "both":
        ncolors = nlevs + 2
    levels = np.concatenate((np.array([0.0, 1e-8]),
                             np.linspace(start=vmin,
                                         stop=vmax,
                                         num=nlevs)[1:]))
    colors = cmap(np.linspace(start=0.0, stop=1.0, num=ncolors))
    cmap, norm = from_levels_and_colors(levels, colors, extend=extend)
    return((cmap, norm))

def main():
    bathymetry_path = ""
    topo_path = "/RECH2/huziy/coupling/coupled-GL-NEMO1h_30min/geophys_452x260_directions_new_452x260_GL+NENA_0.1deg_SAND_CLAY_LDPT_DPTH.fst"




    plot_utils.apply_plot_params()

    with RPN(topo_path) as r:
        assert isinstance(r, RPN)
        topo = r.get_first_record_for_name("ME")
        lons, lats = r.get_longitudes_and_latitudes_for_the_last_read_rec()

        print(lons.shape)

        prj_params = r.get_proj_parameters_for_the_last_read_rec()
        rll = RotatedLatLon(**prj_params)
        bmap = rll.get_basemap_object_for_lons_lats(lons2d=lons, lats2d=lats, resolution="i")


    xx, yy = bmap(lons, lats)

    plt.figure()
    ax = plt.gca()

    lons1 = np.where(lons <= 180, lons, lons - 360)
    topo = maskoceans(lons1, lats, topo)


    topo_clevs = [0, 100, 200, 300, 400, 500, 600, 800, 1000, 1200]
    # bn = BoundaryNorm(topo_clevs, len(topo_clevs) - 1)
    cmap = cm.get_cmap("terrain")

    ocean_color = cmap(0.18)




    cmap, norm = colors.from_levels_and_colors(topo_clevs, cmap(np.linspace(0.3, 1, len(topo_clevs) - 1)))


    add_rectangle(ax, xx, yy, margin=20, edge_style="solid")
    add_rectangle(ax, xx, yy, margin=10, edge_style="dashed")



    im = bmap.pcolormesh(xx, yy, topo, cmap=cmap, norm=norm)
    bmap.colorbar(im, ticks=topo_clevs)
    bmap.drawcoastlines(linewidth=0.3)
    bmap.drawmapboundary(fill_color=ocean_color)
    bmap.drawparallels(np.arange(-90, 90, 10), labels=[1, 0, 0, 1], color="0.3")
    bmap.drawmeridians(np.arange(-180, 190, 10), labels=[1, 0, 0, 1], color="0.3")
    plt.savefig("GL_452x260_0.1deg_domain.png", dpi=300, bbox_inches="tight")

def get_colormap_norm(cmap_type, levels):
    """Create a custom colormap."""
    if cmap_type == "rain":
        cmap, norm = from_levels_and_colors(levels, sns.color_palette("Blues", n_colors=len(levels)),
                                                    extend='max')
    elif cmap_type == "snow":
        cmap, norm = from_levels_and_colors(levels, sns.color_palette("PuRd", n_colors=len(levels)),
                                                    extend='max')
    elif cmap_type == "snow_discrete":    
        colors = ["#DBF069","#5AE463","#E3BE45","#65F8CA","#32B8EB",
                    "#1D64DE","#E97BE4","#F4F476","#E78340","#D73782","#702072"]
        cmap, norm = from_levels_and_colors(levels, colors, extend='max')
    elif cmap_type == "rain_acc":    
        cmap, norm = from_levels_and_colors(levels, sns.color_palette('gist_stern_r', n_colors=len(levels)),
                         extend='max')
    elif cmap_type == "rain_new":
        colors_tuple = pd.read_csv('/home/mpim/m300382/icon_forecasts/cmap_prec.rgba').values    
        cmap, norm = from_levels_and_colors(levels, sns.color_palette(colors_tuple, n_colors=len(levels)),
                         extend='max')

    return(cmap, norm)

def test_cmap_and_norm_from_levels_and_colors():
    data = np.linspace(-2, 4, 49).reshape(7, 7)
    levels = [-1, 2, 2.5, 3]
    colors = ['red', 'green', 'blue', 'yellow', 'black']
    extend = 'both'
    cmap, norm = mcolors.from_levels_and_colors(levels, colors, extend=extend)

    ax = plt.axes()
    m = plt.pcolormesh(data, cmap=cmap, norm=norm)
    plt.colorbar(m)

    # Hide the axes labels (but not the colorbar ones, as they are useful)
    ax.tick_params(labelleft=False, labelbottom=False)

def test_cmap_and_norm_from_levels_and_colors():
    data = np.linspace(-2, 4, 49).reshape(7, 7)
    levels = [-1, 2, 2.5, 3]
    colors = ['red', 'green', 'blue', 'yellow', 'black']
    extend = 'both'
    cmap, norm = mcolors.from_levels_and_colors(levels, colors, extend=extend)

    ax = plt.axes()
    m = plt.pcolormesh(data, cmap=cmap, norm=norm)
    plt.colorbar(m)

    # Hide the axes labels (but not the colorbar ones, as they are useful)
    for lab in ax.get_xticklabels() + ax.get_yticklabels():
        lab.set_visible(False)

def colorbar(nvalues,nticks=None,cmap='gnuplot', start=0.1, stop=0.8,strmap = '%.3f', shrink=1, aspect=20):
    ncolors = len(nvalues)
    if nticks==None:
        nticks = ncolors
    colors = mpltools.color.colors_from_cmap(ncolors, cmap=cmap, start=start, stop=stop)
    cmap, norm = mcolors.from_levels_and_colors(range(ncolors + 1), colors)
    mappable = matplotlib.cm.ScalarMappable(cmap=cmap)
    mappable.set_array([])
    mappable.set_clim(-0.5, ncolors+0.5)
    colorbar = plt.colorbar(mappable,shrink=shrink,aspect=aspect)
    tickpos = np.linspace(0, ncolors,nticks)
    tickvalues = np.linspace(min(nvalues), max(nvalues), nticks)
    tickvalues = [strmap%f for f in tickvalues]
    colorbar.set_ticks(tickpos)
    colorbar.set_ticklabels(tickvalues)
    return colorbar

def setup_colormap(vmin, vmax, nlevs=5,
                   cmap=plt.get_cmap('Greens'),
                   extend='both'):
    """
    create a discrete colormap from a continuous colormap

    Returns the N-level colormap and a normalizer to plot arbitrary
    data using the colormap.  The normalizing places the N levels at
    constant intervals between the vmin and vmax.

    ARGS:
        vmin (float): value to map to the lowest color level
        vmax (float): value to map to the highest color level
        nlevs (integer): number of levels for the colormap (default is 5)
        cmap (:class:`matplotlib.colors.Colormap` instance): colormap to use for the N-level colormap
        extend (string): ({'both'} | 'min' | 'max' | 'neither'}) should the top or bottom of the colormap use an arrow to indicate "and larger" or "and smaller"

    RETURNS:
        tuple (cmap, norm) containing a
        :class:`matplotlib.colors.Colormap` object and a
        :class:`matplotlib.colors.Normalize` object.

    EXAMPLES:
    >>> import matplotlib.pyplot as plt
    >>> import numpy as np
    >>> from timutils.colormap_nlevs import setup_colormap
    >>> data = np.random.rand(100, 100)  * 1000
    >>> mycmap, mynorm = setup_colormap(vmin=200.0, vmax=800.0, extend='max')
    >>> fig, ax = plt.subplots()
    >>> cm = ax.pcolormesh(data, cmap=mycmap, norm=mynorm)
    >>> plt.colorbar(cm)
    >>> plt.title(('setup_colormap example\\n data values 0 to 1000; nlevs=5, vmin=200, vmax=800'))
    >>> plt.show()

    """
    print('setting up colormaps')
    # Pick some of the nicer colors from the palette...
    if extend is "neither":
        ncolors = nlevs - 1
    elif (extend is "min") or (extend is "max"):
        ncolors = nlevs
    elif extend is "both":
        ncolors = nlevs + 1
    levels = np.linspace(start=vmin, stop=vmax, num=nlevs)
    colors = cmap(np.linspace(start=0.0, stop=1.0, num=ncolors))
    cmap, norm = from_levels_and_colors(levels, colors, extend=extend)
    return((cmap, norm))