Шейп-файл и векторный файл не точно перекрываются при построении

Я хочу построить шейп-файл поверх растрового файла в R, но я не могу сделать так, чтобы они идеально перекрывались: растр кажется повернутым на несколько градусов против часовой стрелки. Проблема с проекцией?

введите здесь описание изображения

Пожалуйста, рассмотрите следующие MWE

library(raster)
library(rgdal) 

# Download from http://biogeo.ucdavis.edu/data/gadm2/shp/ITA_adm.zip
shape_file = "ITA_adm1.shp"
# Download from http://sedac.ciesin.columbia.edu/data/set/gpw-v3-population-density/data-download
# Setting Geography: Country, Italy; Data Attributes: Grid
pop_density_file ="w001001.adf"

italy_map <- readOGR(dsn = shape_file, layer = "ITA_adm1")
italy_map_dens <- raster(pop_density_file)

colPal <- colorRampPalette(c("white", "red"))( 500 )

par(mar=c(0,0,0,0))
plot(italy_map_dens, xlim = c(6.70, 18.32), ylim = c(35.2, 47.6), axes=FALSE, box=FALSE, legend=FALSE, col=colPal)
plot(italy_map, xlim = c(6.70, 18.32), ylim = c(35.2, 47.6), border="grey", add=TRUE)

r raster rgdal
person CptNemo    schedule 21.11.2013    source источник
comment
Данные о плотности населения требуют регистрации и авторизации. Фу.   -  person Spacedman    schedule 21.11.2013
comment
Вы можете скачать без авторизации. Под данными для входа есть кнопка.   -  person CptNemo    schedule 22.11.2013

Ответы (2)


arrow_upward
2
arrow_downward

По-видимому, в пакете rgdal была ошибка. Моя проблема решилась после того, как я обновил его до версии 0.8-12.

person CptNemo    schedule 21.11.2013

arrow_upward
1
arrow_downward

Похоже проблема с проекцией. Вам нужно узнать точную проекцию обоих наборов данных и преобразовать один из наборов данных в проекцию другого набора данных.

person Paul Hiemstra    schedule 21.11.2013
comment
Как узнать их проекцию? Если я открою их в QGIS, они прекрасно перекрываются... - person CptNemo; 21.11.2013
comment
Может в свойствах QGis что-то сказано по поводу проекции? В противном случае вам необходимо получить эту информацию от поставщика данных. - person Paul Hiemstra; 21.11.2013
comment
Судя по свойствам слоя в QGIS, они, кажется, имеют один и тот же CRS EPSG: 4326 - WGS 84. Я также пытался экспортировать их, указав CRS, но все же не перекрывая. - person CptNemo; 21.11.2013
comment
Для меня они идеально перекрываются - обратите внимание, что их легче увидеть, если вы нарисуете italy_map_dens > 0, чтобы получить двухцветный растр. Я не уверен, почему вы рисуете с xlim и ylim в обоих случаях - нет смысла иметь xlim, ylim и add=TRUE. Убедитесь, что у вас есть последняя версия rgdal, а также сообщите нам, какие из растровых файлов вы использовали — в наборе данных их несколько. - person Spacedman; 21.11.2013