[R] Enumerating unit cell traversals Re: Interpolate? a line
Michael Bedward
michael.bedward at gmail.com
Wed Sep 15 14:46:02 CEST 2010
Nice posts David. I like the way that 'simple' problems such as this
one give rise to an interesting assortment of approaches and gotchas
:)
Michael
On 15 September 2010 22:41, David Winsemius <dwinsemius at comcast.net> wrote:
>
> On Sep 15, 2010, at 7:55 AM, David Winsemius wrote:
>
>>
>> On Sep 15, 2010, at 7:24 AM, David Winsemius wrote:
>>
>>> Replacing context:
>>>
>>>>> Hello everyone.
>>>>> I have created a 100*100 matrix in R.
>>>>> Let's now say that I have a line that starts from (2,3) point and ends
>>>>> to the
>>>>> (62,34) point. In other words this line starts at cell (2,3) and ends
>>>>> at cell
>>>>> (62,34).
>>>>>
>>>>> Is it possible to get by some R function all the matrix's cells that
>>>>> this line
>>>>> transverses?
>>>>>
>>>>> I would like to thank you for your feedback.
>>>>>
>>>>> Best Regards
>>>>> Alex
>>>
>>> On Sep 15, 2010, at 6:52 AM, Michael Bedward wrote:
>>>
>>>> Hello Alex,
>>>>
>>>> Here is one way to do it. It works but it's not pretty :)
>>>
>>> If you want an alternative, consider that produces the Y cell indices
>>> (since the x cell indices are already 2:62):
>>>
>>> > linefn <- function(x) 3+((34-3)/(62-2)) *(x-2)
>>> > findInterval(linefn(2:62), 3:34)
>>> [1] 1 1 2 2 3 3 4 4 5 5 6 6 7 7 8 8 9 9 10 10 11 11 12
>>> 12 13 13 14
>>> [28] 14 15 15 16 17 17 18 18 19 19 20 20 21 21 22 22 23 23 24 24 25 25 26
>>> 26 27 27 28
>>> [55] 28 29 29 30 30 31 32
>>> # that seems "off" by two
>>> > linefn(62)
>>> [1] 34
>>> > linefn(2)
>>> [1] 3 # but that checks out and I realized those were just indices for
>>> the 3:34 findInterval vector
>>>
>>> > (3:34)[findInterval(linefn(2:62), 3:34)]
>>> [1] 3 3 4 4 5 5 6 6 7 7 8 8 9 9 10 10 11 11 12 12 13 13 14
>>> 14 15 15 16
>>> [28] 16 17 17 18 19 19 20 20 21 21 22 22 23 23 24 24 25 25 26 26 27 27 28
>>> 28 29 29 30
>>> [55] 30 31 31 32 32 33 34
>>>
>>> ( no rounding and I think the logic is clearer.)
>>
>> But I also realized it didn't enumerate all the the cells were crossed
>> either, only indicating which cell was associated with an integer value of
>> x. Also would have even more serious problems if the slope were greater than
>> unity. To enumerate the cell indices that were crossed, try:
>>
>> unique( floor( cbind( seq(2,62, by=0.1), linefn(seq(2,62, by=0.1)) ) ) )
>> [,1] [,2]
>> [1,] 2 3
>> [2,] 3 3
>> [3,] 4 4
>> [4,] 5 4
>> [5,] 5 5
>> [6,] 6 5
>> [7,] 7 5
>> [8,] 7 6
>> snipping interior results
>> [83,] 58 32
>> [84,] 59 32
>> [85,] 60 32
>> [86,] 60 33
>> [87,] 61 33
>> [88,] 62 34
>>
>> That could probably be passed to rect() to illustrate (and check logic):
>>
>> rect(cellidxs[,1], cellidxs[,2], cellidxs[,1]+1, cellidxs[,2]+1,
>> col="red")
>>
>> #redraw line :
>> lines(2:62, 3+(34-3)/(62-2)*(0:60))
>
> And then I got to wondering if every 0.1 was sufficient to catch all the
> corners and discovered I could identify 3 more cell traversals with by=0.01
>
>> cellid2 <-unique( floor(cbind(seq(2,62, by=0.01), linefn(seq(2,62,
>> by=0.01) )) ) )
>> NROW(cellid2) # 91 cells
>> rect(cellid2[,1], cellid2[,2], cellid2[,1]+1, cellid2[,2]+1, col="blue")
>> rect(cellidxs[,1], cellidxs[,2], cellidxs[,1]+1, cellidxs[,2]+1,
>> col="red")
>> lines(2:62, 3+(34-3)/(62-2)*(0:60))
>
> (Trying by=0.001 did not change the count, but did take longer)
>
> --
> David.
>
>>> -
>>> David.
>>>
>>>>
>>>> interp <- approx(c(2, 62), c(3, 34), method="linear", xout=2:62)
>>>> m <- matrix(c(interp$x, round(interp$y)), ncol=2)
>>>> tie <- m[,2] == c(-Inf, m[-nrow(m),2])
>>>> m <- m[ !tie, ]
>>>>
>>>> You might want to examine the result like this...
>>>>
>>>> plot(m) # plots points
>>>> lines(c(2,26), c(3, 34)) # overlay line for comparison
>>>
>>> you can add a grid with
>>> abline(v=2:62, h=3:34)
>>>>
>>>> Michael
>>
>
>
> David Winsemius, MD
> West Hartford, CT
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
More information about the R-help
mailing list