Of course! Showing a GeoTIFF in Python is a common task for geospatial analysis. The best tool for the job depends on whether you want a simple static plot or an interactive map.

Here’s a comprehensive guide covering the most popular libraries: rasterio (for static plots) and folium (for interactive maps).

Prerequisites

First, you need to install the necessary libraries. rasterio is essential for reading the GeoTIFF file, and matplotlib is used for plotting.

pip install rasterio matplotlib folium

You'll also need a sample GeoTIFF file to work with. For this example, we'll use a public Digital Elevation Model (DEM) of the Mount St. Helens area. You can download it from here:

• Download Link: mt_st_helens_cog.tif
• Place this file in the same directory as your Python script or notebook.

Method 1: Static Plot with rasterio and matplotlib

This is the most common and straightforward method for creating publication-quality static images. It's great for reports, papers, and presentations.

Step-by-Step Code

import rasterio
from rasterio.plot import show
import matplotlib.pyplot as plt
# --- 1. Open the GeoTIFF file ---
# The 'with' statement ensures the file is automatically closed
try:
    with rasterio.open('mt_st_helens_cog.tif') as src:
        # --- 2. Read the first band ---
        # Most single-band GeoTIFFs (like DEMs or satellite images) store data in the first band
        dem = src.read(1)
        # --- 3. Plot the data using rasterio's built-in show function ---
        # This is a convenient wrapper around matplotlib
        show(dem, title='Mount St. Helens Elevation')
        # --- 4. Display the plot ---
        plt.show()
except FileNotFoundError:
    print("Error: 'mt_st_helens_cog.tif' not found.")
    print("Please download the file and place it in the same directory as your script.")
except Exception as e:
    print(f"An error occurred: {e}")

Explanation

1. rasterio.open('...'): Opens the GeoTIFF file. This gives you a DatasetReader object, which contains all the metadata (like coordinate system, bounds, etc.) and access to the raster data.
2. src.read(1): Reads the data from the first band of the raster. For a color image, you would read multiple bands (e.g., src.read([1, 2, 3]) for RGB).
3. rasterio.plot.show(): This is the key function. It intelligently uses the raster's metadata (specifically the transform) to correctly place the raster in its geographic coordinate space on the plot. It automatically creates a color ramp based on the data's min/max values.
4. plt.show(): Displays the plot window generated by Matplotlib.

Enhancing the Static Plot

The basic plot is good, but we can make it much better by adding a color bar, a custom colormap, and proper axis labels.

import rasterio
from rasterio.plot import show
import matplotlib.pyplot as plt
# Use a context manager to ensure the file is closed
with rasterio.open('mt_st_helens_cog.tif') as src:
    dem = src.read(1)
    # Create a figure and axes for more control
    fig, ax = plt.subplots(figsize=(10, 8))
    # Define a colormap (e.g., 'terrain' is great for elevation data)
    cmap = plt.cm.terrain
    # Plot the data with a custom colormap
    # The 'ax' parameter tells matplotlib to plot on our specific axes
    im = show(dem, ax=ax, cmap=cmap, title='Mount St. Helens Elevation')
    # Add a colorbar
    # The im object is the "image" created by the show function
    fig.colorbar(im, ax=ax, label='Elevation (meters)')
    # Remove axis ticks and labels for a cleaner map look
    ax.set_xticks([])
    ax.set_yticks([])
    ax.set_xlabel('')
    ax.set_ylabel('')
    # Display the plot
    plt.show()

Method 2: Interactive Map with folium

For web-based applications or when you want to pan, zoom, and get information on hover, an interactive map is perfect. folium integrates beautifully with GeoTIFF data.

Step-by-Step Code

import rasterio
import folium
import numpy as np
from rasterio.plot import show_hist
# --- 1. Open the GeoTIFF and get metadata ---
with rasterio.open('mt_st_helens_cog.tif') as src:
    # Read the first band
    dem = src.read(1)
    # Get the bounds of the raster
    bounds = src.bounds
    # Get the centroid for the initial map center
    # lon, lat order for folium
    center_lon = (bounds.left + bounds.right) / 2
    center_lat = (bounds.bottom + bounds.top) / 2
# --- 2. Create a Folium map centered on the raster ---
# The zoom_start level can be adjusted
m = folium.Map(
    location=[center_lat, center_lon],
    zoom_start=10,
    tiles='OpenStreetMap' # You can use other tile layers like 'Stamen Terrain'
)
# --- 3. Add the GeoTIFF as an overlay ---
# We need to convert the data to a format folium can use (PNG)
# We'll normalize the data to 0-255 for the image
def normalize(array):
    """Normalizes a numpy array to a 0-255 range."""
    min_val = np.nanmin(array)
    max_val = np.nanmax(array)
    return ((array - min_val) / (max_val - min_val) * 255).astype(np.uint8)
normalized_dem = normalize(dem)
# Add the image to the map
# The 'opacity' parameter makes the map semi-transparent over the basemap
folium.raster_layers.ImageOverlay(
    image=normalized_dem,
    bounds=[[bounds.bottom, bounds.left], [bounds.top, bounds.right]],
    colormap=plt.cm.terrain, # Use the same colormap as before
    opacity=0.6,
    name='Elevation Data'
).add_to(m)
# Add a layer control to toggle the overlay
folium.LayerControl().add_to(m)
# --- 4. Display the map ---
# In a Jupyter Notebook, this will display the map directly.
# If running as a .py file, you may need to save it to an HTML file.
m.save('map.html') # Save to an HTML file to view in a browser
# m # In a Jupyter Notebook, just typing 'm' will display it

Explanation

1. rasterio.open() and Metadata: We open the file and extract the raster data (dem) and, crucially, its bounds. The bounds tell folium exactly where to place the image on the map.
2. folium.Map(): Creates a basic interactive map. We set the initial center point using the centroid of our raster's bounds.
3. normalize() function: folium's ImageOverlay expects image data in a standard format (like 0-255 for grayscale). Our elevation data is in meters, so we must scale it.
4. folium.raster_layers.ImageOverlay(): This is the function that overlays our raster data onto the map.
   • image: The normalized NumPy array.
   • bounds: A list of [southwest, northeast] corner coordinates. This is how folium georeferences the image.
   • colormap: Applies a color gradient to the single-band image, making it look like a proper elevation map.
   • opacity: Makes the overlay semi-transparent, so you can see the underlying map tiles.
5. folium.LayerControl(): Adds a small widget in the top right corner to toggle the visibility of your overlay.
6. m.save('map.html'): Saves the interactive map to an HTML file, which you can open in any web browser.

Summary: Which Method to Use?