Ordenamiento y reducción de arrays

El ordenamiento de arrays es útil para obtener mosaicos de calidad personalizados que implican reducir un subconjunto de bandas de imagen según los valores de una banda diferente. En el siguiente ejemplo, se ordena por NDVI y, luego, se obtiene el promedio de un subconjunto de observaciones de la colección con los valores de NDVI más altos:

Editor de código (JavaScript)

// Define a function that scales and masks Landsat 8 surface reflectance images
// and adds an NDVI band.
function prepSrL8(image) {
  // Develop masks for unwanted pixels (fill, cloud, cloud shadow).
  var qaMask = image.select('QA_PIXEL').bitwiseAnd(parseInt('11111', 2)).eq(0);
  var saturationMask = image.select('QA_RADSAT').eq(0);

  // Apply the scaling factors to the appropriate bands.
  var opticalBands = image.select('SR_B.').multiply(0.0000275).add(-0.2);
  var thermalBands = image.select('ST_B.*').multiply(0.00341802).add(149.0);

  // Calculate NDVI.
  var ndvi = opticalBands.normalizedDifference(['SR_B5', 'SR_B4'])
      .rename('NDVI');

  // Replace original bands with scaled bands, add NDVI band, and apply masks.
  return image.addBands(opticalBands, null, true)
      .addBands(thermalBands, null, true)
      .addBands(ndvi)
      .updateMask(qaMask)
      .updateMask(saturationMask);
}

// Define an arbitrary region of interest as a point.
var roi = ee.Geometry.Point(-122.26032, 37.87187);

// Load a Landsat 8 surface reflectance collection.
var collection = ee.ImageCollection('LANDSAT/LC08/C02/T1_L2')
  // Filter to get only imagery at a point of interest.
  .filterBounds(roi)
  // Filter to get only six months of data.
  .filterDate('2021-01-01', '2021-07-01')
  // Prepare images by mapping the prepSrL8 function over the collection.
  .map(prepSrL8)
  // Select the bands of interest to avoid taking up unneeded memory.
  .select('SR_B.|NDVI');

// Convert the collection to an array.
var array = collection.toArray();

// Label of the axes.
var imageAxis = 0;
var bandAxis = 1;

// Get the NDVI slice and the bands of interest.
var bandNames = collection.first().bandNames();
var bands = array.arraySlice(bandAxis, 0, bandNames.length());
var ndvi = array.arraySlice(bandAxis, -1);

// Sort by descending NDVI.
var sorted = bands.arraySort(ndvi.multiply(-1));

// Get the highest 20% NDVI observations per pixel.
var numImages = sorted.arrayLength(imageAxis).multiply(0.2).int();
var highestNdvi = sorted.arraySlice(imageAxis, 0, numImages);

// Get the mean of the highest 20% NDVI observations by reducing
// along the image axis.
var mean = highestNdvi.arrayReduce({
  reducer: ee.Reducer.mean(),
  axes: [imageAxis]
});

// Turn the reduced array image into a multi-band image for display.
var meanImage = mean.arrayProject([bandAxis]).arrayFlatten([bandNames]);
Map.centerObject(roi, 12);
Map.addLayer(meanImage, {bands: ['SR_B6', 'SR_B5', 'SR_B4'], min: 0, max: 0.4});

Configuración de Python

Consulta la página Entorno de Python para obtener información sobre la API de Python y el uso de geemap para el desarrollo interactivo.

import ee
import geemap.core as geemap

Colab (Python)

# Define a function that scales and masks Landsat 8 surface reflectance images
# and adds an NDVI band.
def prep_sr_l8(image):
  # Develop masks for unwanted pixels (fill, cloud, cloud shadow).
  qa_mask = image.select('QA_PIXEL').bitwiseAnd(int('11111', 2)).eq(0)
  saturation_mask = image.select('QA_RADSAT').eq(0)

  # Apply the scaling factors to the appropriate bands.
  optical_bands = image.select('SR_B.').multiply(0.0000275).add(-0.2)
  thermal_bands = image.select('ST_B.*').multiply(0.00341802).add(149.0)

  # Calculate NDVI.
  ndvi = optical_bands.normalizedDifference(['SR_B5', 'SR_B4']).rename('NDVI')

  # Replace the original bands with the scaled ones and apply the masks.
  return (
      image.addBands(optical_bands, None, True)
      .addBands(thermal_bands, None, True)
      .addBands(ndvi)
      .updateMask(qa_mask)
      .updateMask(saturation_mask)
  )


# Define an arbitrary region of interest as a point.
roi = ee.Geometry.Point(-122.26032, 37.87187)

# Load a Landsat 8 surface reflectance collection.
collection = (
    ee.ImageCollection('LANDSAT/LC08/C02/T1_L2')
    # Filter to get only imagery at a point of interest.
    .filterBounds(roi)
    # Filter to get only six months of data.
    .filterDate('2021-01-01', '2021-07-01')
    # Prepare images by mapping the prep_sr_l8 function over the collection.
    .map(prep_sr_l8)
    # Select the bands of interest to avoid taking up unneeded memory.
    .select('SR_B.|NDVI')
)

# Convert the collection to an array.
array = collection.toArray()

# Label of the axes.
image_axis = 0
band_axis = 1

# Get the NDVI slice and the bands of interest.
band_names = collection.first().bandNames()
bands = array.arraySlice(band_axis, 0, band_names.length())
ndvi = array.arraySlice(band_axis, -1)

# Sort by descending NDVI.
sorted = bands.arraySort(ndvi.multiply(-1))

# Get the highest 20% NDVI observations per pixel.
num_images = sorted.arrayLength(image_axis).multiply(0.2).int()
highest_ndvi = sorted.arraySlice(image_axis, 0, num_images)

# Get the mean of the highest 20% NDVI observations by reducing
# along the image axis.
mean = highest_ndvi.arrayReduce(reducer=ee.Reducer.mean(), axes=[image_axis])

# Turn the reduced array image into a multi-band image for display.
mean_image = mean.arrayProject([band_axis]).arrayFlatten([band_names])
m = geemap.Map()
m.center_object(roi, 12)
m.add_layer(
    mean_image, {'bands': ['SR_B6', 'SR_B5', 'SR_B4'], 'min': 0, 'max': 0.4}
)
m

Al igual que en el ejemplo de modelado lineal, separa las bandas de interés del índice de ordenamiento (NDVI) con arraySlice() a lo largo del eje de la banda. Luego, ordena las bandas de interés por índice de ordenamiento con arraySort(). Después de que los píxeles se hayan ordenado por NDVI descendente, usa arraySlice() junto con imageAxis para obtener el 20% de los píxeles de NDVI más altos. Por último, aplica arrayReduce() a lo largo de imageAxis con un reductor medio para obtener el promedio de los píxeles de NDVI más altos. El paso final convierte la imagen del array en una imagen multibanda para su visualización.