为什么sf::st_transform()返回的对象与调用中使用的投影不同?

evrscar2  于 2023-10-13  发布在  其他
关注(0)|答案(1)|浏览(158)

我想使用sf::st_transform()重新投影我的sf对象,但是转换后的对象的投影与我在转换调用中指定的投影不一样。为什么?为什么?

  1. library(sf)
  2. #> Linking to GEOS 3.6.1, GDAL 2.2.3, PROJ 4.9.3
  3. # target proj4string
  4. my_crs <- "+proj=lcc +lat_1=40.66666666666666 +lat_2=41.03333333333333 +lat_0=40.16666666666666 +lon_0=-74 +x_0=300000.0000000001 +y_0=0 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=us-ft +no_defs"
  5. # create test sfc
  6. p1 <- st_sfc(st_point(c(1,2)), crs = "+init=epsg:3857")
  7. # re-project p1 to `my_crs`
  8. p2 <- st_transform(p1, crs = my_crs)
  9. all.equal(my_crs, st_crs(p2)$proj4string)
  10. #> [1] "1 string mismatch"
  11. st_crs(p2)
  12. #> Coordinate Reference System:
  13. #> No EPSG code
  14. #> proj4string: "+proj=lcc +lat_1=40.66666666666666 +lat_2=41.03333333333333 +lat_0=40.16666666666666 +lon_0=-74 +x_0=300000 +y_0=0 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=us-ft +no_defs"

这两个投影之间唯一的区别是proj 4字符串中的+x_0元素。在p2的投影中,尾部的1 e-10已经被删除。

muk1a3rh

muk1a3rh1#

sf::st_transform()使用GDAL进行投影。帮助页面声明:
GDAL不支持某些PROJ.4预测,例如"+proj=wintri",因为它没有逆投影。投影到不支持的投影可以通过lwgeom包的一部分st_transform_proj完成。请注意,不支持的proj4string不能作为参数传递给st_crs,但必须作为字符串给出。
如果您在示例中使用PROJ.4进行投影,则它可以正常工作:

  1. library(sf)
  2. #> Linking to GEOS 3.6.1, GDAL 2.2.3, PROJ 4.9.3
  3. library(lwgeom)
  4. #> Linking to liblwgeom 2.5.0dev r16016, GEOS 3.6.1, PROJ 4.9.3
  5. my_crs <- "+proj=lcc +lat_1=40.66666666666666 +lat_2=41.03333333333333 +lat_0=40.16666666666666 +lon_0=-74 +x_0=300000.0000000001 +y_0=0 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=us-ft +no_defs"
  6. p1 <- st_sfc(st_point(c(1,2)), crs = "+init=epsg:3857")
  7. p2 <- lwgeom::st_transform_proj(p1, crs = my_crs)
  8. all.equal(my_crs, st_crs(p2)$proj4string)
  9. #> [1] TRUE

我对制图学的了解不多,不知道这种差异是由于帮助页面中提到的逆投影问题还是其他原因造成的。
相关信息:
https://stackoverflow.com/a/51663647/1707525
https://github.com/r-spatial/sf/issues/810
https://github.com/r-spatial/sf/issues/1019

展开查看全部

相关问题