Skip to content

find cells of a raster lying on a ring

1 message · Barry Rowlingson

#
On Fri, Jul 13, 2012 at 1:25 AM, Lorenzo Cattarino
<l.cattarino at uq.edu.au> wrote:
stop there!

 You might want to use the raster package, especially if you are using
geo-referenced data:

 > library(raster);library(sp)
 > z = raster(matrix(runif(10000),100,100))
 > projection(z)="+init=epsg:4400"
 > pt = cbind(.2,.6)
 > zd = distanceFromPoints(z,pt)
 > circ = zd>.4 & zd<.41
 > plot(circ)

 Note the circle is here computed as all grid cells within .4 and
.4+cell size distance. Then you can easily get the values in z that
correspond to circ==TRUE...

 > summary(z[circ])
    Min.  1st Qu.   Median     Mean  3rd Qu.     Max.
0.002594 0.230600 0.479600 0.475500 0.690800 0.985400

 Another way would be to rasterize a circular line:

 > cline = cbind(.2+.4*sin(seq(0,2*pi,len=100)),.6+.4*cos(seq(0,2*pi,len=100)))
 > L = SpatialLines(list(Lines(Line(cline),"C")))
 > plot(z)
 > plot(L,add=TRUE)
 > LR = rasterize(L,z)
 > plot(LR)
 > plot(L,add=TRUE)

 but this is a bit slower, and slightly different:

 > summary(z[LR])
    Min.  1st Qu.   Median     Mean  3rd Qu.     Max.
0.003472 0.239600 0.474300 0.476800 0.692900 0.997500

Barry