我有一个带有几个多边形的形状文件,我需要在一个不同的文件或列表文件中的每个多边形中获得mercator坐标。这里的,包含来自我的shapefile中的每个多边形的范围的列表
a<-extent(659302.4, 659802.4, 9860325, 9860825)
b<- extent(647663, 648163, 10033987, 10034487)
c<-extent(609015.1, 609515.1, 9871028, 9871528)
ab<-c(a,b,c)我需要为列表或中的每个区段生成一个墨卡托坐标列表,这个文件包含每个文件的坐标。到目前为止,我唯一的选择是对每个区域进行栅格化,分配CRS坐标,然后转换到mercator,如下所示
bb<-raster(b)
crs(bb)<- "+proj=utm +zone=17 +south +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"
mbb<-projectRaster(bb, crs="+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0
+lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +no_defs")mbb级: RasterLayer尺寸: 12,12,144 (nrow,ncol,ncell)分辨率: 50,50.3 (x,y)范围:-8869218,-8868618,34170.81,34774.41 (xmin,xmax,ymin,ymax)。参考:+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +no_defs
这个方法唯一的问题是非常慢,特别是当我需要得到数百个多边形的墨卡托坐标时。
编辑:我想要的是将投影从UTM更改为Mercator,我考虑的函数允许为区段列表执行此操作,而不是像上面的示例所示那样单独执行。
发布于 2018-02-14 14:20:46
如果要更改范围的投影(而不是栅格),可以使用:
library(raster)
library(dplyr)
library(purrr)
a<-extent(659302.4, 659802.4, 9860325, 9860825)
b<- extent(647663, 648163, 10033987, 10034487)
c<-extent(609015.1, 609515.1, 9871028, 9871528)
ab<-c(a,b,c)
ab <- map(ab, function(x) {as(x,'SpatialPolygons')}) %>%
map(.,function(x){x@proj4string <- CRS("+proj=utm +zone=17 +south +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0");
spTransform(x,CRS("+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +no_defs"))})
# back to extent
new_extent <- ab %>% map(extent)https://stackoverflow.com/questions/48787500
复制相似问题