使用R将shapefile的值提取到点

Question

我想在特定位置提取一个shapefile值。我使用的shapefile可以在这里找到，通过点击Download IHO Sea Areas进行下载。该shapefile包含所有可能的海洋。

我可以使用以下代码读取并绘制它：

require("maptools")
require(rgdal)
require(sp)

ogrListLayers("World_Seas.shp")
shape <- readOGR("World_Seas.shp", layer="World_Seas") 

然而，我想提取特定位置的海洋值，比如说

p <- c(-20, 40)

- user3910073

---

你在那个点上期望什么样的“值”？高程吗？但根据定义，那应该是零。还有其他值吗？ - Felix

Answer 1

可能有更简单的方法，但这是我的做法

require("maptools")
require(rgdal)
require(sp)
library(plyr)
library(dplyr)

setwd("/Users/drisk/Downloads/seas")
ogrListLayers("World_Seas.shp")
shape=readOGR("World_Seas.shp", layer="World_Seas") 

datapol <- data.frame(shape)
pointtoplot <- data.frame(x=-20, y=40)
coordinates(pointtoplot) <- ~ x + y 
proj4string(pointtoplot) <- CRS("+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0")
#
#function over from package sp
test <- data.frame(xx=over(shape, pointtoplot))
combine <- cbind(test, datapol)
combine <- na.omit(combine) #only one point left

输入x=-20，y=40的结果如下：

   xx                 NAME ID Gazetteer_ id
35  1 North Atlantic Ocean 23       1912 35

- MLavoie

Answer 2

您可以使用sp包中的over函数：

library(rgdal)
library(sp)
library(raster)

shape <- shapefile("~/tmp/World_Seas.shp")
head(shape)

plot(shape[shape$ID == 35, ], axes = TRUE)
points(pts)

pts <- SpatialPoints(cbind(-20, 40), 
                     proj4string = CRS(proj4string(shape)))
over(pts, shape)

甚至更简洁：

pts %over% shape

- johannes

---

这是最好的答案，100%确定。 - colin