如何在R中将UTM坐标转换为经纬度

Question

如何在R中将UTM坐标转换为经纬度

rlatitude-longitudedata-conversionutm

13

我正在尝试运行一个物种分布模型，并需要创建背景点来运行我的逻辑回归模型。我刚刚创建了500个随机点，但它们是在UTM坐标系下的，我需要的是经纬度。是否有一种方法可以在R中将它们转换为经纬度？如果可以，请与我分享代码。我对R还比较新。谢谢！

- Caroline4

1个回答

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Robert Hijmans · Accepted Answer

如果你需要经纬度，最好使用该坐标参考系统生成随机点。但如果不需要，你可以按照以下方式操作。

首先，我会生成一个示例坐标集（points）。

 set.seed(1)
 x <- runif(10, -10000, 10000)   
 y <- runif(10, -10000, 10000)   
 points <- cbind(x, y)

现在，假设您知道您的点的坐标参考系统（CRS），您可以创建一个空间对象并分配已知的CRS。在我的示例中，这些点位于UTM，zone 10，datum=WGS84投影中。

library(terra)
v <- vect(points, crs="+proj=utm +zone=10 +datum=WGS84  +units=m")
v
# class       : SpatVector 
# geometry    : points 
# dimensions  : 10, 0  (geometries, attributes)
# extent      : -8764.275, 8893.505, -6468.865, 9838.122  (xmin, xmax, ymin, ymax)
# coord. ref. : +proj=utm +zone=10 +datum=WGS84 +units=m +no_defs

现在我们可以将这些转换为另一个坐标参考系统，例如经度/纬度。

y <- project(v, "+proj=longlat +datum=WGS84")
y
# class       : SpatVector 
# geometry    : points 
# dimensions  : 10, 0  (geometries, attributes)
# extent      : -127.5673, -127.4091, -0.05834327, 0.08873723  (xmin, xmax, ymin, ymax)
# coord. ref. : +proj=longlat +datum=WGS84 +no_defs

你可以这样提取坐标

lonlat <- geom(y)[, c("x", "y")]
head(lonlat, 3)
#             x             y
#[1,] -127.5308 -0.0530354276
#[2,] -127.5117 -0.0583432750
#[3,] -127.4757  0.0337371933

当然你也可以做相反的操作。

back <- project(y, "+proj=utm +zone=10 +datum=WGS84 +units=m")

同样的事情可以使用sf包或旧的sp包来完成。使用sp，创建一个SpatialPoints对象并使用spTransform函数。

library(sp)
sputm <- SpatialPoints(points, proj4string=CRS("+proj=utm +zone=10 +datum=WGS84")) 
spgeo <- spTransform(sputm, CRS("+proj=longlat +datum=WGS84"))
lnlt <- coordinates(spgeo)

我在示例中使用了UTM 10区。但请注意，有60个UTM区域，你需要选择一个。每个区域覆盖6度的地带（360/60=）。如果你的数据跨越很多经度或穿越UTM区域，你不应该使用UTM。对于经度在[-180, 180)之间的情况，你可以按照以下方式计算所需的区域。

utm_zone <- function(longitude) { trunc((180 + longitude) / 6) + 1 }

longs <- c(-122,-119, -118)
utm_zone(min(longs))
# [1] 10
utm_zone(max(longs))
# [1] 11

utm_zone(max(longs))

或者看看像这个这样的地图

UTM在南半球使用“虚假北移”时容易引起混淆，以避免使用负坐标。这是通过在y坐标上加上10000000来实现的，如下所示，使用+south元素。

s <- vect(cbind(174, -44), crs="+proj=longlat +datum=WGS84")
geom(project(s, "+proj=utm +zone=59"))[, c("x", "y")]
#       x          y 
#740526.3 -4876249.1 

geom(project(s, "+proj=utm +south +zone=59"))[, c("x", "y")]
#       x         y 
#740526.3 5123750.9

请注意，我使用“PROJ4”符号来定义坐标参考系统（CRS）。如果使用的基准面（基准面是地球表面形状的模型）是WGS84或NAD83，那么这种方法是有效的。如果不是这种情况，您需要使用“EPSG”代码或“WKT2”描述来定义您的坐标参考系统（CRS）。