Read GRASS 7 raster files from GRASS 7 into R SpatialGridDataFrame objects, and write single columns of R SpatialGridDataFrame objects to GRASS 7. readRAST and writeRAST use temporary binary files and r.out.bin and r.in.bin for speed reasons.

readRAST(vname, cat=NULL, ignore.stderr = get.ignore.stderrOption(),
 NODATA=NULL, plugin=get.pluginOption(), mapset=NULL,
 useGDAL=get.useGDALOption(), close_OK=TRUE, drivername="GTiff",
 driverFileExt=NULL, return_SGDF=TRUE)
writeRAST(x, vname, zcol = 1, NODATA=NULL,
 ignore.stderr = get.ignore.stderrOption(), useGDAL=get.useGDALOption(),
 overwrite=FALSE, flags=NULL, drivername="GTiff")

Arguments

vname

A vector of GRASS 7 raster file names

cat

default NULL; if not NULL, must be a logical vector matching vname, stating which (CELL) rasters to return as factor

ignore.stderr

default taking the value set by set.ignore.stderrOption; can be set to TRUE to silence system() output to standard error; does not apply on Windows platforms

plugin

default taking the value set by set.pluginOption; NULL does auto-detection, changes to FALSE if vname is longer than 1, and a sanity check will be run on raster and current region, and the function will revert to FALSE if mismatch is found; if TRUE, the plugin is available and the raster should be read in its original region and resolution; if the plugin is used, no further arguments other than mapset are respected

mapset

default NULL, if plugin is TRUE, the mapset of the file to be imported will be autodetected; if not NULL and if plugin is TRUE, a character string overriding the autodetected mapset, otherwise ignored

useGDAL

(effectively defunct, only applies to use of plugin) default taking the value set by set.useGDALOption; use plugin and readGDAL if autodetected or plugin=TRUE; or for writing writeGDAL, GTiff, and r.in.gdal, if FALSE using r.out.bin or r.in.bin

close_OK

default TRUE - clean up possible open connections used for reading metadata; may be set to FALSE to avoid the side-effect of other user-opened connections being broken

drivername

default "GTiff"; a valid GDAL writable driver name to define the file format for intermediate files

driverFileExt

default NULL; otherwise string value of required driver file name extension

return_SGDF

default TRUE returning a SpatialGridDataFrame object, if FALSE, return a list with a GridTopology object, a list of bands, and a proj4string; see example below

x

A SpatialGridDataFrame object for export to GRASS as a raster layer

zcol

Attribute column number or name

NODATA

by default NULL, in which case it is set to one less than floor()of the data values, otherwise an integer NODATA value (required to be integer by GRASS r.out.bin)

overwrite

default FALSE, if TRUE inserts "overwrite" into the value of the flags argument if not already there to allow existing GRASS rasters to be overwritten

flags

default NULL, character vector, for example "overwrite"

Value

readRAST returns a SpatialGridDataFrame objects with an data.frame in the data slots, and with the projection argument set. Note that the projection argument set is the the GRASS rendering of proj4, and will differ from the WKT/ESRI rendering returned by readVECT in form but not meaning. They are exchangeable but not textually identical, usually with the +ellps= term replaced by ellipsoid parameters verbatim. If return_SGDF is FALSE, a list with a GridTopology object, a list of bands, and a proj4string is returned, with an S3 class attribute of “gridList”.

Author

Roger S. Bivand, e-mail: Roger.Bivand@nhh.no

Examples

use_sp() run <- FALSE if (nchar(Sys.getenv("GISRC")) > 0 && read.dcf(Sys.getenv("GISRC"))[1,"LOCATION_NAME"] == "nc_basic_spm_grass7") run <- TRUE GV <- Sys.getenv("GRASS_VERBOSE") Sys.setenv("GRASS_VERBOSE"=0) # require(rgdal) ois <- get.ignore.stderrOption() set.ignore.stderrOption(TRUE)
#> [1] FALSE
#> [1] TRUE
if (run) { nc_basic <- readRAST(c("geology", "elevation"), cat=c(TRUE, FALSE), useGDAL=FALSE) nc_basic <- readRAST(c("geology", "elevation"), cat=c(TRUE, FALSE), useGDAL=TRUE) print(table(nc_basic$geology)) }
#> Warning: Discarded datum NAD83 (High Accuracy Reference Network) in CRS definition
#> Warning: non-unique category labels; category number appended
#> Warning: Discarded datum NAD83 (High Accuracy Reference Network) in CRS definition
#> Warning: non-unique category labels; category number appended
#> #> CZfg_217 CZlg_262 CZig_270 CZbg_405 CZve_583 CZam_720 CZg_766 CZam_862 #> 725562 198684 689373 253710 21609 4824 7074 61722 #> CZbg_910 Km_921 CZbg_945 CZam_946 CZam_948 #> 44964 12528 9 4068 873
if (run) { execGRASS("r.stats", flags=c("c", "l", "quiet"), input="geology") }
#> 217 CZfg 725562 #> 262 CZlg 198684 #> 270 CZig 689373 #> 405 CZbg 253710 #> 583 CZve 21609 #> 720 CZam 4824 #> 766 CZg 7074 #> 862 CZam 61722 #> 910 CZbg 44964 #> 921 Km 12528 #> 945 CZbg 9 #> 946 CZam 4068 #> 948 CZam 873
if (run) { boxplot(nc_basic$elevation ~ nc_basic$geology) }
if (run) { nc_basic$sqdem <- sqrt(nc_basic$elevation) } if (run) { if ("GRASS" %in% rgdal::gdalDrivers()$name) { execGRASS("g.region", raster="elevation") dem1 <- readRAST("elevation", plugin=TRUE, mapset="PERMANENT") print(summary(dem1)) execGRASS("g.region", raster="elevation") } } if (run) { writeRAST(nc_basic, "sqdemSP", zcol="sqdem", flags=c("quiet", "overwrite")) execGRASS("r.info", map="sqdemSP") }
#> +----------------------------------------------------------------------------+ #> | Map: sqdemSP Date: Wed Jan 6 13:39:48 2021 | #> | Mapset: rsb Login of Creator: rsb | #> | Location: nc_basic_spm_grass7 | #> | DataBase: /home/rsb/topics/grassdata | #> | Title: | #> | Timestamp: none | #> |----------------------------------------------------------------------------| #> | | #> | Type of Map: raster Number of Categories: 0 | #> | Data Type: DCELL | #> | Rows: 1350 | #> | Columns: 1500 | #> | Total Cells: 2025000 | #> | Projection: Lambert Conformal Conic | #> | N: 228500 S: 215000 Res: 10 | #> | E: 645000 W: 630000 Res: 10 | #> | Range of data: min = 7.45511854848878 max = 12.5031941719687 | #> | | #> | Data Description: | #> | generated by r.in.bin | #> | | #> | Comments: | #> | r.in.bin --overwrite --quiet -d input="/home/rsb/topics/grassdata/nc\ | #> | _basic_spm_grass7/rsb/.tmp/localhost.localdomain/X601" output="sqdem\ | #> | SP" bytes=8 header=0 bands=1 order="native" north=228500 south=21500\ | #> | 0 east=645000 west=630000 rows=1350 cols=1500 anull=6 | #> | | #> +----------------------------------------------------------------------------+ #>
if (run) { execGRASS("g.remove", flags="f", name="sqdemSP", type="raster") } if (run) { writeRAST(nc_basic, "sqdemSP", zcol="sqdem", useGDAL=TRUE, flags=c("quiet", "overwrite")) execGRASS("r.info", map="sqdemSP") }
#> +----------------------------------------------------------------------------+ #> | Map: sqdemSP Date: Wed Jan 6 13:39:48 2021 | #> | Mapset: rsb Login of Creator: rsb | #> | Location: nc_basic_spm_grass7 | #> | DataBase: /home/rsb/topics/grassdata | #> | Title: | #> | Timestamp: none | #> |----------------------------------------------------------------------------| #> | | #> | Type of Map: raster Number of Categories: 0 | #> | Data Type: DCELL | #> | Rows: 1350 | #> | Columns: 1500 | #> | Total Cells: 2025000 | #> | Projection: Lambert Conformal Conic | #> | N: 228500 S: 215000 Res: 10 | #> | E: 645000 W: 630000 Res: 10 | #> | Range of data: min = 7.45511854848878 max = 12.5031941719687 | #> | | #> | Data Description: | #> | generated by r.in.bin | #> | | #> | Comments: | #> | r.in.bin --overwrite --quiet -d input="/home/rsb/topics/grassdata/nc\ | #> | _basic_spm_grass7/rsb/.tmp/localhost.localdomain/X158" output="sqdem\ | #> | SP" bytes=8 header=0 bands=1 order="native" north=228500 south=21500\ | #> | 0 east=645000 west=630000 rows=1350 cols=1500 anull=6 | #> | | #> +----------------------------------------------------------------------------+ #>
if (run) { print(system.time(sqdemSP <- readRAST(c("sqdemSP", "elevation"), useGDAL=TRUE, return_SGDF=FALSE))) }
#> Warning: Discarded datum NAD83 (High Accuracy Reference Network) in CRS definition
#> user system elapsed #> 0.549 0.065 0.719
if (run) { print(system.time(sqdemSP <- readRAST(c("sqdemSP", "elevation"), useGDAL=TRUE, return_SGDF=TRUE))) }
#> Warning: Discarded datum NAD83 (High Accuracy Reference Network) in CRS definition
#> user system elapsed #> 0.539 0.083 0.740
if (run) { print(system.time(sqdemSP <- readRAST(c("sqdemSP", "elevation"), useGDAL=FALSE, return_SGDF=TRUE))) }
#> Warning: Discarded datum NAD83 (High Accuracy Reference Network) in CRS definition
#> user system elapsed #> 0.527 0.087 0.721
if (run) { print(system.time(sqdemSP <- readRAST(c("sqdemSP", "elevation"), useGDAL=FALSE, return_SGDF=FALSE))) }
#> Warning: Discarded datum NAD83 (High Accuracy Reference Network) in CRS definition
#> user system elapsed #> 0.523 0.090 0.724
if (run) { str(sqdemSP) mat <- do.call("cbind", sqdemSP$dataList) str(mat) }
#> List of 3 #> $ grid :Formal class 'GridTopology' [package "sp"] with 3 slots #> .. ..@ cellcentre.offset: num [1:2] 630005 215005 #> .. ..@ cellsize : num [1:2] 10 10 #> .. ..@ cells.dim : int [1:2] 1500 1350 #> $ dataList :List of 2 #> ..$ sqdemSP : num [1:2025000] 11.9 11.9 11.9 11.9 12 ... #> ..$ elevation: num [1:2025000] 142 141 141 142 143 ... #> $ proj4string:Formal class 'CRS' [package "sp"] with 1 slot #> .. ..@ projargs: chr "+proj=lcc +lat_0=33.75 +lon_0=-79 +lat_1=36.1666666666667 +lat_2=34.3333333333333 +x_0=609601.22 +y_0=0 +ellps="| __truncated__ #> - attr(*, "class")= chr "gridList" #> num [1:2025000, 1:2] 11.9 11.9 11.9 11.9 12 ... #> - attr(*, "dimnames")=List of 2 #> ..$ : NULL #> ..$ : chr [1:2] "sqdemSP" "elevation"
if (run) { print(system.time(SGDF <- sp::SpatialGridDataFrame(grid=sqdemSP$grid, proj4string=sqdemSP$proj4string, data=as.data.frame(sqdemSP$dataList)))) }
#> user system elapsed #> 0.002 0.000 0.002
if (run) { summary(SGDF) }
#> Length Class Mode #> 2025000 SpatialGridDataFrame S4
if (run) { execGRASS("g.remove", flags="f", name="sqdemSP", type="raster") execGRASS("r.mapcalc", expression="basins0 = basins - 1") execGRASS("r.stats", flags="c", input="basins0") }
#> 1 116975 #> 3 75480 #> 5 1137 #> 7 80506 #> 9 7472 #> 11 348209 #> 13 51456 #> 15 81959 #> 17 29652 #> 19 267883 #> 21 89465 #> 23 63433 #> 25 91354 #> 27 190239 #> 29 59727 #> * 470053
if (run) { basins0 <- readRAST("basins0") print(table(basins0$basins0)) }
#> Warning: Discarded datum NAD83 (High Accuracy Reference Network) in CRS definition
#> #> 1 3 5 7 9 11 13 15 17 19 21 #> 116975 75480 1137 80506 7472 348209 51456 81959 29652 267883 89465 #> 23 25 27 29 #> 63433 91354 190239 59727
if (run) { basins0 <- readRAST("basins0", plugin=FALSE) print(table(basins0$basins0)) }
#> Warning: Discarded datum NAD83 (High Accuracy Reference Network) in CRS definition
#> #> 1 3 5 7 9 11 13 15 17 19 21 #> 116975 75480 1137 80506 7472 348209 51456 81959 29652 267883 89465 #> 23 25 27 29 #> 63433 91354 190239 59727
if (run) { execGRASS("g.remove", flags="f", name="basins0", type="raster") } Sys.setenv("GRASS_VERBOSE"=GV) set.ignore.stderrOption(ois)
#> [1] TRUE