# Problem points into a shape polygon

I want to determine which points from a set of hundreds of geographic coordinates are located into a shape polygon. I was able to convert all the coordinates (latitude and longitude ) into a SpatialPointsDataFrame and I can read and plot the shapefile in R as well. A graph of my shp file together with the points is:

From similar questions I was trying to use the function point.in.polygon in the sp package as well as the the over() function.

``````Polygon: polygon1,
Points: spdf

pt.in.poly <- sp::over( spdf, polygon1, fn = NULL)
``````

or:

``````new_shape <- point.in.poly(spdf, polygon1)
``````

but in both cases the result I get is:

``````Error: identicalCRS(x, y) is not TRUE
``````

The SpatialPointsDataFrame was created with the coordinates: CRS("+proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0") and the shapefile have the CRS arguments: +proj=cea +lon_0=0 +lat_ts=30 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0 , so they seem to be the same, I am not sure, Im not an expert on this.

If its a problem of coordinates system how can I fix it? or what am I doing wrong?

You have to transform your polygon (shapefile) Coordinate Reference System (CRS) with the function `spTransform` from the `sp` package. You can simply do this:

``````polygon.t <- spTransform(polygon, CRSobj = "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0")
``````

I also give you a complete reproducible example below:

# Get example data:

``````# Load libraries
library('rgdal')
library('raster')

# Get SpatialPolygonsDataFrame object example
polygon <- getData('GADM', country = 'URY', level = 0)

# Plot data
plot(polygon)

# Make sample points
sample.spdf <- spsample(polygon, n = 200, type = "random")
sample.spdf <- spsample(sample.spdf, n = 200, type = "random") # get points outside polygon
sample.spdf <- SpatialPointsDataFrame(sample.spdf, data = data.frame("ID" = 1:200, "V2" = sample(x = 1:1000, size = 200))) # add random data
``````

# Plot data:

``````# Plot points and polygon
plot(polygon)
plot(sample.spdf, add = TRUE, pch = 19, cex = 1)
box()
``````

``````# Interactive map with Mapview
library(mapview)
mapview(polygon) + mapview(sample.spdf)
``````

# Spatial Query: points in polygon with different CRS:

``````# Apply transformation to sample.spdf
sample.spdf.t <- spTransform(sample.spdf, CRSobj = "+proj=cea +lon_0=0 +lat_ts=30 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0")

projection(sample.spdf.t) # [1] "+proj=cea +lon_0=0 +lat_ts=30 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"
projection(polygon) # [1] "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"

pt.in.poly.sp <- sp::over(sample.spdf.t, polygon, fn = NULL)
# Error: identicalCRS(x, y) is not TRUE
pt.in.poly.rgeos <- rgeos::gIntersection(sample.spdf.t, polygon)
# spgeom1 and spgeom2 have different proj4 strings

# Plot points transformed + polygon
plot(polygon)
plot(sampleSP.t, add = TRUE, pch = 19, cex = 0.5)
``````

# Spatial Query: points in polygon with same CRS:

``````# Ask layers CRS
projection(sample.spdf) # [1] [1] "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"
projection(polygon) # [1] "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"

pt.in.poly.sp <- sp::over(sample.spdf, polygon, fn = NULL)
pt.in.poly.rgeos <- rgeos::gIntersection(sample.spdf, polygon)

# Plot points in polygon
plot(polygon)
plot(sample.spdf, pch = 19, cex = 1, add = TRUE)
plot(pt.in.poly.rgeos, pch = 19, cex = 0.5, col = 'green', add = TRUE)
box()
``````

``````# Interactive map with Mapview
mapview(polygon) + mapview(sample.spdf) + mapview(pt.in.poly.rgeos)
``````

• Thanks @Guzmán, it was a very clear and great answer, I was able to solve it! Commented Feb 16, 2017 at 14:46