dimanche 3 janvier 2021

Detecting upper edge from two overlayed rasters, focal matrix and conditional (ifelse) statements

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