Generating GPS Coordinates / creating a raster of evenly distributed points 1km apart

Question

i would like to cover an area defined by a bbox lat long coordinates with a raster of gps points 1 km apart. Currently i generate 2000 points for a bboxbbox=8.9771580802,47.2703623267,13.8350427083,50.5644529365 the following way:

as.data.frame(cbind(runif(2000,8.9771580802 ,13.8350427083),runif(2000,47.2703623267,50.5644529365)))

Since runif is a normal distribution, i think i just have to increase the amount of points to cover the whole area the way i need it. Is there a more clever way to do it? How many points would i need?

UPDATE

I thought i maybe can use the package sp to do the job but im still not realy familiar the with settings:

  longitudes <- c(8.9771580802, 13.8350427083)
  latitudes <- c(47.2703623267, 50.5644529365)
  bounding_box <- matrix(c(longitudes, latitudes), nrow = 2, byrow = TRUE, dimnames = list(NULL, c("min", "max")))
  projection <- "+proj=longlat" 
  sp_box<-Spatial(bbox = bounding_box, proj4string = CRS(projection))
  p_sample<-spsample(sp_box, 10000, type="regular")

If i understand correctly this will give me a number of points evenly distributed within my coordinates. spsample has an option for cell size but i dont grasp it yet. BR Andreas

runif will give a random distribution, ie not evenly spaced. do you need points on a 1km raster? — Waldi, Jul 22 '20 at 08:35
Could you specify, what you mean by "1km apart" and "raster"? Are we talking about a grid of squares, where the distance between one point and its direct neighbour is exactly 1km (which also implies that the distance to a point diagonal to the current point is sqrt(2)km)? Sorry, maybe I'm just not enough into spatial analysis to grasp your question correctly. — mabreitling, Jul 22 '20 at 09:01
Hi @mabreitling, thank you for you reply. Yes its exactly what i mean : "we are talking about a grid of squares, where the distance between one point and its direct neighbour is exactly 1km (which also implies that the distance to a point diagonal to the current point is sqrt(2)km)" — Andreas, Jul 22 '20 at 09:10
It looks as if the `destPoint` function in the `geosphere` [package](https://cran.r-project.org/web/packages/geosphere/index.html) should give you what you want: Function destPoint returns the location of point given a point of origin, and a distance and bearing. — Limey, Jul 22 '20 at 09:58

mabreitling · Accepted Answer · 2020-07-22T21:56:02.140

I'm not too much into spatial data and analysis but maybe it helps as a first step (I took different coordinates to get a reproducible example, which fits the dimensions of Germany where I have some feelings for the dimensions). I'm sure that there is a more elegant way but it should give you what you need. geosphere::destPoint() is ued to compute the points for a given distance and direction, geosphere::distGeo() computes the north-south/west-east distance of the given box to compute how many points we need for each direction. expand_grid() is then used to compute every combination for the computed border points.

Be also aware that I changed the distance between the points to 10,000 meters or 10 km to get fewer points and a nicer plot. You would have to change the numbers accordingly

nw <- c(5.8 55) 
se <- c(15.1, 47)

lon1 <- nw[1]
lat1 <- nw[2]
lon2 <- se[1]
lat2 <- se[2]


#(1) compute the border points in y direction, going south from the nw-point 
# while keeping lon constant

lat <- geosphere::destPoint(nw, 180, 1:floor(geosphere::distGeo(c(lon1,lat1), 
                                                                 c(lon1,lat2))/10000)*10000) 
lat <- as_tibble(lat) 

#(2) compute the border point in x direction (analog to above)
lon <- geosphere::destPoint(nw, 90, 1:floor(geosphere::distGeo(c(lon1,lat1), 
                                                                c(lon2,lat1))/10000)*10000)
lon <- as_tibble(lon)

# use expand_grid() to compute all combinations 
grid <- tidyr::expand_grid(lat$lat, lon$lon) 
names(grid) <- c("lat", "lon") #nicer names

### for visualizing what we've done, map germany with a grid overlay
germany  <- rnaturalearth::ne_countries(type =  "countries", 
            country = "germany", returnclass = "sf")

ggplot2::ggplot(data = germany)+
  ggplot2::geom_sf()+
  ggplot2::geom_point(data = grid, mapping = aes(x = lon, y = lat), size = 0.01)

Thank you very much mabreitling! One additional question: the initial coordinates nw,se: where are they come from? If I would like to scale this approach for each Bundesland for example, how would i go about it? — Andreas, Jul 22 '20 at 16:37
Well, those were just approximate values for the bounding box of Germany (which is not totally correct as you can see). nw is the combination the latitude value from the northernmost (somewhere at the border to Denmark) and the longitude of the westernmost point. For each Bundesland those values could be computed from the respective shapefile. — mabreitling, Jul 22 '20 at 16:40

Generating GPS Coordinates / creating a raster of evenly distributed points 1km apart

1 Answers1