使用 tidyverse 在 UTM 和十进制度之间进行转换-解网

问：

我正在使用 and 包将坐标从 UTM 转换为十进制度。我可以进行转换，但是，使用该函数后我丢失了一列。我需要保留这些列，以便可以将其连接到另一个数据帧;有没有办法在不使用该函数的情况下提取纬度/纬度？tidyversesfuniqueIDst_coordinates()uniqueIDst_coordinates()

例

library(tidyverse)
library(sf)

dat <- data.frame(
  uniqueID = seq(1,5,1),
  UTM_X = c(305334,NA,302685,300026,298030),
  UTM_Y = c(5320733,5320926,NA,5320882,5321002)
)

dat_sf <- dat %>%
  na.omit() %>%  # removes rows with missing values
  st_as_sf(coords = c('UTM_X','UTM_Y'), crs = 32616) %>% # let's R know the UTM is in zone 16
  st_transform(crs = 4326) %>% # converts UTM to decimal degree but is now stored as an sf_POINT
  st_coordinates() # extracts the lat/lon into two columns BUT drops the uniqueID

dat_sf
#>              X        Y
#> [1,] -89.61023 48.01022
#> [2,] -89.68139 48.00992
#> [3,] -89.70818 48.01037

^{创建于 2023-10-31 使用 reprex v2.0.2}

理想输出

dat_sf

  uniqueID      lat       lon
1        1 48.01022 -89.61023
2        4 48.00992 -89.68139
3        5 48.01037 -89.70818

R dplyr 坐标空间坐标变换

dat_sf <- dat %>%
  na.omit() %>% 
  st_as_sf(coords = c('UTM_X','UTM_Y'), crs = 32616) %>% 
  st_transform(crs = 4326) %>% 
  mutate(lat = st_coordinates(.)[, 2],
         lon = st_coordinates(.)[, 1]) %>%
  data.frame() %>%
  select(!geometry)

dat_sf
#>   uniqueID      lat       lon
#> 1        1 48.01022 -89.61023
#> 4        4 48.00992 -89.68139
#> 5        5 48.01037 -89.70818

^{创建于 2023-10-31 使用 reprex v2.0.2}

1赞 margusl 11/1/2023 #2

变化略有不同：

dat %>%
  na.omit() %>%
  st_as_sf(coords = c('UTM_X','UTM_Y'), crs = 32616) %>%
  st_transform(crs = 4326) %>%
  { bind_cols(st_drop_geometry(.), st_coordinates(.)) } %>% 
  rename(lon = X, lat = Y)
#>   uniqueID       lon      lat
#> 1        1 -89.61023 48.01022
#> 2        4 -89.68139 48.00992
#> 3        5 -89.70818 48.01037

1赞 Jonathan V. Solórzano 11/1/2023 #3

有几种方法可以达到所需的输出。这是另一个使用，和 .drop_naunnestst_drop_geometry

dat |>
  # Drop rows with NA in columns that start with UTM
  drop_na(starts_with("UTM")) |> 
  # Create sf object using UTM_X and UTM_Y as coordinates
  st_as_sf(coords = c("UTM_X","UTM_Y"), 
           crs = 32616) |> 
  # Transform to lat lon WGS84
  st_transform(crs = 4326) |>
  # Get coordinates as a data.frame
  mutate(coords = as.data.frame(st_coordinates(geometry))) |>
  # Unnest the data frame to get data at the same level as uniqueID
  unnest(coords) |>
  # Drop geometry column
  st_drop_geometry() |>
  # Rename coordinates columns
  rename("lon" = "X",
         "lat" = "Y")

# A tibble: 3 × 3
#  uniqueID   lon   lat
#*    <dbl> <dbl> <dbl>
#1        1 -89.6  48.0
#2        4 -89.7  48.0
#3        5 -89.7  48.0

1赞 Spacedman 11/2/2023 #4

你的问题是试图将所有这些强行放入一个管道中。没有理由不分阶段执行此操作。

例如，如果您这样做：

pts = st_transform(st_as_sf(na.omit(dat), coords=c("UTM_X","UTM_Y"), crs=32616), 4326)

然后你就得到了你需要的大部分位：

> head(pts)
Simple feature collection with 3 features and 1 field
Geometry type: POINT
Dimension:     XY
Bounding box:  xmin: -89.70818 ymin: 48.00992 xmax: -89.61023 ymax: 48.01037
Geodetic CRS:  WGS 84
  uniqueID                   geometry
1        1 POINT (-89.61023 48.01022)
4        4 POINT (-89.68139 48.00992)
5        5 POINT (-89.70818 48.01037)

然后，您可以获取 X-Y 坐标并将其重新组合在一起：

> pts = cbind(pts$uniqueID, st_coordinates(pts))
> head(pts)
               X        Y
[1,] 1 -89.61023 48.01022
[2,] 4 -89.68139 48.00992
[3,] 5 -89.70818 48.01037

然而，重要的是将所有这些放在一个函数中：

UTMto4326 = function(dat,X="UTM_X", Y="UTM_Y"){
 pts = st_transform(st_as_sf(na.omit(dat), coords=c(X,Y), crs=32616), 4326)
 cbind(st_drop_geometry(pts), st_coordinates(pts))
}

现在我想想，你可以让它变得更有用和通用：

transform_coords = function(dat, X, Y, crs_from, crs_to){
    pts = st_transform(st_as_sf(na.omit(dat), coords=c(X,Y), crs=crs_from), crs_to)
    cbind(st_drop_geometry(pts), st_coordinates(pts))
}

每当您想将两个坐标列从一个系统转换为另一个系统时，请使用：

> transform_coords(dat, "UTM_X","UTM_Y", 32616, 4326)
  uniqueID         X        Y
1        1 -89.61023 48.01022
4        4 -89.68139 48.00992
5        5 -89.70818 48.01037
> transform_coords(dat, "UTM_X","UTM_Y", 32616, 3875)
  uniqueID        X        Y
1        1 16803910 11960957
4        4 16806944 11967048
5        5 16808160 11969303

这还允许您检查反向转换是否使您回到开始的位置，减去具有两个调用的 NA 项：

> transform_coords(transform_coords(dat, "UTM_X","UTM_Y", 32616, 4326), "X","Y", 4326, 32616)
  uniqueID      X       Y
1        1 305334 5320733
4        4 300026 5320882
5        5 298030 5321002
> na.omit(dat)
  uniqueID  UTM_X   UTM_Y
1        1 305334 5320733
4        4 300026 5320882
5        5 298030 5321002

这也只使用基础 R 和包，这是此转换真正需要的全部内容，因此保持简单是最佳做法。sf

使用 tidyverse 在 UTM 和十进制度之间进行转换

Using tidyverse to convert between UTM and decimal degree

评论

评论