计算从多边形到空间点的最小距离

问题描述 投票:0回答:1

我无法找到多边形和空间点之间的距离。我试图测量区域(多边形)和以前的锡矿(空间点)之间的距离。我怀疑我可能在投影方面做错了。我似乎能够计算出各区之间的距离,但不能计算区和矿的距离。

#Open file
mapping<- readOGR(dsn="C:/Users/noble/OneDrive - London School of Economics/LSE/EC465/Extended Essay/Stack Exchange Question/gadm2.Peninsula.dbf")


#Define the districts
Kinta <- mapping[mapping@data$NAME_1 == "Perak" &
             mapping@data$NAME_2 == "Kinta", ]
Gombak <- mapping[mapping@data$NAME_1 == "Selangor" &
             mapping@data$NAME_2 == "Gombak", ]

#Set a projection
EPSG.3375<-"+proj=omerc +lat_0=4 +lonc=102.25 +alpha=323.0257964666666 +k=0.99984 +x_0=804671 +y_0=0 +ellps=GRS80 +units=m +no_defs"
Gombak.km<-spTransform(Gombak,CRS(EPSG.3375))
Kinta.km<-spTransform(Kinta,CRS(EPSG.3375))

#Calculate the distance
gDistance(Kinta.km, Gombak.km, byid=TRUE)

#Open data on tin mines and define as spatial points
tin.01<-read.csv("C:/Users/noble/OneDrive - London School of Economics/LSE/EC465/Extended Essay/Stack Exchange Question/Tin Mine Location_01.csv")
coordinates(tin.01)<-8:9

proj4string(tin.01) <- CRS("+proj=omerc +lat_0=4 +lonc=102.25 +alpha=323.0257964666666 +k=0.99984 +x_0=804671 +y_0=0 +ellps=GRS80 +units=m +no_defs")

#Find distance between district and mines
gDistance(Kinta.km,tin.01,byid=TRUE)

多边形之间距离的输出似乎是正确的:

> gDistance(Kinta.km, Gombak.km, byid=TRUE)
   59
71 100676

但我对地区和矿山之间距离的输出肯定是错误的:

> gDistance(Kinta.km,tin.01,byid=TRUE)
    59
1 661153.5
2 661152.6
3 661153.0
4 661152.7
5 661151.8
6 661152.9
7 661153.1
8 661153.3
9 661153.2

我打算做的是:1)计算马来西亚所有地区与所有以前的锡矿之间的距离; 2)提取最近的锡矿到每个区的距离。

这是我正在使用的数据的link。我正在使用区域多边形的GADM数据,并且我亲自编写了历史锡矿的位置。非常感谢所有帮助。谢谢。

r gis spatial shapefile sp
1个回答
0
投票

您当前的方法存在两个问题。 1.)coordinates任务不太正确。它应该是long / lat,但是你将它指定为lat / long。 2.)以当前方式直接设置CRS实际上不会以必要的方式改变点。您需要首先分配适当的长/纬度CRS,然后执行spTransform操作。

#Open data on tin mines and define as spatial points
tin.01 <- read.csv("random/Tin Mine Location_01.csv")
coordinates(tin.01) <- c("Longitude","Latitude")  ## Should be `9:8` if you wanted to stick with indexes, but using the names here is generally lower risk.
proj4string(tin.01) <- CRS(proj4string(Kinta)) ## Setting the initial projection to the same one the polygons are using. You should change this if your original data source uses some other known long/lat projection.

tin.km <- spTransform(tin.01,CRS(EPSG.3375)) ## Creating a transformed set of points for the distance calculation.

#Find distance between district and mines
gDistance(Kinta.km,tin.km,byid=TRUE) ## '0' distance means the mine is inside the district.

          59
1 194384.372
2 223773.999
3      0.000
4  36649.914
5 102944.361
6      0.000
7      0.000
8   6246.066
9      0.000
© www.soinside.com 2019 - 2024. All rights reserved.