[R-sig-Geo] calculating shortest distance from a set of SpatialPoints to a polygon

[Malcolm Fairbrother]

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