Преобразование шейп-файла в растровый

Question

У меня проблема с растеризацией шейп-файла для создания точек в сетке 0,5 * 0,5. Файл формы представляет собой классификацию уровней риска (Low-0, Medium-100, High-1000, Very High-1500) глобальных коралловых рифов для интегрированных угроз.

Я вытащил код из другого примера, который отлично работает, но когда я пытаюсь использовать его для своих данных, я ничего не получаю от функции построения графика. Ниже ссылки на шейп и мой код:

Reefs At Risk: Global Integreated Threats

# Read shapefile into R 
library(rgdal) 
library(raster)  

int.threat.2030 <- readOGR(dsn = "Global_Threats/Integrated_Future", 
          layer = "rf_int_2030_poly") 

## Set up a raster "template" for a 0.5 degree grid 
ext <- extent(-110, -50, 0, 35) 
gridsize <- 0.5 
r <- raster(ext, res=gridsize) 

## Rasterize the shapefile 
rr <- rasterize(int.threat.2030, r) 

## Plot raster 
plot(rr) 

Любые идеи, где я мог бы быть неправильно? Это проблема с самим шейп-файлом?

Просьба и спасибо!

Answer 1

Вы предположили, что многоугольники были в долготу/Lat координаты, но они не являются:

library(raster) 
library(rgdal) 
p <- shapefile('Global_Threats/Integrated_Future/rf_int_2030_poly.shp') 
p 

#class  : SpatialPolygonsDataFrame 
#features : 63628 
#extent  : -18663508, 14601492, -3365385, 3410115 (xmin, xmax, ymin, ymax) 
#coord. ref. : +proj=cea +lon_0=-160 +lat_ts=0 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0 
#variables : 3 
#names  : ID, THREAT, THREAT_TXT 
#min values :  1,  0, Critical 
#max values : 63628, 2000, Very High 

Вы можете изменить проекцию

pgeo <- spTransform(p, CRS('+proj=longlat +datum=WGS84')) 

, а затем сделать что-то вроде:

ext <- floor(extent(pgeo)) 
rr <- raster(ext, res=0.5) 
rr <- rasterize(pgeo, rr, field=1) 

Или сохранить оригинальную CRS и сделать что-то вроде:

ext <- extent(p) 
r <- raster(ext, res=50000) 
r <- rasterize(p, r, field=1) 
plot(r) 

Обратите внимание, что вы растеризуете очень маленькие полигоны в большие растровые ячейки. Полигон считается «внутри», если он покрывает центр ячейки (т. Е. Предполагая случай, когда многоугольники покрывают несколько ячеек). Поэтому для этих данных вам нужно будет использовать гораздо более высокое разрешение (а затем, возможно, суммировать результаты). В качестве альтернативы вы можете растрировать центроиды полигонов.

Но ни одно из вышесказанных не имеет отношения к действительности, так как вы делаете это все назад. Многоугольники, очевидно, получены из растра (посмотрите, насколько они являются блочными), а растровый доступен в наборе данных, на который вы указываете!

Таким образом, вместо растеризации, сделайте следующее:

x <- raster('Global_Threats/Integrated_Future/rf_int_2030') 
x 
#class  : RasterLayer 
#dimensions : 25456, 80150, 2040298400 (nrow, ncol, ncell) 
#resolution : 500, 500 (x, y) 
#extent  : -20037508, 20037492, -6363885, 6364115 (xmin, xmax, ymin, ymax) 
#coord. ref. : NA 
#data source : C:\temp\Global_Threats\Integrated_Future\rf_int_2030 
#names  : rf_int_2030 
#values  : 0, 2000 (min, max) 
#attributes : 
# ID COUNT THREAT_TXT 
# 0 80971  Low 
# 100 343535  Medium 
# 1000 322231  High 
# 1500 168518 Very High 
# 2000 83598 Critical 

Здесь черчения часть Палаван:

e <- extent(c(-8990636, -8929268, 1182946, 1256938)) 
plot(x, ext=e) 
plot(p, add=TRUE) 

Если вам нужна более низкое разрешение см raster::aggregate. Для другой системы координат координат см. raster::projectRaster.