数组排序对于获取自定义质量拼接图非常有用，这涉及根据其他波段中的值缩减部分图像波段。以下示例按 NDVI 排序，然后获取集合中 NDVI 值最高的部分观测值的平均值：
// 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});
import ee import geemap.core as geemap
# 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
与线性建模示例中一样，沿波段轴使用 arraySlice() 将感兴趣的波段与排序索引 (NDVI) 分隔开来。然后，使用 arraySort() 按排序索引对感兴趣的频段进行排序。将像素按 NDVI 值从高到低排序后，沿 imageAxis 使用 arraySlice() 即可获取 NDVI 值最高的 20% 像素。最后，使用平均值 reducer 沿 imageAxis 应用 arrayReduce()，以获取最高 NDVI 像素的平均值。最后一步是将阵列图像转换回多波段图像以进行显示。