[R-sig-Geo] Randomly moving a locality (within set limits)

Barry Rowlingson b.rowlingson at lancaster.ac.uk
Fri Aug 22 18:59:24 CEST 2014


On Fri, Aug 22, 2014 at 4:58 PM, MacQueen, Don <macqueen1 at llnl.gov> wrote:
> My first thought would be to use spTransform() to convert to a (local)
> cartesian coordinate system, do the shift [perhaps using
> maptools::elide()], and then convert back.
>
> Identifying appropriate (local) coordinate systems for each locality could
> be a chore, but other than that it looks to me like a conceptually simple
> process.

This SO question:
http://stackoverflow.com/questions/9186496/determining-utm-zone-to-convert-from-longitude-latitude

has a simple formula for computing the UTM longitude zone:

long2UTM <- function(long) {
    (floor((long + 180)/6) %% 60) + 1
}

SFlong <- -122.4192
long2UTM(SFlong)
# [1] 10

So San Francisco (at that longitude) is in UTM 10 zone, which is epsg
code 32600 + 10 (for WGS84 spheroids). So the CRS string is
"+init=epsg:32610". You can construct this in R with a bit of
paste-fu.

For things in the southern hemisphere its 32700+zone number. See here:

http://spatialreference.org/ref/epsg/32610/
http://spatialreference.org/ref/epsg/32710/

That gets you a coordinate system in metres, don't forget to multiply
by 1000 to get km!

Here's a complete and slightly untested function (surely this exists
somewhere..):

 UTMzone = function(lat,long){
z = (floor((long + 180)/6) %% 60) + 1
 ns = ifelse(lat>=0,32600,32700)
 paste("+init=epsg:",ns+z,sep="")
 }

Here's the UTM zones it thinks are for a bunch of +10degree latitude
points across the globe:

> UTMzone(10,seq(-180,180,len=20))
 [1] "+init=epsg:32601" "+init=epsg:32604" "+init=epsg:32607" "+init=epsg:32610"
 [5] "+init=epsg:32613" "+init=epsg:32616" "+init=epsg:32619" "+init=epsg:32623"
 [9] "+init=epsg:32626" "+init=epsg:32629" "+init=epsg:32632" "+init=epsg:32635"
[13] "+init=epsg:32638" "+init=epsg:32642" "+init=epsg:32645" "+init=epsg:32648"
[17] "+init=epsg:32651" "+init=epsg:32654" "+init=epsg:32657" "+init=epsg:32601"

they look plausible, but may or may not be right, I've not checked...

Barry



More information about the R-sig-Geo mailing list