[R-sig-Geo] Including SAGA grid as GDAL-supported format

[Alexander Brenning]

Hi,

a GDAL interface to SAGA grids would of course be very welcome. Right 
now the easiest (but indirect) way of importing SAGA grids is with

   read.sgrd("mysagagrid")

which uses rsaga.sgrd.to.esri (to a tempfile()) and then read.ascii.grid

In the next RSAGA release I will consider using readBin() in read.sgrd() 
  in order to read the SAGA grid file directly as you suggested.

Cheers
  Alex




Tomislav Hengl wrote:
> Dear r-sig-geo,
> 
> I would like to initiate the processes of registering the SAGA grid  
> format under GDAL (the SAGA 2.0.4 source code is available for  
> download here  
> http://sourceforge.net/projects/saga-gis/files/SAGA%20-%202.0/SAGA%202.0.4/saga_2.0.4_src.zip/download).
> 
> Do I have to follow some formal procedure, or do I have to prepare the  
> GDAL driver myself (as explained at  
> http://www.gdal.org/gdal_drivertut.html)? Any help/suggestions are  
> welcome (apparently it should not be too complicated).
> 
> SAGA grid consists of tree types of files:
> 1. "*.sgrd" - the header file with name, data format, XLL, YLL, rows  
> columns, cell size, z-factor and no data value;
> 2. "*.sdat" - the raw data file;
> 3. "*.hgrd" - the history file;
> 
> Here are some examples hot to read and write the SAGA grids to R:
> 
>> library(gstat)
>> library(RSAGA)
>> library(spatstat)
> 
>> data(meuse.grid)
>> coordinates(meuse.grid) <- ~x+y
>> gridded(meuse.grid) <- TRUE
>> proj4string(meuse.grid) = CRS("+init=epsg:28992")
> 
> # write to SAGA grid format;
>> write.asciigrid(meuse.grid["soil"], "meuse_soil.asc")
>> rsaga.esri.to.sgrd(in.grids="meuse_soil.asc",  
>> out.sgrd="meuse_soil.sgrd", in.path=getwd())
> 
> # read SAGA grid format:
>> sgrd <- matrix((unlist(strsplit(readLines(file("meuse_soil.sgrd")),  
>> split="\t= "))), ncol=2, byrow=T)
>> sgrd
>        [,1]              [,2]
>   [1,] "NAME"            "meuse_soil"
>   [2,] "DESCRIPTION"     "UNIT"
>   [3,] "DATAFILE_OFFSET" "0"
>   [4,] "DATAFORMAT"      "FLOAT"
>   [5,] "BYTEORDER_BIG"   "FALSE"
>   [6,] "POSITION_XMIN"   "178460.0000000000"
>   [7,] "POSITION_YMIN"   "329620.0000000000"
>   [8,] "CELLCOUNT_X"     "78"
>   [9,] "CELLCOUNT_Y"     "104"
> [10,] "CELLSIZE"        "40.0000000000"
> [11,] "Z_FACTOR"        "1.000000"
> [12,] "NODATA_VALUE"    "-9999.000000"
> [13,] "TOPTOBOTTOM"     "FALSE"
> 
> # read the raw data: 4bit, numeric (FLOAT), byte order small;
>> sdat <- readBin("meuse_soil.sdat", what="numeric", size=4,  
>> n=as.integer(sgrd[8,2])*as.integer(sgrd[9,2]))
>> sdat.sp <- as.im(list(x=seq(from=as.integer(sgrd[6,2]),  
>> length.out=as.integer(sgrd[8,2]), by=as.integer(sgrd[10,2])),  
>> y=seq(from=as.integer(sgrd[7,2]), length.out=as.integer(sgrd[9,2]),  
>> by=as.integer(sgrd[10,2])), z=matrix(sdat,  
>> nrow=as.integer(sgrd[8,2]), ncol=as.integer(sgrd[9,2]))))
>> sdat.sp <- as(sdat.sp, "SpatialGridDataFrame")
> # replace the mask value with NA's:
>> sdat.sp at data[[1]] <-  
>> ifelse(sdat.sp at data[[1]]==as.integer(sgrd[12,2]), NA,  
>> sdat.sp at data[[1]])
>> spplot(sdat.sp)
> 
> 
> Of course, it would be much easier to have this in a single line:
> 
>> meuse.grid <- readGDAL("meuse_soil.sgrd")
> 
> or
> 
>> writeGDAL(meuse.grid["soil"], "meuse_soil.sgrd", "SAGA")
> 
> 
> 
> PS: A new version of SAGA has just been released few days ago.
> 
> 
> thank you,
> 
> T. Hengl
> http://home.medewerker.uva.nl/t.hengl/
> 
> Connected discussion:  
> https://stat.ethz.ch/pipermail/r-sig-geo/attachments/20080130/cd5c3748/attachment.pl
> 
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at stat.math.ethz.ch
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> 

-- 
Alexander Brenning
brenning at uwaterloo.ca - T +1-519-888-4567 ext 35783
Department of Geography and Environmental Management
University of Waterloo
200 University Ave. W - Waterloo, ON - Canada N2L 3G1
http://www.fes.uwaterloo.ca/geography/faculty/brenning/