使用R进行Proj4投影重投影

Question

我正在尝试将WGS84坐标重新投影到MGA Zone 53，这是基于GDA94基准的UTM投影。但我的结果是无穷大，这绝对是不正确的。我使用R的proj4包进行操作：

> 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

> 

这里出了什么问题？

- fmark

---

你尝试过使用 silent=F 参数来查看是否有任何警告吗？ - VitoshKa

我已经看到了，它只是警告Z值，而这里我并不关心。 - fmark

Answer 1

问题在于ptransform需要弧度而不是角度。函数proj4:::project默认为角度。如果将角度转换为弧度，则使用ptransform的结果相同。

- mdsumner

---

如果你看了帮助页面，这个想法本来就会出现在你的脑海里。当情况紧急时，这是一个好主意；） - mdsumner

Answer 2

proj4包从哪里获取？

如果您可以安装rgdal，请尝试：

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函数期望一个矩阵，假定源数据是经纬度/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

这些是一组数字，可能是坐标或其他数据。

- mdsumner

---

啊，我看到CRAN上有proj4，抱歉，这里可能出了些问题 - 我会检查一下。 - mdsumner

proj4 在 Windows 的 64 位版本上不可用，这就是我最初认为它不在 CRAN 上的原因。 - mdsumner