[R] spTransform longlat to utm

[MacQueen, Don]

r-sig-geo is a better place to ask this question.

There doesn't appear to be anything wrong with how you are using
spTransform().


  "x" should be longitude
  "y" should be latitude

More precisely, the first column of your matrix
  matrix(c(x,y), ncol=2)
should be longitude, the second column latitude. So you had it right the
first time.

The most likely explanation is that you have your CRS string wrong. That
is, it would be something other than
  "+proj=merc +zone=32u +datum=WGS84"

When I search the proj.4 codes in QGIS, I don't find any with "32u". Try
UTM zone 32 (32N):


> sp2 <- spTransform(sp1,CRS(" +proj=utm +zone=32 +ellps=WGS84
>+datum=WGS84 +units=m +no_defs +towgs84=0,0,0"))
> sp2
SpatialPoints:
     coords.x1 coords.x2
[1,]  497734.1   5478009
Coordinate Reference System (CRS) arguments: +proj=utm +zone=32
+ellps=WGS84 +datum=WGS84 +units=m +no_defs +towgs84=0,0,0

I think this at least shows that you had your x,y correct the first time.

-Don


-- 
Don MacQueen

Lawrence Livermore National Laboratory
7000 East Ave., L-627
Livermore, CA 94550
925-423-1062





On 10/18/12 8:32 AM, "Alexander Belousov" <mynewsmailbox at gmx.de> wrote:

>Dear all,
>
>I am trying to project my LongLat-maps to a plane.
>The ultimate purpose is to do a search of points in vicinity of other
>points using overlay-commands (sp) with radius in km.
>
>I am applying spTransform (package rgdal) and it gives my some curious
>results.
>
>An example.
>Let's take a point lying somewhere in Germany, zone=32U
>
>x <- 8.968735
>y <- 49.454735
>
>After conversion I sould get something like
>
>Easting:  426858 (km)
>Northing: 5427937 (km)
>
>sp1 <- SpatialPoints(matrix(c(x,y), ncol=2), proj4string =
>CRS("+proj=longlat +datum=WGS84"))
>
>sp1Transformed <-  spTransform(sp1, CRS("+proj=merc +zone=32u
>+datum=WGS84"))
>coordinates(sp1Transformed)
>
>       coords.x1   coords.x2
>[1,] 998395.0133 6319888.068
>
>The result is an obvious nonsense.
>
>Well, after some deliberation I swapped the original coordinate columns:
>
>x <- 49.454735
>y <- 8.968735
>
>sp2 <- SpatialPoints(matrix(c(x,y), ncol=2), proj4string =
>CRS("+proj=longlat +datum=WGS84"))
>
>sp2Transformed <-  spTransform(sp2, CRS("+proj=merc +zone=32u
>+datum=WGS84"))
>coordinates(sp2Transformed)
>
>      coords.x1  coords.x2
>[1,] 5505275.918 995840.692
>
>Now the northing (which comes first after swapping) looks better, at
>least it falls within Germany, which is however not to say about the
>easting.
>
>Apparently I am missing an important point by the transformation.
>Any help will be greatly appreciated!
>
>Regards
>
>Alex
>
>______________________________________________
>R-help at r-project.org mailing list
>https://stat.ethz.ch/mailman/listinfo/r-help
>PLEASE do read the posting guide
>http://www.R-project.org/posting-guide.html
>and provide commented, minimal, self-contained, reproducible code.