Zooming in with Satellites

This lesson shows how to use Earth Engine and geemap to find and visualize the clearest available satellite image over a site near a specified date. By adjusting filters and using sorting, you can display the most usable image (with minimal clouds obscuring the image).

import folium
import ee
import geemap
from datetime import datetime, timedelta

# Authenticate your Google account with Earth Engine
ee.Authenticate()

# Write your project ID here, in quotes
ee.Initialize(project = "emerge-lessons")

# Define your point of interest
latitude = 27.536873
longitude = -81.469549

map = folium.Map(location=[latitude, longitude], tiles="Cartodb dark_matter", zoom_start=9)

Note: Above, we used the folium Python library to make the map. Another option is to use the following:

Map = geemap.Map(center=(lat, lon), zoom=10)

This would have a similar result. We used folium to make the map instead of geemap because folium allows the map to remain displayed even after the notebook is uploaded online, whereas the geemap.Map() method may not be as easily formatted.

Like in Chapter 1, we need to define a function (from this tutorial) to add Google Earth Engine data to a map in a way that allows it to be interactively displayed.

def add_ee_layer(self, ee_image_object, vis_params, name):
    """Adds a method for displaying Earth Engine image tiles to folium map."""
    map_id_dict = ee.Image(ee_image_object).getMapId(vis_params)
    folium.raster_layers.TileLayer(
        tiles=map_id_dict['tile_fetcher'].url_format,
        attr='Map Data &copy; <a href="https://earthengine.google.com/">Google Earth Engine</a>',
        name=name,
        overlay=True,
        control=True
    ).add_to(self)

folium.Map.add_ee_layer = add_ee_layer

# Create a point at the longitude and latitude
point = ee.Geometry.Point([longitude, latitude])

# Choose your target date (when you'd like to be "closest" to)
target_date_str = '2024-01-01'
target_date = datetime.strptime(target_date_str, '%Y-%m-%d')

# Define a time window (e.g. 60 days before and after)
# You can make this wider or narrower depending on how many images you want to consider
window_days = 60
start_date = (target_date - timedelta(days=window_days)).strftime('%Y-%m-%d')
end_date = (target_date + timedelta(days=window_days)).strftime('%Y-%m-%d')

# Build a Sentinel-2 ImageCollection:
#  - Only images over your point
#  - Within your date window
#  - Sorted so the least cloudy image comes first
collection = (
    ee.ImageCollection("COPERNICUS/S2_SR_HARMONIZED")
    .filterBounds(point)                        # Covers the point
    .filterDate(start_date, end_date)           # Within date window
    .sort('CLOUDY_PIXEL_PERCENTAGE')            # Sort least cloudy first
)

# Grabs the least cloudy image in your window
image = collection.first()

# Only add the image to the map if there is one found in the window
if image.getInfo():  # Will be False if there are no images found at all
    # Get the image's capture date
    image_date = ee.Date(image.get("system:time_start")).format("YYYY-MM-dd").getInfo()
    # Get the cloudiness percentage
    cloud_pct = image.get('CLOUDY_PIXEL_PERCENTAGE').getInfo()
    # Visualization parameters: true color, stretch from 0 to 3000 as usual for S2
    vis_params = {'bands': ['B4', 'B3', 'B2'], 'min': 0, 'max': 3000}
    # Add the least cloudy image, layer name shows its date and cloud %
    map.add_ee_layer(image, vis_params, f'Least Cloudy: {image_date} ({cloud_pct:.2f}% cloud)')

# Show the map
display(map)

Make this Notebook Trusted to load map: File -> Trust Notebook