0

I'm trying to randomly generate a point within a 10 kilometer radius of a point on the globe in R, and cannot figure out how to in a direct way. The goal is to have an arbitrary point on the globe with coordinates given in WGS84 coordinates, pick a random direction, then a random distance 0-10km and place a new point there. My current solution works, and relies on picking a random distance using runif, then using st_buffer to create a ring of that radius, subtracting the intersection of the ring with a slightly smaller ring to approximate a circle, and then using st_sample to generate a point on the ring. This works, but is incredibly slow, and I need to perform this operation millions of times in my program, so I need to find a more efficient method if possible. The overall goal is to replicate the DHS anonymization algorithm in R on a data set.

I have set up a working means of turning the generated random distances and angles into x and y coordinates in meters that can be added to a point's coordinates to displace it the correct way, but my sample is an sf data frame using WGS84, and the coordinates are therefore not in meters.here are some example coords

I'm not sure how to properly convert the points or translations without having issues involving the curvature of the Earth messing with the distances and angles.

  • 5
    Welcome to SO, PatrickGlenn! (1) Questions on SO (especially in R) do much better if they are reproducible and self-contained, see https://stackoverflow.com/q/5963269 , [mcve], and https://stackoverflow.com/tags/r/info. (2) Please don't post only images of data, see https://meta.stackoverflow.com/a/285557 (and https://xkcd.com/2116/). (3) Try `geosphere::destPoint(cbind(-74.00684, 40.70878), 130, 622)` finds a location 130 degrees (southeast), 622 meters away from StackExchange HQ (roughly). – r2evans Jul 31 '23 at 19:56
  • How accurate do you need to be? If you are willing to let the Earth be a perfect sphere, start with https://en.wikipedia.org/wiki/Spherical_circle , then look at https://en.wikipedia.org/wiki/Circular_sector to generate the 2-D coordinates you need to establish the Spherical CIrcle of interest. – Carl Witthoft Jul 31 '23 at 20:36

1 Answers1

3

We can use geosphere::destPoint to find the coordinates along some bearing (degrees) and range (meters) from a known point. For an example, the Stack Exchange office is advertised as being in NYC, so we can do

stackexchange <- cbind(-74.00684444589426, 40.708763941777285)
geosphere::destPoint(p=stackexchange, b=145, d=600)
#            lon      lat
# [1,] -74.00277 40.70434

enter image description here

This is vectorized as well. For instance,

  • p (the origin coordinates) can be a matrix of 1 or more rows, and if b and d are each length 1 then the output is the same dimensions as the p input, just offset;
  • b and d can be length greater than 1 (same lengths), so if p is just 1 row (one coordinate), then it will give you multiple points from that single origin p;
  • if all three p, b, and d are the same length (number rows of p, length of b and d), then you get independent bearings/ranges from the points.
r2evans
  • 141,215
  • 6
  • 77
  • 149