Image Metadata for Timeseries

Hi @bluesky314 ,

Since you are using mosaicking, the output images could be a mosaic from multiple tiles if the area of interest is not fully covered by a single tile.

If you are interested in of which tiles the mosaic is made, I suggest using Mosaicking.TILE and doing the mosaicking with evalcscript. In this case the metadata of the tiles used to create the output mosaic can be extracted by the updateOutputMetadata function.

Here’s an example evalcsript to do the mosaicking by selecting the tile with least cloud coverage:

//VERSION=3
var scenes_list = [];
function setup() {
    return {
        input: [{
            // Add the dataMask band for data filtering
            bands: ["B02", "B03", "B04", "dataMask"]
        }],
        output: {
            bands: 3
        },
        // Use TILE instead of SIMPLE so the metadata can be accessed by `scenes` object
        mosaicking: Mosaicking.TILE
    };
}

// Extract the metadata of the tiles
function updateOutputMetadata(scenes, inputMetadata, outputMetadata) {
    let metadata = [];
    for (i=0; i<scenes_list.length; i++) {
      let acquisition_object = {
        "date": {},
        "tile_id":{},
        "cloud_coverage":{}
      };
      Object.assign(acquisition_object, {"date": scenes.tiles[i].date});
      Object.assign(acquisition_object, {"tile_id": scenes.tiles[i].tileOriginalId});
      Object.assign(acquisition_object, {"cloud_coverage": scenes.tiles[i].cloudCoverage});
      metadata.push(acquisition_object)
    }
    outputMetadata.userData = {
      "metadata": metadata
    }
}

function evaluatePixel(samples) {
    let index_for_mosaicking = get_index_for_mosaicking(samples);
    if (!scenes_list.includes(index_for_mosaicking)) {
      scenes_list.push(index_for_mosaicking)
    }
    return [samples[index_for_mosaicking].B04, samples[index_for_mosaicking].B03, samples[index_for_mosaicking].B02];
}

// Define a function to get the index of sample with the least cloud coverage
function get_index_for_mosaicking(samples) {
    let sample_with_data = [];
    for (i=0; i < samples.length; i++) {
      // Select the index only if there is data
      if (samples[i].dataMask == 1) {
      sample_with_data.push(i);
      }
    }
    // Return the first one which has the least cloud coverage
    return sample_with_data[0];
}

With the above evalscript, you can have the image and the metadata of the tiles used to generate the output image in a JSON file by editing the python script as following.

To display the images:

def get_true_color_request(time_interval):
    return SentinelHubRequest(
        evalscript=evalscript_true_color,
        input_data=[
            SentinelHubRequest.input_data(
                data_collection=DataCollection.SENTINEL2_L1C,
                time_interval=time_interval,
                mosaicking_order='leastCC'
            )
        ],
        responses=[
            SentinelHubRequest.output_response('default', MimeType.PNG),
            SentinelHubRequest.output_response('userdata', MimeType.JSON)
        ],
        bbox=betsiboka_bbox,
        size=betsiboka_size,
        config=config
    )

list_of_requests = [get_true_color_request(slot) for slot in slots]
list_of_requests = [request.download_list[0] for request in list_of_requests]

data = SentinelHubDownloadClient(config=config).download(list_of_requests, max_threads=5)

ncols = 4
nrows = 3
aspect_ratio = betsiboka_size[0] / betsiboka_size[1]
subplot_kw = {'xticks': [], 'yticks': [], 'frame_on': False}

fig, axs = plt.subplots(ncols=ncols, nrows=nrows, figsize=(5 * ncols * aspect_ratio, 5 * nrows),
                        subplot_kw=subplot_kw)

for idx, monthly_data in enumerate(data):
    image = monthly_data['default.png']
    ax = axs[idx // ncols][idx % ncols]
    ax.imshow(np.clip(image * 2.5/255, 0, 1))
    ax.set_title(f'{slots[idx][0]}  -  {slots[idx][1]}', fontsize=10)

plt.tight_layout()

To have the metadata of the tiles used to generate the output images:

metadata = {}
for idx, monthly_data in enumerate(data):
    acquisition_dates = monthly_data['userdata.json']
    metadata[f'{slots[idx][0]}-{slots[idx][1]}'] = acquisition_dates['metadata']
print(metadata)

Note: the cloud coverage in the metadata is the cloud coverage of that specific tile, not the cloud coverage of the output image.

To get the cloud coverage of the output image, you may request for SCL or CLM band and compute the cloud coverage from the output data.