How to use the earthpy.plot function in earthpy

# Landsat 8 red band is band 4 at [3]
# Landsat 8 near-infrared band is band 5 at [4]
ndvi = es.normalized_diff(arr_st[4], arr_st[3])


###############################################################################
# Plot NDVI With Colorbar Legend of Continuous Values
# ----------------------------------------------------

# You can plot NDVI with a colorbar legend of continuous values using the
# ``plot_bands`` function from the ``earthpy.plot`` module.

titles = ["Landsat 8 - Normalized Difference Vegetation Index (NDVI)"]

# Turn off bytescale scaling due to float values for NDVI
ep.plot_bands(
    ndvi, cmap="RdYlGn", cols=1, title=titles, scale=False, vmin=-1, vmax=1
)


###############################################################################
# Classify NDVI
# -------------

# Next, you can classify the NDVI to categorize the results into useful classes.
# Values under 0 will be classified together as no vegetation. Additional classes
# will be created for bare area # and low, moderate, and high vegetation areas.

# Create classes and apply to NDVI results
ndvi_class_bins = [-np.inf, 0, 0.1, 0.25, 0.4, np.inf]
ndvi_landsat_class = np.digitize(ndvi, ndvi_class_bins)

# ---------------
# To begin, open your DEM layer as a numpy array using Rasterio. Below you set all
# terrain values < 0 to ``nan``. Next, plot the data using ``ep.plot_bands()``.

# Set the home directory and get the data for the exercise
os.chdir(os.path.join(et.io.HOME, "earth-analytics"))
dtm = "data/vignette-elevation/pre_DTM.tif"

# Open the DEM with Rasterio
with rio.open(dtm) as src:
    elevation = src.read(1)
    # Set masked values to np.nan
    elevation[elevation < 0] = np.nan

# Plot the data
ep.plot_bands(
    elevation,
    scale=False,
    cmap="gist_earth",
    title="DTM Without Hillshade",
    figsize=(10, 6),
)
plt.show()

####################################################################################
# Create the Hillshade
# --------------------
# Once the DEM is read in, call ``es.hillshade()`` to create the hillshade.

# Create and plot the hillshade with earthpy
hillshade = es.hillshade(elevation)

)
plt.show()

####################################################################################
# Change the Angle Altitude of the Sun
# -------------------------------------
# Another variable you can adjust for hillshade is what angle of the sun.
# The ``angle_altitude`` parameter values range from 0 to 90. 90 represents the sun
# shining from directly above the scene. The default value for ``angle_altitude`` in
# ``es.hillshade()`` is 30 degrees.

# Adjust the azimuth value
hillshade_angle_10 = es.hillshade(elevation, altitude=10)

# Plot the hillshade layer with the modified angle altitude
ep.plot_bands(
    hillshade_angle_10,
    scale=False,
    cbar=False,
    title="Hillshade with Angle Altitude set to 10 Degrees",
    figsize=(10, 6),
)
plt.show()

####################################################################################
# Overlay a DEM on top of the Hillshade
# -------------------------------------
# A hillshade can be used to visually enhance a DEM.
# To overlay the data, use the ``ep.plot_bands()`` function in EarthPy combined with
# ``ax.imshow()``. The alpha setting sets the tranparency value for the hillshade layer.

# Plot the DEM and hillshade at the same time

# Landsat 8 red band is band 4 at [3]
# Landsat 8 near-infrared band is band 5 at [4]
ndvi = es.normalized_diff(arr_st[4], arr_st[3])


###############################################################################
# Plot NDVI With Colorbar Legend of Continuous Values
# ----------------------------------------------------
#
# You can plot NDVI with a colorbar legend of continuous values using the
# ``plot_bands`` function from the ``earthpy.plot`` module.

titles = ["Landsat 8 - Normalized Difference Vegetation Index (NDVI)"]

# Turn off bytescale scaling due to float values for NDVI
ep.plot_bands(
    ndvi, cmap="RdYlGn", cols=1, title=titles, scale=False, vmin=-1, vmax=1
)


###############################################################################
# Classify NDVI
# -------------
#
# Next, you can categorize (or classify) the NDVI results into useful classes.
# Values under 0 will be classified together as no vegetation. Additional classes
# will be created for bare area and low, moderate, and high vegetation areas.

# Create classes and apply to NDVI results
ndvi_class_bins = [-np.inf, 0, 0.1, 0.25, 0.4, np.inf]
ndvi_landsat_class = np.digitize(ndvi, ndvi_class_bins)

# plot NDVI with scaling turned off in order for the proper range of values
# (-1 to 1) to be displayed. You can use the ``cmap=`` parameter to adjust
# the colormap for the plot

NDVI = es.normalized_diff(array_stack[4], array_stack[3])
ep.plot_bands(NDVI, scale=False, cmap="RdYlGn")
plt.show()

##################################################################################
# Adjust the Number of Columns for a Multi Band Plot
# ---------------------------------------------------
#
# The number of columns used while plotting multiple bands can be changed in order
# to change the arrangement of the images overall.

ep.plot_bands(array_stack, cols=2)
plt.show()