我有一大堆包裹的多边形形状文件(50个县),我想循环遍历每个包裹并计算从每个包裹的中心到最近的铁路的距离,全国存储在线形文件中 .
我提出了两种不同的措施:来自geospheres包的point2Line和来自rgeos包的gDistance . 每个人都有不同的问题 .
对于point2line,我的代码看起来像这样:
rail_projection <-spTransform(rail_file,CRS(proj4string(parcels)) #project to same CRS
rail_crop<-crop(rail_projection,extent(parcels)) #crop to extent of parcels
centroids<-gCentroid(parcels,byid=T) # take centroid of parcel polygons.
m<-gDistance(rail_crop,centroids,byid=T)
但是这会发出以下警告:
Warning messages:
1: In RGEOSDistanceFunc(rail_crop, centroids, byid=T, "rgeos_distance") :
Spatial object 1 is not projected; GEOS expects planar coordinates
2: In RGEOSDistanceFunc(rail_crop, centroids, byid=T, "rgeos_distance") :
Spatial object 2 is not projected; GEOS expects planar coordinates
它还给了我一个nparcels X(我认为行数矩阵)的数字,我不太确定测量 . 我想用KM或里程或其他东西进行测量 .
使用dist2line的第二种方法如下所示:dist <-dist2Line(centroids,rail_crop)它工作正常,但唯一的问题是需要很长时间!要循环到60个左右的县,可能需要几天时间 . 另外,如果最近的铁路在另一个县,我宁愿不必裁剪铁路形状文件 . 使用整个shapefile运行代码需要更长的时间 .
所以我基本上寻找的是一种方法,以一种相当有效的方式从点到线形文件获得最近距离,以公里或英里的距离吐出距离,或类似的可解释的东西 . 如果这可以通过gDistance以没有这些错误的方式完成那么好!如果有办法从dist2Line或其他方法加速,那也很棒 . 如果您对如何使用缓冲区裁剪轨道shapefile有任何想法,我也很感激 . (由于包裹数据是专有的,我无法发布代码)
我仍然是一个有空间的初学者,如果我遇到麻烦或者之前已经回答了这个问题,我很抱歉 . 我环顾四周,找不到能够在这里工作的解决方案 .
谢谢!
使用数据编辑
所以我想我想出了我的警告信息来自哪里this . 我使用的是三维坐标系而不是平面坐标系 . 我运行了一些代码,其中包含我从互联网上为威斯康星州的一个县下载的shapefile,并且gDistance吐出一个矩阵,其中包裹的数量x从每个地块到每条铁路线的距离(我假设)的行数 . 有谁知道这些距离的测量结果是什么?如何将它们转换为KM或英里?我也应该使用"preferred"平面投影 . 此外,任何关于如何使用缓冲区裁剪shapefile的想法都将非常感激 . 文件/代码可以找到here
1 回答
我建议调查Universal Transverse Mercator system . 我不是最好的投影坐标系,但它是我在场上最常使用的坐标系(只要研究区域不是太大) .
找到特定于您所在地区的UTM区域后,您可以使用
sp::spTransform()
重新投影您的shapefile,以便gDistance()
可以以公里或米为单位返回您的距离 .这是覆盖肯尼亚的UTM区域的proj4string:
+proj=utm +zone=37 +south +ellps=WGS84 +datum=WGS84 +units=m +no_defs
. 如您所见,此处的单位为米(+units=m
) .不要太粗鲁,但似乎你对空间分析有些新意 . 您可能值得花时间研究相对较新的
sf
软件包,它简化并标准化了空间数据的功能和类 . 以下是您可以查看的简要介绍:http://strimas.com/r/tidy-sf/ . 该软件包中的类似函数是st_transform()
和st_distance()
.