Search code examples
rgisrgdal

Extract shapefile value to point with R


I want to extract a shapefile value at a particular location. The shapefile I am using is can be found here, downloaded by clicking Download IHO Sea Areas. The shape file contains all the possible seas.

I can read it and plot it using:

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

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

However, I would like to extract the sea value for a particular location, say

p <- c(-20, 40)

Solution

  • there is a probably an easier way but here is my take

    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
    

    output for your point x=-20, y=40

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