1

I want get 4 objects from histogram of a dataset in 136 line of code. getting one object is ok but for getting 4 objects (NDVI_before_histogram', 'NBR_before_histogram', 'gNDVI_before_histogram','MIRBI_before_histogram) it seem need a function. How can I do this?

the code:

// Select rangeland class from ESRI 2021
var remapped = image.select(['remapped']);
var rangeland = remapped.addBands(ee.Image(1).updateMask(remapped.gt(1)).rename('range'));
print('Rangeland:', rangeland);

//clip rangeland with boundry
var helen= rangeland.select(1).clip(table).reduceToVectors({
geometry: table,
scale: 10,
 })

var before = ee.Image('COPERNICUS/S2_SR_HARMONIZED/20210611T071621_20210611T072134_T39SVR')
  .clip(helen);

var after = ee.Image('COPERNICUS/S2_SR_HARMONIZED/20210626T071619_20210626T072134_T39SVR')
  .clip(helen);


// Optional viewing of cropped images
Map.centerObject(helen, 12);
Map.addLayer(before, {bands: ['B4', 'B3', 'B2'], min: 0, max: 3000}, 'Before fire');
Map.addLayer(after, {bands: ['B4', 'B3', 'B2'], min: 0, max: 3000}, 'After Fire');

// Calculate indices for the before image
var gNDVIbefore= before.normalizedDifference(['B8', 'B3']).rename('gNDVI_before');
print(gNDVIbefore, 'gNDVIbefore')
var NDVIbefore = before.normalizedDifference(['B8', 'B4']).rename('NDVI_before');
print(NDVIbefore, 'NDVIbefore')
var NBRbefore = before.normalizedDifference(['B8', 'B12']).rename('NBR_before');
print(NBRbefore, 'NBRbefore')
var MIRBIbefore = before.expression(
  '10 * B12 - 9.8 * B11 + 2',
  {
    'B12': before.select('B12'),
    'B11': before.select('B11')
  }
).rename('MIRBI_before');
print(MIRBIbefore, 'MIRBIbefore')

// Calculate indices for the after image
var gNDVIafter= after.normalizedDifference(['B8', 'B3']).rename('gNDVI_after');
print(gNDVIafter,'gNDVIafter')
var NDVIafter = after.normalizedDifference(['B8', 'B4']).rename('NDVI_after');
var NBRafter = after.normalizedDifference(['B8', 'B12']).rename('NBR_after');
var MIRBIafter = after.expression(
  '10 * B12 - 9.8 * B11 + 2',
  {
    'B12': after.select('B12'),
    'B11': after.select('B11')
  }
).rename('MIRBI_after');

// Calculate differences in indices
var dgNDVI = gNDVIbefore.subtract(gNDVIafter);
var dNDVI = NDVIbefore.subtract(NDVIafter);
var dNBR = NBRbefore.subtract(NBRafter);
var dMIRBI = MIRBIbefore.subtract(MIRBIafter);

// Calculate dataset before the event
var multiindex_before = NDVIbefore.addBands(NBRbefore).addBands(gNDVIbefore).addBands(MIRBIbefore)
print(multiindex_before, 'multiindex_before')
Map.addLayer(multiindex_before, {}, 'multiindex_before')

// Calculate multi-index after the event
var multiindex_after = NDVIafter.addBands(NBRafter).addBands(gNDVIafter).addBands(MIRBIafter)
print(multiindex_after, 'multiindex_after')
Map.addLayer(multiindex_after, {}, 'multiindex_after')

// Calculate the difference in multi-index
var dmultiindex = multiindex_before.subtract(multiindex_after)
//.rename(['NDVI_before', 'NBR_before', 'gNDVI_before', 'MIRBI_before'])
print(dmultiindex.bandNames(), '.bandNames')

// Add the difference in multi-index to the map
Map.addLayer(dmultiindex, {}, 'Difference in Multi-Index');



////////////////////////////////////////////////////////
// Calculate the histogram of the dmultiindex band, including mean and variance for additional information.
////////////////////////////////////////////////////////
var histogram = dmultiindex.select(['NDVI_before', 'NBR_before', 'gNDVI_before', 'MIRBI_before']).reduceRegion({
  reducer: ee.Reducer.histogram({maxBuckets: 255})
    .combine(ee.Reducer.mean(), '', true)
    .combine(ee.Reducer.variance(), '', true), 
  geometry: helen, 
  scale: 10,
  bestEffort: true
});

// Print histogram, mean, and variance
print('Histogram, Mean, and Variance:', histogram);

// Chart the histogram
print(ui.Chart.image.histogram(dmultiindex.select(['NDVI_before', 'NBR_before', 'gNDVI_before', 'MIRBI_before']), helen, 30)
  .setOptions({
    title: 'Histogram of dmultiindex Band',
    hAxis: {title: 'dmultiindex Value'},
    vAxis: {title: 'Frequency'},
    legend: 'none',
    lineWidth: 1,
    colors: ['orange']
  }));
  
// Return the DN that maximizes interclass variance in dmultiindex (in the region).
var otsu = function(histogram) {
  var counts = ee.Array(ee.Dictionary(histogram).get('histogram'));
  var means = ee.Array(ee.Dictionary(histogram).get('bucketMeans'));
  var size = means.length().get([0]);
  var total = counts.reduce(ee.Reducer.sum(), [0]).get([0]);
  var sum = means.multiply(counts).reduce(ee.Reducer.sum(), [0]).get([0]);
  var mean = sum.divide(total);
  
  var indices = ee.List.sequence(1, size);
  
  // Compute between sum of squares, where each mean partitions the data.
  var bss = indices.map(function(i) {
    var aCounts = counts.slice(0, 0, i);
    var aCount = aCounts.reduce(ee.Reducer.sum(), [0]).get([0]);
    var aMeans = means.slice(0, 0, i);
    var aMean = aMeans.multiply(aCounts)
        .reduce(ee.Reducer.sum(), [0]).get([0])
        .divide(aCount);
    var bCount = total.subtract(aCount);
    var bMean = sum.subtract(aCount.multiply(aMean)).divide(bCount);
    return aCount.multiply(aMean.subtract(mean).pow(2)).add(
           bCount.multiply(bMean.subtract(mean).pow(2)));
  });
  
  print(ui.Chart.array.values(ee.Array(bss), 0, means));
  
  // Return the mean value corresponding to the maximum BSS.
  return means.sort(bss).get([-1]);
};

var threshold = otsu(histogram.get('NDVI_before_histogram', 'NBR_before_histogram', 'gNDVI_before_histogram','MIRBI_before_histogram'));
print('threshold', threshold);

var thrdmultiindex = dmultiindex.select(['NDVI_before', 'NBR_before', 'gNDVI_before', 'MIRBI_before']).gt(threshold);
Map.addLayer(thrdmultiindex)

Map.addLayer(thrdmultiindex.selfMask(), {palette: 'orange'}, 'thrdmultiindex');

//////////////////////////////////// otsu for after dataset///////////////////////////////////////

var histogram = multiindex_after.select('multiindex_after').reduceRegion({
  reducer: ee.Reducer.histogram({maxBuckets: 255})
    .combine(ee.Reducer.mean(), '', true)
    .combine(ee.Reducer.variance(), '', true), 
  geometry: helen, 
  scale: 10,
  bestEffort: true
});

// Print histogram, mean, and variance
print('Histogram, Mean, and Variance:', histogram);

// Chart the histogram
print(ui.Chart.image.histogram(multiindex_after.select('multiindex_after'), helen, 30)
  .setOptions({
    title: 'Histogram of multiindex_after Band',
    hAxis: {title: 'multiindex_afte Value'},
    vAxis: {title: 'Frequency'},
    legend: 'none',
    lineWidth: 1,
    colors: ['orange']
  }));
  
// Return the DN that maximizes interclass variance in dmultiindex_after (in the region).
var otsu = function(histogram) {
  var counts = ee.Array(ee.Dictionary(histogram).get('histogram'));
  var means = ee.Array(ee.Dictionary(histogram).get('bucketMeans'));
  var size = means.length().get([0]);
  var total = counts.reduce(ee.Reducer.sum(), [0]).get([0]);
  var sum = means.multiply(counts).reduce(ee.Reducer.sum(), [0]).get([0]);
  var mean = sum.divide(total);
  
  var indices = ee.List.sequence(1, size);
  
  // Compute between sum of squares, where each mean partitions the data.
  var bss = indices.map(function(i) {
    var aCounts = counts.slice(0, 0, i);
    var aCount = aCounts.reduce(ee.Reducer.sum(), [0]).get([0]);
    var aMeans = means.slice(0, 0, i);
    var aMean = aMeans.multiply(aCounts)
        .reduce(ee.Reducer.sum(), [0]).get([0])
        .divide(aCount);
    var bCount = total.subtract(aCount);
    var bMean = sum.subtract(aCount.multiply(aMean)).divide(bCount);
    return aCount.multiply(aMean.subtract(mean).pow(2)).add(
           bCount.multiply(bMean.subtract(mean).pow(2)));
  });
  
  print(ui.Chart.array.values(ee.Array(bss), 0, means));
  
  // Return the mean value corresponding to the maximum BSS.
  return means.sort(bss).get([-1]);
};

var threshold = otsu(histogram.get('multiindex_after_histogram'));
print('threshold of multiindex_after', threshold);

var thrmultiindex_after = multiindex_after.select('multiindex_after').gt(threshold);
Map.addLayer(thrmultiindex_after)

Map.addLayer(thrmultiindex_after.selfMask(), {palette: 'orange'}, 'multiindex_after');

https://code.earthengine.google.com/8dbb18ab5764f9a0f85e12c2449b3e30

1
  • Please Edit the Question to contain a minimal reproducible code sample, or this is likely to be closed.
    – Vince
    Commented Jul 9 at 18:43

2 Answers 2

0

The primary issue with your function call is that the get method only accepts a single argument, but you are attempting to pass four. Here’s the line causing the issue:

var threshold = otsu(histogram.get('NDVI_before_histogram', 'NBR_before_histogram', 'gNDVI_before_histogram','MIRBI_before_histogram'));

To fix this, you need to retrieve each histogram individually and then pass the necessary components to your otsu function:

var otsu = function(params) {
  var counts = params.counts
  var means = params.means
  ...
};

For example, if you're processing the NDVI_before_histogram:

var NDVI_before_histogram = ee.Dictionary(histogram.get('NDVI_before_histogram'));
var threshold = otsu({
  counts: ee.List(NDVI_before_histogram.get('histogram')),
  means: ee.List(NDVI_before_histogram.get('bucketMeans'))
});

Addressing the error in line 136 is a good start, but there seem to be numerous other issues throughout your code that suggest a fundamental misunderstanding of how certain GEE functions operate and the arguments they require. I recommend conducting a thorough review of your entire script to ensure a better understanding of these functions and their usage. It may be helpful to consult the GEE documentation or seek tutorials that clarify these concepts before revisiting the forum for further assistance.

1
  • thanks Nicolas. retrieving each histogram individually get 4 thresholds. I need a threshold for a dataset (dmultiindex) in line 70. is there any way to get just one threshold? Commented Jul 10 at 14:05
0

The issues are produced because is mandatory to have 4 thresholds and 4 dmultiindex variables (or 8 because there are before and after events). I fixed almost all the code but I left a small part for you to do on your own.

A sample of main modifications looks as follows:

/*
var threshold = otsu(histogram.get('NDVI_before_histogram', 
                                   'NBR_before_histogram', 
                                   'gNDVI_before_histogram',
                                   'MIRBI_before_histogram'));
*/

var threshold_NDVI = otsu(histogram.get('NDVI_before_histogram'));

print('threshold_NDVI', threshold_NDVI);

var threshold_NBR = otsu(histogram.get('NBR_before_histogram'));

print('threshold_NBR', threshold_NBR);

var threshold_gNDVI = otsu(histogram.get('gNDVI_before_histogram'));

print('threshold_gNDVI', threshold_gNDVI);

var threshold_MIRBI = otsu(histogram.get('MIRBI_before_histogram'));

print('threshold_MIRBI', threshold_MIRBI);



var thrdmultiindex_NDVI = dmultiindex.select(['NDVI_before', 
                                         'NBR_before', 
                                         'gNDVI_before', 
                                         'MIRBI_before']).gt(threshold_NDVI);

Map.addLayer(thrdmultiindex_NDVI);

var thrdmultiindex_NBR = dmultiindex.select(['NDVI_before', 
                                         'NBR_before', 
                                         'gNDVI_before', 
                                         'MIRBI_before']).gt(threshold_NBR);

Map.addLayer(thrdmultiindex_NBR);

var thrdmultiindex_MIRBI = dmultiindex.select(['NDVI_before', 
                                         'NBR_before', 
                                         'gNDVI_before', 
                                         'MIRBI_before']).gt(threshold_MIRBI);

Map.addLayer(thrdmultiindex_MIRBI);



Map.addLayer(thrdmultiindex_MIRBI.select('MIRBI_before')
                                 .selfMask(), {palette: 'orange'}, 'thrdmultiindex_MIRBI');

Complete code can be accessed here.

After running it in the GEE code editor, results can be observed in following picture. Orange layer (the last one) is masked for the very negative threshold MIRBI value.

enter image description here

1
  • thank you very much. in this code you edited, there are 4 threshold that made for each bands. But I want just one threshold that should be for my dataset. it means a threshold for all bands together not separate band. Commented Jul 10 at 20:05

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.