> extract does not sort the values according to the SpatialLine direction
but sorts it somehow on its extraction order. > I think it is a bug or at
least a not ideal behavior in extract() on spatial lines.

This is by design because it is simpler (faster) than extracting "along the
line".
I have added an option "along" in version 2.1-23 of raster (now on R-Forge)
to get the values in the expected order (along=TRUE).

Robert


On Fri, Apr 5, 2013 at 5:24 AM, Matteo Mattiuzzi <
matteo.mattiuzzi@boku.ac.at> wrote:

> Hi, my doubt is that extract does not sort the values according to the
> SpatialLine direction but sorts it somehow on its extraction order. I think
> it is a bug or at least a not ideal behavior in extract() on spatial lines.
> I think Robert can clarify this?
>
>
> library(raster)
> # raster example in extract:
>
>
>  r <- raster(ncol=36, nrow=18)
>  r[] <- 1:ncell(r)
>  cds1 <- rbind(c(-50,0), c(0,60), c(40,5), c(15,-45), c(-10,-25))
>  cds2 <- rbind(c(80,20), c(140,60), c(160,0), c(140,-55))
>  lines <- SpatialLines(list(Lines(list(Line(cds1)), "1"),
> Lines(list(Line(cds2)), "2") ))[1] # only 1 line!
>
>
> plot(r)
> lines(lines)
> x11()
> plot(extract(r, lines)[[1]]) # following the line the values should be
> sorted differently.
>
>
>
>
> #
> rasGround  <- raster("rasGroundCrop.tif", crs = "+proj=tmerc +lat_0=0
> +lon_0=15 +k=0.9999 +x_0=500000 +y_0=-5000000 +ellps=bessel
> +towgs84=653.0,-212.0,449.0 +units=m +no_defs")
> Road34 <- shapefile("Profil")
> projection(Road34) <- projection(rasGround)
>
>
>
> plot(rasGround)
> plot(Road34, add = TRUE, col = "black", lwd = 3)
> plot(extract(rasGround, Road34)[[1]], type = "l")
>
>
>
> Matteo
>
> >>> Omphalodes Verna  04/05/13 1:24 PM >>>
> Dear Matteo!
>
> Thanks
> for suggestion. I tried your code, but it seems that it is not working.
> Please
> see my code and in attached files are files.
>
> library(sp)
> library(raster)
> library(rgdal)
>
>
> rasGround  <-
> raster("rasGroundCrop.tif", crs = "+proj=tmerc +lat_0=0
> +lon_0=15 +k=0.9999 +x_0=500000 +y_0=-5000000 +ellps=bessel
> +towgs84=653.0,-212.0,449.0 +units=m +no_defs")
> Road34 <- readOGR(".",
> "Profil",
> p4s = "+proj=tmerc +lat_0=0 +lon_0=15 +k=0.9999
> +x_0=500000 +y_0=-5000000 +ellps=bessel +towgs84=653.0,-212.0,449.0
> +units=m
> +no_defs")
>
>
> plot(rasGround)
> plot(Road34, add = TRUE, col = "black",
> lwd = 3)
> plot(extract(rasGround, Road34)[[1]], type =
> "l") #Not works
>
> Maybe this could be a solution but
> I lose direction:
>
> pointsXY <- rasterToPoints(rasGround)
> rasX <- rasGround
> rasY <- rasGround
> rasX[] <- pointsXY[, 1]
> rasY[] <- pointsXY[, 2]
>
> rasXYZ <- stack(rasGround, rasX, rasY)
> segment <- extract(rasXYZ, Road34)[[1]]
> IDX <- order(segment[, 2], segment[, 3])
> segment <- segment[IDX, ]
> plot(segment[,1], type = "l")
>
> Thanks, OV
>
>
>
>
> ----- Original Message -----
> From: Matteo Mattiuzzi
> To: r-sig-geo@r-project.org; omphalodes.verna@yahoo.com
> Cc:
> Sent: Thursday, March 21, 2013 7:52 PM
> Subject: Re: [R-sig-Geo] vertical profile of raster DEM along road
> (package raster)
>
> Dear Omphalodes,
>
>
> library(raster)
> library(sp)
> set.seed(2)
> data(volcano)
> r <- raster(volcano)
> plot(r)
> l1 = cbind(c(seq(0,1 , by = 0.1)), runif(11, 0, 1))
> S1 = Lines(list(Line(l1)), ID = "a")
> lines(S1, asp = 1, col = "red")
>
>
> S1<-SpatialLines(list(S1)) # I'm not practic in sp, but this at least
> works!
> extract(r,S1)
>
>
>
> Matteo
>
> >>> Omphalodes Verna  03/21/13 6:06 PM >>>
> Dear list.
>
> I am asking for advice for getting verticale profile of DEM along line
> object. I konw, it is easy for one straight line (distance and extract
> functions) also with using a combination of GRASS and R (
> http://casoilresource.lawr.ucdavis.edu/drupal/node/375 ). Is in R
> possible to do this.
>
> Thanks to all. OV
>
> Here is a sample code:
>
> library(raster)
> library(sp)
> set.seed(2)
> data(volcano)
> r <- raster(volcano)
> plot(r)
> l1 = cbind(c(seq(0,1 , by = 0.1)), runif(11, 0, 1))
> S1 = Lines(list(Line(l1)), ID = "a")
> lines(S1, asp = 1, col = "red")
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo@r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
>

	[[alternative HTML version deleted]]

