[R-sig-Geo] Step characteristics on raster maps

Tomislav Hengl hengl at spatial-analyst.net
Tue Oct 5 09:27:32 CEST 2010


Hi Basile,

If you work with large data, then you should definitively consider using 
SAGA GIS. Here are some examples:

--------------

# download gridded data:
 > 
download.file("http://geomorphometry.org/sites/default/files/volcano_maungawhau.zip", 
destfile=paste(getwd(), "volcano_maungawhau.zip", sep="/"))
 > unzip("volcano_maungawhau.zip")
 > volcano <- readGDAL("volcano_maungawhau.asc")
 > writeGDAL(volcano, "volcano_dem.sdat", "SAGA")

# create a shape file:
 > tmp <- 
list(Lines(list(Line(data.frame(x=c(2667700,2667700),y=c(6479520,6478720)))), 
"1"))
 > trans <- SpatialLinesDataFrame(SpatialLines(tmp), 
data.frame(ID=c("1")), match.ID=FALSE)
 > writeOGR(trans, "trans.shp", "trans", "ESRI Shapefile")

# overlay shapes over grids:
 > rsaga.get.usage("shapes_grid", 1)
SAGA CMD 2.0.4
library path:   C:/PROGRA~1/R/R-211~1.1/library/RSAGA/saga_vc/modules
library name:   shapes_grid
module name :   Get Grid Data for Shapes
Usage: 1 -SHAPES <str> -GRIDS <str> -RESULT <str> [-INTERPOL <num>]
   -SHAPES:<str>         Shapes
         Shapes (input)
   -GRIDS:<str>          Grids
         Grid list (input)
   -RESULT:<str>         Shapes (Grid Information)
         Shapes (output)
   -INTERPOL:<num>       Interpolation
         Choice
         Available Choices:
         [0] Nearest Neighbor
         [1] Bilinear Interpolation
         [2] Inverse Distance Interpolation
         [3] Bicubic Spline Interpolation
 > rsaga.geoprocessor("shapes_grid", 1, param=list(SHAPES="trans.shp", 
GRIDS="volcano_dem.sgrd", RESULT="trans_ov.shp", INTERPOL=0))
# get's an error message; looks like a bug :(

# convert segments to points:
 > rsaga.geoprocessor("shapes_points", 5, 
param=list(POINTS="trans_pt.shp", LINES="trans.shp", ADD=TRUE, DIST=20))
# now overlay points and lines:
 > trans.pt <- readOGR("trans_pt.shp", "trans_pt")
 > ov.trans <- overlay(volcano, trans.pt)
 > summary(ov.trans)
Object of class SpatialPointsDataFrame
Coordinates:
               min     max
coords.x1 2667700 2667700
coords.x2 6478740 6479520
Is projected: NA
proj4string : [NA]
Number of points: 41
Data attributes:
    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.
   109.0   134.0   158.0   150.7   166.0   195.0

--------------


T. Hengl
[http://www.wewur.wur.nl/popups/vcard.aspx?id=HENGL001]


Op 30-9-2010 20:57, Mathieu Basille schreef:
> Dear Robert,
>
> I just understood the interest of 'crosstab' with 'mask', this is pretty
> neat! Thanks for the suggestion.
>
> However, I can see some potential drawbacks with this approach: as my
> objective is to describe each step (each segment), I should first cut
> each 'SpatialLines' into a list of 'Lines', each 'Lines' being a segment
> (and not the whole trajectory). I apology if I'm totally wrong, I'm not
> that familiar with SpatialLines...
>
> But then, what about crossing/overlaping segments, i.e. several segments
> that all fall in the same pixel (which usually happens a lot)? With
> linesToRaster, only the last segment would be kept. So that I should run
> independently the approach for each 'Line' (each segment) just to avoid
> these crossings, and that would result in as many new rasters as the
> number of segments (and I'm talking about hundreds of thousands here,
> over large rasters).
>
> Basically, to rephrase my problem, here is what I want to achieve. Given
> a set of points (say a SpatialPointsDataFrame, or SpatialLinesDataFrame,
> or ltraj), I'd like to be able to compute new variables that give for
> each point summaries or characteristics (e.g. proportion of each type of
> the raster if it is categorical, or mean if is continuous) of the
> segment from that point to the next (or previous, it doesn't matter).
>
> It seems to be a simple problem in terms of low-level functions
> (segments over raster), but a complex one in terms of data structure (an
> single object with a set of individual trajectories which are themselves
> sets of segments). And I have to admit that I have some troubles going
> from the former to the latter.
>
> Thanks again for your time,
> Mathieu.
>
>
> Le 2010-09-30 11:28, Robert J. Hijmans a écrit :
>> Perhaps you can do something like
>>
>> r is a Raster* object
>> line is a SpatialLines* object
>>
>> library(raster)
>> rr<- linesToRaster(line, r)
>> rm<- mask(r, rr)
>> crosstab(rm, rr)
>>
>> Robert
>>
>> On Thu, Sep 30, 2010 at 5:46 AM, Mathieu Basille
>> <basille at ase-research.org>  wrote:
>>> Dear list members,
>>>
>>> I'm trying to compute characteristics along steps (i.e. segments between two
>>> points), based on underlying raster maps. The steps originally come from
>>> radiotracking data, converted to ltraj objects (adehabitat). The idea is to
>>> compute (for example) the habitat composition corresponding to each step
>>> instead of each points, as we are interesting in the movement path.
>>>
>>> I tried different solutions, as I would like to do it with R. I did not find
>>> any solution using adehabitat (or the development versions adehabitatMA/LS);
>>> 'join' only works for point characteristics, not steps. I also tried using
>>> sp and SpatialLinesDataFrame, but overlay does not seem to work with
>>> SpatialLines(DataFrame) and SpatialPixelsDataFrame:
>>>
>>> Error in function (classes, fdef, mtable)  :
>>>   unable to find an inherited method for function "overlay", for signature
>>> "SpatialPixelsDataFrame", "SpatialLinesDataFrame"
>>>
>>> I also investigated packages raster, trip, and rgeos, without success. Maybe
>>> the low level functions of rgeos could be used, but it seems a bit out of my
>>> skills (and time available) at the moment.
>>>
>>> Another solution might be to use spgrass6 in conjunction with GRASS, but I'm
>>> not familiar enough with GRASS to judge if it is a viable alternative...
>>>
>>> I'd welcome any hints/thoughts on this question.
>>> Sincerely,
>>> Mathieu Basille.
>>>
>>>
>>> --
>>>
>>> ~$ whoami
>>> Mathieu Basille, Post-Doc
>>>
>>> ~$ locate
>>> Laboratoire d'Écologie Comportementale et de Conservation de la Faune
>>> + Centre d'Étude de la Forêt
>>> Département de Biologie
>>> Université Laval, Québec
>>>
>>> ~$ info
>>> http://ase-research.org/basille
>>>
>>> ~$ fortune
>>> ``If you can't win by reason, go for volume.''
>>> Calvin, by Bill Watterson.
>>>
>>> _______________________________________________
>>> R-sig-Geo mailing list
>>> R-sig-Geo at stat.math.ethz.ch
>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>
>



More information about the R-sig-Geo mailing list