Dear list, I have a large set of locations scattered over a land area covering about 200,000 km?. That area borders the sea, and I'd like to find the shortest distance from each point location to the sea/coast. Point locations are in lat/long. I've turned the locations into a SpatialPoints object, with a proj4string. And I have a SpatialPolygons object consisting of a single polygon, defining the sea areas (same proj4string as for the SpatialPoints object). My problem is I'm not clear on the next step. I tried using "gDistance" from the "rgeos" package: dists <- gDistance(locations, sea) but got two warnings: Warning messages: 1: In RGEOSDistanceFunc(spgeom1, spgeom2, byid, "rgeos_distance") : Spatial object 1 is not projected; GEOS expects planar coordinates 2: In RGEOSDistanceFunc(spgeom1, spgeom2, byid, "rgeos_distance") : Spatial object 2 is not projected; GEOS expects planar coordinates So should I set the proj4string to NA? And/or is gDistance even the right tool for this job? Maybe it's designed for instances where distances are small, not in the hundreds of kilometres? Any advice on how to do this appropriately would be much appreciated. - Malcolm
calculating shortest distance from a set of SpatialPoints to a polygon
5 messages · O'Hanlon, Simon J, Roger Bivand, Malcolm Fairbrother
Hi Malcom, I think you can try spatstat::nncross. You will need to turn your sp objects into spatstat objects. So coerce your locations into a point pattern of class 'ppp', which I believe is done easily enough with maptools::as.ppp and your polygons into a spatial lines pattern of class 'psp' using maptools::as.psp. Then it should be as simple as: dists <-nncross(locations.ppp,sea.psp) HTH. Simon -----Original Message----- From: r-sig-geo-bounces at r-project.org [mailto:r-sig-geo-bounces at r-project.org] On Behalf Of Malcolm Fairbrother Sent: 13 June 2012 13:46 To: r-sig-geo at r-project.org Subject: [R-sig-Geo] calculating shortest distance from a set of SpatialPoints to a polygon Dear list, I have a large set of locations scattered over a land area covering about 200,000 km?. That area borders the sea, and I'd like to find the shortest distance from each point location to the sea/coast. Point locations are in lat/long. I've turned the locations into a SpatialPoints object, with a proj4string. And I have a SpatialPolygons object consisting of a single polygon, defining the sea areas (same proj4string as for the SpatialPoints object). My problem is I'm not clear on the next step. I tried using "gDistance" from the "rgeos" package: dists <- gDistance(locations, sea) but got two warnings: Warning messages: 1: In RGEOSDistanceFunc(spgeom1, spgeom2, byid, "rgeos_distance") : Spatial object 1 is not projected; GEOS expects planar coordinates 2: In RGEOSDistanceFunc(spgeom1, spgeom2, byid, "rgeos_distance") : Spatial object 2 is not projected; GEOS expects planar coordinates So should I set the proj4string to NA? And/or is gDistance even the right tool for this job? Maybe it's designed for instances where distances are small, not in the hundreds of kilometres? Any advice on how to do this appropriately would be much appreciated. - Malcolm _______________________________________________ R-sig-Geo mailing list R-sig-Geo at r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Hi Malcom, I forgot one step (and as.ppp is from spatstat not maptools package). as.psp is expecting a SpatialLines object not a SpatialPolygons object. You can change your SpatialPolygons to SpatialLines easily enough, so try this: library(spatstat) sea.l <- as(sea, "SpatialLines") sea.psp <- as.psp(sea.l) locations.ppp <- as.ppp(locations) dists <- nncross(locations.ppp, sea.psp) HTH, Simon -----Original Message----- From: r-sig-geo-bounces at r-project.org [mailto:r-sig-geo-bounces at r-project.org] On Behalf Of Malcolm Fairbrother Sent: 13 June 2012 13:46 To: r-sig-geo at r-project.org Subject: [R-sig-Geo] calculating shortest distance from a set of SpatialPoints to a polygon Dear list, I have a large set of locations scattered over a land area covering about 200,000 km?. That area borders the sea, and I'd like to find the shortest distance from each point location to the sea/coast. Point locations are in lat/long. I've turned the locations into a SpatialPoints object, with a proj4string. And I have a SpatialPolygons object consisting of a single polygon, defining the sea areas (same proj4string as for the SpatialPoints object). My problem is I'm not clear on the next step. I tried using "gDistance" from the "rgeos" package: dists <- gDistance(locations, sea) but got two warnings: Warning messages: 1: In RGEOSDistanceFunc(spgeom1, spgeom2, byid, "rgeos_distance") : Spatial object 1 is not projected; GEOS expects planar coordinates 2: In RGEOSDistanceFunc(spgeom1, spgeom2, byid, "rgeos_distance") : Spatial object 2 is not projected; GEOS expects planar coordinates So should I set the proj4string to NA? And/or is gDistance even the right tool for this job? Maybe it's designed for instances where distances are small, not in the hundreds of kilometres? Any advice on how to do this appropriately would be much appreciated. - Malcolm _______________________________________________ R-sig-Geo mailing list R-sig-Geo at r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
On Wed, 13 Jun 2012, Malcolm Fairbrother wrote:
Dear list, I have a large set of locations scattered over a land area covering about 200,000 km?. That area borders the sea, and I'd like to find the shortest distance from each point location to the sea/coast. Point locations are in lat/long. I've turned the locations into a SpatialPoints object, with a proj4string. And I have a SpatialPolygons object consisting of a single polygon, defining the sea areas (same proj4string as for the SpatialPoints object). My problem is I'm not clear on the next step. I tried using "gDistance" from the "rgeos" package: dists <- gDistance(locations, sea) but got two warnings: Warning messages: 1: In RGEOSDistanceFunc(spgeom1, spgeom2, byid, "rgeos_distance") : Spatial object 1 is not projected; GEOS expects planar coordinates 2: In RGEOSDistanceFunc(spgeom1, spgeom2, byid, "rgeos_distance") : Spatial object 2 is not projected; GEOS expects planar coordinates
Wait a moment! The warnings are there for a purpose (and would also affect the use of functions in spatstat, which also assumes planar coordinates). Never delete important data (in this case the metadata of your points, saying which spatial reference system they possess). If your data are in geographical coordinates, you need to state that specifically, and look for suitable functions. One such is perhaps dist2Line() in geosphere. The alternative (and I advise doing both) is to project your points and coastine, and use gDistance(locations, sea, byid=TRUE), which will return distances in the metric of the projection (probably metres).
So should I set the proj4string to NA? And/or is gDistance even the right tool for this job? Maybe it's designed for instances where distances are small, not in the hundreds of kilometres?
Never set proj4string to NA if it is known (NA means not available, it doesn't mean ignore). The function works on planar geometries, and projection distortions for well-chosed projections at that scale are not large. Usually, wrongly defined datum metadata are more of a problem. Roger
Any advice on how to do this appropriately would be much appreciated. - Malcolm
_______________________________________________ R-sig-Geo mailing list R-sig-Geo at r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Roger Bivand Department of Economics, NHH Norwegian School of Economics, Helleveien 30, N-5045 Bergen, Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43 e-mail: Roger.Bivand at nhh.no
Dear Roger and Simon, Thanks very much for the suggestions. I got both the gDistance and nncross options to work (after learning rather more about projections and their use in R), and they're producing identical, sensible results. Interestingly, nncross is faster than gDistance, and both are way faster than dist2Line (which produces very similar though not identical results). I think my main problem was understanding projections, how to use proj4string, and cluing in that "+proj=utm" is projected whereas "+proj=longlat" isn't! Rather new material for me. Cheers, Malcolm
On 13 Jun 2012, at 14:22, Roger Bivand wrote:
On Wed, 13 Jun 2012, Malcolm Fairbrother wrote:
Dear list, I have a large set of locations scattered over a land area covering about 200,000 km?. That area borders the sea, and I'd like to find the shortest distance from each point location to the sea/coast. Point locations are in lat/long. I've turned the locations into a SpatialPoints object, with a proj4string. And I have a SpatialPolygons object consisting of a single polygon, defining the sea areas (same proj4string as for the SpatialPoints object). My problem is I'm not clear on the next step. I tried using "gDistance" from the "rgeos" package: dists <- gDistance(locations, sea) but got two warnings: Warning messages: 1: In RGEOSDistanceFunc(spgeom1, spgeom2, byid, "rgeos_distance") : Spatial object 1 is not projected; GEOS expects planar coordinates 2: In RGEOSDistanceFunc(spgeom1, spgeom2, byid, "rgeos_distance") : Spatial object 2 is not projected; GEOS expects planar coordinates
Wait a moment! The warnings are there for a purpose (and would also affect the use of functions in spatstat, which also assumes planar coordinates). Never delete important data (in this case the metadata of your points, saying which spatial reference system they possess). If your data are in geographical coordinates, you need to state that specifically, and look for suitable functions. One such is perhaps dist2Line() in geosphere. The alternative (and I advise doing both) is to project your points and coastine, and use gDistance(locations, sea, byid=TRUE), which will return distances in the metric of the projection (probably metres).
So should I set the proj4string to NA? And/or is gDistance even the right tool for this job? Maybe it's designed for instances where distances are small, not in the hundreds of kilometres?
Never set proj4string to NA if it is known (NA means not available, it doesn't mean ignore). The function works on planar geometries, and projection distortions for well-chosed projections at that scale are not large. Usually, wrongly defined datum metadata are more of a problem. Roger
Any advice on how to do this appropriately would be much appreciated. - Malcolm
_______________________________________________ R-sig-Geo mailing list R-sig-Geo at r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
-- Roger Bivand Department of Economics, NHH Norwegian School of Economics, Helleveien 30, N-5045 Bergen, Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43 e-mail: Roger.Bivand at nhh.no