fix issue where readRAST fails with grass >=7.6

[ltalluto]

The most recent version changes to reading WKT projections from grass by default with grass >= 7.6, which is not working for rasters. I was getting no arguments in initialization list whenever I tried to read a raster.

This seemed like a simple fix, so I just did it myself rather than creating an issue, but happy for feedback if this is somehow breaking something.

[rsbivand]

No, I must have a reporoducible example and full analysis. For rgrass7 0.2-5, GRASS 7.8.5 on Fedora 33 (PROJ, GDAL, GRASS etc. installed from source), and the small NC standard data set, and all the examples work as expected. WKT2 is always to be prefered for GRASS >= 7.6. I suspect that your installation(s) of sp and/or rgdal are faulty, probably because the PROJ data directory is not present. You do not say which OS you are using, or how you installed GRASS, GDAL or PROJ.

[ltalluto]

Indeed, sorry for the lack of info. I've tried 2 test systems now. First: macOS 10.15.7, R 3.6.3, GRASS 7.6.1 (using the mac binary installer), CRAN binary installations of rgrass7 0.2-5, sp 1.4-5, rgdal 1.4-4, raster 3.3-11. PROJ is included with the binary of rgdal and the the directory is where it is supposed to be. Manually loading rgdal prints a message with the proj directory. Both proj and sp are working normally.

Second system: macOS 10.15.7, R 4.0.1, GRASS 7.8.5, CRAN binary installs of rgrass7 0.2-5, sp 1.4-5. rgdal 1.5-23, raster 3.4-5.

The following fails reliably with the first system, succeeds on the second. I was also getting the error on a Windows system last week, but it's not available at the moment for testing, so I will have to come back to it later. If I get a moment, I can try up- and downgrading rgdal on the two macs to see if it is reproducible across systems just by changing rgdal versions. The first system was also working fine on rgrass7 0.2-3.

library(rgrass7)
library(raster)
use_sp()
gisBase = "/Applications/GRASS-7.8.app/Contents/Resources"
ra = raster(system.file("external/test.grd", package = "raster"))
ext = as.character(as.vector(extent(ra)))
nr=nrow(ra)
nc = ncol(ra)
res = as.character(res(ra))
ra = as(ra, "SpatialGridDataFrame")
initGRASS(gisBase, home=tempdir(), gisDbase = tempdir(), location = "test", mapset = "PERMANENT", 
    override = TRUE)
execGRASS("g.proj", flags = "c", proj4 = sp::proj4string(ra), intern=TRUE)
execGRASS("g.region", n = ext[4], s = ext[3], e = ext[2], w = ext[1], rows=nr, 
    cols=nc, nsres = res[2], ewres = res[1])

writeRAST(ra, "test")
out = readRAST("test") ## fails

[rsbivand]

rgdal < 1.5 are built against PROJ < 6. So although GRASS >= 7.6 do provide WKT2, sp/rgdal cannot use WKT2 to create a CRS object. I'll add a version check for rgdal in that condition.

[ltalluto]

Ahh, good to know, thanks! Time to update my test machine I think.

[rsbivand]

Please try remotes::install_github("rsbivand/rgrass7"), to try out https://github.com/rsbivand/rgrass7/commit/d8a6705d6f8da38bb660850927a307fcd5f82a47 to check that old rgdal is properly trapped.

[ltalluto]

Now I get:

initGRASS(gisBase, home=tempdir(), gisDbase = tempdir(), location = "test", mapset = "PERMANENT", 
     override = TRUE)
Error in if (substr(x$proj4, 1, 1) == "+") cat("projection ", paste(strwrap(x$proj4),  : 
  missing value where TRUE/FALSE needed

Grass still starts, and I am able to read a raster now.

[rsbivand]

I tried to repair the test criteria, but would be interested in seeing this:

x <- initGRASS(gisBase, home=tempdir(), gisDbase = tempdir(), location = "test", mapset = "PERMANENT", override = TRUE)
str(x)
str(x$proj4)

[ltalluto]

Interesting, that error message does not print if I assign the result of initGRASS to x. It also does not print if I run the next two lines (g.proj and g.region), and then re-run initGRASS.

str(x)
List of 25
 $ GISDBASE     : chr "/var/folders/p5/l0zg51f17nb8_lfv6888gxh40000gn/T//Rtmpq44ufD"
 $ LOCATION_NAME: chr "test"
 $ MAPSET       : chr "PERMANENT"
 $ GRASS_GUI    : chr "text"
 $ projection   : chr "0"
 $ zone         : chr "0"
 $ n            : num 1
 $ s            : num 0
 $ w            : num 0
 $ e            : num 1
 $ t            : num 1
 $ b            : num 0
 $ nsres        : num 1
 $ nsres3       : num 1
 $ ewres        : num 1
 $ ewres3       : num 1
 $ tbres        : num 1
 $ rows         : int 1
 $ rows3        : int 1
 $ cols         : int 1
 $ cols3        : int 1
 $ depths       : int 1
 $ cells        : chr "1"
 $ cells3       : chr "1"
 $ proj4        : chr NA
 - attr(*, "class")= chr "gmeta"

str(x$proj4)
 chr NA

[ltalluto]

Whoops, sorry, didn't see the new commit. Working fine now!