\[TALLER\ INTERPOLACION\]

install.packages( c( “sf”, “stars”,“leaflet”, “gstat”,“automap”,“raster”, “RColorBrewer”))

library(sf) library(stars) library(leaflet) library(gstat) library(automap) library(raster) library(RColorBrewer)

getwd()

setwd(“C:/Users/asus_/OneDrive/Escritorio/INTERPOLACION”)

ruta_aoi=“./municipios_muestreo.shp” ruta_puntos=“./puntos_muestreo.shp” ruta_raster=“./dem_srtm_9377.tif”

raster_dem=read_stars(ruta_raster, RasterIO = list(bands = 1))
raster_dem
plot(raster_dem)
crs_raster=st_crs(raster_dem)
crs_raster

puntos = st_read(ruta_puntos)
crs_puntos=st_crs(puntos)
crs_puntos


if (identical(crs_raster, crs_puntos)) {
  cat("El sistema de referencia de coordenadas del raster y de los puntos es el mismo:\n")
  print(crs_raster)
} else {
  cat("Los sistemas de referencia de coordenadas son diferentes.\n")
  cat("Raster:\n")
  print(crs_raster)
  cat("Puntos:\n")
  print(crs_puntos)
}

puntos_reproj <-  st_transform(puntos, crs_raster)

(puntos_reproj)


var="pH" 
puntos[var] 



st_crs(puntos_reproj)
crs(raster_dem)


  
  g = gstat(formula = pH ~ 1, data = puntos_reproj,set=list(idp=2))


  z = predict(g, raster_dem)

print(z)
z = z["var1.pred",,]


names(z) = "pH"

b = seq(4, 7.5, 0.1)


plot(z, breaks = b, col = hcl.colors(length(b)-1, "Spectral"), reset = FALSE)


plot(st_geometry(puntos_reproj), pch = 3, add = TRUE)


contour(z, breaks = b, add = TRUE)


interpolacion_idw = "./ph_idw.tif"

(formato geotiff) write_stars(z, dsn = interpolacion_idw)

```
v_emp_ok = variogram(pH ~ 1, puntos_reproj)

plot(v_emp_ok)

v_mod_ok = autofitVariogram(pH ~ 1, as(puntos_reproj, “Spatial”))

g = gstat(formula = pH ~ 1, model = v_mod_ok$var_model, data = puntos)

z = predict(g, raster_dem)

print(z)

prediccion = z[“var1.pred”,,] names(prediccion) = “pH”

b_predict = seq(4, 7, 0.1)

plot(prediccion, breaks = b_predict, col = hcl.colors(length(b_predict)-1, “Spectral”), reset = FALSE)

varianza = z[“var1.var”,,] names(varianza) = “varianza pH”

b_var = seq(0.1, 1.4, 0.1)

plot(varianza, breaks = b_var, col = hcl.colors(length(b_var)-1, “Spectral”), reset = FALSE)