重新投影 SpatialPointsDataFrame 不会改变范围答案

【问题标题】：Re-projecting a SpatialPointsDataFrame does not change the extent重新投影 SpatialPointsDataFrame 不会改变范围
【发布时间】：2018-03-24 20:10:25
【问题描述】：

我正在尝试将与蝙蝠位置 (SpatialPointsDataFrame) 相关的数据叠加到科罗拉多州 (SpatialPolygonsDataFrame)。两个对象的CRS不同：

crs(colorado)
#CRS arguments:
# +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0 
crs(BatRange_data)
#CRS arguments:
# +proj=utm +zone=13 +datum=NAD83 +units=m +no_defs +ellps=GRS80 +towgs84=0,0,0

当我重新投影蝙蝠数据时，CRS 确实发生了变化，但范围不受影响。

之前：

BatRange_data
#class       : SpatialPointsDataFrame 
#features    : 2456 
#extent      : 139996.3, 748812, 4094998, 4535103  (xmin, xmax, ymin, ymax)
#coord. ref. : +proj=utm +zone=13 +datum=NAD83 +units=m +no_defs +ellps=GRS80 +towgs84=0,0,0 
#variables   : 17

之后：

geo_proj = "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"
crs(BatRange_data) = geo_proj

BatRange_trnsfrmd = spTransform(BatRange_data,geo_proj)
BatRange_trnsfrmd
#class       : SpatialPointsDataFrame 
#features    : 2456 
#extent      : 139996.3, 748812, 4094998, 4535103  (xmin, xmax, ymin, ymax)
#coord. ref. : +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0 
#variables   : 17

我无法在多边形上绘制点，因为我的多边形对象的范围与我的点对象的范围不同：

colorado
#class       : SpatialPolygonsDataFrame 
#features    : 1 
#extent      : -109.0608, -102.042, 36.99223, 41.00561  (xmin, xmax, ymin, ymax)
#coord. ref. : +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0 
#variables   : 13

为什么 SpatialPointsDataFrame 的范围在重新投影和转换时没有变化？

可重现的示例（我自己对此进行了测试，并将其发送给了朋友；提供的内容可重现）：

#Libraries
x = c('raster','sp','rgdal')

# If you don't know if you've installed the packages for some of these 
# libraries, run this:
# install.packages(x)

lapply(x, library, character.only=TRUE)
rm(x)

# Read in and format data -------------------------------------------------
BatRange_data = new("SpatialPoints", coords = structure(c(179935.719907205, 179935.938979813, 179935.938979813, 176598.335967664, 176598.335967664, 4499963.43180688, 4499963.30060606, 4499963.30060606, 4489332.4211975, 4489332.4211975), .Dim = c(5L, 2L), .Dimnames = list(NULL, c("coords.x1", "coords.x2"))), bbox = structure(c(176598.335967664, 4489332.4211975, 179935.938979813, 4499963.43180688), .Dim = c(2L, 2L), .Dimnames = list(c("coords.x1", "coords.x2"), c("min", "max"))) , proj4string = new("CRS"    , projargs = NA_character_))

#Use state bounds from gadm website:
us = getData("GADM", country="USA", level=1)

#Extract states (need to uppercase everything)
co = "Colorado" 

colorado = us[match(toupper(co),toupper(us$NAME_1)),]

#Re-project bat data to 'longlat'
geo_proj = 
  "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"
crs(BatRange_data) = geo_proj

BatRange_trnsfrmd =
  spTransform(BatRange_data,geo_proj)

plot(BatRange_trnsfrmd)
plot(colorado,add=TRUE)

【问题讨论】：

这个例子对我们来说是不可重现的，因为它依赖于“Armstrong_Edited”；并且不必要地复杂地加载大约 30 个库，而三个就足够了。（或两个：library(raster); library(rgdal)）
"Armstrong_Edited" 用于制作BatRange_data，为此我在末尾添加了dput 输出。所以，这对你来说是可重现的。如果您不想加载所有库，请删除您不想要的库。加载需要几秒钟，所以它并不完全“复杂”。

标签： r r-raster sp rgdal

【解决方案1】：

在您的代码中：

geo_proj = "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"
crs(BatRange_data) = geo_proj

BatRange_trnsfrmd = spTransform(BatRange_data,geo_proj)
BatRange_trnsfrmd

在尝试重新投影之前，您将 BatRange_data 的原始投影覆盖为纬度/经度。 spTransform 因此什么也不做，因为它试图从 epsg:4326 重新投影到 epsg:4326。这就是为什么您的范围没有改变的原因。

因此，您可能应该删除这一行：

crs(BatRange_data) = geo_proj

一切都应该工作。

HTH。

【讨论】：