[R] Enumerating unit cell traversals Re: Interpolate? a line

David Winsemius dwinsemius at comcast.net
Wed Sep 15 15:23:56 CEST 2010 

On Sep 15, 2010, at 8:58 AM, Rainer M Krug wrote:

> -----BEGIN PGP SIGNED MESSAGE-----
> Hash: SHA1
>
> On 15/09/10 14:46, Michael Bedward wrote:
>> 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?
>
> This is a known problem in computer graphics and the algorythm used to
> identify these cells is the "Bresenham's Line Drawing Algorithm"
>
> You can check out one version at
>
> http://www.falloutsoftware.com/tutorials/dd/dd4.htm
>
> or
>
> http://en.wikipedia.org/wiki/Bresenham's_line_algorithm

After viewing the Wikipedia description and the illustration I think  
that may not be the same problem. Note that the cells traversed only  
at their very outer corners do not get chosen. The goal is to  
"approximate a line" rather than to enumerate traversals.

-- 
David.


>
> Cheers,
>
> Rainer
>
>>>>>>>
>>>>>>> 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.
>>>
>>
>> ______________________________________________
>> 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.
>
>
> - --
> Rainer M. Krug, PhD (Conservation Ecology, SUN), MSc (Conservation
> Biology, UCT), Dipl. Phys. (Germany)
>
> Centre of Excellence for Invasion Biology
> Natural Sciences Building
> Office Suite 2039
> Stellenbosch University
> Main Campus, Merriman Avenue
> Stellenbosch
> South Africa
>
> Tel:        +33 - (0)9 53 10 27 44
> Cell:       +27 - (0)8 39 47 90 42
> Fax (SA):   +27 - (0)8 65 16 27 82
> Fax (D) :   +49 - (0)3 21 21 25 22 44
> Fax (FR):   +33 - (0)9 58 10 27 44
> email:      Rainer at krugs.de
>
> Skype:      RMkrug
> -----BEGIN PGP SIGNATURE-----
> Version: GnuPG v1.4.10 (GNU/Linux)
> Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org/
>
> iEYEARECAAYFAkyQwwUACgkQoYgNqgF2egrDSACeJPgbW5tQtqh7RyUMIezeNNL9
> 1+QAn2Vgow6lgeCuzZ2N2ysgxVTxmNm4
> =NFGN
> -----END PGP SIGNATURE-----

David Winsemius, MD
West Hartford, CT

More information about the R-help mailing list