0

I decided to find out how the FWI index works for operational monitoring of the fire situation in the region. This is my first project on Google Engine and my programming knowledge leaves much to be desired, so please do not criticize. The index is calculated based on a set of indicators, one of which is FFMC (Fine Fuel Moisture Code). The value for calculating the indicator is 80 (start FFMC) on the very first day of monitoring. Further calculation is based on data for the previous day. To obtain information for a specific day, it is necessary to calculate all previous days, since the indicator for a certain date depends on the indicator of the previous one.

How can one iterate through these equations, saving the final indicator value for the previous date and using it for calculations at the current moment?

Roughly speaking, to calculate the F2 indicator for 2105-04-04, you need to know F1, F3 - F2, and so on.

var aoi = ee.Geometry.Polygon(
        [[[29.269278980991366, 51.95773487804608],
          [29.269278980991366, 51.22745421938314],
          [30.670035816928866, 51.22745421938314],
          [30.670035816928866, 51.95773487804608]]], null, false);

var RGBPalette = ['#8B0000', '#FFFF00', '#FFFACD', '#F0E68C', '#ADFF2F', '#7FFF00', '#7CFC00'];
Map.setCenter (30.0149, 51.6353, 10);
var startDate = ee.Date('2015-04-03');
var prevDate = (startDate.advance(-1, 'day'));
print (startDate);
print (prevDate);

var t = ee.ImageCollection("ECMWF/ERA5_LAND/DAILY_AGGR").filterDate(startDate).filterBounds(aoi).median();

var r = ee.ImageCollection("ECMWF/ERA5_LAND/DAILY_AGGR").filterDate(prevDate).filterBounds(aoi).median();

var w = ee.ImageCollection('NASA/GLDAS/V021/NOAH/G025/T3H').filterDate(startDate).filterBounds(aoi).median();

var temper = t.select('temperature_2m');
var precip_today = t.select('total_precipitation_sum');
var precip_prev = r.select('total_precipitation_sum');
var windsp = w.select('Wind_f_inst');
var dewpo = t.select ('dewpoint_temperature_2m');

var T = temper.subtract(ee.Number(273.15));   
var R1 = precip_today.multiply(ee.Number(1000));
var R0 = precip_prev.multiply(ee.Number(1000));

var W = windsp.multiply(ee.Number(1));
var dewpoi = dewpo.subtract(ee.Number(273.15));   
var H =  temper.expression('100 - 5 * (T - Td)',{'T': T.select('temperature_2m'),'Td': dewpoi.select('dewpoint_temperature_2m')});

Map.addLayer(T.clip(aoi), tempParams, 'Температура', false);
Map.addLayer(R1.clip(aoi), precParams, 'Осадки сегодня', false);
Map.addLayer(R0.clip(aoi), precParams, 'Осадки вчера', false);
Map.addLayer(W.clip(aoi), windParams, 'Скорость ветра', false);
Map.addLayer(H.clip(aoi), humidParams, 'Относительная влажность', false);

//FFMC
var startFFMC = ee.Number(80)
var Mo = R0.expression('147.2 * (101 - startFFMC)/(59.5 + startFFMC)',{'startFFMC': startFFMC});
Map.addLayer(Mo.clip(aoi), RGBPalette, 'Mo', false);

var Rf = R0.expression('Ro - 0.5',{'Ro': R0});
Map.addLayer(Rf.clip(aoi), RGBPalette, 'Rf', false);

var Mr = Mo<=ee.Number(150.0) ? T.expression ('Mo + 42.5 * Rf * (2.718**(-100/(251-Mo)))*(1-(2.718**(-6.93/Rf)))', {'Mo': Mo, 'Rf': Rf}) : 
T.expression ('Mo + 42.5 * Rf * (2.718**(-100/(251-Mo)))*(1-(2.718**(-6.93/Rf)))+0.0015*((Mo-150)**2)*(Rf**(0.5))', {'Mo': Mo, 'Rf': Rf});
Map.addLayer(Mr.clip(aoi), RGBPalette, 'Mr', false);

var Mrn = Mo.expression ('Rf > 0.5 ? Mr : Mo', {'Rf': Rf, 'Mo': Mo, 'Mr': Mr});
Map.addLayer(Mrn.clip(aoi), RGBPalette, 'Mrn', false);

var Ed = Rf.expression('(0.942 * H**0.679) + (11*(2.718**(H-100)/10))+0.18*(21.1 - T)*(1-(2.718**-(0.115*H)))',{'H': H,'T': T}).rename('EMC for drying');
Map.addLayer(Ed.clip(aoi), RGBPalette, 'Ed', false);

var Ew = Rf.expression('(0.618 * H**0.753) + (10*(2.718**(H-100)/10))+0.18*(21.1 - T)*(1-(2.718**(-0.115*H)))',{'H': H,'T': T}).rename('EMC for wetting');
Map.addLayer(Ew.clip(aoi), RGBPalette, 'Ew', false);

var Ko = Rf.expression('0.424 * (1-((H/100)**1.7)) + (0.0694 * (W**0.5)) * (1-((H/100)**8))',{'W': W,'H': H}).rename('Ko');
Map.addLayer(Ko.clip(aoi), RGBPalette, 'Ko', false);

var Kd = Ko.expression('1*Ko * 0.581 * (2.718**0.0365*T)',{'Ko': Ko,'T': T}).rename('Kd');
Map.addLayer(Kd.clip(aoi), RGBPalette, 'Kd', false);

var Kl = H.expression('(0.424 * (1-((100 - H)/100)**1.7) + (0.0694 * (W**0.5)) * (1-(100-H)/100)**8)',{'W': W,'H': H}).rename('Kl');
Map.addLayer(Kl.clip(aoi), RGBPalette, 'Kl', false);

var Kw = Kl.expression('1*Kl * 0.581 * (2.718**0.0365*T)',{'Kl': Kl,'T': T}).rename('Kw');
Map.addLayer(Kw.clip(aoi), RGBPalette, 'Kw', false);

var M1 = Ed.expression('Ed + (Mr - Ed) * 10**-Kd',{'Ed': Ed,'Mr': Mr,'Kd' : Kd}).rename('M1');
Map.addLayer(M1.clip(aoi), RGBPalette, 'M1', false);

var M2 = Ed.expression('Ew - (Ew - Mr) * 10**-Kd', {'Ew': Ew, 'Mr': Mr, 'Kd' : Kd}).rename('M2');
Map.addLayer(M2.clip(aoi), RGBPalette, 'M2', false);

var Mn = Mrn.expression ('Mrn>Ed ? M1 : Mrn<Ew && Mrn < Ed ? M2 : Ew<=Mrn<=Ed ? Mrn: 0', {'Mrn':Mrn,'Ed':Ed,'M1':M1,'Ew':Ew,'M2':M2});
Map.addLayer(Mn.clip(aoi), RGBPalette, 'Mn', false);

var F = Mn.expression ('59.5 * (250 - Mn)/(147.2 + Mn)', {'Mn':Mn});
Map.addLayer(F.clip(aoi), RGBPalette, 'F', false);

var day1 = ee.Date('2015-04-04');
var t1 = ee.ImageCollection("ECMWF/ERA5_LAND/DAILY_AGGR").filterDate(day1).filterBounds(aoi).median();
var r1 = ee.ImageCollection("ECMWF/ERA5_LAND/DAILY_AGGR").filterDate(day1).filterBounds(aoi).median();
var w1 = ee.ImageCollection('NASA/GLDAS/V021/NOAH/G025/T3H').filterDate(day1).filterBounds(aoi).median();
var temper1 = t1.select('temperature_2m');
var precip_today1 = t1.select('total_precipitation_sum');
var precip_prev1 = r1.select('total_precipitation_sum');
var windsp1 = w1.select('Wind_f_inst');
var dewpo1 = t1.select ('dewpoint_temperature_2m');
var T1 = temper1.subtract(ee.Number(273.15));   
var R11 = precip_today1.multiply(ee.Number(1000));
var R01 = precip_prev1.multiply(ee.Number(1000));
var W1 = windsp1.multiply(ee.Number(1));
var dewpoi1 = dewpo1.subtract(ee.Number(273.15));   
var H1 =  temper1.expression('100 - 5 * (T1 - Td1)',{'T1': T1.select('temperature_2m'),'Td1': dewpoi1.select('dewpoint_temperature_2m')});
var Mo1 = R0.expression('147.2 * (101 - F)/(59.5 + F)',{'F': F});
var Rf1 = R01.expression('Ro1 - 0.5',{'Ro1': R01});
var Mr1 = Mo1<=ee.Number(150.0) ? T.expression ('Mo1 + 42.5 * Rf1 * (2.718**(-100/(251-Mo1)))*(1-(2.718**(-6.93/Rf1)))', {'Mo1': Mo1, 'Rf1': Rf1}) : 
T1.expression ('Mo1 + 42.5 * Rf1 * (2.718**(-100/(251-Mo1)))*(1-(2.718**(-6.93/Rf1)))+0.0015*((Mo1-150)**2)*(Rf1**(0.5))', {'Mo1': Mo1, 'Rf1': Rf1});
var Mrn1 = Mo1.expression ('Rf1 > 0.5 ? Mr1 : Mo1', {'Rf1': Rf1, 'Mo1': Mo1, 'Mr1': Mr1});
var Ed1 = Rf1.expression('(0.942 * H1**0.679) + (11*(2.718**(H1-100)/10))+0.18*(21.1 - T1)*(1-(2.718**-(0.115*H1)))',{'H1': H1,'T1': T1}).rename('EMC for drying');
var Ew1 = Rf1.expression('(0.618 * H1**0.753) + (10*(2.718**(H1-100)/10))+0.18*(21.1 - T1)*(1-(2.718**(-0.115*H1)))',{'H1': H1,'T1': T1}).rename('EMC for wetting');
var Ko1 = Rf1.expression('0.424 * (1-((H1/100)**1.7)) + (0.0694 * (W1**0.5)) * (1-((H1/100)**8))',{'W1': W1,'H1': H1}).rename('Ko1');
var Kd1 = Ko1.expression('1*Ko1 * 0.581 * (2.718**0.0365*T1)',{'Ko1': Ko1,'T1': T1}).rename('Kd1');
var Kl1 = H1.expression('(0.424 * (1-((100 - H1)/100)**1.7) + (0.0694 * (W1**0.5)) * (1-(100-H1)/100)**8)',{'W1': W1,'H1': H1}).rename('Kl1');
var Kw1 = Kl1.expression('1*Kl1 * 0.581 * (2.718**0.0365*T1)',{'Kl1': Kl1,'T1': T1}).rename('Kw1');
var M11 = Ed1.expression('Ed1 + (Mr1 - Ed1) * 10**-Kd1',{'Ed1': Ed1,'Mr1': Mr1,'Kd1' : Kd1}).rename('M11');
var M21 = Ed1.expression('Ew1 - (Ew1 - Mr1) * 10**-Kd1', {'Ew1': Ew1, 'Mr1': Mr1, 'Kd1' : Kd1}).rename('M21');
var Mn1 = Mrn1.expression ('Mrn1>Ed1 ? M11 : Mrn1<Ew1 && Mrn1 < Ed1 ? M21 : Ew1<=Mrn1<=Ed1 ? Mrn1: 0', {'Mrn1':Mrn1,'Ed1':Ed1,'M11':M11,'Ew1':Ew1,'M21':M21});
var F1 = Mn1.expression ('59.5 * (250 - Mn1)/(147.2 + Mn1)', {'Mn1':Mn1});
Map.addLayer(F1.clip(aoi), RGBPalette, 'F1', false);

1 Answer 1

0

Without going into the specifics of the algorithm itself, here's how you can structure your code to perform your calculation:

var aoi = ee.Geometry.Polygon([[
  [29.269278980991366, 51.95773487804608],
  [29.269278980991366, 51.22745421938314],
  [30.670035816928866, 51.22745421938314],
  [30.670035816928866, 51.95773487804608]
]])

var startDate = ee.Date('2024-01-01')
var endDate = ee.Date('2024-02-01') // Exclusive

var collection = ee.ImageCollection("ECMWF/ERA5_LAND/DAILY_AGGR")
  .filterDate(startDate, endDate)
  .filterBounds(aoi)

var numberOfDays = endDate.difference(startDate, 'days')
var dayOffsets = ee.List.sequence(0, numberOfDays.subtract(1))

// The initial condition of your algorithm - pick something appropriate for your case
var initialImage = ee.Image(0)

// This calls the function calculateSomething() one time for each dayOffset.
// First time around, the second argument, prevImage, will be initialImage.
var result = ee.Image(
  dayOffsets.iterate(
    calculateSomething,
    initialImage
  )
)

Map.addLayer(result, {bands: 'temperature_2m', min: 270, max: 274})
Map.centerObject(aoi)


function calculateSomething(dayOffset, prevImage) {
  // Pick up the date for current day offset
  var date = startDate.advance(ee.Number(dayOffset), 'days')

  // Filter down collection to imagery from the current date
  var image = collection
    .filterDate(date.getRange('day'))
    .median()
    
  // Implement this yourself. This currently just takes mean between previous and current image
  // 
  // When iterating over this function, the image being returned here
  // will be prevImage next time this function is called
  return ee.Image(prevImage).add(image).divide(2) 
}

https://code.earthengine.google.com/bf7c5d644285a9d85ef43fa5369e9064

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.