Downloading multiple datasets using the same reduceRegion in Earth Engine

Question

I'm trying to get the average value within a week timeframe, for every band, for every dataset of Sentinel 5P on Google Earth Engine, for a 10km radius around a point.

The part I'm having trouble with is efficiently joining (currently with an inner join?) the various Sentinel 5P datasets together.

This is currently resulting in an empty dictionary as the result of .getInfo(), where as if I use reduceRegion on a single image (i.e. the mean of a single dataset filtered by date) I get a reasonable response. I've done this before by looping through the databases, but I'm assuming there is a much more efficient way.

How can I efficiently apply a reduceRegion across multiple earth engine datasets at the same time?

start_date = '01-01-2019'
next_date = '01-08-2019'
geometry = ee.Geometry.Point(lon=lon, lat=lat).buffer(1e5) #10 km around an arbitrary point

EARTHENGINE_DBS = [
    "COPERNICUS/S5P/OFFL/L3_NO2",
    "COPERNICUS/S5P/OFFL/L3_SO2",
    "COPERNICUS/S5P/OFFL/L3_CO",
    "COPERNICUS/S5P/OFFL/L3_HCHO",
    "COPERNICUS/S5P/OFFL/L3_O3",
]

# -- Make a date filter to get images in this date range. --
dateFilter = ee.Filter.date(start_date, next_date)

# --- Select imaged collection ---
ic = ee.ImageCollection(EARTHENGINE_DBS[0])
for db in EARTHENGINE_DBS[1:]:

    # --- Define the new data ---
    _new = ee.ImageCollection(db).filter(dateFilter)

    # --- Define an inner join. ---
    innerJoin = ee.Join.inner()

    # --- Specify an equals filter for image timestamps. ---
    filterTimeEq = ee.Filter.equals(** {
            'leftField': 'system:time_start',
            'rightField': 'system:time_start'
    })

    # --- Apply the join. ---
    ic = innerJoin.apply(ic, _new, filterTimeEq)

    # --- flatten ---
    ic = ic.map(lambda feature: ee.Image.cat(feature.get('primary'), feature.get('secondary')))

# --- aggregate ic ---
agg_ic = ic.mean()

mean_dict = agg_ic.reduceRegion(**{
    'reducer': ee.Reducer.mean(),
    'geometry':geometry,
    'scale': 1000,
    'bestEffort':True,
    'tileScale':16
})

mean_dict.getInfo()

Answer 1

Instead of joining, you could perhaps merge all the image collections and calculate your aggregate from that. Each image in that collection must (?) have the same bands. You can map over each collection and fill in missing bands with masked out images. This might get you started:

import ee

ee.Initialize()

lon = 0
lat = 10
start_date = '2019-01-01'
next_date = '2019-01-08'
geometry = ee.Geometry.Point(lon=lon, lat=lat).buffer(1e4)

EARTHENGINE_DBS = [
    "COPERNICUS/S5P/OFFL/L3_NO2",
    "COPERNICUS/S5P/OFFL/L3_SO2",
    "COPERNICUS/S5P/OFFL/L3_CO",
    "COPERNICUS/S5P/OFFL/L3_HCHO",
    "COPERNICUS/S5P/OFFL/L3_O3",
]

collections = ee.List([
    ee.ImageCollection(ee.String(name))
            .filterDate(start_date, next_date)
            .filterBounds(geometry)   
            .map(lambda image: image.float())
    for name in EARTHENGINE_DBS
])

band_names = ee.List(collections.iterate(
    lambda c, acc: ee.List(acc).cat(ee.ImageCollection(c).first().bandNames()),
    ee.List([]))
).distinct()

def add_additional_bands(image):
    empty_bands = ee.List.repeat(ee.Image(), band_names.size())
    empty_image = ee.ImageCollection(empty_bands).toBands().rename(band_names)
    return image.float().addBands(empty_image).select(band_names).float()

collection = ee.ImageCollection(collections \
    .iterate(
        lambda collection, acc: ee.ImageCollection(acc)
            .merge(ee.ImageCollection(collection).map(add_additional_bands)),
        ee.ImageCollection([])
    )
)

agg_ic = collection.mean()
mean_dict = agg_ic.reduceRegion(
    reducer=ee.Reducer.mean(),
    geometry=geometry,
    scale=1000,
    bestEffort=True
)

mean_dict.getInfo()