Band math question

[gagecarto]

I am working to quantify a change between a current and future dataset. I am already using TiTiler for tiles and some basic stats.

For example purposes, imagine a 2 band, categorical raster which has integer values 1-5. Band one represents current time and band two is future. Below is a tabular representation of 15 cells (b1 and b2).

What I'd like to summarize is how many of each value changed to another value. The table below should illustrate this concept.

Any ways to accomplish something like this in TiTiler?

[zacdezgeo]

Hello!

TiTiler is an excellent tool for tiling and visualizing geospatial raster data, but there might be better fits for directly calculating these types of statistics. Some functionality for summary statistics is available, but I don't see a direct path that would cover your use case, given you are interested in the counts from one category to another. There may be some complicated way of changing the values of your first band and applying some band calculations, like in this example. I would not necessarily recommend the approach.

Instead, you could use Rasterio to read your COG and calculate your needed statistics. Since I don't have access to your data, I can't provide a tailored solution, but I've drafted a Python script below that demonstrates a general approach to achieve this:

import rasterio
import numpy as np
import csv

def calculate_change_matrix(current_band, future_band):
    change_matrix = np.zeros((5, 5), dtype=int)
    for current_value, future_value in zip(np.nditer(current_band), np.nditer(future_band)):
        if 1 <= current_value <= 5 and 1 <= future_value <= 5:
            change_matrix[current_value-1, future_value-1] += 1
    return change_matrix

def write_change_matrix_to_csv(change_matrix, csv_file_path):
    with open(csv_file_path, mode='w', newline='') as file:
        writer = csv.writer(file)
        writer.writerow(['From', 'To', 'Count'])
        for i in range(change_matrix.shape[0]):
            for j in range(change_matrix.shape[1]):
                writer.writerow([i + 1, j + 1, change_matrix[i, j]])

# Paths to your COGs
cog_current_path = 'path_to_cog.tif'

# Read the first and second band
with rasterio.open(cog_current_path) as src:
    current_band = src_current.read(1)  # Read the first band 
    future_band = src_future.read(2)  # Read the second band

# Calculate the change matrix
change_matrix = calculate_change_matrix(current_band, future_band)

# Write the change matrix to a CSV file
csv_file_path = 'data/change_matrix.csv'
write_change_matrix_to_csv(change_matrix, csv_file_path)

That being said, TiTiler could be very useful for visualizing where these changes have occurred on a map. It allows you to overlay and compare your two images, offering a visual insight into the changes. For an example of how TiTiler can be used for such visual comparisons, you might find this notebook interesting. It compares satellite imagery before and after an earthquake in Turkey.

I hope this helps!

[vincentsarago]

As @zacharyDez mentioned, this is not really a task for TiTiler, but you could totally add a new endpoint which will do the processing and show the results.

That's said, TiTiler is best when reading part of dataset. What you're trying to achieve is to work on a full dataset (at full resolution), which will mean that you'll put a lot in memory and might have to let the user wait a bit before having the result ready.

[gagecarto]

This is all really helpful.. I was thinking there might be some band math magic but @zacharyDez's Rasterio solution could work really well. Thanks as always for such a great set of tools!