Calculating Mean of Divided LANDSAT C2 ARD Earth Explorer Tiles

I need to calculate the mean of LANDSAT C2 ARD Earth Explorer image data such as evapotranspiration. The tif tiles overlap each other and cover more than 1 degree latitude/1 degree longitude, so I need to divide the area of interest I am assigned with into blocks smaller than that and calculate the mean of the image data for each block.

My area of interest covers 7 degrees longitude and 4 degrees latitude. That means I can divide it into 28 or more blocks. My blocks can be 0.5x0.5 degrees, 1x1 degrees, etc.; it doesn't matter as long as they are divided evenly.

I am very much a beginner at Matlab. I just started programming this year and hadn't used Matlab before I got this job with my university. These are the instructions from my professor:

"Matrix 1 (size: N x 3) is for image data

column 1: latitude

column 2: longitude

column 3: image data (e.g., LST)

Matrix 2 is for block location (latitude)

Assuming there are N x M blocks after you divide the entire area evenly.

matrix size: N x M

The value of (i,j)-th element is the latitude of the top-left corner.

Matrix 3 is the same as Matrix 2, but for longitude.

Matrix 4 is similar to Matrix 2, but for the width of the block along latitude.

Matrix 5 is the same as Matrix 4, but for the length of the block along longitude.

With the above 5 matrices, you should be able to calculate any statistics (of LST) within every blocks. Then you will produce a final N x M matrix of the statistics."

I have the corners of my area of interest saved as parameters. I think I understand what to do for matricies 2-5, but I'm not sure what to do for matrix 1. My tiles are 5000x5000 matricies, so how would I condense/move that into column 3 of matrix 1? Columns 1 and 2 of that matrix will have the latitude and longitude of that particular tile/image. I eventually need to make a line of code that divides the image into blocks and calculates the mean of the image data for each block automatically, all in one. aaaaah

Answer 1

[Abhishek Chakram]

Hi Eric Williams.

It is my understanding that you want to calculate the mean of divided LANDSAT C2 ARD Earth Explorer Tiles. Here is a sample code for the same:

% Define random values for the matrices
numBlocksLatitude = 4;
numBlocksLongitude = 7;

blockLatitude = [35, 35.5, 36, 36.5; 37, 37.5, 38, 38.5; 39, 39.5, 40, 40.5; 41, 41.5, 42, 42.5];
blockLongitude = [-120, -119.5, -119, -118.5, -118, -117.5, -117; -120, -119.5, -119, -118.5, -118, -117.5, -117; -120, -119.5, -119, -118.5, -118, -117.5, -117; -120, -119.5, -119, -118.5, -118, -117.5, -117];
blockWidth = [0.5, 0.5, 0.5, 0.5; 0.5, 0.5, 0.5, 0.5; 0.5, 0.5, 0.5, 0.5; 0.5, 0.5, 0.5, 0.5];
blockLength = [1, 1, 1, 1, 1, 1, 1; 1, 1, 1, 1, 1, 1, 1; 1, 1, 1, 1, 1, 1, 1; 1, 1, 1, 1, 1, 1, 1];

% Generate random image data
imageData = [rand(5000, 1) * 10 + 30, rand(5000, 1) * 10 - 120, rand(5000, 1)];

% Calculate the mean within each block
meanData = zeros(numBlocksLatitude, numBlocksLongitude);

for i = 1:numBlocksLatitude
    for j = 1:numBlocksLongitude
        % Get the coordinates of the current block
        lat = blockLatitude(i, j);
        lon = blockLongitude(i, j);
        width = blockWidth(i, j);
        length = blockLength(i, j);
        
        % Determine the indices of the current block within the image data
        latIndices = imageData(:, 1) >= lat & imageData(:, 1) < lat + width;
        lonIndices = imageData(:, 2) >= lon & imageData(:, 2) < lon + length;
        blockIndices = latIndices & lonIndices;
        
        % Extract the image data within the current block
        blockData = imageData(blockIndices, 3);
        
        % Calculate the mean of the image data within the current block
        meanData(i, j) = mean(blockData);
    end
end

% Display the mean data
disp(meanData);

In this example, “numBlocksLatitude” and “numBlocksLongitude” is traverserd loop and mean is calculated which is further stored in the array “meanData”.

Hope this answers your question.

Best Regards,

Abhishek Chakram