Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Efficient way to calculate grid quadrants a line passes through

I have a 2-dimmensional unit grid, and a bunch of line segments that start and end at any rational number. I need an efficient way to calculate which grid cells the line passes through. For example, the line:

From (2.1, 3.9) to (3.8, 4.8) passes through the grid cells with lower left points (2, 3), (2, 4), and (3, 4).

Is there a quick, efficient way to calculate these quadrants from the line's endpoints?

I'll be working in R, but an answer in Python or pseudocode would work too. Thanks!

like image 280
Bonnie Scott Avatar asked Dec 11 '12 01:12

Bonnie Scott


1 Answers

Folks who work with spatial data deal with this kind of question all the time, so it may be worth piggy-backing on their efforts. Here's a solution that uses R's raster package (and functions from the sp package on which it depends):

library(raster)

## Create a SpatialLines object
a <- c(2.1, 3.9) 
b <- c(3.8, 4.8)
## Method #1 -- Uses functions from the sp package.
SL <- SpatialLines(list(Lines(list(Line(rbind(a,b))), "ab")))
## Method #2 -- Uses readWKT() from the rgeos package. Easier to read.
# library(rgeos)
# string <- paste0("LINESTRING(", paste(a, b, collapse=", "), ")")
# SL <- readWKT(string)

## Create a raster object
m <- 10
n <- 10
mat <- matrix(seq_len(m*n), nrow = m, ncol = n)
r <- raster(mat, xmn = 0, xmx = n, ymn = 0, ymx = m) 

## Find which cells are intersected & get coordinates of their lower-left corners
ii <- extract(r, SL, cellnumbers=TRUE)[[1]][, "cell"]
floor(xyFromCell(r, ii))
#      x y
# [1,] 2 4
# [2,] 3 4
# [3,] 2 3

## Confirm that this is correct with a plot
image(r)
plot(as(rasterize(SL, r), "SpatialPolygons"), 
     border = "darkgrey", lwd = 2, add = TRUE)
lines(SL)

enter image description here

like image 122
Josh O'Brien Avatar answered Oct 11 '22 15:10

Josh O'Brien