Skip to content

Match coordinates to regional polygon

2 messages · Miluji Sb, Roger Bivand

Original
Miluji Sb
# 
Dear all,

Re-sending because the message was too big due to the attachment,

I am trying to match a large number of coordinates (attached) sub-national
regions using GADM shapefile.

Coordinates:
https://drive.google.com/file/d/1PUsi4d0wP7hB6Aps6UmpXsnIPSD3I1sT/view?usp=sharing

Shapefile:
https://drive.google.com/drive/folders/1DxyygjQNeg2GIS9doaFM8yM2LvBYZi8E?usp=sharing


This is the code I am using;

###
library(data.table)
library(ncdf4)
library(rgdal)
library(raster)
library(tidyverse)

## Read coordinates
coords <- read.csv("coords.csv")

## Read shapefile
world <-readOGR(dsn = "/gadm36_levels_shp", layer = "gadm36_2")

## Convert to spatial points
pts <- coords
coordinates(pts) <- ~lon+lat

## Provide the same CRS to the coordinates as the shapefile
proj4string(pts) <- proj4string(world)

## Use over function to match coordinates to polygons
pts_district <- over(pts, world)$NAME_1
###

However, less than 20% of the coordinates are matched to a polygon. What am
I doing wrong? Any help will be greatly appreciated.

Regards,

Milu
Roger Bivand
#
On Tue, 10 Mar 2020, Miluji Sb wrote:

            

      
      
      
    

        
Apart from rgdal, the remainer are pure distraction. Parsimony is 
important.

            

      
      
      
    

        
This is the full GADM world map of admin units at level 2:

            

      
      
      
    

        
[1] 45962    13

            

      
      
      
    

        
This is a 0.5 degree world grid

            

      
      
      
    

        
min    max
lon -179.75 179.75
lat  -89.75  89.75

            

      
      
      
    

        
cellcentre.offset cellsize cells.dim
lon           -179.75      0.5       720
lat            -89.75      0.5       360

            

      
      
      
    

        
min    max
lon -179.75 179.75
lat  -89.75  89.75

            

      
      
      
    

        
plot(pts[!is.na(pts_district$GID_0),], axes=TRUE)

looks very much like world land mass to me. Oceans are big, your polygons 
omit Antarctica, and Plate Carree area (0.5x0.5 degree squares) is 
misleading.

Hope this clarifies,

Roger