search

rsbivand / rgrass

Interpreted interface between the GRASS geographical information system and R
https://rsbivand.github.io/rgrass/
26 stars 9 forks source link

Issues with projection when using RGRASS7 #32

Closed manomuthus closed 3 years ago

manomuthus commented 3 years ago

I have just started using rgrass7 and getting error messages when I try to use readRAST. Please see below

initGRASS(gisBase="C:/Program Files/GRASS GIS 7.8", mapset = "PERMANENT", override = T)

rast_sgdf <- as(rast, "SpatialGridDataFrame")
writeRAST(rast_sgdf, vname="rast_grass", flags = "overwrite")

trying to read rast_grass 

readRAST("rast_grass")
Error in sp::CRS(getLocationProj()) : 
  PROJ4 argument-value pairs must begin with +: XY location (unprojected)

I read somewhere that setting the projection using proj4 will solve this issue, so tried the following 

execGRASS("g.proj", flags = "c", proj4 = "+proj=utm +zone=45 +datum=WGS84 +units=m +no_defs")

but now I get the following error message

readRAST("rast_grass")
Error in sp::CRS(getLocationProj()) : 
  PROJ4 argument-value pairs must begin with +: PROJCS["unknown",
    GEOGCS["wgs84",
        DATUM["WGS_1984",
            SPHEROID["WGS_1984",6378137,298.257223563]],
        PRIMEM["Greenwich",0,
            AUTHORITY["EPSG","8901"]],
        UNIT["degree",0.0174532925199433,
            AUTHORITY["EPSG","9122"]]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",87],
    PARAMETER["scale_factor",0.9996],
    PARAMETER["false_easting",500000],
    PARAMETER["false_northing",0],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]],
    AXIS["Easting",EAST],
    AXIS["Northing",NORTH]]

What am I doing wrong?

data

dput(rast)
new("RasterLayer", file = new(".RasterFile", name = "", datanotation = "FLT4S", 
    byteorder = "little", nodatavalue = -3.4e+38, NAchanged = FALSE, 
    nbands = 1L, bandorder = "BIL", offset = 0L, toptobottom = TRUE, 
    blockrows = 16L, blockcols = 128L, driver = "", open = FALSE), 
    data = new(".SingleLayerData", values = c(NA, NA, NA, NA, 
    NA, 1208.45349121094, 1207.15441894531, 1207.19006347656, 
    1207.94274902344, 1207.89782714844, 1207.94274902344, 1209.11303710938, 
    1210.1640625, 1208.50305175781, 1209.26330566406, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 1208.56213378906, 
    1207.99792480469, 1207.19848632812, 1207.94274902344, 1207.94274902344, 
    1208.01794433594, 1208.43347167969, 1208.52490234375, 1210.61877441406, 
    1209.56994628906, 1212.35241699219, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, 1208.55212402344, 1208.19763183594, 
    1207.61181640625, 1207.9072265625, 1208.43347167969, 1208.53259277344, 
    1208.56506347656, 1208.56506347656, 1209.91125488281, 1209.72192382812, 
    1217.07543945312, 1219.02136230469, 1219.72485351562, 1217.35620117188, 
    1214.68359375, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    1208.24816894531, 1208.21142578125, 1207.71325683594, 1207.91076660156, 
    1208.19030761719, 1208.56506347656, 1208.53442382812, 1208.56506347656, 
    1211.30285644531, 1215.82885742188, 1214.61682128906, 1221.67639160156, 
    1220.3544921875, 1219.9892578125, 1217.69470214844, 1221.14685058594, 
    1219.81823730469, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    1208.55249023438, 1208.02697753906, 1208.017578125, 1207.400390625, 
    1208.65942382812, 1208.56506347656, 1208.53625488281, 1208.76940917969, 
    1208.76940917969, 1208.76940917969, 1208.76940917969, 1209.5849609375, 
    1217.29479980469, 1222.43395996094, 1220.53955078125, 1222.69262695312, 
    1221.01965332031, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    1208.05090332031, 1207.9970703125, 1207.96997070312, 1208.05041503906, 
    1208.76940917969, 1208.76940917969, 1208.76940917969, 1208.76940917969, 
    1208.76940917969, 1209.5849609375, 1209.5849609375, 1209.5849609375, 
    1222.4541015625, 1220.85107421875, 1221.82141113281, 1216.05773925781, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 1208.44177246094, 
    1208.52172851562, 1208.26806640625, 1208.46484375, 1208.78295898438, 
    1207.81494140625, 1208.16723632812, 1208.69140625, 1209.5849609375, 
    1209.5849609375, 1209.5849609375, 1209.720703125, 1213.89770507812, 
    1214.76220703125, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, 1209.41796875, 1208.56091308594, 1208.10803222656, 
    1208.53002929688, 1208.42626953125, 1207.88256835938, 1208.81066894531, 
    1208.83459472656, 1209.5849609375, 1209.52807617188, 1209.81579589844, 
    1214.59741210938, 1223.80065917969, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 1208.49841308594, 
    1209.14697265625, 1208.54370117188, 1208.89050292969, 1209.81579589844, 
    1209.81579589844, 1214.78759765625, 1221.84899902344, 1223.90710449219, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, 1208.6904296875, 1208.94995117188, 1209.23315429688, 
    1209.67565917969, 1209.67565917969, 1214.1171875, 1215.07580566406, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, 1208.83178710938, 1208.76831054688, 1209.67565917969, 
    1209.67565917969, 1211.46594238281, 1214.94274902344, 1215.39904785156, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, 1209.25732421875, 1208.86218261719, 1209.048828125, 
    1209.67565917969, 1209.67565917969, 1211.82763671875, 1215.14599609375, 
    1215.67370605469, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, 1209.06079101562, 1209.17163085938, 
    1208.79455566406, 1209.67565917969, 1209.67565917969, 1213.59326171875, 
    1220.58471679688, 1217.48803710938, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 1210.30834960938, 
    1209.05725097656, 1210.30346679688, 1209.67565917969, 1209.67565917969, 
    1215.83471679688, 1222.40283203125, 1218.2578125, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, 1209.61926269531, 1210.16052246094, 1210.56823730469, 
    1210.79309082031, 1210.46276855469, 1214.73291015625, 1216.96398925781, 
    1219.07653808594, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, 1209.46716308594, 1209.22875976562, 
    1210.47521972656, 1210.0283203125, 1210.0283203125, 1212.83312988281, 
    1216.43481445312, 1219.17358398438, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 1210.00036621094, 
    1209.9521484375, 1209.96520996094, 1210.0283203125, 1210.0283203125, 
    1213.79821777344, 1217.90661621094, 1219.22473144531, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, 1209.54418945312, 1210.26818847656, 1210.0283203125, 
    1210.03039550781, 1210.0283203125, 1212.66833496094, 1216.18420410156, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, 1209.45166015625, 1209.70434570312, 1209.33239746094, 
    1210.0283203125, 1210.0283203125, 1211.15905761719, 1212.17797851562, 
    1214.48999023438, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    NA, NA, NA, NA, NA, NA, NA, NA, 1209.94604492188, 1210.36206054688, 
    1210.103515625, 1210.0283203125, 1210.0283203125, 1211.44226074219, 
    1214.18041992188, NA, NA), offset = 0, gain = 1, inmemory = TRUE, 
        fromdisk = FALSE, isfactor = FALSE, attributes = list(), 
        haveminmax = TRUE, min = 1207.15441894531, max = 1223.90710449219, 
        band = 1L, unit = "", names = "nakkhu_hosp_10m"), legend = new(".RasterLegend", 
        type = character(0), values = logical(0), color = logical(0), 
        names = logical(0), colortable = logical(0)), title = character(0), 
    extent = new("Extent", xmin = 333003.9801, xmax = 333403.9801, 
        ymin = 3060402.4038, ymax = 3060602.4038), rotated = FALSE, 
    rotation = new(".Rotation", geotrans = numeric(0), transfun = function () 
    NULL), ncols = 40L, nrows = 20L, crs = new("CRS", projargs = "+proj=utm +zone=45 +datum=WGS84 +units=m +no_defs"), 
    history = list(), z = list())
rsbivand commented 3 years ago

Please see #27 and #31, and try to install from source from github, or from the link to a Windows binary package https://github.com/rsbivand/rgrass7/files/6285516/rgrass7_0.2-6.zip built some time ago. Please report back. Note that you should set the projection immediately after initialising GRASS, as the location you create should have the same projection as the files being added. If you can confirm tha this removes your problem, I'll submit to CRAN soon.

manomuthus commented 3 years ago

I unistalled previous veriosn and reinstalled rgrass7 (0.2-6) from your link. But I still get the following error.

Error in sp::CRS(getLocationProj()) : 
  PROJ4 argument-value pairs must begin with +: PROJCS["unknown",
    GEOGCS["wgs84",
        DATUM["WGS_1984",
            SPHEROID["WGS_1984",6378137,298.257223563]],
        PRIMEM["Greenwich",0,
            AUTHORITY["EPSG","8901"]],
        UNIT["degree",0.0174532925199433,
            AUTHORITY["EPSG","9122"]]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",87],
    PARAMETER["scale_factor",0.9996],
    PARAMETER["false_easting",500000],
    PARAMETER["false_northing",0],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]],
    AXIS["Easting",EAST],
    AXIS["Northing",NORTH]]

I set the projection immediatly after initialising GRASS using execGRASS("g.proj", flags = "c", proj4 = "+proj=utm +zone=45 +datum=WGS84 +units=m +no_defs"). The proj4 string is derived using crs(rast). I hope I am doing everything right.

rsbivand commented 3 years ago

I'll re-try on Windows (I'm in Linux now), the binary package is some months old. Please also add all the code you used to apply the projection (immediately after instantiating the GRASS session), it may be your use of "PERMANENT", but without the code it is hard to say.

manomuthus commented 3 years ago

Thank you. The complete code is provided in the initial comment.

rsbivand commented 3 years ago

Windows, R 4.1.1, GRASS 7.8.5, CRAN sp and rgrass7, no g.proj: 

> library(raster)
Loading required package: sp
> load("rast.rda")
> library(rgrass7)
Loading required package: XML
GRASS GIS interface loaded with GRASS version: (GRASS not running)
> initGRASS(gisBase="C:/Program Files/GRASS GIS 7.8", mapset = "PERMANENT", override=TRUE)
gisdbase    C:/Users/RB/AppData/Local/Temp/RtmpOUxb7k
location    file225037446716
mapset      PERMANENT
rows        1
columns     1
north       1
south       0
west        0
east        1
nsres       1
ewres       1
projection:
 XY location (unprojected)
> rast_sgdf <- as(rast, "SpatialGridDataFrame")
> use_sp()
> writeRAST(rast_sgdf, vname="rast_grass", flags = "overwrite")
 100%
> rg <- readRAST("rast_grass")
Error in sp::CRS(getLocationProj()) : NA
> packageVersion("sp")
[1] '1.4.5'
> packageVersion("rgrass7")
[1] '0.2.5'
>

Old updated rgrass7, no g.proj:

> install.packages("rgrass7_0.2-6.zip")
Installing package into 'C:/Users/RB/Documents/R/win-library/4.1'
(as 'lib' is unspecified)
inferring 'repos = NULL' from 'pkgs'
package 'rgrass7' successfully unpacked and MD5 sums checked
> library(raster)
Loading required package: sp
> library(rgrass7)
Loading required package: XML
GRASS GIS interface loaded with GRASS version: (GRASS not running)
> packageVersion("sp")
[1] '1.4.5'
> packageVersion("rgrass7")
[1] '0.2.6'
> load("rast.rda")
> initGRASS(gisBase="C:/Program Files/GRASS GIS 7.8", mapset = "PERMANENT", override=TRUE)
gisdbase    C:/Users/RB/AppData/Local/Temp/Rtmp0EUcCj
location    file1c05a116da1
mapset      PERMANENT
rows        1
columns     1
north       1
south       0
west        0
east        1
nsres       1
ewres       1
> rast_sgdf <- as(rast, "SpatialGridDataFrame")
> use_sp()
> writeRAST(rast_sgdf, vname="rast_grass", flags = "overwrite")
 100%
> rg <- readRAST("rast_grass")
Error in sp::CRS(getLocationProj()) : NA

Setting g.proj with old update rgrass7 shows a further problem:

>  library(raster)
Loading required package: sp
> library(rgrass7)
Loading required package: XML
GRASS GIS interface loaded with GRASS version: (GRASS not running)
> packageVersion("sp")
[1] '1.4.5'
> packageVersion("rgrass7")
[1] '0.2.6'
>  load("rast.rda")
> initGRASS(gisBase="C:/Program Files/GRASS GIS 7.8", mapset = "PERMANENT", override=TRUE)
gisdbase    C:/Users/RB/AppData/Local/Temp/Rtmpwx11mg
location    file18546c4872e2
mapset      PERMANENT
rows        1
columns     1
north       1
south       0
west        0
east        1
nsres       1
ewres       1
> execGRASS("g.proj", flags = "c", proj4 = "+proj=utm +zone=45 +datum=WGS84 +units=m +no_defs")
Default region was updated to the new projection, but if you have multiple
mapsets `g.region -d` should be run in each to update the region from the
default
Projection information updated
> getLocationProj()
ERROR 1: PROJ: proj_as_wkt: SQLite error on SELECT name, ellipsoid_auth_name, ellipsoid_code, prime_meridian_auth_name, prime_meridian_code, area_of_use_auth_name, area_of_use_code, publication_date, deprecated FROM geodetic_datum WHERE auth_name = ? AND code = ?: no such column: area_of_use_auth_name
[1] ""
> Sys.getenv("PROJ_LIB")
[1] "C:/Users/RB/Documents/R/win-library/4.1/rgdal/proj"
> rgdal::rgdal_extSoftVersion()
          GDAL GDAL_with_GEOS           PROJ             sp           EPSG
       "3.2.1"         "TRUE"        "7.2.1"        "1.4-5"      "v10.008"

This is not the same problem, but is because rgdal on load (used by sp::CRS()) on Windows sets PROJ_LIB, leading g.proj to fail in looking up the WKT2 version as it is not the same EPSG version (and proj.db was re-organized between EPSG 9 and 10):

> packageVersion("rgdal")
[1] '1.5.27'
> library(RSQLite)
> db <- dbConnect(SQLite(), dbname=file.path("C:/Program Files/GRASS GIS 7.8/share/proj", "proj.db"))
> (metadata <- dbReadTable(db, "metadata"))
           key                                                    value
1 EPSG.VERSION                                                   v9.8.6
2    EPSG.DATE                                               2020-01-22
3 ESRI.VERSION                                            ArcMap 10.8.0
4    ESRI.DATE                                               2019-12-01
5  IGNF.SOURCE https://geodesie.ign.fr/contenu/fichiers/IGNF.v3.1.0.xml
6 IGNF.VERSION                                                    3.1.0
7    IGNF.DATE                                               2019-05-24
> dbDisconnect(db)

Knowing this, a further work-around has been added to getLocationProj(), and a new Windows binary created: (rgrass7_0.2-6.zip. With this (and correcting an omission of correctly setting the GRASS region) and rast.rda rast.rda.zip:

> install.packages("rgrass7_0.2-6.zip")
Installing package into 'C:/Users/RB/Documents/R/win-library/4.1'
(as 'lib' is unspecified)
inferring 'repos = NULL' from 'pkgs'
package 'rgrass7' successfully unpacked and MD5 sums checked
> library(rgrass7)
Loading required package: XML
GRASS GIS interface loaded with GRASS version: (GRASS not running)
>  library(raster)
Loading required package: sp
> load("rast.rda")
>  initGRASS(gisBase="C:/Program Files/GRASS GIS 7.8", mapset = "PERMANENT", override=TRUE)
gisdbase    C:/Users/RB/AppData/Local/Temp/RtmpUzAsbn
location    file2f0444971c68
mapset      PERMANENT
rows        1
columns     1
north       1
south       0
west        0
east        1
nsres       1
ewres       1
> execGRASS("g.proj", flags = "c", proj4 = "+proj=utm +zone=45 +datum=WGS84 +units=m +no_defs")
Default region was updated to the new projection, but if you have multiple
mapsets `g.region -d` should be run in each to update the region from the
default
Projection information updated
> getLocationProj()
ERROR 1: PROJ: proj_as_wkt: SQLite error on SELECT name, ellipsoid_auth_name, ellipsoid_code, prime_meridian_auth_name, prime_meridian_code, area_of_use_auth_name, area_of_use_code, publication_date, deprecated FROM geodetic_datum WHERE auth_name = ? AND code = ?: no such column: area_of_use_auth_name
[1] "+proj=utm +no_defs +zone=45 +a=6378137 +rf=298.257223563 +towgs84=0.000,0.000,0.000 +type=crs  +to_meter=1"
> rast_sgdf <- as(rast, "SpatialGridDataFrame")
>  use_sp()
> writeRAST(rast_sgdf, vname="rast_grass", flags = "overwrite")
 100%
>  rg <- readRAST("rast_grass")
ERROR 1: PROJ: proj_as_wkt: SQLite error on SELECT name, ellipsoid_auth_name, ellipsoid_code, prime_meridian_auth_name, prime_meridian_code, area_of_use_auth_name, area_of_use_code, publication_date, deprecated FROM geodetic_datum WHERE auth_name = ? AND code = ?: no such column: area_of_use_auth_name
Creating BIL support files...
Exporting raster as floating values (bytes=8)
 100%
Warning message:
In showSRID(uprojargs, format = "PROJ", multiline = "NO", prefer_proj = prefer_proj) :
  Discarded datum unknown in Proj4 definition
> summary(rg)
Object of class SpatialGridDataFrame
Coordinates:
     min max
[1,]   0   1
[2,]   0   1
Is projected: TRUE
proj4string :
[+proj=utm +zone=45 +ellps=WGS84 +units=m +no_defs]
Grid attributes:
  cellcentre.offset cellsize cells.dim
1               0.5        1         1
2               0.5        1         1
Data attributes:
   rast_grass
 Min.   : NA
 1st Qu.: NA
 Median : NA
 Mean   :NaN
 3rd Qu.: NA
 Max.   : NA
 NA's   :1
> rast
class      : RasterLayer
dimensions : 20, 40, 800  (nrow, ncol, ncell)
resolution : 10, 10  (x, y)
extent     : 333004, 333404, 3060402, 3060602  (xmin, xmax, ymin, ymax)
crs        : +proj=utm +zone=45 +datum=WGS84 +units=m +no_defs
source     : memory
names      : nakkhu_hosp_10m
values     : 1207.154, 1223.907  (min, max)

> summary(rast_sgdf)
Object of class SpatialGridDataFrame
Coordinates:
       min     max
s1  333004  333404
s2 3060402 3060602
Is projected: TRUE
proj4string :
[+proj=utm +zone=45 +datum=WGS84 +units=m +no_defs]
Grid attributes:
   cellcentre.offset cellsize cells.dim
s1            333009       10        40
s2           3060407       10        20
Data attributes:
 nakkhu_hosp_10m
 Min.   :1207
 1st Qu.:1209
 Median :1210
 Mean   :1211
 3rd Qu.:1213
 Max.   :1224
 NA's   :594
Warning message:
In wkt(obj) : CRS object has no comment
> execGRASS("r.info", map="rast_grass")
 +----------------------------------------------------------------------------+
 | Map:      rast_grass                     Date: Wed Sep 29 10:25:32 2021    |
 | Mapset:   PERMANENT                      Login of Creator: RB              |
 | Location: file2f0444971c68                                                 |
 | DataBase: C:/Users/RB/AppData/Local/Temp/RtmpUzAsbn                        |
 | Title:                                                                     |
 | Timestamp: none                                                            |
 |----------------------------------------------------------------------------|
 |                                                                            |
 |   Type of Map:  raster               Number of Categories: 0               |
 |   Data Type:    DCELL                                                      |
 |   Rows:         20                                                         |
 |   Columns:      40                                                         |
 |   Total Cells:  800                                                        |
 |        Projection: UTM (zone 45)                                           |
 |            N: 3060602.4038    S: 3060402.4038   Res:    10                 |
 |            E: 333403.9801    W: 333003.9801   Res:    10                   |
 |   Range of data:    min = 1207.15441894531  max = 1223.90710449219         |
 |                                                                            |
 |   Data Description:                                                        |
 |    generated by RINBIN~1                                                   |
 |                                                                            |
 |   Comments:                                                                |
 |    RINBIN~1 --overwrite -d input="C:/Users/RB/AppData/Local/Temp/RtmpUz\   |
 |    Asbn/file2f0444971c68/PERMANENT/.tmp/unknown/X916" output="rast_gras\   |
 |    s" bytes=8 header=0 bands=1 order="native" north=3060602.4038 south=\   |
 |    3060402.4038 east=333403.9801 west=333003.9801 rows=20 cols=40 anull\   |
 |    =1206                                                                   |
 |                                                                            |
 +----------------------------------------------------------------------------+
> execGRASS("r.stats", flags="A", input="rast_grass")
 100%
*
> execGRASS("g.region", flags="p")
projection: 1 (UTM)
zone:       45
datum:      wgs84
ellipsoid:  wgs84
north:      1
south:      0
west:       0
east:       1
nsres:      1
ewres:      1
rows:       1
cols:       1
cells:      1
> execGRASS("g.region", raster="rast_grass")
> execGRASS("g.region", flags="p")
projection: 1 (UTM)
zone:       45
datum:      wgs84
ellipsoid:  wgs84
north:      3060602.4038
south:      3060402.4038
west:       333003.9801
east:       333403.9801
nsres:      10
ewres:      10
rows:       20
cols:       40
cells:      800
>  rg <- readRAST("rast_grass")
ERROR 1: PROJ: proj_as_wkt: SQLite error on SELECT name, ellipsoid_auth_name, ellipsoid_code, prime_meridian_auth_name, prime_meridian_code, area_of_use_auth_name, area_of_use_code, publication_date, deprecated FROM geodetic_datum WHERE auth_name = ? AND code = ?: no such column: area_of_use_auth_name
Creating BIL support files...
Exporting raster as floating values (bytes=8)
 100%
Warning message:
In showSRID(uprojargs, format = "PROJ", multiline = "NO", prefer_proj = prefer_proj) :
  Discarded datum unknown in Proj4 definition
> summary(rg)
Object of class SpatialGridDataFrame
Coordinates:
         min     max
[1,]  333004  333404
[2,] 3060402 3060602
Is projected: TRUE
proj4string :
[+proj=utm +zone=45 +ellps=WGS84 +units=m +no_defs]
Grid attributes:
  cellcentre.offset cellsize cells.dim
1            333009       10        40
2           3060407       10        20
Data attributes:
   rast_grass
 Min.   :1207
 1st Qu.:1209
 Median :1210
 Mean   :1211
 3rd Qu.:1213
 Max.   :1224
 NA's   :594

In summary, the further problem beyond those resolved in April, was that recent CRAN rgdal binaries bundle a more recent and incompatible version of proj.db than GRASS Windows stand-alone 7.8.5. The error message comes from PROJ called by g.proj and is not captured by stderr, so will be shown when a version mismatch is present.

library(raster)
load("rast.rda")
rast
rast_sgdf <- as(rast, "SpatialGridDataFrame")
td <- tempdir()
SG <- as(rast_sgdf, "SpatialGrid") #set region on initiation
library(rgrass7)
initGRASS(gisBase="C:/Program Files/GRASS GIS 7.8", home=td, SG=SG, mapset = "PERMANENT", override=TRUE)
execGRASS("g.proj", flags = "c", proj4 = slot(crs(rast), "projargs"))
writeRAST(rast_sgdf, vname="rast_grass", flags = "overwrite")
execGRASS("g.region", raster="rast_grass")
rg <- readRAST("rast_grass")
(rrg <- as(rg, "RasterLayer"))
all.equal(rast, rrg)
manomuthus commented 3 years ago

Many thanks. It works now with rgdal 1.5-27.

rsbivand commented 3 years ago

source 0.2-6 released on CRAN, macOS and Windows binaries to follow.