Skip to content

C API for sp classes?

9 messages · Roger Bivand, wouter buytaert, Barry Rowlingson +2 more

Original
wouter buytaert
# 
Hi,

I want to write some R functions for calculating hydrological 
catchment characteristics from a DEM (topographic index, overland flow 
delay function etc.).

Looks like using a DEM with SpatialGridDataFrame class as an input is a 
good idea.

Then what is the best way to extract the data matrix from this class to 
pass it on to a .C() function which does the actual calculations?

cheers,

wouter
Roger Bivand
#
On Tue, 7 Mar 2006, Edzer J. Pebesma wrote:

            

      
      
      
    

        
Is the overhead of learning and using the .Call() interface too great? As 
Edzer says, since the AttributeList can contain multiple vectors of 
possibly different types, it may be helpful to get a better grip on them 
that way. The geometry slots are more complicated, and would need more 
thought, though, but would be needed to make a raster, as Edzer shows.

Interesting issue. What would be returned? I guess in fact that a .C() 
interface with a more complicated R function to pick everything apart and 
check types is the way to go for prototyping. What would you use as 
baseline for comparison to see if your output is plausible?

Roger

            

      
      
      
    

        

      
      
    

            

      
      
      
    

        

      
      
    

        

  
    

      
      
    

  


  


      

    

    

      
      

        


  

        

      Edzer Pebesma
    

    

      
      #
    

  


  

      
Wouter, if you want to use the .C interface, you need to pass
the maps as simple vectors. I would pass one with the map
values, and one with the topology. Below is some example
code how to get them. An R function should be used to do
the coercion to double arrays, and wrap the .C call.

Please note that a SpatialGridDataFrame can actually hold
more than one attribute (map).
--
Edzer

For getting the values:

 > library(sp)
 > data(meuse.grid)
 > gridded(meuse.grid)=~x+y
 > fullgrid(meuse.grid)=T
 > meuse.grid at data@att[[1]] # first attribute
 [1] NA NA NA NA NA NA NA NA NA NA ...
 >

Getting the topology:

 > meuse.grid at grid
                       x      y
cellcentre.offset 178460 329620
cellsize              40     40
cells.dim             78    104
 > class(meuse.grid at grid)
[1] "GridTopology"
attr(,"package")
[1] "sp"
 > unlist(gridtopology(meuse.grid))
Error in unlist(gridtopology(meuse.grid)) :
        couldn't find function "gridtopology"
 > unlist(gridparameters(meuse.grid))
cellcentre.offset1 cellcentre.offset2          cellsize1          cellsize2
            178460             329620                 40                 40
        cells.dim1         cells.dim2
                78                104
 > as.double(unlist(gridparameters(meuse.grid)))
[1] 178460 329620     40     40     78    104
 >

        
Wouter Buytaert wrote:

        

            

      
      
      
    

        

      
      
    

  
  


  


      

    

    

      
      

        


  

        

      wouter buytaert
    

    

      
      #
    

  


  

      
Thanks for the info!

Learning the .Call interface is certainly no problem of course. I just 
thought about going with .C() because:

- the "writing R extensions" manual says you should consider .C before 
using .Call
- it is easy to pass the map values as a double vector, and only a few 
additional parameters have to be passed on (resolution, number of 
rows/cols etc...), so using .C() is quite straightforward.
- I prefer doing a datacheck in R before passing the data to .C()

The resulting R function doesn't seem very complicated. I'm a newbie 
though, so all recommendations welcome:

topidx <- function(map) {

if(!(class(map) == "SpatialGridDataFrame"))
 	printf("Map should be of SpatialGridDataFrame class")

values <- map at data@att[[1]]
ew_res <- map at grid@cellsize[1]
ns_res <- map at grid@cellsize[2]
cols <- map at grid@cells.dim[1]
rows <- map at grid@cells.dim[2]

# here comes al sorts of data check and NA handling...

   .C("topidx",
 	PACKAGE = "mypackage",
 	as.double(values),
 	as.integer(rows),
 	as.integer(cols),
 	as.double(ew_res),
 	as.double(ns_res),
 	topidxmap = double(rows*cols))$topidxmap

# return is a double vector which can be postprocessed here...
}

cheers
wouter

        
On Tue, 7 Mar 2006, Roger Bivand wrote:

        

            

      
      
      
    

        

      
      
    

            

      
      
      
    

        

      
      
    

            

      
      
      
    

  
  


  


      

    

    

      
      

        


  

        

      Barry Rowlingson
    

    

      
      #
    

  


  

      
Wouter Buytaert wrote:

        

            

      
      
      
    

        
Recommendation number 1 - read all about Object-Oriented Programming!

  You shouldn't have a function that tests the class of an object, you 
should write a generic function and a method for your class.

  With S3 methods this was easy, you wrote a generic function called 
'topidx' and then a method for SpatialWhateverObject called 
topidx.SpatialWhateverObject. Then your topidx.SpatialWhateverObject() 
function could be pretty sure it was getting something that was a 
SpatialWhateverObject. I've never got my head round the new S4 class 
stuff, but sp uses it widely so check out Edzer's code, he's a master :)

  The reason explicit checks for class(foo) are wrong is because your 
object might not be of that class, but might inherit from that class, ie 
be a superclass of SpatialGridDataFrame. Add to that the fact that 
class() can be a vector:

  > class(glm(y~x))
  [1] "glm" "lm"

  If you really have to check to see if the object you are being passed 
is of the right sort, use 'inherits':

  > inherits(glm(y~x),'lm')
  [1] TRUE
  > inherits(glm(y~x),'glm')
  [1] TRUE

This should be sufficient to guarantee that the object implements the 
functionality you need to do your work.

Otherwise, I vote for .Call() :)

  Barry

  
  


  


      

    

    

      
      

        


  

        

      Tim Keitt
    

    

      
      #
    

  


  

      
On Tue, 2006-03-07 at 16:24 +0000, Barry Rowlingson wrote:

            

      
      
      
    

        
That's a rather odd aspect of S methods -- even if you only want to call
the function for a single class, you have to create a generic function.
In cases where I'm not wanting polymorphic behavior (the function really
only makes sense for a single class), I tend to write ordinary functions
and either test for class membership or just document that the passed
object should be of the appropriate class. When the function will be
called on multiple classes, then it is worth defining a generic
function. Making everything generic is too much overhead.

THK

  
    

      
      
    

  


  


      

    

    

      
      

        


  

        

      Barry Rowlingson
    

    

      
      #
    

  


  

      
Tim Keitt wrote:

        

            

      
      
      
    

        
one of many :)

            

      
      
      
    

        
I think this stems from the S3 way that what other languages call 
"methods of objects", R has as global namespace functions. In Python for 
example, the method name is tied to the object, so you do: foo.bar() 
instead of bar(foo). S3 shoe-horned the OO paradigm into the existing 
functional syntax instead of extending the syntax such that methods were 
true properties of objects. S4 methods continue this (IMHO) poor design 
choice, and then tries to tidy up by having 'namespaces' declared in 
packages - but at least it has a more formal method for declaring 
methods instead of S3's reliance on <generic>.<class> naming scheme!

  If I could transparently bind all the statistical and graphical 
functionality of R with the Python language I'd be happy :)[1]

Barry

[1] Yeah, I know there are R-Python links but they're not quite 
transparent enough....

  
  


  


      

    

    

      
      

        


  

        

      Tim Keitt
    

    

      
      #
    

  


  

      
On Tue, 2006-03-07 at 17:39 +0000, Barry Rowlingson wrote:

            

      
      
      
    

        
I couldn't agree more. It seems namespaces in R confuse binding names to
namespaces with encapsulation (hiding symbols within packages). For
efficiency reasons, your forced to export all your public symbols into
the global namespace anyway. Ugh!

As you point out, it would be nice if classes defined namespaces so that
methods could be attached to the class as in C++ and Python. Of course
programming S/R would not be nearly so much fun if it didn't have so
many quirks... :-)

            

      
      
      
    

        
I do prefer S4 methods. On the other hand, I found they really slowed my
rate of programming. Several of my packaging projects (DBI.PgSQL) got
derailed while I came up to speed with S4 methods.

THK

            

      
      
      
    

  
    

      
      
    

  


  


      

    

    

      
      

        


  

        

      Edzer Pebesma
    

    

      
      #
    

  


  

      
Barry Rowlingson wrote:

        

            

      
      
      
    

        
Wouter, I have never completely understood Barry's reservations against
the S4 object orientation framework, which he never managed (or wanted?) to
get a real grip on; please don't let it intimidate you. There are many 
ways of
checking for classes; I would prefer is(map, "SpatialGridDataFrame") because
it will then work as well for classes that are derived from 
SpatialGridDataFrame.

Another issue in this thread: R namespaces and S4 classes are two rather
different issues; at some areas they even seem not to work together too 
well.
Ask Roger (in private, please!) for details, or search through r-devel.

            

      
      
      
    

        
I managed to have package gstat using exclusively .Call() at some
stage, but had to work into poorly documented areas of this interface:
R has quite a different way of dealing with strings than that documented
in the green book (Programming with data, J. Chambers).

The final pleasure was mostly that of aesthetics, and had nothing
to do with "good code".
--
Edzer