Proj4 reprojection using R
-
09-10-2019 - |
Question
I am trying to reproject coordinates from WGS84 to MGA Zone 53, a UTM projection based on the GDA94 datum. I get infinity as my result, which is definitely incorrect. I am using R's proj4
package like so:
> library(proj4)
> df <- data.frame("x" = c(131.1, 131.102, 131.1106, 133.34), "y" = c(-13.23, -13.243, -13.22, -22.66))
> df
x y
1 131.1000 -13.230
2 131.1020 -13.243
3 131.1106 -13.220
4 133.3400 -22.660
> ptransform(data = df, src.proj = "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs", dst.proj = "+proj=utm +zone=53 +south +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs")
$x
[1] Inf Inf Inf Inf
$y
[1] Inf Inf Inf Inf
$z
[1] 0 0 0 0
>
What is going wrong here?
Solution
The problem is that ptransform expects radians, not degrees. The function proj4:::project defaults to degrees. The results are the same with ptransform if you convert to radians.
OTHER TIPS
Where is the proj4 package obtained from?
Try rgdal if you can install it:
df <- data.frame("x" = c(131.1, 131.102, 131.1106, 133.34), "y" = c(-13.23, -13.243, -13.22, -22.66))
library(rgdal)
## project expects a matrix, assumes source is longlat/WGS84
project(as.matrix(df), "+proj=utm +zone=53 +south +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs")
[,1] [,2]
[1,] 77177.18 8534132
[2,] 77416.79 8532695
[3,] 78310.75 8535258
[4,] 329440.68 7493165