Я строю карту для северо-восточной части США. Фон карты должен быть либо картой высоты, либо картой средней годовой температуры. У меня есть два растра с Worldclim.org, которые дают мне эти переменные, но мне нужно обрезать их в соответствии с интересующими меня состояниями. Любые предложения о том, как это сделать. Это то, что я до сих пор:
#load libraries
library (sp)
library (rgdal)
library (raster)
library (maps)
library (mapproj)
#load data
state<- data (stateMapEnv)
elevation<-raster("alt.bil")
meantemp<-raster ("bio_1.asc")
#build the raw map
nestates<- c("maine", "vermont", "massachusetts", "new hampshire" ,"connecticut",
"rhode island","new york","pennsylvania", "new jersey",
"maryland", "delaware", "virginia", "west virginia")
map(database="state", regions = nestates, interior=T, lwd=2)
map.axes()
#add site localities
sites<-read.csv("sites.csv", header=T)
lat<-sites$Latitude
lon<-sites$Longitude
map(database="state", regions = nestates, interior=T, lwd=2)
points (x=lon, y=lat, pch=17, cex=1.5, col="black")
map.axes()
library(maps) #Add axes to main map
map.scale(x=-73,y=38, relwidth=0.15, metric=T, ratio=F)
#create an inset map
# Next, we create a new graphics space in the lower-right hand corner. The numbers are proportional distances within the graphics window (xmin,xmax,ymin,ymax) on a scale of 0 to 1.
# "plt" is the key parameter to adjust
par(plt = c(0.1, 0.53, 0.57, 0.90), new = TRUE)
# I think this is the key command from http://www.stat.auckland.ac.nz/~paul/RGraphics/examples-map.R
plot.window(xlim=c(-127, -66),ylim=c(23,53))
# fill the box with white
polygon(c(0,360,360,0),c(0,0,90,90),col="white")
# draw the map
map(database="state", interior=T, add=TRUE, fill=FALSE)
map(database="state", regions=nestates, interior=TRUE, add=TRUE, fill=TRUE, col="grey")
Объекты высот и средних значений - это те, которые должны быть обрезаны до экстента области вложенных объектов. Любой вклад поможет
Ответы:
Я бы отказался от использования
maps
пакета и нашел бы шейп-файл состояния. Затем загрузите это в R, используяrgdal
, и затем сделайте некоторую работу по наложению многоугольника.Как это? Шейп-файл gadm довольно подробный, вы можете найти более обобщенный.
источник
Вот подход с использованием
extract()
изraster
пакета. Я протестировал его с данными о высоте и средней температуре с веб-сайта WorldClim (я ограничиваю этот пример только высотой, температура работает аналогично), и соответствующий шейп-файл с границами штатов США находится здесь . Просто загрузите данные .zip и распакуйте их в свой рабочий каталог.Вам нужно загрузить
rgdal
иraster
библиотеки, чтобы продолжить.Давайте импортируем американский шейп-файл, используя сейчас
readOGR()
. После настройки CRS шейп-файла я создаю подмножество, содержащее желаемые состояния. Обратите внимание на использование заглавных и маленьких начальных букв!Затем импортируйте растровые данные, используя их,
raster()
и обрежьте их по экстенту ранее созданного подмножества состояний.В качестве последнего шага вам необходимо определить те пиксели растра высот, которые лежат в границах заданных полигонов состояния. Для этого используйте функцию «маска».
Вот очень простой график результатов:
Ура,
Флориан
источник