r-projectRaster在使用自己的crs进行投影时会产生不同的输出



对于一个函数,我想使用用户设置的CRS重新投影光栅输入——它本身就是使用裁剪和掩码后的输出。我认为用现有的crs投影不会有任何作用,只会返回输入光栅。令我惊讶的是,事实并非如此。下面是一个可复制的示例

创建伪光栅:

library(raster)
# Download country polygin, in this case Malawi
mwi <- raster::getData("GADM", country = "MWI", level = 0)
# Create dummy raster
grid <- raster::raster() # 1 degree raster
grid <- raster::disaggregate(grid, fact = 12) # 5 arcmin
grid <- raster::crop(grid, mwi)
values(grid) <- rep(1, ncell(grid)) # Set a value
# The input raster with dimensions 94,39,3666
grid <- raster::mask(grid, mwi)
plot(grid)
grid
# Reproject the raster using its own crs. I use ngb as it is a categorical variable.
# This raster has dimensions 102, 47, 4794 so it seems a lot of white space (NA) is added.
own_crs <- crs(grid)
grid_reproj <- raster::projectRaster(grid, crs = own_crs, method = "ngb")
plot(grid_reproj)
grid_reproj
# To remove the white space I use trim
# This results in a waster with dimensions 93, 39, 3627
grid_trim <- raster::trim(grid_reproj)
plot(grid_trim)
grid_trim
# I also decided to compare the maps visually with mapview
library(mapview)
# There seems to be a trim function in mapview which I set to FALSE
# Also use the browser for easy viewing
options(viewer = NULL)
mapviewOptions(trim = FALSE)
mapviewGetOption("trim")
mapview(grid, col.regions = "green", legend = F) +
mapview(grid_reproj, col.regions = "red", legend = F) +
mapview(grid_trim, col.regions = "blue", legend = F) 

比较地图,我观察到两件事:

(1( 除了单个网格单元外,grid和gridtrim几乎相同顶部。也许这是由于四舍五入,

(2( grid_proj具有更大的尺寸和不同的程度。它还与其他地图相比,该地图似乎发生了轻微的偏移。这通过修剪进行校正,所以我认为这些实际上是NA细胞,差异可能是有关mapview如何显示地图。

因此,我的主要问题是,当rasterProject使用同样的程度?为什么修剪后的结果会有所不同?

要以右侧的方式投影光栅数据,必须提供另一个raster*对象,而不是唯一的crs。

如果您提供crs,则会计算出一个新的范围,该范围比避免丢失数据所需的范围大一点。当然,在这种crs实际上相同的奇怪情况下,返回输入不变是有意义的。

最新更新