我正在尝试使用 R 将包含北向和东向数据的向量转换为十进制度。到目前为止,我已经能够使用 rgdal 包来执行投影,但我遇到了一个问题。UTM 区域的参数必须是单个字符串,但我的数据中有多个 UTM 区域。有没有办法使用 rgdal::SpatialPoints(( 函数来设置多个 UTM 区域?我尝试在 CRS(( 函数中输入设置字符串并使用向量,但我收到一条警告消息,指出仅使用了第一个区域。
df<- data.frame(X = c(774869, 771437, 1051883, 524468),
Y = c(414498, 403790, 184967, 779682),
WGS.1984.UTM.Zone = c("57N", "57N", "54N", "59N"))
library(rgdal)
sputm <- SpatialPoints(df[c("X", "Y")], proj4string=CRS("+proj=utm +zone=57N +datum=WGS84"))
spgeo <- spTransform(sputm, CRS("+proj=longlat +datum=WGS84"))
#resultant coordinates are all transformed relative to UTM zone 57N
df$projargs<- paste("+proj=utm +zone=", df$WGS.1984.UTM.Zone, " +datum=WGS84", sep = "")
#trying again with set projargs character strings
spUTM<- SpatialPoints(df[c("X", "Y")], proj4string=CRS(df$projargs))
警告消息:
1: In if (!is.na(projargs(( { : 条件的长度> 1,并且仅使用第一个元素
2: In if (!is.na(projargs(( { : 条件的长度> 1,并且仅使用第一个元素
3: In if (is.na(projargs(( uprojargs <- projargselse uprojargs <- paste(unique(unlist(strsplit(projargs, : 条件的长度> 1,并且仅使用第一个元素
这是sf
包的另一种解决方案。
使用st_as_sf()
创建 sf 对象。 请注意新的geometry
列。
library(sf)
library(dplyr)
df_sf <- st_as_sf(df,
coords = c("X", "Y"),
remove = FALSE,
crs = projargs)
df_sf
Simple feature collection with 4 features and 4 fields
geometry type: POINT
dimension: XY
bbox: xmin: 524468 ymin: 184967 xmax: 1051883 ymax: 779682
epsg (SRID): 32657
proj4string: +proj=utm +zone=57 +datum=WGS84 +units=m +no_defs
X Y WGS.1984.UTM.Zone projargs geometry
1 774869 414498 57N +proj=utm +zone=57N +datum=WGS84 POINT (774869 414498)
2 771437 403790 57N +proj=utm +zone=57N +datum=WGS84 POINT (771437 403790)
3 1051883 184967 54N +proj=utm +zone=54N +datum=WGS84 POINT (1051883 184967)
4 524468 779682 59N +proj=utm +zone=59N +datum=WGS84 POINT (524468 779682)
使用更新几何列的st_transform()
重新投影。 我们可以用st_coordinates()
将它们拉到新列中,它为每一行返回一个两列矩阵。
df_sf_longLat <- st_transform(df_sf, crs = "+proj=longlat +datum=WGS84") %>%
mutate(long_degrees = st_coordinates(.)[,1],
lat_degrees = st_coordinates(.)[,2])
df_sf_longLat
Simple feature collection with 4 features and 6 fields
geometry type: POINT
dimension: XY
bbox: xmin: 159.2216 ymin: 1.667156 xmax: 163.9555 ymax: 7.053618
epsg (SRID): 4326
proj4string: +proj=longlat +datum=WGS84 +no_defs
X Y WGS.1984.UTM.Zone projargs geometry long_degrees lat_degrees
1 774869 414498 57N +proj=utm +zone=57N +datum=WGS84 POINT (161.4747 3.746529) 161.4747 3.746529
2 771437 403790 57N +proj=utm +zone=57N +datum=WGS84 POINT (161.4435 3.64983) 161.4435 3.649830
3 1051883 184967 54N +proj=utm +zone=54N +datum=WGS84 POINT (163.9555 1.667156) 163.9555 1.667156
4 524468 779682 59N +proj=utm +zone=59N +datum=WGS84 POINT (159.2216 7.053618) 159.2216 7.053618
函数未矢量化。 所以我们需要循环
library(rgdal)
spUTM <- lapply(df$projargs, function(x)
SpatialPoints(df[c("X", "Y")], proj4string=CRS(x)))
str(spUTM)
#List of 4
# $ :Formal class 'SpatialPoints' [package "sp"] with 3 slots
# .. ..@ coords : num [1:4, 1:2] 774869 771437 1051883 524468 414498 ...
# .. .. ..- attr(*, "dimnames")=List of 2
# .. .. .. ..$ : NULL
# .. .. .. ..$ : chr [1:2] "X" "Y"
# .. ..@ bbox : num [1:2, 1:2] 524468 184967 1051883 779682
# .. .. ..- attr(*, "dimnames")=List of 2
# .. .. .. ..$ : chr [1:2] "X" "Y"
# ...
或者使用purrr
中的map
library(dplyr)
library(purrr)
df %>%
mutate(spUTM = map(projargs, ~
SpatialPoints(df[c('X', 'Y')], proj4string = CRS(.x))))