У меня есть шейп-файл с высоким разрешением от NOAA, который охватывает большую часть земного шара. Моя цель:
1 - обрезать этот набор данных в меньшую область (уже достигнута ОК);
2 - Наложить ограниченную область на базовый график R (все еще возникают проблемы).
Вот мой код:
library(rgeos)
library(rgdal)
download.file(url = 'https://dnc.nga.mil/DNCSiteContent/Shapefiles/shapefile01.zip', destfile = '~/Desktop/shapefile01.zip')
unzip('~/Desktop/shapefile01.zip', exdir = '~/Desktop/shapefile01')
closeAllConnections()
system('rm -rf ~/Desktop/shapefile01.zip')
targetfolder <- dir('~/Desktop', pattern = 'shapefile01', full.names = TRUE) # need full directory path for rgdal::readOGR
shape <- rgdal::readOGR(dsn = targetfolder, layer = 'cd01')
bbx <- rgeos::readWKT('POLYGON((-42.05 -23.03, -41.95 -23.03, -41.95 -22.93, -42.05 -22.93, -42.05 -23.03))')
proj4string(bbx) <- proj4string(shape)
mysite <- rgeos::gIntersection(shape, bbx)
plot(mysite)
В идеале я хочу раскрасить местность на карте так, чтобы она выглядела так (здесь отредактировано с помощью Photoshop):
Но mysite, объект класса SpatialLines с 28 отрезками, которые должны быть единым целым, но на самом деле не связаны между собой. Поэтому, когда я рисую весь объект, контур линии этой прибрежной зоны выглядит просто отлично. Но поскольку они представляют собой несколько меньших отрезков, R не знает, как их соединить, чтобы я мог заполнить местность. Есть ли способ сказать R, чтобы цвет заполнил континент? Может быть, объединить все эти сегменты в один непрерывный континентальный контур?
Кстати, в идеале я хотел бы сделать это в базе R, потому что это будет позже добавлено (т.е. наложено) на другой график с заливкой градиентом цвета для океана вокруг той же географической области (сгенерированной из другого пространственного файла) .
UPDATE
Используя пакет spatstat, я теперь могу очень легко визуализировать все отрезки, поэтому я подумал, что было бы неплохо обновить это здесь, чтобы помочь визуализировать проблему. Вот код:
library(spatstat)
mysite2 <- spatstat::as.psp(mysite)
plot(mysite2)
Имеет смысл? Я думаю, что было бы неплохо иметь постоянное решение, потому что такова природа этих высококачественных шейп-файлов береговой линии, предоставленных NOAA.