裁剪光栅到最小程度不在R中工作

问题描述 投票:0回答:1

我正在编写一个脚本,它将采用任意三个栅格,并尽可能地将它们裁剪掉。所有三个栅格都具有相同的分辨率和投影。但是,裁剪到最小程度不会改变三个栅格的范围。我尝试过setExtent,同样的事情发生了。如果有人能提出建议,我会非常感激。这是示例代码:

 library(raster)

#Projection of all three rasters
newproj<- "+proj=utm +zone=4 +datum=WGS84 +units=m +no_defs +ellps=WGS84
+towgs84=0,0,0" 

#Create three rasters with varying extents
raster1p<- raster(crs = newproj)
extent(raster1p)<- c(531247, 691847, 7856684, 7987884) 
res(raster1p)<- c(100, 100)
values(raster1p)<- NA

raster2p<- raster(crs = newproj)
extent(raster2p)<- c(533550.8, 646550.8, 7881307, 7973807) 
res(raster2p)<- c(100, 100)
values(raster2p)<- NA

raster3p<- raster(crs = newproj)
extent(raster3p)<- c(525739, 689839, 7857305, 7996505) 
res(raster3p)<- c(100, 100)
values(raster3p)<- NA


#Find minimum extent
xmin1<- c(xmin(extent(raster1p)), xmin(extent(raster2p)), xmin(extent(raster3p)))
xmax1<- c(xmax(extent(raster1p)), xmax(extent(raster2p)), xmax(extent(raster3p)))
ymin1<- c(ymin(extent(raster1p)), ymin(extent(raster2p)), ymin(extent(raster3p)))
ymax1<- c(ymax(extent(raster1p)), ymax(extent(raster2p)), ymax(extent(raster3p)))
xmin_new<- min(xmin1)
xmax_new<- min(xmax1)
ymin_new<- min(ymin1)
ymax_new<- min(ymax1)
newextent=c(xmin_new, xmax_new, ymin_new, ymax_new)

#Crop rasters to minimum extent
crop(raster1p, newextent)
crop(raster2p, newextent)
crop(raster3p, newextent)

#Compare extents
extent_check<- c(extent(raster1p), extent(raster2p), extent(raster3p))

但是,当我查看extent_check以查看扩展区现在是否匹配时,我发现扩展区根本没有改变:

> extent_check
[[1]]
class       : Extent 
xmin        : 531247 
xmax        : 691847 
ymin        : 7856684 
ymax        : 7987884 

[[2]]
class       : Extent 
xmin        : 533550.8 
xmax        : 646550.8 
ymin        : 7881307 
ymax        : 7973807 

[[3]]
class       : Extent 
xmin        : 525739 
xmax        : 689839 
ymin        : 7857305 
ymax        : 7996505 

知道我可能做错了吗?感谢您的时间

r crop raster
1个回答
0
投票

我认为这不是一件坏事,而是一个误解(尽管你的代码中存在错误)。

示例数据

library(raster)  
prj <- "+proj=utm +zone=4 +datum=WGS84"
r1 <- raster(res=100, ext=extent(c(531247, 691847, 7856684, 7987884)), crs=prj, vals=NA)
r2 <- raster(res=100, ext=extent(c(533550.8, 646550.8, 7881307, 7973807)), crs=prj, vals=NA)
r3 <- raster(res=100, ext=extent(c(525739, 689839, 7857305, 7996505)), crs=prj, vals=NA)

找到“最小范围”

e <- intersect(intersect(extent(r1), extent(r2)), extent(r3))

请注意,结果与您的结果不同,因为您使用

xmin_new <- min(xmin1)ymin_new <- min(ymin1)

它应该在哪里

xmin_new <- max(xmin1)ymin_new <- max(ymin1)

现在裁剪

r1e <- crop(r1, e)
r2e <- crop(r2, e)
r3e <- crop(r3, e)

检查生成的范围

t(sapply(c(r1e, r2e, r3e), function(i) as.vector(extent(i))))
#         [,1]     [,2]    [,3]    [,4]
#[1,] 533547.0 646547.0 7881284 7973784
#[2,] 533550.8 646550.8 7881307 7973807
#[3,] 533539.0 646539.0 7881305 7973805

它们并不完全相同,因为这是不可能的,因为栅格不对齐。他们的“起源”是不同的

t(sapply(c(r1e, r2e, r3e), origin))
#      [,1] [,2]
#[1,]  47.0  -16
#[2,] -49.2    7
#[3,]  39.0    5

为了使它们对齐,你需要做这样的事情

r1e <- crop(r1, e)
r2e <- resample(r2, r1e)
r3e <- resample(r3, r1e)

t(sapply(c(r1e, r2e, r3e), function(i) as.vector(extent(i))))
#       [,1]   [,2]    [,3]    [,4]
#[1,] 533547 646547 7881284 7973784
#[2,] 533547 646547 7881284 7973784
#[3,] 533547 646547 7881284 7973784
© www.soinside.com 2019 - 2024. All rights reserved.