为什么在带有 sf 列(POINT)的 df 上应用时 sf::st_transform returns 空几何?
Why sf::st_transform returns empty geometries when applied on a df with sf column (POINT)?
我正在尝试 st_transform 数据帧中的 sf 点(CRS:4326,基于度)到局部度量投影(EPSG:23867,DGN95 / UTM 区域 47N)。
所以我有一个数据框 "waypoints",其中包含我的所有属性以及经/纬度列。我使用 "lat" 和 "lon" 作为输入 (waypoints_sf) 在数据框中创建了一个 sf 点几何。 (另外,我在两个新列中将坐标保存为北距和东距)
waypoints_sf <- st_as_sf(waypoints,
coords = c("lat", "lon"),
crs = 4326)
coordinates_tmp <- st_coordinates(waypoints_sf)
colnames(coordinates_tmp) <- c("E","N")
waypoints_sf <- cbind(waypoints_sf,coordinates_tmp) # save coordinates as column E and N
class(waypoints_sf)
[1] "sf" "grouped_df" "tbl_df" "tbl" "data.frame"
class(waypoints_sf$geometry)
[1] "sfc_POINT" "sfc"
st_crs(waypoints_sf$geometry)
[1] Coordinate Reference System:
EPSG: 4326
proj4string: "+proj=longlat +datum=WGS84 +no_defs"
现在,当我尝试 st_transform 这个数据框时,它不会抛出错误或警告,但它的 sf 点是空几何 c(NaN, NaN) 不管我怎么努力。
waypoints_sf <- st_transform(waypoints_sf, 23867)
waypoints_sf$geometry[1]
[1] Geometry set for 1 feature (with 1 geometry empty)
geometry type: POINT
dimension: XY
bbox: xmin: NA ymin: NA xmax: NA ymax: NA
epsg (SRID): 23867
proj4string: +proj=utm +zone=47 +ellps=WGS84 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
POINT EMPTY
所以我做错了什么,也许我真的错过了一些简单的东西,我真的不知道为什么它不起作用?
如果它不适用于 st_transform 可能使用 N 和 E raw,那么简单的解决方法是什么?
目前对我有用的解决方法("waypoints" 是简单的 df):
waypoints_sf <- waypoints
waypoints_sf <- waypoints_sf %>% mutate(longitude=lon,latitude=lat) # save if later needed
coordinates(waypoints_sf) <- c("lon","lat")
proj4string(waypoints_sf) <- CRS("+proj=longlat +datum=WGS84") # set crs
waypoints_sf <- spTransform(test, CRS("+init=epsg:23867")) # transform
waypoints_sf <- st_as_sf(waypoints_sf,crs = 23867)
我怀疑您调换了纬度和经度,这导致在 EPSG:23867
的上下文中出现非法组合
考虑这段代码(注意我已经交换了你的经纬度坐标!)
library(sf)
library(dplyr)
library(leaflet)
waypoints_sf <- readr::read_csv("sample.csv") %>%
st_as_sf(coords = c("lon", "lat"), # note the order!
crs = 4326)
# this works...
waypoints_trans <- st_transform(waypoints_sf, crs = 23867)
# sanity check - a leaflet plot of the original sf data frame
leaflet(data = waypoints_sf) %>%
addTiles() %>%
addCircleMarkers(fillColor = "red",
stroke = FALSE)
作为完整性检查和最终确认,传单中的地图;这看起来可以吗?我是,然后你的坐标被翻转了。
我正在尝试 st_transform 数据帧中的 sf 点(CRS:4326,基于度)到局部度量投影(EPSG:23867,DGN95 / UTM 区域 47N)。
所以我有一个数据框 "waypoints",其中包含我的所有属性以及经/纬度列。我使用 "lat" 和 "lon" 作为输入 (waypoints_sf) 在数据框中创建了一个 sf 点几何。 (另外,我在两个新列中将坐标保存为北距和东距)
waypoints_sf <- st_as_sf(waypoints,
coords = c("lat", "lon"),
crs = 4326)
coordinates_tmp <- st_coordinates(waypoints_sf)
colnames(coordinates_tmp) <- c("E","N")
waypoints_sf <- cbind(waypoints_sf,coordinates_tmp) # save coordinates as column E and N
class(waypoints_sf)
[1] "sf" "grouped_df" "tbl_df" "tbl" "data.frame"
class(waypoints_sf$geometry)
[1] "sfc_POINT" "sfc"
st_crs(waypoints_sf$geometry)
[1] Coordinate Reference System:
EPSG: 4326
proj4string: "+proj=longlat +datum=WGS84 +no_defs"
现在,当我尝试 st_transform 这个数据框时,它不会抛出错误或警告,但它的 sf 点是空几何 c(NaN, NaN) 不管我怎么努力。
waypoints_sf <- st_transform(waypoints_sf, 23867)
waypoints_sf$geometry[1]
[1] Geometry set for 1 feature (with 1 geometry empty)
geometry type: POINT
dimension: XY
bbox: xmin: NA ymin: NA xmax: NA ymax: NA
epsg (SRID): 23867
proj4string: +proj=utm +zone=47 +ellps=WGS84 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
POINT EMPTY
所以我做错了什么,也许我真的错过了一些简单的东西,我真的不知道为什么它不起作用? 如果它不适用于 st_transform 可能使用 N 和 E raw,那么简单的解决方法是什么?
目前对我有用的解决方法("waypoints" 是简单的 df):
waypoints_sf <- waypoints
waypoints_sf <- waypoints_sf %>% mutate(longitude=lon,latitude=lat) # save if later needed
coordinates(waypoints_sf) <- c("lon","lat")
proj4string(waypoints_sf) <- CRS("+proj=longlat +datum=WGS84") # set crs
waypoints_sf <- spTransform(test, CRS("+init=epsg:23867")) # transform
waypoints_sf <- st_as_sf(waypoints_sf,crs = 23867)
我怀疑您调换了纬度和经度,这导致在 EPSG:23867
的上下文中出现非法组合考虑这段代码(注意我已经交换了你的经纬度坐标!)
library(sf)
library(dplyr)
library(leaflet)
waypoints_sf <- readr::read_csv("sample.csv") %>%
st_as_sf(coords = c("lon", "lat"), # note the order!
crs = 4326)
# this works...
waypoints_trans <- st_transform(waypoints_sf, crs = 23867)
# sanity check - a leaflet plot of the original sf data frame
leaflet(data = waypoints_sf) %>%
addTiles() %>%
addCircleMarkers(fillColor = "red",
stroke = FALSE)
作为完整性检查和最终确认,传单中的地图;这看起来可以吗?我是,然后你的坐标被翻转了。