I am a beginner user of R as GIS. My main aim is detect only the upper (highest) edge between two classes (forest, non forest) trought of an elevational gradient.
To do the edge detection I used two raster layers, the data:
Veg_DEM_Test <- "https://ift.tt/38VZCfl"
Veg_Class_Test <- "https://ift.tt/3rOSlGN"
1.-Veg_DEM_Test 1km resolution (r1).
2.-Veg_Class_Test 1km resolution (r2). A Raster layer with the classes, where the class 1 is non-forest and the class 2 is forest.
Then, I followed the next steps:
First: I created a focal matrix to run through all overlayed raster. Second: I maked a function with conditional statements for detection of upper edge between forest non-forest.
Where:
x=Veg_DEM_Test and y=Veg_Class_Test.
The conditional statements to detect the highest limit are: When the focal cell of the matrix has the vegetation class as forest==2, ej.(y[5]==2) and at same time this central cell are below (elevation) of at least one from the 4 cells neighbours. This upper neighbour cell is class vegetation non forest==1 ej. (x[5]<x[2]&x[2]==1). The result, should be a contour line that outline the upper edges of forest with the non forest. Like this image.
Expected result https://www.dropbox.com/s/sg7t3spm0f6vjy1/Veg_upper_edge.jpg?dl=0
But I am stuck here I have a issues to do it: The wrote code looks right, apparently runs well but doesn't produce any raster result.
Please, hoping someone can point me in the right direction.
Pd: I vectorized the overlayed raster, because at the first attempt R show me a error message: Error in (function (x, fun, filename = "", recycle = TRUE, forcefun = FALSE, : cannot use this formula, probably because it is not vectorized
Find below the code.
R version: 3.6.3 R Studio Version: 1.3.959 Raster package: 3.4-5
library(raster)
#URL data
Veg_DEM_Test <- "https://www.dropbox.com/s/tklr5roxhd6v1n6/Veg_DEM_Test.tif?dl=0"
Veg_Class_Test <- "https://www.dropbox.com/s/t1s18rump8zfreh/Veg_Class_Test.tif?dl=0"
#Download data
download.file(Veg_DEM_Test, destfile=paste0(getwd(),"/", "Veg_DEM_Test.tif"), method="auto", mode="wb", timeout="6000")
download.file(Veg_Class_Test, destfile=paste0(getwd(),"/", "Veg_Class_Test.tif"), method="auto", mode="wb", timeout="6000")
#Load Rasters
r1 <- raster("...\\Veg_DEM_Test.tif")
r1 #View the characteristics of the layer
extent(r1) #view their characteristics of the Tile
r2 <- raster("...\\Veg_Class_Test.tif")
r2 #View the characteristics of the layer
extent(r2) #view their characteristics of the Tile
extent (r1)== extent (r2)
#Creating the Window Matrix Move
f = matrix(c(0,1,0,
1,0,1,
0,1,0), nrow=3, byrow = T)
#Creating function and upper edge
rast_edge_TEST <-focal(na.rm = T, w=f, pad= T, (overlay(r1, r2,fun=Vectorize(function(x,y){ifelse(y[5]==2 & x[5]<x[2] & y[2]==1,-1,
ifelse(y[5]==2 & x[5]<x[4] & y[4]==1,-1,
ifelse(y[5]==2 & x[5]<x[6] & y[6]==1,-1,
ifelse(y[5]==2 & x[5]<x[8] & y[8]==1,-1,NA))))}))))
rast_edge_TEST
#Plotting
x11()
plot(rast_edge_TEST)
Aucun commentaire:
Enregistrer un commentaire