Обрезать растровый слой с помощью SpatialPolygonDataFrame

Question

У меня есть растровая сетка, которую я хочу обрезать в соответствии с границами суши на карте мира, предоставленной данными пакета 'maptools'. Сделав некоторые исследования, я обнаружил, что я должен использовать crop() функция, а затем mask() функция, но я получаю сообщение об ошибке.

Вот мой код:

# load the worldmap SpatialPolygonDataFrame
library(maptools)
data(wrld_simpl)
ll=CRS("+init=epsg:4326")
world<-spTransform(wrld_simpl, ll)
ext<-extent(-10.417,31.917,34.083,71.083)
# get only region covering europe 
world@bbox<-as.matrix(ext)


# create the origin in WGS84 CRS
ll = CRS("+init=epsg:4326") 
origin = SpatialPoints(cbind(7,40),ll) 
# create the raster grid 
grid = raster() 
# define extent and resolution
e = extent(5,18,40,50)
extent(grid)<-e
res(grid)=c(0.08,0.08) 
# fill the raster values
grid[] = runif(1:20250)

# crop the grid with worldmap to erase sea areas
test<-crop(grid,world)
rasterize(world,test)
Found 246 region(s) and 3768 polygon(s)
Error in if (x1a > rxmx) { : argument is of length zero

Answer 1

Урожай не обрезает до полигонов страны, как вы думаете. На самом деле, crop(grid, world) на самом деле ничего не обрезает, так как размер мира больше, чем размер сетки. Вы можете проверить это, построив график test, Это так же, как grid,

Вместо этого после grid[] = runif(1:ncell(grid)) Вы могли бы сделать:

# we're going to make 3 plots
par(mfrow=c(1,3))

# crop the Worldmap to grid 
worldcrop<-crop(world, grid)
plot(worldcrop)

# rasterize output, give cells value of NAME(seas are NA)
worldcropr = rasterize(worldcrop,grid, field='NAME', fun='first')
plot(worldcropr)

# mask random grid by worldcropr
gridmask = mask(x=grid, mask=worldcropr)
plot(gridmask)

# number of land pixels:
sum(!is.na(gridmask@data@values))

# number of sea pixels:
sum(is.na(gridmask@data@values))

Я верю, что это то, что вы после.