0
votes

I am having trouble finding the distance between polygons and spatial points. I am trying to measure the distance between districts (polygons) and former tin mines (spatial points). I suspect I may have done something wrong with the projection. I seem to be able to calculate the distance between the districts, but not for districts and mines.

#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)

My output for the distance between polygons appear to be correct:

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

But my output for distance between the districts and mines are definitely wrong:

> 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

What I intend to do is to: 1) calculate the distance between all districts in Malaysia with all the former tin mines; and 2) extract the distance of the closest tin mine to each district.

Here's a link to the data I'm using. I'm using GADM data for the district polygons and I've hand coded the location of historical tin mines myself. Would appreciate all help given. Thank you.

1

1 Answers

0
votes

There are two issues with your current approach. 1.) The coordinates assignment isn't quite correct. It should be long/lat, but you're assigning it as lat/long. 2.) Directly setting the CRS in the current way will not actually alter the points in the necessary way. You need to assign an appropriate long/lat CRS first, and then perform a spTransform action.

#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