Растр и многоугольник не выравниваются, даже если размеры и поляризации совпадают - PullRequest
0 голосов
/ 22 февраля 2019

Я пытаюсь наложить пространственный многоугольник (полученный из сетки, созданной icosa) на глобальную растровую карту.Я удостоверился, что ограничивающие рамки / экстенты объектов равны, и что они оба проецируются на один и тот же CRS, но полигон не совпадает с экстентом растра правильно.Есть идеи?(базовое изображение доступно здесь: https://imgur.com/a/c1oVm6y) misaligned raster and spatial polygon

# required packages
library(icosa)
library(raster)
library(rgdal)
library(rgeos)

# define projection string
equi <- CRS("+proj=longlat +datum=WGS84")

# load data and set extent
lad <- raster("Data/Maps/Ladinian_grey.jpg", crs = equi)
extent(lad) = c(-180, 180, -90, 90)

# generate grid and convert to spatial polygon
large_grid <- hexagrid(c(2, 3))
large_grid <- newsp(large_grid)
map_grid <- SpPolygons(large_grid, res = 50)
map_grid <- spTransform(map_grid, equi)

# plot
plot(lad)
plot(map_grid, col = NA, add = TRUE)

1 Ответ

0 голосов
/ 22 февраля 2019

Трудно увидеть, что происходит без данных, и странные вещи, которые вы делаете, изменяя экстенты объектов.В приведенном ниже примере показано, что он должен работать:

library(raster)
library(icosa)

r <- raster()
values(r) <- 1:ncell(r)
large_grid <- hexagrid(c(2, 3))
large_grid <- newsp(large_grid)
map_grid <- SpPolygons(large_grid, res = 50)
plot(r)
lines(map_grid)

enter image description here

Плохое выравнивание, которое вы показываете , может быть вызвано изменением размераГрафик после рисования.

Он также отлично работает с вашими данными и (упрощенным) кодом:

library(icosa)
library(raster)

lad <- raster("f91xQIt.jpg")
extent(lad) = c(-180, 180, -90, 90)
crs(lad) <- "+proj=longlat +datum=WGS84"

large_grid <- hexagrid(c(2, 3))
large_grid <- newsp(large_grid)
map_grid <- SpPolygons(large_grid, res=50)

# I use image so that it is clearer where all the lines come from
image(lad)
plot(map_grid, col = NA, add = TRUE)

enter image description here

И степенькак и ожидалось

map_grid
#class       : SpatialPolygons 
#features    : 362 
#extent      : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
#crs         : +proj=longlat +a=6371007 +b=6371007 
...