del 1 al 566 no tocar!!!

Índice

2 Calculo de Areas 2018-01-01


library(rgee)
ee_Initialize()
## -- rgee 1.0.9 --------------------------------------- earthengine-api 0.1.259 -- 
##  v email: not_defined
##  v Initializing Google Earth Engine:
 v Initializing Google Earth Engine:  DONE!
## 
 v Earth Engine user: users/tarredwall 
## --------------------------------------------------------------------------------

1 Calculo de Areas “2019-01-01” “2019-04-01”

(volver al índice)

Extraemos la imagen satelital Sentinel filtrada por la region de los lagos

mask <- st_read("region_los_lagos.shp", 
        quiet = TRUE) %>% 
        sf_as_ee()

# convertimos el shp en geometria:
region <- mask$geometry()

start <- ee$Date("2019-01-01")
finish <- ee$Date("2019-04-01")
cc <- 20

sentinel1 = ee$ImageCollection('COPERNICUS/S2')$filterDate(start, finish)$filterBounds(region)$
  filter(ee$Filter$lt("CLOUDY_PIXEL_PERCENTAGE", cc))

first <- sentinel1$median()

# definimos los parametros de visualizacion:
vizParams <- list(
  bands = c("B8","B5" , "B3"),
  #  bands = c("B2", "B3"),
  min = 100,
  max = 1000,
  gamma = 2
)
Map$setCenter(lon = -73.079, lat = -42.611, zoom = 6)
Map$addLayer(first, vizParams, "Landsat 8 image")

Creamos un cuadrado arbitrario

box_2 <- ee$Geometry$Rectangle(coords = c(-73,-43,-72.2,-42),
                             proj = "EPSG:4326",
                             geodesic = FALSE)
Map$setCenter(lon = -73.079, lat = -42.611, zoom = 7)
Map$addLayer(box_2)
# obtenenos informacion de las medias de una imagen
# mean <- first$reduceRegion(
#   reducer = ee$Reducer$mean(),
#   geometry = box_2,
#   scale = 30,
#   bestEffort= TRUE
# )
# print(mean$getInfo())

construimos un raster con NDGI mayor a 0

Lo que aca hacemos es construir una imagen raster con una columna llamada nd con dos categorias para cada pixel: 0 si es rojo, 1 si es verde

getNDGI <- function(image) 
{
  #image$normalizedDifference(c("B3", "B4")) 
     image$normalizedDifference(c("B3", "B4")) > 0
}

ndgi_01 <- getNDGI(first)
## Warning: Ops.ee.image.Image will be deprecated in rgee v.1.1.0. Please install
## rgeeExtra (https://github.com/r-earthengine/rgeeExtra). Deeply sorry for the
## inconveniences.
ndgiParams <- list(palette = c(
  "#d73027", "#f46d43", "#fdae61",
  "#fee08b", "#d9ef8b", "#a6d96a",
  "#66bd63", "#1a9850"
))

Map$setCenter(lon = -73.079, lat = -42.611, zoom = 6)
Map$addLayer(ndgi_01, ndgiParams, "NDGI")

Intersectamos con los límites políticos de los glaciares

mask_0 <- st_read("Lim_glaciares_v2_Simplify.shp",
        quiet = TRUE) %>%
        sf_as_ee()

# convertimos el shp en geometria:
region_0 <- mask_0$geometry()

sale_4 <- ndgi_01$clip(region_0)
Map$setCenter(lon = -73.079, lat = -42.611, zoom = 7)
Map$addLayer(sale_4, ndgiParams, "NDGI")

Los rojos interiores son glaciares 0 = rojo

muestras_in <- sale_4$sampleRegions(
    collection = ee$Feature(region_0),
    scale = 1500,
    tileScale = 16,
    geometries = TRUE
  )

muestras_in_gee_rojas = muestras_in$filter(ee$Filter$eq('nd', 0))
#muestras_in_gee_verdes = muestras_in$filter(ee$Filter$eq('nd', 1))

#Obtencion de puntos glaciares:
Map$setCenter(-72.7, -43.4, 6)
Map$addLayer(
  eeObject =  muestras_in_gee_rojas,
  visParams = {},
  name = "puntos glaciares"
)

Obtencion de puntos no glaciares:

box_2 <- ee$Geometry$Rectangle(coords = c(-73,-43,-72.2,-42),
                             proj = "EPSG:4326",
                             geodesic = FALSE)
Map$setCenter(lon = -73.079, lat = -42.611, zoom = 7)
Map$addLayer(box_2)
diff1 <- box_2$difference(region_0, ee$ErrorMargin(1))
Map$addLayer(diff1, list(color = "FFFF00"), "diff1")
sale_3 <- ndgi_01$clip(diff1)
Map$setCenter(lon = -73.079, lat = -42.611, zoom = 7)
Map$addLayer(sale_3, ndgiParams, "NDGI")

MUESTRAS

muestras_out <- sale_3$sampleRegions(
    collection = ee$Feature(diff1),
    scale = 1500,
    tileScale = 16,
    geometries = TRUE
  )

muestras_in_gee_verdes = muestras_in$filter(ee$Filter$eq('nd', 1))
#muestras_in_gee_rojas = muestras_in$filter(ee$Filter$eq('nd', 0))

#Obtencion de puntos glaciares:
Map$setCenter(-72.7, -43.4, 6)
Map$addLayer(
  eeObject =  muestras_in_gee_verdes,
  visParams = {},
  name = "puntos glaciares"
)
union_total <- muestras_in_gee_verdes$merge(muestras_in_gee_rojas)

3 Aplicamos Random Forest

(volver al índice)

mask <- st_read("region_los_lagos.shp", 
        quiet = TRUE) %>% 
        sf_as_ee()

# convertimos el shp en geometria:
region <- mask$geometry()

start <- ee$Date("2019-01-01")
finish <- ee$Date("2019-04-01")
cc <- 20

sentinel1 = ee$ImageCollection('COPERNICUS/S2')$filterDate(start, finish)$filterBounds(region)$
  filter(ee$Filter$lt("CLOUDY_PIXEL_PERCENTAGE", cc))

first <- sentinel1$median()

# Use these bands for prediction.
bands <- c("B2", "B3", "B4", "B8")

# This property of the table stores the land cover labels.
label <- "nd"

# Overlay the points on the imagery to get training.
training <- first$select(bands)$sampleRegions(
  collection = union_total,
  properties = list(label),
  scale = 10
)

# Train a CART classifier with default parameters.
trained <- ee$Classifier$smileRandomForest(10)$train(training, label, bands)

# Classify the image with the same bands used for training.
classified <- first$select(bands)$classify(trained)

# Viz parameters.
# viz_img <- list(bands = c("B4", "B3", "B2"), max = 0.4)
vizParams <- list(
  bands = c("B8","B5" , "B3"),
  #  bands = c("B2", "B3"),
  min = 100,
  max = 1000,
  gamma = 2
)

viz_class <- list(palette = c("red", "green", "blue"), min = 0, max = 2)

Map$addLayer(first, vizParams, name = "image") +
Map$addLayer(classified, viz_class, name = "classification")
verdes_in <- classified$select('classification')$eq(0)

# Calcula el área de píxeles en kilómetros cuadrados 
verdes_in_km = verdes_in$multiply(ee$Image$pixelArea())$divide(1000*1000)

verdes_in_km_area = verdes_in_km$reduceRegion(
    reducer= ee$Reducer$sum(),
    geometry= box_2,
    scale= 100  ,
    maxPixels= 1e9
  )

verdes_in_km_area$getInfo()
## $classification
## [1] 346.3606

2 Calculo de Areas “2018-01-01” “2018-04-01”

(volver al índice)

Extraemos la imagen satelital Sentinel filtrada por la region de los lagos

mask <- st_read("region_los_lagos.shp", 
        quiet = TRUE) %>% 
        sf_as_ee()

# convertimos el shp en geometria:
region <- mask$geometry()

start <- ee$Date("2018-01-01")
finish <- ee$Date("2018-04-01")
cc <- 20

sentinel1 = ee$ImageCollection('COPERNICUS/S2')$filterDate(start, finish)$filterBounds(region)$
  filter(ee$Filter$lt("CLOUDY_PIXEL_PERCENTAGE", cc))

first <- sentinel1$median()

# definimos los parametros de visualizacion:
vizParams <- list(
  bands = c("B8","B5" , "B3"),
  #  bands = c("B2", "B3"),
  min = 100,
  max = 1000,
  gamma = 2
)
Map$setCenter(lon = -73.079, lat = -42.611, zoom = 6)
Map$addLayer(first, vizParams, "Landsat 8 image")
box_2 <- ee$Geometry$Rectangle(coords = c(-73,-43,-72.2,-42),
                             proj = "EPSG:4326",
                             geodesic = FALSE)
Map$setCenter(lon = -73.079, lat = -42.611, zoom = 7)
Map$addLayer(box_2)
# obtenenos informacion de las medias de una imagen
# mean <- first$reduceRegion(
#   reducer = ee$Reducer$mean(),
#   geometry = box_2,
#   scale = 30,
#   bestEffort= TRUE
# )
# print(mean$getInfo())

construimos un raster con NDGI mayor a 0

Lo que aca hacemos es construir una imagen raster con una columnallamada nd con dos categorias para cada pixel: 0 si es rojo, 1 si es verde

getNDGI <- function(image) 
{
  #image$normalizedDifference(c("B3", "B4")) 
     image$normalizedDifference(c("B3", "B4")) > 0
}

ndgi_01 <- getNDGI(first)
## Warning: Ops.ee.image.Image will be deprecated in rgee v.1.1.0. Please install
## rgeeExtra (https://github.com/r-earthengine/rgeeExtra). Deeply sorry for the
## inconveniences.
ndgiParams <- list(palette = c(
  "#d73027", "#f46d43", "#fdae61",
  "#fee08b", "#d9ef8b", "#a6d96a",
  "#66bd63", "#1a9850"
))

Map$setCenter(lon = -73.079, lat = -42.611, zoom = 6)
Map$addLayer(ndgi_01, ndgiParams, "NDGI")

Intersectamos con los límites políticos de los glaciares

mask_0 <- st_read("Lim_glaciares_v2_Simplify.shp",
        quiet = TRUE) %>%
        sf_as_ee()

# convertimos el shp en geometria:
region_0 <- mask_0$geometry()

sale_4 <- ndgi_01$clip(region_0)
Map$setCenter(lon = -73.079, lat = -42.611, zoom = 7)
Map$addLayer(sale_4, ndgiParams, "NDGI")

Los pixeles rojos son glaciares

muestras_in <- sale_4$sampleRegions(
    collection = ee$Feature(region_0),
    scale = 1500,
    tileScale = 16,
    geometries = TRUE
  )

muestras_in_gee_rojas = muestras_in$filter(ee$Filter$eq('nd', 0))
#muestras_in_gee_verdes = muestras_in$filter(ee$Filter$eq('nd', 1))

#Obtencion de puntos glaciares:
Map$setCenter(-72.7, -43.4, 6)
Map$addLayer(
  eeObject =  muestras_in_gee_rojas,
  visParams = {},
  name = "puntos glaciares"
)

Obtencion de puntos no glaciares:

box_2 <- ee$Geometry$Rectangle(coords = c(-73,-43,-72.2,-42),
                             proj = "EPSG:4326",
                             geodesic = FALSE)
Map$setCenter(lon = -73.079, lat = -42.611, zoom = 7)
Map$addLayer(box_2)
diff1 <- box_2$difference(region_0, ee$ErrorMargin(1))
Map$addLayer(diff1, list(color = "FFFF00"), "diff1")
sale_3 <- ndgi_01$clip(diff1)
Map$setCenter(lon = -73.079, lat = -42.611, zoom = 7)
Map$addLayer(sale_3, ndgiParams, "NDGI")
muestras_out <- sale_3$sampleRegions(
    collection = ee$Feature(diff1),
    scale = 1500,
    tileScale = 16,
    geometries = TRUE
  )

muestras_in_gee_verdes = muestras_in$filter(ee$Filter$eq('nd', 1))
#muestras_in_gee_rojas = muestras_in$filter(ee$Filter$eq('nd', 0))

#Obtencion de puntos glaciares:
Map$setCenter(-72.7, -43.4, 6)
Map$addLayer(
  eeObject =  muestras_in_gee_rojas,
  visParams = {},
  name = "puntos glaciares"
)
union_total <- muestras_in_gee_verdes$merge(muestras_in_gee_rojas)

3 Aplicamos Random Forest

(volver al índice)

mask <- st_read("region_los_lagos.shp", 
        quiet = TRUE) %>% 
        sf_as_ee()

# convertimos el shp en geometria:
region <- mask$geometry()

start <- ee$Date("2018-01-01")
finish <- ee$Date("2018-04-01")
cc <- 20

sentinel1 = ee$ImageCollection('COPERNICUS/S2')$filterDate(start, finish)$filterBounds(region)$
  filter(ee$Filter$lt("CLOUDY_PIXEL_PERCENTAGE", cc))

first <- sentinel1$median()

# Use these bands for prediction.
bands <- c("B2", "B3", "B4", "B8")

# This property of the table stores the land cover labels.
label <- "nd"

# Overlay the points on the imagery to get training.
training <- first$select(bands)$sampleRegions(
  collection = union_total,
  properties = list(label),
  scale = 10
)

# Train a CART classifier with default parameters.
trained <- ee$Classifier$smileRandomForest(10)$train(training, label, bands)

# Classify the image with the same bands used for training.
classified <- first$select(bands)$classify(trained)

# Viz parameters.
# viz_img <- list(bands = c("B4", "B3", "B2"), max = 0.4)
vizParams <- list(
  bands = c("B8","B5" , "B3"),
  #  bands = c("B2", "B3"),
  min = 100,
  max = 1000,
  gamma = 2
)

viz_class <- list(palette = c("red", "green", "blue"), min = 0, max = 2)

Map$addLayer(first, vizParams, name = "image") +
Map$addLayer(classified, viz_class, name = "classification")
verdes_in <- classified$select('classification')$eq(0)

# Calcula el área de píxeles en kilómetros cuadrados 
verdes_in_km = verdes_in$multiply(ee$Image$pixelArea())$divide(1000*1000)

verdes_in_km_area = verdes_in_km$reduceRegion(
    reducer= ee$Reducer$sum(),
    geometry= box_2,
    scale= 100  ,
    maxPixels= 1e9
  )

verdes_in_km_area$getInfo()
## $classification
## [1] 368.7491
III IMPORTANTE LO ULTIMO

Ahora lo que necesitamos es la obtencion clasificada de glaciar por glaciar:

Lo primero que podriamos hacer es obtener muestras dentro de los limites glciares que sean rojos en la region de arica y tarapaca

mask <- st_read("regiones_separadas/Region_15.shp", 
        quiet = TRUE) %>% 
        sf_as_ee()

# convertimos el shp en geometria:
region <- mask$geometry()
Map$setCenter(lon = -69.0055, lat = -18.5774, zoom = 7)
Map$addLayer(region)
start <- ee$Date("2019-03-01")
finish <- ee$Date("2020-03-01")
cc <- 20

sentinel1 = ee$ImageCollection('COPERNICUS/S2')$filterDate(start, finish)$filterBounds(region)$
  filter(ee$Filter$lt("CLOUDY_PIXEL_PERCENTAGE", cc))

first <- sentinel1$median()

# definimos los parametros de visualizacion:
vizParams <- list(
  bands = c("B8","B5" , "B3"),
  #  bands = c("B2", "B3"),
  min = 100,
  max = 1000,
  gamma = 2
)
Map$setCenter(lon = -69.0055, lat = -18.5774, zoom = 7)
Map$addLayer(first,vizParams, "Landsat 8 image")
getNDGI <- function(image) 
{
  #image$normalizedDifference(c("B3", "B4")) 
     image$normalizedDifference(c("B3", "B4")) < 0
}

ndgi_01 <- getNDGI(first)
## Warning: Ops.ee.image.Image will be deprecated in rgee v.1.1.0. Please install
## rgeeExtra (https://github.com/r-earthengine/rgeeExtra). Deeply sorry for the
## inconveniences.
ndgiParams <- list(palette = c(
  "#d73027", "#f46d43", "#fdae61",
  "#fee08b", "#d9ef8b", "#a6d96a",
  "#66bd63", "#1a9850"
))

Map$setCenter(lon = -69.0055, lat = -18.5774, zoom = 7)
Map$addLayer(ndgi_01, ndgiParams, "NDGI")

Intersectamos con los límites políticos de los glaciares

glaciar_01 <- st_read("GlaciaresxRegion/Simple_ARICA_Y_PARINACOTA.shp",
        quiet = TRUE) %>%
        sf_as_ee()

# convertimos el shp en geometria:
region_0 <- glaciar_01$geometry()

sale_4 <- ndgi_01$clip(region_0)
Map$setCenter(lon = -69.0055, lat = -18.5774, zoom = 7)
Map$addLayer(sale_4, ndgiParams, "NDGI")

muestras de areas rojas de los graciares

Los rojos interiores son glaciares 0 = rojo

muestras_in <- sale_4$sampleRegions(
    collection = ee$Feature(region_0),
    scale = 100,
    tileScale = 16,
    geometries = TRUE
  )


###Los pixels rojos son los 0


muestras_in_gee_rojas = muestras_in$filter(ee$Filter$eq('nd', 0))

#Obtencion de puntos glaciares:
Map$setCenter(lon = -69.0055, lat = -18.5774, zoom = 7)
Map$addLayer(
  eeObject =  muestras_in_gee_rojas,
  visParams = {},
  name = "puntos glaciares"
)

debemos obtener los puntos verdes de afuera:

mask <- st_read("regiones_separadas/Region_15.shp", 
        quiet = TRUE) %>% 
        sf_as_ee()

# convertimos el shp en geometria:
region <- mask$geometry()
Map$setCenter(lon = -69.0055, lat = -18.5774, zoom = 7)
Map$addLayer(region)
diff1 <- region$difference(region_0, ee$ErrorMargin(1))
Map$addLayer(diff1, list(color = "FFFF00"), "diff1")
sale_3 <- ndgi_01$clip(diff1)
Map$setCenter(lon = -69.0055, lat = -18.5774, zoom = 7)
Map$addLayer(sale_3, ndgiParams, "NDGI")

Ahora extraemos las muestras de color verde, esto es igual a 1 en la zonas de fuera los limites graciares.

muestras_in <- sale_3$sampleRegions(
    collection = ee$Feature(diff1),
    scale = 1000,
    tileScale = 16,
    geometries = TRUE
  )

muestras_in_gee_verdes = muestras_in$filter(ee$Filter$eq('nd', 1))
#muestras_in_gee_verdes = muestras_in$filter(ee$Filter$eq('nd', 1))

#Obtencion de puntos glaciares:
Map$setCenter(lon = -69.0055, lat = -18.5774, zoom = 7)
Map$addLayer(eeObject =  muestras_in_gee_verdes,visParams = {}, name = "puntos glaciares")
union_total <- muestras_in_gee_verdes$merge(muestras_in_gee_rojas)

3 Aplicamos Random Forest

(volver al índice)

mask <- st_read("regiones_separadas/Region_15.shp", 
        quiet = TRUE) %>% 
        sf_as_ee()

# convertimos el shp en geometria:
region <- mask$geometry()

start <- ee$Date("2018-01-01")
finish <- ee$Date("2018-04-01")
cc <- 20

sentinel1 = ee$ImageCollection('COPERNICUS/S2')$filterDate(start, finish)$filterBounds(region)$
  filter(ee$Filter$lt("CLOUDY_PIXEL_PERCENTAGE", cc))

first <- sentinel1$median()

# Use these bands for prediction.
bands <- c("B2", "B3", "B4", "B8")

# This property of the table stores the land cover labels.
label <- "nd"

# Overlay the points on the imagery to get training.
training <- first$select(bands)$sampleRegions(
  collection = union_total,
  properties = list(label),
  scale = 10
)

# Train a CART classifier with default parameters.
trained <- ee$Classifier$smileRandomForest(10)$train(training, label, bands)

# Classify the image with the same bands used for training.
classified <- first$select(bands)$classify(trained)

# Viz parameters.
# viz_img <- list(bands = c("B4", "B3", "B2"), max = 0.4)
vizParams <- list(
  bands = c("B8","B5" , "B3"),
  #  bands = c("B2", "B3"),
  min = 100,
  max = 1000,
  gamma = 2
)

viz_class <- list(palette = c("red", "green", "blue"), min = 0, max = 2)

Map$addLayer(first, vizParams, name = "image") +
Map$addLayer(classified, viz_class, name = "classification")

Este paso en muy, muy importante, y es el area poligonal oficial donde nosotros estudiamos la variacion glaciar. Ahora debemos cortar por los limites graciares y obtener los shps. No nos interesan las zonas no glaciares!!!

mask_0 <- st_read("GlaciaresxRegion/Simple_ARICA_Y_PARINACOTA.shp")
## Reading layer `Simple_ARICA_Y_PARINACOTA' from data source `C:\Users\usuario\Desktop\ds\ds_rgee\cart\GlaciaresxRegion\Simple_ARICA_Y_PARINACOTA.shp' using driver `ESRI Shapefile'
## Simple feature collection with 327 features and 55 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: -69.81277 ymin: -19.03571 xmax: -68.91922 ymax: -17.64658
## Geodetic CRS:  WGS 84
reg_15_2 <-  mask_0 %>% 
filter(COD_GLA =="CL101010001")
varia <- reg_15_2 %>%
        sf_as_ee()
# convertimos el shp en geometria:
region_0 <- varia$geometry()

sale_4 <- classified$clip(region_0)
Map$setCenter(lon = -69.0055, lat = -18.5774, zoom = 7)
Map$addLayer(sale_4, viz_class, "NDGI")
fff <- sale_4$select('classification')$eq(0)

class_areas = fff$reduceRegion(
reducer= ee$Reducer$sum(),
geometry= region_0,
scale= 30 ,
maxPixels= 1e9
)

class_areas$getInfo()
## $classification
## [1] 282.3686

Anexo.

glaciar_001 <- read_sf("GlaciaresxRegion/Simple_ARICA_Y_PARINACOTA.shp")
 glaciar_001
## Simple feature collection with 327 features and 55 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: -69.81277 ymin: -19.03571 xmax: -68.91922 ymax: -17.64658
## Geodetic CRS:  WGS 84
## Warning: `...` is not empty.
## 
## We detected these problematic arguments:
## * `needs_dots`
## 
## These dots only exist to allow future extensions and should be empty.
## Did you misspecify an argument?
## # A tibble: 327 x 56
##       Id COD_GLA NOMBRE CLASIFICA REGION COMUNA DATUM HUSO    ESTE  NORTE
##    <int> <chr>   <chr>  <chr>     <chr>  <chr>  <chr> <chr>  <int>  <int>
##  1    11 CL1010~ S/N    GLACIAR ~ ARICA~ PUTRE  WGS ~ 19SUR 490679 8.01e6
##  2     1 CL1010~ VOLCA~ GLACIAR ~ ARICA~ PUTRE  WGS ~ 19SUR 484875 7.99e6
##  3     8 CL1010~ S/N    GLACIAR ~ ARICA~ PUTRE  WGS ~ 19SUR 494376 7.97e6
##  4    11 CL1010~ S/N    GLACIAR ~ ARICA~ PUTRE  WGS ~ 19SUR 494388 7.97e6
##  5    14 CL1010~ S/N    GLACIAR ~ ARICA~ PUTRE  WGS ~ 19SUR 494057 7.97e6
##  6    18 CL1010~ S/N    GLACIARE~ ARICA~ PUTRE  WGS ~ 19SUR 494908 7.97e6
##  7     2 CL1010~ S/N    GLACIARE~ ARICA~ PUTRE  WGS ~ 19SUR 493272 7.98e6
##  8    37 CL1010~ VOLCA~ GLACIAR ~ ARICA~ PUTRE  WGS ~ 19SUR 484474 7.99e6
##  9    36 CL1010~ VOLCA~ GLACIAR ~ ARICA~ PUTRE  WGS ~ 19SUR 484792 7.99e6
## 10    12 CL1010~ S/N    GLACIARE~ ARICA~ PUTRE  WGS ~ 19SUR 490628 8.01e6
## # ... with 317 more rows, and 46 more variables: FUENTE_DIG <chr>,
## #   FUENTE_FEC <chr>, INVENT_FEC <chr>, COD_BNA <chr>, AREA_Km2 <dbl>,
## #   COD_CUEN <chr>, NOMB_CUEN <chr>, COD_REGION <chr>, COD_SCUEN <chr>,
## #   COD_SSCUEN <chr>, N_CUENDRE <chr>, PRECM <chr>, LARGO_PROM <dbl>,
## #   LMAXTOTAL <dbl>, L_MAX_ACUM <dbl>, L_MAX_EXP <dbl>, LMAXABLAC <dbl>,
## #   ANCHO_PROM <dbl>, ESP_MED <dbl>, AREA_EXP <dbl>, AREA_CUB <dbl>,
## #   AREA_ABLAC <dbl>, AREA_ACUM <dbl>, VOL_M3 <dbl>, HMAX <dbl>, HMEDIA <dbl>,
## #   HMINTOTAL <dbl>, HMEDABLAC <dbl>, HMINEXP <dbl>, LATITUD <dbl>,
## #   LONGITUD <dbl>, ORIENACUM <chr>, ORIENABLAC <chr>, ORIENTA <chr>,
## #   CLAS_WGI <chr>, CLAS_2_CUB <chr>, ZONA_GLACI <chr>, HMEDIANA <dbl>,
## #   PENDIENTE <dbl>, ERROR_KM2 <dbl>, ERROR_PORC <dbl>, InPoly_FID <dbl>,
## #   SimPgnFlag <int>, MaxSimpTol <dbl>, MinSimpTol <dbl>,
## #   geometry <MULTIPOLYGON [°]>