Skip to content

'Customizing' spsample

2 messages · Jean-Paul Kibambe Lubamba, Paul Hiemstra

Original
Jean-Paul Kibambe Lubamba
# 
Hello everybody,

I am using 'spsample' to derive points along a line and it works perfectly
by doing the following :

xx <- readShapeLines("line2",proj4string=CRS("+proj=utm +zone=34
+datum=WGS84"))
sppts <- spsample(xx, n=50, type="random")
plot(sppts)

However, I was wondering if the user can somehow 'customize' the sampling
process, as for example, by separating sampling points by a given
distance. And this distance being in the same units as the line. The 'n'
argument in spsample will therefore vary according with the length of each
line.

Any help is welcome!


Jean-Paul
2 days later
Paul Hiemstra
#
Jean-Paul Kibambe Lubamba wrote:
Hi Jean-Paul,

My suggestion would be to take the function sample.Line (from sp) as a 
starting point and modify it to your specifications. You can find the 
source in the /R directory (spsample.R) of the source package, which you 
can download from CRAN.

cheers,
Paul

sample.Line
function (x, n, type, offset = runif(1), proj4string = 
CRS(as.character(NA)),
    ...)
{
    offset = offset[1]
    if (missing(n))
        n <- as.integer(NA)
    if (!is.finite(n) || n < 1)
        return(NULL)
    cc = coordinates(x)
    lengths = LineLength(cc, longlat = FALSE, sum = FALSE)
    if (any(abs(lengths) < .Machine$double.eps)) {
        wl <- which(abs(lengths) < .Machine$double.eps)
        cc <- cc[-(wl), ]
        lengths <- lengths[-(wl)]
    }
    csl = c(0, cumsum(lengths))
    maxl = csl[length(csl)]
    if (type == "random")
        pts = runif(n) * maxl
    else if (type == "stratified")
        pts = ((1:n) - runif(n))/n * maxl
    else if (type == "regular")
        pts = ((1:n) - (1 - offset))/n * maxl
    else stop(paste("type", type, "not available for Line"))
    int = findInterval(pts, csl, all.inside = TRUE)
    where = (pts - csl[int])/diff(csl)[int]
    xy = cc[int, , drop = FALSE] + where * (cc[int + 1, , drop = FALSE] -
        cc[int, , drop = FALSE])
    if (nrow(xy) < 1)
        return(NULL)
    SpatialPoints(xy, proj4string)
}
<environment: namespace:sp>