Search code examples
rggplot2gisspatial

How to create a matrix of evenly-spaced points within an angled polygon, given the corner coordinates [R]


Given some example random data, with UTM coordinates for each corner:

       test<-structure(list(name = c("P11C1", "P11C2", "P11C3", "P11C4"), 
    east = c(6404807.016, 6404808.797, 6404786.695, 6404784.761
    ), north = c(497179.4834, 497159.1862, 497156.6599, 497176.4444
    ), plot_num = c(11, 11, 11, 11)), row.names = c(NA, -4L), class = c("tbl_df", 
"tbl", "data.frame"))

If we plot this as a polygon. we can see a tilted rectangle (this is because this shape is generated using real differential-GPS captured coordinates on the ground):

library(ggplot2)
ggplot(test) + geom_polygon(aes(east, north))

random rectangle

  • My question is, how can I generate points among custom dimensions that are evenly spaced within this polygon? For instance, if I want to generate a grid of evenly spaced 10x11 points within this grid. Can anyone suggest a neat to do this, given the corner points? I have hundreds of discrete plots for which I would then like to loop/map a solution over. I assume this involves some simple geometry, but with the added confusion of a tilted plot, I've gotten really confused and could not find a similar solution here on SO or elsewhere! FYI in this instance I am not expecting projection to be an issue since it is UTM coordinates, but a spatial solution that accounts for global projections would be cool to see, too!

Solution

  • You could use this little function:

    gridify <- function(x, y, grid_x = 10, grid_y = 10) {
      x <- sort(x)
      y <- sort(y)
      xvals <- do.call(rbind, Map(function(a, b) seq(b, a, length = grid_x),
          a = seq(x[1], x[3], length = grid_y), 
          b = seq(x[2], x[4], length = grid_y)))
      yvals <- do.call(rbind, Map(function(a, b) seq(a, b, length = grid_y),
                                  a = seq(y[1], y[3], length = grid_x), 
                                  b = seq(y[2], y[4], length = grid_x)))
      as.data.frame(cbind(x = c(xvals), y = c(t(yvals))))
    }
    

    So for example, to plot a 10 by 11 grid, we would do:

    ggplot(test) + 
      geom_polygon(aes(east, north)) +
      geom_point(data = gridify(x = test$east, y = test$north, grid_x = 11),
                 aes(x, y), color = 'red') +
      coord_equal()
    

    enter image description here

    And we can scale to arbitrary numbers of points:

    library(ggplot2)
    
    ggplot(test) + 
      geom_polygon(aes(east, north)) +
      geom_point(data = gridify(x = test$east, y = test$north, 50, 50),
                 aes(x, y), color = 'red') +
      coord_equal()
    

    enter image description here