rparallel-processingterraclamp# R: Use terra to rescale large raster files and run in parallel

I have some very large raster files and want to rescale the values of these rasters to 0-1. Here are the information of my large rasters:

```
class : SpatRaster
dimensions : 161239, 132862, 1 (nrow, ncol, nlyr)
resolution : 25, 25 (x, y)
extent : 2634975, 5956525, 1385025, 5416000 (xmin, xmax, ymin, ymax)
coord. ref. : +proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +ellps=GRS80 +units=m +no_defs
source : mean_temp.tif
name : mean
min value : -9.1
max value : 15.4
```

To rescale the raster values to 0~1, I separate the raster values to three parts:

- For pixels with values <= 5th quantile, all of these pixels get values 0.
- Pixel values that are greater or equal to 95% will get 1.
- Pixel values between the 5th quantile and 95th quantile with be

normalized to 0~1.

I use `terra`

and `sf`

R packages for rescaling and make smaller tiles of the large raster.

Here is how I do it:

```
library(terra)
library(sf)
min_offset <- rast("test.tif")
quantiles_min <- spatSample(min_offset, 10^4, "regular")
q2 <- as.data.frame(quantile(quantiles_min, probs=c(0.05, 0.95), na.rm=TRUE))
print(q2)
#Read in data in sf format.
c_shape <- read_sf("test.shp")
# create an initial grid for centroid determination
c_grid <- st_make_grid(c_shape, cellsize = c(5e5, 5e5), square = TRUE) |>
st_as_sf()
c_grid_spat <- vect(c_grid)
## make tiles of large rasters:
filename <- paste0("/tile/tile", "_.tif")
tile <- makeTiles(min_offset,c_grid_spat, filename, overwrite = TRUE)
tile_rasters <- c(lapply(tile, rast))
Rescaling <- function(x){
##clamp to get the range of values that are need to be rescaled.
middle_min <- clamp(x,q2[1,1], q2[2,1], values = FALSE)
#use normalization to rescale the middle range values.
normalized_mid <- ((middle_min-q2[1,1])/(q2[2,1]- q2[1,1]))|> round(digits = 1)
## Values that are lower or equal to 5th quantile became 0
min_q5 <- clamp(x, upper = q2[1,1], value=FALSE)|> round(digits = 1)
min_repl5 <- init(min_q5, fun=0) |> mask(min_q5)
##values that are greater or equal to 95% get 1.
min_q95 <- clamp(x,q2[2,1],Inf, value = FALSE)|> round(digits = 1)
min_repl95 <- init(min_q95, fun=1) |> mask(min_q95)
##save raster
final_min<- merge(normalized_mid,min_repl5,min_repl95)
return(final_min)
}
x <- rast(tile_rasters[2]);Rescaling(x) #this works
MIs <- function(w){
r <- rast(w)
y <- Rescaling(r)
wrap(y,proxy=TRUE)
} #This function work. when > i=1; r <- tile_rasters[i] ;MI_tile <- MIs(r)
###Test function###
test <- MIs(tile_rasters[2])#this work with next line.
test <- terra::unwrap(test)#this work
```

When I tried to apply the two functions above in `foreach`

parallel, I encountered some errors:

```
library(parallelly)
library(doParallel)
library(foreach)
cl <- makeClusterPSOCK(length(tile_rasters), autoStop = TRUE)
registerDoParallel(cl)
test <- foreach(x =1:length(tile_rasters),
.packages = "terra") %dopar% {
resu <- MIs(tile_rasters[i])
return(resu)
}
stopCluster(cl)
```

Error message:

```
Error in { : task 1 failed - "NULL value passed as symbol address"
```

Does anyone has clues that how I can solve this problem to have my function parallelized on each raster tile?

For the test.shp file, you can find it from here, for the test.tif file, download from here.

To merge all the tile results into a final raster, I use `vrt`

function, but it seems to still take quite some memory space.
Here is how I did it:

```
t.lst <- list.files("/AllTiles/", pattern=".tif", full.names=TRUE)
MIs_final <- function(t.lst, fout="") {
r <- vrt(t.lst)
if (fout != "") {
writeRaster(r, fout, overwrite=TRUE)
fout
} else {
wrap(r)
}
}
```

Solution

You wrap the SpatRaster that is send back from a node, but you also need to wrap the input SpatRaster. But I would send the filenames instead. Otherwise you would probably run out of memory.

Something like this

Set up

```
library(terra)
f <- system.file("ex/elev.tif", package="terra")
r <- rast(f)
qnt <- global(r, quantile, probs=c(0.05, 0.95), na.rm=TRUE, maxcell=10^4) |> unlist()
qrng <- diff(qnt)
tiles <- makeTiles(r, 20, overwrite = TRUE)
```

A function that takes the file (this function should be much faster than yours; perhaps reducing the benefit of parallelization)

```
rescaler1 <- function(f){
x <- rast(f)
norm <- (x-qnt[1])/qrng
rcl <- rbind(c(-Inf, qnt[1], 0), c(qnt[2], Inf, 1))
cls <- classify(x, rcl, other=NA)
cover(cls, norm) |> round(1)
}
```

or easier to understand (but probably slower)

```
rescaler2 <- function(f) {
x <- rast(f)
x <- ifel(x < qnt[1], 0,
ifel(x > qnt[2], 1, (x-qnt[1]) / qrng))
round(x, 1)
}
```

And now you can do

```
a <- rescaler1(tiles[2])
b <- rescaler2(tiles[2])
```

And a function like this

```
MIs <- function(f, fout="") {
r <- rescaler1(f)
if (fout != "") {
writeRaster(r, fout, overwrite=TRUE)
fout
} else {
wrap(r)
}
}
```

And you need to export `qnt`

and `qrng`

to the nodes

Here is an alternative approach with multiple cores

```
fresc <- function(x){
x <- ifelse(x < qnt[1], 0,
ifelse(x > qnt[2], 1, (x-qnt[1]) / qrng))
round(x, 1)
}
r <- rast(f)
a <- app(r, fresc)
library(parallel)
cl <- makeCluster(4)
clusterExport(cl, c("qnt", "qrng"))
b <- app(r, fresc, cores=cl)
```

- Installing R on Linux: configure: error: libcurl >= 7.28.0 library and headers are required with support for https
- How to do ensembles with time series using AICc?
- planes3d expands and draws the area based on the sphere's radius
- How to extract tag code itself using R, rvest
- How to Display or Print Contents of Environment in R
- How to use Windows user credentials for proxy authentication in R/RStudio
- R reticulate specifying python executable to use
- Replace multiple Instances of a variable name in an R function and save the modified function
- Standardizing address formatting in R
- How to fix "failed to load cairo DLL" in R?
- Using grepl to filter columns names in specific range of columns
- changing the legends in ggplot2 to have groups of similar labels
- How to keep only unique rows but ignore a column?
- convert string date to R Date FAST for all dates
- Add subgroup text to plotly pie chart
- R Shiny : adjust height of DT datatable when fillContainer=TRUE,
- Why do R external pointers' "unusual copying semantics" mean they should not be used stand-alone?
- How to extract somo character after a string with a number of word which can change in R
- What does `se` stand for in geom_smooth(..., se = FALSE)?
- How to find number of rows greater than any values in R
- Align text and reduce space between text and parentheses in plotly hover info box
- Remove outer box of geom_bar plot with broken y-axis
- How to use lag/lead in mutate with an initial value?
- Is it possible to have a Shiny ConditionalPanel whose condition is a global variable?
- counting elements in one list in another list
- How to vectorize nested loops in R?
- Replace NA values with an incrementing sequence starting from the previous non-NA value
- How can I calculate the number of uniques in a row within a species matrix?
- How to perform operations on pairs of rows, based on a "distinguishing" column's values
- Mutate variable based on previous observations