Using R map data to determine associated state for a coordinate?

Werner Wernersen wrote:
Hi!

I have no idea if this is maybe an easy task utilizing
R since I read there is 
geographical map data in some package:

I have a huge number of geographical points with their
coordinates in Germany. 
Now I want to determine for each point in which
"Bundesland" = state it is located.

Can anybody tell me if this is doable relatively easy
in R and if so give me 
some hints or links how to do it?

Thanks a million,
   Werner
Hello Werner,

two building blocks, but don't know if the precision meets your needs.

1. Do you have a good map of Germany *with* Federal States?

* If YES and if it's free:
==> I would be interested! Please post it's source.

* If NO:
==> Downloadable map data are available on:
     http://www.vdstech.com/map_data.htm

2. The following approach reads and converts a shapefile with functions 
from maptools and then follows the example of inside.owin() from the 
spatstat package.
The example code will work very well, but since yesterday, when we
released a new version of maptools depending on the sp package, it can
look like:
library(maptools)
Loading required package: foreign
Loading required package: sp
nc <- readShapePoly(system.file("shapes/sids.shp", package="maptools")[1])
plot(nc, lwd=2, border="grey")
bbox(nc)
min       max
r1 -84.32385 -75.45698
r2  33.88199  36.58965
x <- runif(1000, min=-84.32385, max=-75.45698)
y <- runif(1000, min=33.88199, max=36.58965)
xypts <- SpatialPoints(cbind(x, y))
plot(xypts, add=TRUE, pch=19, cex=0.2)
where_am_i <- overlay(xypts, nc)
plot(xypts[is.na(where_am_i),], add=TRUE, pch=19, cex=0.2, col="grey80")  
summary(where_am_i)
Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's 
   1.00   29.00   51.00   53.11   77.00  100.00  450.00 

and in this case the points would be read into the SpatialPoints object 
directly. Should overlay have trouble with the "huge" number of points, 
you could take them in smaller batches, storing the intermediate results. 
As Thomas said, you need the boundaries of the "Bundesland" first, and the 
accuracy of your results will depend on the degree of detail of the 
boundary polygons.

Roger Bivand
Hope that helps

Thomas Petzoldt

##########################################################
library(maptools)
library(spatstat)

ger <- read.shape("germany.shp")
plot(ger)

pger <- Map2poly(ger)
sx<- pger[[13]]
lines(sx, type="l", col="red") # Saxony ;-)

## Create an owin (observation window) object
# direction of coordinates must be reversed, in some cases
# if error message: remove rev()'s
saxony <- owin(poly=list(x=rev(sx[,1]), y=rev(sx[,2])))

# random points in rectangle
x <- runif(1000, min= 6, max=15)
y <- runif(1000, min=46, max=56)

ok <- inside.owin(x, y, saxony)

points(x[ok], y[ok])
points(x[!ok], y[!ok], pch=".")
##########################################################

______________________________________________
R-help at stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Helleveien 30, N-5045 Bergen,
Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: Roger.Bivand at nhh.no

Using R map data to determine associated state for a coordinate?

Thread (3 messages)