Using NDVI thresholds as a condition for successive elaborations in Google Earth Engine

Question

From Sentinel-1 and Sentinel-2 it was created an image collection of 255 images, each with 5 bands:

ndvi,
VV_corr (VV corrected for the incidence angle),
VV_diff (VV_corr difference between two successive acquisitions),
VV_diff_abs (absolute value of VV_corr difference between two successive acquisitions),
meanNdvi (mean ndvi between two successive acquisitions)

Now I want to subtract at a pixel level each VV_corr with the minimum VV_corr in the image collection based on the ndvi value. I should calculate the minimum VV_corr for each pixel with ndvi comprised between 0-0.1, 0.1-0.2, 0.2-0.3, 0.3-0.4, 0.4-0.5 etc. Than evaluate each pixel: if its ndvi is comprises between 0 and 0.1 it should be subtracted with the corrisponding minimum VV_corr, else if it is comprises between 0.1 and 0.2 it should be subtracted with the corrisponding minimum VV_corr etc.

I don't how to approach the problem without using if/else statements. At the moment I kind of did it but with a lot of copy and paste.

var geometry = 
    /* color: #d63000 */
    /* shown: false */
    ee.Geometry.Point([12.198218946865174, 44.47799625529383]),
    area = 
    /* color: #98ff00 */
    /* shown: false */
    ee.Geometry.Polygon(
        [[[10.423257686605991, 43.3608941765534],
          [13.213785030355991, 43.735146720411635],
          [12.444742061605991, 45.21664268167095],
          [10.064828731527866, 44.90526427338962]]]);
          
          
//Corine land cover
var corine = ee.Image('COPERNICUS/CORINE/V20/100m/2012')
//aree agricole = 1; altre coperture di suolo = 0
var areeAgricole = corine.eq(211)     //seminativi non irrigui
                  .or(corine.eq(212)) //seminativi irrigui
                  .or(corine.eq(213)) //risaie
                  .or(corine.eq(231)) //prati stabili
                  .or(corine.eq(241)) //colture temporanee
                  .or(corine.eq(242)) //sistemi colturali complessi
                  .or(corine.eq(243)) //aree agrarie e spazi naturali
                  .or(corine.eq(321)) //aree a pascolo e praterie
                  .or(corine.eq(333)) //aree con vegetazione rada
                  .or(corine.eq(322));//brughiere

//function to correct theta_i
var toGammaVV = function (image) {
      return image.addBands(image.select('VV').subtract(image.select('angle')
        .multiply(Math.PI/180.0).cos().log10().multiply(10.0)).rename('VV_corr'))
}


//collection sentinel-1               
var s1 = ee.ImageCollection("COPERNICUS/S1_GRD") 
              .filterBounds(geometry) 
              .filter(ee.Filter.eq('orbitProperties_pass','ASCENDING')) 
              .filter(ee.Filter.eq('instrumentMode','IW')) 
              .filter(ee.Filter.eq('relativeOrbitNumber_start',117)) 
              .filterDate('2015-08-01','2020-11-01') 
              .map(toGammaVV)
              .select('VV_corr')
              .map(function(image){return image.updateMask(areeAgricole)})
              
              
//sentinel-2
var s2 = ee.ImageCollection("COPERNICUS/S2") 
  .filterBounds(geometry) 
  .filterDate("2015-08-01","2020-11-01") 
  .filter(ee.Filter.lt("CLOUDY_PIXEL_PERCENTAGE", 30)) 
  .map(function (img) { 
    var ndvi = img.normalizedDifference(["B8","B4"])
    .rename("ndvi");
    return img.addBands(ndvi)})
  .map(function(image){return image.updateMask(areeAgricole)})


//combine s1 and s2     
var maxDiffFilter = ee.Filter.maxDifference({ 
  difference: 10 * 24 * 60 * 60 * 1000,  //differenza massima 10giorni
  leftField: "system:time_start", //operatore 1 
  rightField: "system:time_start" //operatore 2
});

var saveBestJoin = ee.Join.saveBest({ 
  matchKey: "bestImage",
  measureKey: "timeDiff"
});

var join = saveBestJoin.apply(s1, s2, maxDiffFilter); 

var plot = join.map(function(img) {
  var cat = ee.Image.cat([img, img.get("bestImage")]);
  return ee.Image(cat.select(['VV_corr','ndvi']))
});


//make sure the images are ordered
var list = plot.sort('system:time_start').toList(plot.size());

//couple the images 2 by 2
var pair = list.slice(0,-1).zip(list.slice(1));

//get mean NDVI and VV diff
var diff = ee.ImageCollection(pair.map(function(f) { 
  var diff =  ee.Image(ee.List(f).get(1)).select('VV_corr').subtract(ee.Image(ee.List(f).get(0)).select('VV_corr'))
  var ndvi =  (ee.Image(ee.List(f).get(1)).select('ndvi').add(ee.Image(ee.List(f).get(0)).select('ndvi'))).divide(2)
  var diff_abs = diff.abs()
  return ndvi.addBands([
    ee.Image(ee.List(f).get(1)).select('VV_corr'),
    ee.Image(ee.List(f).get(1)).select('ndvi'),
    diff,
    diff_abs]).rename(['meanNdvi','VV_corr','ndvi','VV_diff','VV_diff_abs'])
}));


//check 
// Map.addLayer(ee.Image(diff.get(3)), {},'diff (i,i+1)', false)
// Map.addLayer(ee.Image(list.get(3)),{},'i', false)
// Map.addLayer(ee.Image(list.get(4)),{},'i+1', false)



//sample over the images in the collection
var sampColl = function (img) {
  return ee.Image(img).sampleRegions({
    collection: area, 
    scale: 13000, //min scale
    tileScale: 16,
    })
}
var samp = diff.map(sampColl).flatten();


// Generate chart from sample
var chart = ui.Chart.feature.byFeature(samp,'meanNdvi','VV_diff')
                .setChartType('ScatterChart')
                .setOptions({
                  title: 'Radar backscatter difference as a function of NDVI',
                  hAxis: {
                    title: 'NDVI(mean)',
                    titleTextStyle: {italic: false, bold: true},
                    viewWindow: {min: 0, max: 0.8}
                  },
                  vAxis: {
                    title: 'VV difference(i,i+1)',
                    titleTextStyle: {italic: false, bold: true}
                  },
                  pointSize: 4,
                  dataOpacity: 0.6,
                  colors: ['1d6b99', 'cf513e'],
                });


// Generate chart abs from sample
var chart_abs = ui.Chart.feature.byFeature(samp,'meanNdvi','VV_diff_abs')
                .setChartType('ScatterChart')
                .setOptions({
                  title: 'Radar backscatter difference (absolute value) as a function of NDVI',
                  hAxis: {
                    title: 'NDVI(mean)',
                    titleTextStyle: {bold: true},
                    viewWindow: {min: 0, max: 0.8}
                  },
                  vAxis: {
                    title: 'VV difference(i,i+1)',
                    titleTextStyle: {bold: true}
                  },
                  pointSize: 4,
                  dataOpacity: 0.4,
                  colors: ['1d6b99', 'cf513e'],
                  trendlines: {
                    0: {  
                      type: 'linear',  
                      color: 'red',
                      lineWidth: 2,
                    },
                  },
                });
print('SOIL MOISTURE CHANGE DETECTION WITH NDVI (Gao et al., 2017)')
print("If NDVI increases, the radar backscattering difference \n between two acquisitions decreases.\n Also the sensitivity to soil moisture variations decreases")
print(chart);
print(chart_abs);


//how to loop the following???

var ndvi01 = diff.map(function(img){return img.updateMask(img.select('ndvi').lte(0.1).and(img.select('ndvi').gt(0.0)))}).set('ndvi',0.1);
var ndvi02 = diff.map(function(img){return img.updateMask(img.select('ndvi').lte(0.2).and(img.select('ndvi').gt(0.1)))}).set('ndvi',0.2);
var ndvi03 = diff.map(function(img){return img.updateMask(img.select('ndvi').lte(0.3).and(img.select('ndvi').gt(0.2)))}).set('ndvi',0.3);
var ndvi04 = diff.map(function(img){return img.updateMask(img.select('ndvi').lte(0.4).and(img.select('ndvi').gt(0.3)))}).set('ndvi',0.4);
var ndvi05 = diff.map(function(img){return img.updateMask(img.select('ndvi').lte(0.5).and(img.select('ndvi').gt(0.4)))}).set('ndvi',0.5);
var ndvi06 = diff.map(function(img){return img.updateMask(img.select('ndvi').lte(0.6).and(img.select('ndvi').gt(0.5)))}).set('ndvi',0.6);
var ndvi07 = diff.map(function(img){return img.updateMask(img.select('ndvi').lte(0.7).and(img.select('ndvi').gt(0.6)))}).set('ndvi',0.7);
var ndvi08 = diff.map(function(img){return img.updateMask(img.select('ndvi').lte(0.8).and(img.select('ndvi').gt(0.7)))}).set('ndvi',0.8);

// print(ndvi01)

var maxDiff01 = ndvi01.select('VV_diff_abs').max().rename('VV_diff_max').addBands(ndvi01.select('VV_corr').min().rename('VV_min')).copyProperties(ndvi01)
var maxDiff02 = ndvi02.select('VV_diff_abs').max().rename('VV_diff_max').addBands(ndvi02.select('VV_corr').min().rename('VV_min')).copyProperties(ndvi02)
var maxDiff03 = ndvi03.select('VV_diff_abs').max().rename('VV_diff_max').addBands(ndvi03.select('VV_corr').min().rename('VV_min')).copyProperties(ndvi03)
var maxDiff04 = ndvi04.select('VV_diff_abs').max().rename('VV_diff_max').addBands(ndvi04.select('VV_corr').min().rename('VV_min')).copyProperties(ndvi04)
var maxDiff05 = ndvi05.select('VV_diff_abs').max().rename('VV_diff_max').addBands(ndvi05.select('VV_corr').min().rename('VV_min')).copyProperties(ndvi05)
var maxDiff06 = ndvi06.select('VV_diff_abs').max().rename('VV_diff_max').addBands(ndvi06.select('VV_corr').min().rename('VV_min')).copyProperties(ndvi06)
var maxDiff07 = ndvi07.select('VV_diff_abs').max().rename('VV_diff_max').addBands(ndvi07.select('VV_corr').min().rename('VV_min')).copyProperties(ndvi07)
var maxDiff08 = ndvi08.select('VV_diff_abs').max().rename('VV_diff_max').addBands(ndvi08.select('VV_corr').min().rename('VV_min')).copyProperties(ndvi08)



var coll = ee.ImageCollection([
  maxDiff01, 
  maxDiff02, 
  maxDiff03, 
  maxDiff04, 
  maxDiff05,
  maxDiff06, 
  maxDiff07,
  maxDiff08])

var chart_maxdiff = ui.Chart.image.series({
                  imageCollection: coll.select('VV_diff_max'), 
                  region: area, 
                  reducer: ee.Reducer.median(), 
                  scale:500, 
                  xProperty: 'ndvi'})
                  .setOptions({
                  title: 'Maximum radar backscattering difference (mean value over the image) for differents NDVI classes',
                  hAxis: {
                    title: 'NDVI',
                    titleTextStyle: {bold: true},
                    viewWindow: {min: 0.1, max: 0.8}
                  },
                  vAxis: {
                    title: 'Max VV difference',
                    titleTextStyle: {bold: true},
                    
                  },
                  pointSize: 2,
                  colors: ['red'],
                  })
print('g(NDVI) is the empirical function which describes \n the maximum change in radar signal for each NDVI value')                  
print(chart_maxdiff)           

print('Surface soil moisture at time t2 is given by:')                 
print('SSM(t2) = SSM(t1)+H(diff-sigma(t1,t2))')
print('where:')
print('H(diff-sigma(t1,t2))=(0.15*sigma0(ndvi)-min sigma0(ndvi)) / max diff sigma0(ndvi)')
       


var numFunction = function (img) {
  var ndvi = img.select('ndvi');
  var num01 = ndvi.gte(0.0).and(ndvi.lt(0.1)).select('VV_corr').subtract(ee.Image(maxDiff01).select('VV_min'))
  var num02 = ndvi.gte(0.1).and(ndvi.lt(0.2)).select('VV_corr').subtract(ee.Image(maxDiff02).select('VV_min'))
  var num03 = ndvi.gte(0,2).and(ndvi.lt(0.3)).select('VV_corr').subtract(ee.Image(maxDiff03).select('VV_min'))
  var num04 = ndvi.gte(0.3).and(ndvi.lt(0.4)).select('VV_corr').subtract(ee.Image(maxDiff04).select('VV_min'))
  var num05 = ndvi.gte(0.4).and(ndvi.lt(0.5)).select('VV_corr').subtract(ee.Image(maxDiff05).select('VV_min'))
  var num06 = ndvi.gte(0.5).and(ndvi.lt(0.6)).select('VV_corr').subtract(ee.Image(maxDiff06).select('VV_min'))
  var num07 = ndvi.gte(0.6).and(ndvi.lt(0.7)).select('VV_corr').subtract(ee.Image(maxDiff07).select('VV_min'))
  var num08 = ndvi.gte(0.7).and(ndvi.lt(0.8)).select('VV_corr').subtract(ee.Image(maxDiff08).select('VV_min'))
  
}

//   return num01.rename('VV-VVmin')
// }

// var numColl = diff.map(numFunction)
// print(numColl.first())
// Map.addLayer(ee.Image(diff.first()),{bands:'VV_corr', min:-20, max:0}, 'VV')
// Map.addLayer(ee.Image(maxDiff01),{bands:'VV_min', min:-20, max:0}, 'VV min')
// Map.addLayer(ee.Image(numColl.first()),{bands:'VV_corr', min:-20, max:0}, 'VV - VV min')

This is the link to the code: https://code.earthengine.google.com/56a339b63146288a2ad53edb8403710e

Daniel Wiell · Accepted Answer · 2022-12-20 12:06:47Z

0

It was pretty clean to start with. You can just map over a list of NDVI thresholds.

var step = 0.1
var coll = ee.ImageCollection(
  ee.List.sequence(0, 0.7, step)
    .map(function (min) {
      min = ee.Number(min)
      var max = min.add(step)
      var masked = diff.map(function(img) {
        return img
          .updateMask(
            img.select('ndvi').gt(min).and(img.select('ndvi').lte(max))
          )
      })
      var vvDiffMax = masked.select('VV_diff_abs').max().rename('VV_diff_max')
      var vvMin = masked.select('VV_corr').min().rename('VV_min')
      return vvDiffMax
        .addBands(vvMin)
        .set('ndvi', max)
    })
)

https://code.earthengine.google.com/6e6a07797faa1a5313b897f59416d1d5

answered Dec 20, 2022 at 12:06

Daniel Wiell

14.5k2 gold badges12 silver badges30 bronze badges

Hi Daniel thank you very much! Than, to subtract VV_corr and VV_min based on NDVI, shoud I use a similar approach mapping over a list of ndvi values?
– benedetta
Commented Dec 20, 2022 at 13:19
code.earthengine.google.com/0834e2c2a891f24307779e4618aaba17
– benedetta
Commented Dec 20, 2022 at 13:24

Add a comment |

Stack Exchange Network

Using NDVI thresholds as a condition for successive elaborations in Google Earth Engine

1 Answer 1

Your Answer

Not the answer you're looking for? Browse other questions tagged
google-earth-engine
sentinel-2
sentinel-1
if-else
or ask your own question.

Hot Network Questions

Using NDVI thresholds as a condition for successive elaborations in Google Earth Engine

1 Answer 1

Your Answer

Sign up or log in

Post as a guest

Not the answer you're looking for? Browse other questions tagged google-earth-enginesentinel-2sentinel-1if-else or ask your own question.

Related

Hot Network Questions

Not the answer you're looking for? Browse other questions tagged
google-earth-engine
sentinel-2
sentinel-1
if-else
or ask your own question.