Skip to content
Snippets Groups Projects
Commit 6e366f54 authored by Francesco Sabatini's avatar Francesco Sabatini
Browse files

Fixed bug in A98 optimized

parent 5b0a64e8
No related branches found
No related tags found
No related merge requests found
...@@ -25,11 +25,8 @@ PredExtr <- function(x.shp, myfunction=NA, output=NA, ...@@ -25,11 +25,8 @@ PredExtr <- function(x.shp, myfunction=NA, output=NA,
print(paste("output will be:", output)) print(paste("output will be:", output))
if(typp=="raster"){ mypredictor <- raster(toextract)} else if(typp=="raster"){ mypredictor <- raster(toextract)} else
mypredictor <- readOGR(toextract) mypredictor <- readOGR(toextract)
colnames(mypredictor@data)
header.shp <- readOGR(x.shp) header.shp <- readOGR(x.shp)
crs(mypredictor) <- crs(header.shp) #should be verified beforehand! crs(mypredictor) <- crs(header.shp) #should be verified beforehand!
## Divide in chunks if requested ## Divide in chunks if requested
if(chunkn>1 & !is.na(chunk.i)){ if(chunkn>1 & !is.na(chunk.i)){
...@@ -47,14 +44,16 @@ PredExtr <- function(x.shp, myfunction=NA, output=NA, ...@@ -47,14 +44,16 @@ PredExtr <- function(x.shp, myfunction=NA, output=NA,
cl <- makeForkCluster(ncores, outfile="" ) cl <- makeForkCluster(ncores, outfile="" )
registerDoParallel(cl) registerDoParallel(cl)
print("start main foreach loop")
if(typp=="raster"){ if(typp=="raster"){
out <- foreach(i=1:length(header.shp), .packages=c('raster'), .combine=rbind) %dopar% { print("start main foreach loop")
out <- foreach(i=1:length(header.shp), .packages=c('raster'), .combine=rbind) %dopar% {
tmp <- raster::extract(mypredictor, header.shp[i,], tmp <- raster::extract(mypredictor, header.shp[i,],
buffer=min(10000, max(250, header.shp@data$lc_ncrt[i])), fun=myfunction) buffer=min(10000, max(250, header.shp@data$lc_ncrt[i])), fun=myfunction)
} }
if(!is.na(output)) {write.csv(out, file = output)} else{return(out)} if(!is.na(output)) {write.csv(out, file = output)} else{return(out)}
} else { } else {
print("Match sequentially")
out <- sp::over(x=header.shp, y=mypredictor) out <- sp::over(x=header.shp, y=mypredictor)
i.toassign <- which(is.na(out[,1])) i.toassign <- which(is.na(out[,1]))
toassign <- header.shp[i.toassign,] toassign <- header.shp[i.toassign,]
...@@ -67,7 +66,7 @@ PredExtr <- function(x.shp, myfunction=NA, output=NA, ...@@ -67,7 +66,7 @@ PredExtr <- function(x.shp, myfunction=NA, output=NA,
as_Spatial() as_Spatial()
} else { } else {
mypredictor.expl <- mypredictor} mypredictor.expl <- mypredictor}
print("start main foreach loop")
nearest <- foreach(i=1:length(toassign), .packages=c('raster'), .combine=rbind) %dopar% { nearest <- foreach(i=1:length(toassign), .packages=c('raster'), .combine=rbind) %dopar% {
print(i) print(i)
## create a subset of geoentities based on a 5° buffer radius around each target plot. ## create a subset of geoentities based on a 5° buffer radius around each target plot.
...@@ -84,13 +83,16 @@ PredExtr <- function(x.shp, myfunction=NA, output=NA, ...@@ -84,13 +83,16 @@ PredExtr <- function(x.shp, myfunction=NA, output=NA,
} }
) )
# find nearest neighbour # find nearest neighbour
nearest.tmp <- tryCatch(geosphere::dist2Line(toassign[i,], tmp.mypredictor), nearest.tmp <- tryCatch(tmp.mypredictor@data[geosphere::dist2Line(toassign[i,], tmp.mypredictor)[,"ID"],],
error = function(e){data.frame(i=i, distance=NA, lon=NA, lat=NA,ID=NA)} error = function(e){
ee <- myptedictor@data[1,]
ee[1,] <- rep(NA, ncol(mypredictor))
}
) )
return(nearest.tmp) return(nearest.tmp)
} }
out[i.toassign,] <- tmp.mypredictor@data[nearest[,"ID"],] out[i.toassign,] <- nearest
} }
if(!is.na(output)) {write.csv(out, file = output)} else{return(out)} if(!is.na(output)) {write.csv(out, file = output)} else{return(out)}
} }
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment