Calcular 23 índices de vegetación usados en agricultura con imagenes Landsat y Sentinel

Calcular 23 índices de vegetación usados en agricultura con imagenes Landsat y Sentinel

En esta entrada se va explicar como calcular rápidamente cerca de 23 índices de vegetación (IV) usados comúnmente en aplicaciones agrícolas ya sea para estimar variables como cobertura, índice de área foliar, biomasa, entre otros. Para ello vamos a usar la herramienta VICAL (Jiménez-Jiménez et al, 2022) y Google Earth Engine (GEE). La sección de ayuda de VICAL se encuentra aquí.

VICAL cuenta con tres archivos que puede importar a su Script, estos son:

// archivo para las colecciones de imagenes);
var imp = require('users/InifapCenidRaspa/VICAL:Exportaciones');
// archivo para los indices de vegetacion
var imp2= require('users/InifapCenidRaspa/VICAL:VegetationIndex'); 
// archivo para las visualizaciones
var St= require('users/InifapCenidRaspa/VICAL:Style');

Variables generales

Antes de importar el conjunto de colecciones de imágenes (Landsat o sentinel-2 o ambos) se deben definir ciertas variables que son útiles para filtrar esta colección, estos son: i) un punto o polígono; ii) *intervalo de fechas y iii) valor umbral de nubes en las imágenes*. Para este ejemplo usamos polígonos de unas parcelas agrícolas de un Distrito de Riego en México. co

var fecha = ['2021-01-01', '2022-03-18']; //Fecha inicial y final 
//polígono o punto
var table = ee.FeatureCollection("projects/calcium-verbena-328905/assets/Bate"); 
var p_nubes= 30;//Porcentaje de nubes

Selección de Colección de imagenes

Con estas colecciones de imágenes se pueden calcular series de tiempo de diferentes índices de vegetación.

Landsat

Si desea usar imágenes LandSat (4, 5, 7, 8 y 9) corregidas atmosféricamente libre de nubes se puede usar el siguiente código, donde se crea una función para unir las colecciones de imágenes y ordenarlos por fecha. Para ello se usa el archivo imp descrito al inicio de la entrada. Se realiza el ajuste espectral de los sensores del Landsat, 4, 5 y 7.

function ColeccionImagenSR(fecha, recorte, umbral)
{
  // se importan las colecciones de imagenes usando el archivo "imp"
  var L9sr = imp.ColeccionLandsatSR(fecha, 'LC09', recorte, umbral);
  var L8sr = imp.ColeccionLandsatSR(fecha, 'LC08', recorte, umbral);
  var L7sr = imp.ColeccionLandsatSR(fecha, 'LE07', recorte, umbral);
  var L5sr = imp.ColeccionLandsatSR(fecha, 'LT05', recorte, umbral);
  var L4sr = imp.ColeccionLandsatSR(fecha, 'LT04', recorte, umbral);
  //Los datos de ETM y ETM+ se ajustan espectrales a OLI y OLI-2 
  var L7a = L7sr.map(imp.TMaOLI);
  var L5a = L5sr.map(imp.TMaOLI);
  var L4a = L4sr.map(imp.TMaOLI);
  // Une las tres series de imágenes
  var serieT =L9sr.merge(L8sr).merge(L7a).merge(L5a).merge(L4a).sort('system:time_start');
    return serieT;
}
//Se importa la colección usando la función anterior 
var l8Sergio=ColeccionImagenSR(fecha, table, p_nubes);
//podemos imprimir las imagenes usando la funcion print() para ver 
//si se llevo a cabo el filtrado de la colección de imagenes (Figura 6.1) 
print (l8Sergio);

Sentinel-2

Si se desea usar imágenes Sentinel-2 corregidas atmosféricamente libre de nubes se puede usar el siguiente código.

//Se importa la colección usando el siguiente código
var S2sr = imp.ColeccionImagenSentinelSR(fecha, table, p_nubes);
//podemos imprimir las imagenes usando la funcion print() para ver 
//si se llevo a cabo el filtrado de la colección de imagenes (Figura 6.2)
print (S2sr);

Landsat y Sentinel-2

Si se desea usar imágenes LandSat y Sentinel-2 corregidas atmosféricamente libre de nubes se puede usar el siguiente código, los datos se ajustaron espectralmente a las bandas espectrales de Landsat 8. Se emplean las funciones descritas en la Sección anterior:

function ColeccionImagenAMBOS(fecha, recorte, umbral)
{
  //Se lee la función para imagenes Landsat con ajuste espectral
  var L8Conjunto=ColeccionImagenSR(fecha, recorte, umbral)
  //colección entinel
  var S2sr = imp.ColeccionImagenSentinelSR(fecha, recorte, umbral);
  var S2a = S2sr.map(imp.MSIaOLI); //Se hace el ajuste espectral de sentinel-2 a Landsat
  var serieT = S2a.merge(L8Conjunto).sort('system:time_start');
    return serieT;
}

//Se importa la colección usando el siguiente código
var S2B = ColeccionImagenAMBOS(fecha, table, p_nubes);
//podemos imprimir las imagenes usando la funcion print() para ver 
//si se llevo a cabo el filtrado de la coleccion de imagenes (Figura 6.3)
print (S2sr);

Índices de vegetación

Para usar el código de índices de vegetación de VICAL se tienen que usar el archivo imp2; y se mandan a llamar estos IV usando los nombres de la columna ExpresionGEE que se muestran en la Tabla 1.

Por ejemplo, para calcular el NDVI (Fila 12 de Tabla 1)con la colección de imágenes LandSat y Sentinel-2 de la sección anterior se usaría el siguiente código:

//Normalized Difference Vegetation Index- NDVI
var ivs = ee.ImageCollection(S2B.map(imp2.NDVI));
//para imprimir y ver la banda de NDVI
print (ivs);

Ahora, para visualizar en el mapa el NDVI de la primera imagen de la colección y recortado para la zona quedaría como se muestra a continuación; donde para ver la paleta de colores se utilizaría el archivo st de VICAL.

//Se obtiene el NDVI de la primera imagen en la colección 
var iv = ivs.first();
//Paleta de colores donde se usa el archivo st
var ivVis = {min :0, max : 1, palette : St.paletaIV};
Map.addLayer(iv.clip(table), ivVis,'NDVI'); //Indice
//se centra el mapa a la zona
 Map.centerObject(table, 13);

En la siguiente Figura se muestra el mapa de NDVI para la zona de interés. Figure74.png

Para revisar el código del ejemplo clic aquí

Si se desea visualizar el NDVI de una imagen en particular se debe convertir a lista de después mandar a llamarlo.

Table 1. Nombre de los Índices de vegetación en VICAL

Número Indice AbreviaturaExpresionGEECoeficientes
1 Atmospherically resistant vegetation index ARVI ARVI γ=1.0
1 Atmospherically resistant vegetation index ARVI ARVI γ=1.0
2 Adjusted transformed soil-adjusted vegetation index ATSAVI ATSAVI
3 Difference vegetation index DVI DVI
4 Enhanced vegetation index EVI EVI C1=6.0, C2= 7.5; L=1.0
5 Enhanced vegetation index EVI2 EVI2 C1=2.4
6 Green normalized difference vegetation index GNDVI GNDVI
7 Modified soil adjusted vegetation index MSAVI2 MSAVI2
8 Moisture stress index MSI MSI
9 Modified triangular vegetation index MTVI MTVI
10 Modified triangular vegetation index-2 MTVI2 MTVI2
11 Normalized difference tillage index (NDTI) NDTI NDTI
12 Normalized difference vegetation index NDVI NDVI
13 Normalized difference water index NDWI NDWI
14 Optimized soil adjusted vegetation index OSAVI OSAVI X=0.16
15 Renormalized difference vegetation index RDVI RDVI
16 Redness index RI RI
17 Ratio vegetation index RVI RVI
18 Soil adjusted vegetation index SAVI SAVI L=0.5
19 Triangular vegetation index TVI TVI
20 Transformed soil adjusted vegetation index TSAVI TSAVI a= 1 ; b=0;
21 Visible atmospherically resistant index VARI VARI
22 Vegetation index number or simple ratio VIN VIN
23 Wide dynamic range vegetation index WDRVI WDRVI α=0.2

Did you find this article valuable?

Support SERGIO IVÁN JIMÉNEZ JIMÉNEZ by becoming a sponsor. Any amount is appreciated!