0
votes

I was asked whether R can work with shapefiles - I never worked with shapefiles myself before, but I am sure, others must have come across this kind of question!

I have two shapefiles: a) shapefile 1 (PolygonSamples.shp) contains a list of polygons which are distributed all over Germany (attached is a sample). The polygons might be smaller, equal or larger than the polygon of one postal codes polygon.

b) shapefile 2 lists the german postal codes and can be downloaded from https://blog.oraylis.de/2010/05/german-map-spatial-data-for-plz-postal-code-regions/

The question is now: How to 'match' the two shapefiles to get a dataframe that lists which polygon in shapefile 1 matches which postal codes(s) of shapefile 2. The result ideally looks like

 Polygon ID (shapefile 1)     Postal Code (shapefile 2)
         1                                80995
         2                                80997
         2                                80999
         3                                81247

Nothing of what I found matches really my question. For example From a shapefile with polygons/areas, and points (lat,lon), figure out which polygon/area each point belongs to? In R seems close, but I don't manage to get the desired dataframe (or datatable) output.

library(maps)
library(maptools)

# Polygons
tmp_dir <- "C:/Users/.../"
polygons <- readShapeSpatial(sprintf('%s/polygons.shp', tmp_dir)
plot(polygons)


# Postal codes
dir <- "C:/Users/..../"
postcode <- readShapeSpatial(sprintf('%s/post_pl.shp', dir)
plot(postcode)

The missing codes snipplet would read something like

    result_table <- match(polygons_ID, postcode, 
                      data1= polygon, data2 = postcode, 
                      by = "coordinates in the shapefile"

Sample of polygons in a shapefile (.shp) incl. other spatial files (.dbf,.prj, .qpj,.shx) can be send.

Any help is really VERY much appreciated!

PS: R version 3.2.3, 64 bit, RStudio on Windows 7

2
use function over from library(rgeos)G. Cocca

2 Answers

0
votes

Unfortunately I did not find an answer in R, but I could figure out how to match the two independent shapefiles in QGIS.

The main problem: The custom shapefile uses in the .prj file the geocoding Google Mercator (EPSG = 900913), while the downloaded postal code shapefile uses EPSG 4326.

QGIS does not automatically recognize these .prj files as projection files. One has to set them by hand.

Most importantly: Google Mercator (EPSG = 900913) was changed to EPSG= 3857. So for the custom shapefile I had to set – by hand! – the CRS to WGS 84/Pseudo-Mercator EPSG = 3857.

Now I could right click on the custom shape layer -> save as …. And Change the CRS to EPSG 4326. Thus the new custom shapefile now has the same projection like the downloaded postal code shapefile, and they can be joined by location.

(PS: Although I have a solution to do the conversion by hand, I would love to do this in R, because I need the resulting file for analysis.)

0
votes

Check out: https://gis.stackexchange.com/questions/140504/extracting-intersection-areas-in-r?newreg=033544fa0f5349bcb8167d78867c8073

It gives you which shapefiles in dataset B overlap with a shapefile in dataset A as well as how much area in each of B's shapefiles is present in the target shapefile.