Skip to content
Prev 6452 / 29559 Next

Intersection of two coordinates (triangulation)

On Fri, Sep 11, 2009 at 7:43 PM, Breitbach, Nils <breitbach at uni-mainz.de> wrote:
How far have you got at the moment? What form is your data in?
Shapefiles? CSV files?

 The way to do coordinate transforms in R is with the spTransform
function from the rgdal package. Here I'll make some fake data and
show you how it's done:

 map=data.frame(x=runif(100,0,5),y=runif(100,20,22),s=1:100)
 plot(map$x,map$y)

this is just a plain R data frame. We need to make it spatial, so we
give it coordinates:

 coordinates(map)=~x+y
 plot(map)

 now we see a different view - the aspect ratio is now set to 1.
Currently the points don't have a spatial reference. We'll give them
one, by setting its 'proj4string':

  proj4string(map)=CRS("+init=epsg:4326")

 PROJ4 is a library for doing map transformations, and the sp package
uses it. The "epsg:4326" bit is the epsg code for WGS84 (see
www.epsg.org for more info on epsg codes).

 Now do:

 plot(map);axis(1);axis(2)

 to see the coordinates.

So how do we transform? spTransform to a new CRS!

Look up your EPSG code either from www.epsg.org or maybe here:

http://spatialreference.org/ref/epsg/32632/

 and transform....

 map2 = spTransform(map,CRS("+init=epsg:32632"))
 plot(map2);axis(1);axis(2)

 (Note I'm not sure if the numbers I chose are valid in that zone).

Gauss-Krueger zone 3 seems to be 2166:

http://spatialreference.org/ref/epsg/2166/

Hope this helps!

Barry