rspatial / terra

R package for spatial data handling https://rspatial.github.io/terra/reference/terra-package.html
GNU General Public License v3.0
533 stars 86 forks source link

mean in mosaic breaks the raster #1450

Open JKupzig opened 6 months ago

JKupzig commented 6 months ago

I am trying to merge rasters with same crs and resolution but different extent to one bigger raster. For overlapping areas the function mean should be applied. When I apply mosaic to all my rasters the result is broken. (s. picture below). When I am using merge instead, it works.

Further I have encountered, that using a subset of rasters or doing the mosaic-merge in two steps works. Also, when I apply a different function than mean (e.g, first) it works as expected. Here is a minimal reproducable example

library(terra)

x1 <- rast(xmin=3323433, xmax=3668916, ymin=5774887, ymax=6119500, res=1000, vals=1, crs=terra::crs("EPSG:5677"))
x2 <- rast(xmin=3478201, xmax=3882596, ymin=5757906, ymax=6095527, res=1000, vals=2, crs=terra::crs("EPSG:5677"))
x3 <- rast(xmin=3278500, xmax=3625980, ymin=5410296, ymax=5828826, res=1000, vals=3, crs=terra::crs("EPSG:5677"))
x4 <- rast(xmin=3555086, xmax=3941507, ymin=5550139, ymax=5916728, res=1000, vals=4, crs=terra::crs("EPSG:5677"))
x5 <- rast(xmin=3459229, xmax=3875606, ymin=5228500, ymax=5661015, res=1000, vals=5, crs=terra::crs("EPSG:5677"))
x6 <- rast(xmin=3373358, xmax=3664922, ymin=5244482, ymax=5629050, res=1000, vals=6, crs=terra::crs("EPSG:5677"))
x  <- list(x1, x2, x3, x4, x5, x6)
x_  <- list(x1, x2, x3, x4)

# all together
collection <- terra::sprc(x)
m1 <- mosaic(collection)
plot(m1)

m1

# stepwise approach
collection <- terra::sprc(x_)
m2 <- mosaic(collection)
plot(m2)

m2

m3 <- mosaic(m2, x5, x6)
plot(m3)

m3

# all together merge
collection <- terra::sprc(x)
m4 <- merge(collection)
plot(m4)

m4

I am working on Ubuntu Yammy and this is my SessionInfo

R version 4.1.2 (2021-11-01)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 22.04.3 LTS

Matrix products: default
BLAS:   /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3
LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.20.so

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C               LC_TIME=de_DE.UTF-8        LC_COLLATE=en_US.UTF-8     LC_MONETARY=de_DE.UTF-8    LC_MESSAGES=en_US.UTF-8   
 [7] LC_PAPER=de_DE.UTF-8       LC_NAME=C                  LC_ADDRESS=C               LC_TELEPHONE=C             LC_MEASUREMENT=de_DE.UTF-8 LC_IDENTIFICATION=C       

attached base packages:
[1] parallel  stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] terra_1.7-71      dplyr_1.1.4       doParallel_1.0.17 iterators_1.0.14  foreach_1.5.2    

loaded via a namespace (and not attached):
 [1] Rcpp_1.0.12        pillar_1.9.0       compiler_4.1.2     class_7.3-20       tools_4.1.2        xts_0.13.2         gstat_2.1-1        lifecycle_1.0.4    tibble_3.2.1      
[10] lattice_0.20-45    pkgconfig_2.0.3    rlang_1.1.3        DBI_1.2.2          cli_3.6.2          e1071_1.7-14       generics_0.1.3     vctrs_0.6.5        hms_1.1.3         
[19] classInt_0.4-10    grid_4.1.2         tidyselect_1.2.0   spacetime_1.3-1    glue_1.7.0         sf_1.0-15          R6_2.5.1           fansi_1.0.6        sp_2.1-3          
[28] tzdb_0.4.0         readr_2.1.5        magrittr_2.0.3     intervals_0.15.4   codetools_0.2-18   units_0.8-5        utf8_1.2.4         KernSmooth_2.23-20 proxy_0.4-27      
[37] FNN_1.1.4          zoo_1.8-12 
chris-english commented 6 months ago

Correcting origins

?mosaic
The SpatRasters must have the same origin and spatial resolution.
origin(x1)
[1]  433 -113
origin(x2)
[1] 201 -94
origin(x3)
[1] 500 296
origin(x4)
[1]  86 139
origin(x5)
[1] 229 500
origin(x6)
[1] 358 482

origin(x1) <- c(0,0)
origin(x2) <- c(0,0)
origin(x3) <- c(0,0)
origin(x4) <- c(0,0)
origin(x5) <- c(0,0)
origin(x6) <- c(0,0)
x_lst = list(x1, x2, x3, x4, x5, x6)
collect_orig <- sprc(x_lst)
m1 = mosaic(collect_orig)
plot(m1)

Should resolve with additional data prep.

JKupzig commented 5 months ago

Hi, thanks for the suggestions. Sorry, for the late reply. I've just checked with my original data and unfortunately, my original rasters do have the same origin...

For now, I've found a workaround to get one raster out of my six by having all of my rasters in the same extent and resolution (and all cells outside the original rasters are filled with NA), so I can use simple raster calculation.

dimfalk commented 5 months ago

@rhijmans I'd expect this to be related with #1262, see also my answer on SO here.