Closed olenko closed 2 years ago
rsbivand
commented
2 years ago The environment variables will not help. rgdal and terra for Windows ship with their own copies (identical) of the PROJ shared files. Thw Windows GRASS standalone ships with its copy (which may not be identical), and OSeo4W ships with its copy (which may not be identical). Are you using the GRASS stansalone or OSGeo4W?
olenko
commented
2 years ago I am using the GRASS standalone, but I also installed GRASS with OSGeo4W to check and got the same error. Also, a few students sent me the same error.
rsbivand
commented
2 years ago OK, I will try to replicate later today I hope; this is similar to: https://github.com/rsbivand/rgrass/issues/32, https://github.com/rsbivand/rgrass/issues/31, https://github.com/rsbivand/rgrass/issues/27. Please report the output of print.gmeta(). Which PROJ version is GRASS using?
rsbivand
commented
2 years ago I can replicate the problem, which seems to be the two parallel PROJ installations interfering with each other:
> elev <- readRAST("elev_state_500m")
Error : XML content does not seem to be XML: 'ains DATABASE.LAYOUT.VERSION.MINOR = 0 whereas a number >= 2 is expected. It comes from another PROJ installation.'
In addition: Warning messages:
1: Package rgrass7 is transitioning to package rgrass for GRASS GIS 8.
'readRAST' is deprecated. Use 'read_RAST' instead.
2: In system(cmd0, intern = TRUE) :
running command 'g.proj.exe --interface-description' had status 5
Error in parseGRASS(cmd, legacyExec = legacyExec) : g.proj not parsed
> gmeta()
Error : XML content does not seem to be XML: 'ains DATABASE.LAYOUT.VERSION.MINOR = 0 whereas a number >= 2 is expected. It comes from another PROJ installation.'
In addition: Warning message:
In system(cmd0, intern = TRUE) :
running command 'g.proj.exe --interface-description' had status 5
Error in parseGRASS(cmd, legacyExec = legacyExec) : g.proj not parsed
> elev <- read_RAST("elev_state_500m")
Error : XML content does not seem to be XML: 'ains DATABASE.LAYOUT.VERSION.MINOR = 0 whereas a number >= 2 is expected. It comes from another PROJ installation.'
In addition: Warning message:
In system(cmd0, intern = TRUE) :
running command 'g.proj.exe --interface-description' had status 5
Error in parseGRASS(cmd, legacyExec = legacyExec) : g.proj not parsed
> elev <- read_RAST("elev_state_500m", return_format="terra")
proj_create: C:/Users/RB/Documents/R/win-library/4.1/rgdal/proj\proj.db contains DATABASE.LAYOUT.VERSION.MINOR = 0 whereas a number >= 2 is expected. It comes from another PROJ installation.
proj_create: no database context specified
ERROR 1: PROJ: proj_create_from_database: C:/Users/RB/Documents/R/win-library/4.1/rgdal/proj\proj.db contains DATABASE.LAYOUT.VERSION.MINOR = 0 whereas a number >= 2 is expected. It comes from another PROJ installation.
Checking GDAL data type and nodata value...
Using GDAL data type <Float64>
Exporting raster data to RRASTER format...
ROUTGD~1 complete. File
<C:\Users\RB\AppData\Local\Temp\grass8-RB-5428\RtmpIP2VBg\file4a0673d4b31.grd>
created.
> There is however a workaround, as elev was read but without a coordinate reference system:
> elev
class : SpatRaster
dimensions : 280, 767, 1 (nrow, ncol, nlyr)
resolution : 1000, 1000 (x, y)
extent : 151000, 918000, 27000, 307000 (xmin, xmax, ymin, ymax)
coord. ref. :
source : file4a0673d4b31.grd
name : elev_state_500m
min value : -6.592731
max value : 1899.827
> library(raster)
> elev1 <- raster(elev)
> elev2 <- as(elev1, "SpatialGridDataFrame")
> summary(elev2)
Object of class SpatialGridDataFrame
Coordinates:
min max
s1 151000 918000
s2 27000 307000
Is projected: NA
proj4string : [NA]
Grid attributes:
cellcentre.offset cellsize cells.dim
s1 151500 1000 767
s2 27500 1000 280
Data attributes:
elev_state_500m
Min. : -6.59
1st Qu.: 14.94
Median : 86.69
Mean : 212.95
3rd Qu.: 249.62
Max. :1899.83
NA's :83108
>rgdal and terra set PROJ_LIB to the versions of PROJ metadata files, especially proj.db they ship with, which may not match those GRASS ships with. Since I can replicate the problem, I can try to debug it.
olenko
commented
2 years ago Thank you very much! It works.
rsbivand
commented
2 years ago For this location (looking at PERMANENT/PROJ_EPSG), this restores the CRS:
> crs(elev1) <- "EPSG:3358"
> crs(elev1)
Coordinate Reference System:
Deprecated Proj.4 representation:
+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=GRS80 +units=m
+no_defs
WKT2 2019 representation:
PROJCRS["NAD83(HARN) / North Carolina",
BASEGEOGCRS["NAD83(HARN)",
DATUM["NAD83 (High Accuracy Reference Network)",
ELLIPSOID["GRS 1980",6378137,298.257222101,
LENGTHUNIT["metre",1]]],
PRIMEM["Greenwich",0,
ANGLEUNIT["degree",0.0174532925199433]],
ID["EPSG",4152]],
CONVERSION["SPCS83 North Carolina zone (meters)",
METHOD["Lambert Conic Conformal (2SP)",
ID["EPSG",9802]],
PARAMETER["Latitude of false origin",33.75,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8821]],
PARAMETER["Longitude of false origin",-79,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8822]],
PARAMETER["Latitude of 1st standard parallel",36.1666666666667,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8823]],
PARAMETER["Latitude of 2nd standard parallel",34.3333333333333,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8824]],
PARAMETER["Easting at false origin",609601.22,
LENGTHUNIT["metre",1],
ID["EPSG",8826]],
PARAMETER["Northing at false origin",0,
LENGTHUNIT["metre",1],
ID["EPSG",8827]]],
CS[Cartesian,2],
AXIS["easting (X)",east,
ORDER[1],
LENGTHUNIT["metre",1]],
AXIS["northing (Y)",north,
ORDER[2],
LENGTHUNIT["metre",1]],
USAGE[
SCOPE["Engineering survey, topographic mapping."],
AREA["United States (USA) - North Carolina - counties of Alamance; Alexander; Alleghany; Anson; Ashe; Avery; Beaufort; Bertie; Bladen; Brunswick; Buncombe; Burke; Cabarrus; Caldwell; Camden; Carteret; Caswell; Catawba; Chatham; Cherokee; Chowan; Clay; Cleveland; Columbus; Craven; Cumberland; Currituck; Dare; Davidson; Davie; Duplin; Durham; Edgecombe; Forsyth; Franklin; Gaston; Gates; Graham; Granville; Greene; Guilford; Halifax; Harnett; Haywood; Henderson; Hertford; Hoke; Hyde; Iredell; Jackson; Johnston; Jones; Lee; Lenoir; Lincoln; Macon; Madison; Martin; McDowell; Mecklenburg; Mitchell; Montgomery; Moore; Nash; New Hanover; Northampton; Onslow; Orange; Pamlico; Pasquotank; Pender; Perquimans; Person; Pitt; Polk; Randolph; Richmond; Robeson; Rockingham; Rowan; Rutherford; Sampson; Scotland; Stanly; Stokes; Surry; Swain; Transylvania; Tyrrell; Union; Vance; Wake; Warren; Washington; Watauga; Wayne; Wilkes; Wilson; Yadkin; Yancey."],
BBOX[33.83,-84.33,36.59,-75.38]],
ID["EPSG",3358]]
> elev2 <- as(elev1, "SpatialGridDataFrame")
> cat(wkt(elev2), "\n")
PROJCRS["NAD83(HARN) / North Carolina",
BASEGEOGCRS["NAD83(HARN)",
DATUM["NAD83 (High Accuracy Reference Network)",
ELLIPSOID["GRS 1980",6378137,298.257222101,
LENGTHUNIT["metre",1]]],
PRIMEM["Greenwich",0,
ANGLEUNIT["degree",0.0174532925199433]],
ID["EPSG",4152]],
CONVERSION["SPCS83 North Carolina zone (meters)",
METHOD["Lambert Conic Conformal (2SP)",
ID["EPSG",9802]],
PARAMETER["Latitude of false origin",33.75,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8821]],
PARAMETER["Longitude of false origin",-79,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8822]],
PARAMETER["Latitude of 1st standard parallel",36.1666666666667,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8823]],
PARAMETER["Latitude of 2nd standard parallel",34.3333333333333,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8824]],
PARAMETER["Easting at false origin",609601.22,
LENGTHUNIT["metre",1],
ID["EPSG",8826]],
PARAMETER["Northing at false origin",0,
LENGTHUNIT["metre",1],
ID["EPSG",8827]]],
CS[Cartesian,2],
AXIS["easting (X)",east,
ORDER[1],
LENGTHUNIT["metre",1]],
AXIS["northing (Y)",north,
ORDER[2],
LENGTHUNIT["metre",1]],
USAGE[
SCOPE["Engineering survey, topographic mapping."],
AREA["United States (USA) - North Carolina - counties of Alamance; Alexander; Alleghany; Anson; Ashe; Avery; Beaufort; Bertie; Bladen; Brunswick; Buncombe; Burke; Cabarrus; Caldwell; Camden; Carteret; Caswell; Catawba; Chatham; Cherokee; Chowan; Clay; Cleveland; Columbus; Craven; Cumberland; Currituck; Dare; Davidson; Davie; Duplin; Durham; Edgecombe; Forsyth; Franklin; Gaston; Gates; Graham; Granville; Greene; Guilford; Halifax; Harnett; Haywood; Henderson; Hertford; Hoke; Hyde; Iredell; Jackson; Johnston; Jones; Lee; Lenoir; Lincoln; Macon; Madison; Martin; McDowell; Mecklenburg; Mitchell; Montgomery; Moore; Nash; New Hanover; Northampton; Onslow; Orange; Pamlico; Pasquotank; Pender; Perquimans; Person; Pitt; Polk; Randolph; Richmond; Robeson; Rockingham; Rowan; Rutherford; Sampson; Scotland; Stanly; Stokes; Surry; Swain; Transylvania; Tyrrell; Union; Vance; Wake; Warren; Washington; Watauga; Wayne; Wilkes; Wilson; Yadkin; Yancey."],
BBOX[33.83,-84.33,36.59,-75.38]],
ID["EPSG",3358]]
> Note to self: this is possibly also an rgdal issue when more than one PROJ thread is active. Unknown whether it spills over into terra.
rsbivand
commented
2 years ago And going forward, terra will be an effective resolution:
> library(rgrass7)
Loading required package: XML
GRASS GIS interface loaded with GRASS version: GRASS 8.0.1 (2022)
and location: nc_spm_08_grass7
> library(terra)
terra 1.5.21
> elev <- read_RAST("elev_state_500m", return_format="terra")
Checking GDAL data type and nodata value...
Using GDAL data type <Float64>
Exporting raster data to RRASTER format...
ROUTGD~1 complete. File
<C:\Users\RB\AppData\Local\Temp\grass8-RB-12436\RtmpwJegRI\file206028f634cf.grd>
created.
> inMemory(elev)
[1] FALSE
> elev
class : SpatRaster
dimensions : 280, 767, 1 (nrow, ncol, nlyr)
resolution : 1000, 1000 (x, y)
extent : 151000, 918000, 27000, 307000 (xmin, xmax, ymin, ymax)
coord. ref. : +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=GRS80 +units=m +no_defs
source : file206028f634cf.grd
name : elev_state_500m
min value : -6.592731
max value : 1899.827
> Sys.getenv("PROJ_LIB")
[1] "C:\\Program Files\\GRASS GIS 8.0\\share\\proj"
> elev1 <- raster::raster(elev)
> Sys.getenv("PROJ_LIB")
[1] "C:/Users/RB/Documents/R/win-library/4.1/rgdal/proj"where we see the location of PROJ_LIB change as sp and rgdal are loaded by raster. If you can convert your course materials to use terra rather than sp which uses rgdal to handle CRS, the problem should not occur, because PROJ_LIB is not set to the PROJ files bundled with rgdal or sf. I'll explore whether I can adapt rgdal to avoid resetting PROJ_LIB. @VLucet @florisvdh @hellik : is this an extra reason to prefer terra at least until I fix rgdal (and sf) if I can?
rsbivand
commented
2 years ago @olenko please try install.packages("rgdal", repos="http://R-Forge.R-project.org") to install a trial version of rgdal avoiding overwriting PROJ_LIB, which should restore your previous workflow.
> install.packages("rgdal", repos="http://R-Forge.R-project.org")
Installing package into ‘C:/Users/RB/Documents/R/win-library/4.1’
(as ‘lib’ is unspecified)
trying URL 'http://R-Forge.R-project.org/bin/windows/contrib/4.1/rgdal_1.5-30.zip'
Content type 'application/zip' length 25874825 bytes (24.7 MB)
downloaded 24.7 MB
package ‘rgdal’ successfully unpacked and MD5 sums checked
The downloaded binary packages are in
C:\Users\RB\AppData\Local\Temp\grass8-RB-3176\RtmpEPFCSb\downloaded_packages
> library(rgrass7)
Loading required package: XML
GRASS GIS interface loaded with GRASS version: GRASS 8.0.1 (2022)
and location: nc_spm_08_grass7
> use_sp()
>
> elev <- readRAST("elev_state_500m", ignore.stderr=TRUE)
Warning messages:
1: Package rgrass7 is transitioning to package rgrass for GRASS GIS 8.
'readRAST' is deprecated. Use 'read_RAST' instead.
2: In system(syscmd, intern = intern, ignore.stderr = ignore.stderr, :
running command 'g.proj.exe -w' had status 884
3: In showSRID(uprojargs, format = "PROJ", multiline = "NO", prefer_proj = prefer_proj) :
Discarded datum Unknown based on GRS80 ellipsoid in Proj4 definition
> library(sp)
> cat(wkt(elev), "\n")
PROJCRS["unknown",
BASEGEOGCRS["unknown",
DATUM["Unknown based on GRS80 ellipsoid",
ELLIPSOID["GRS 1980",6378137,298.257222101,
LENGTHUNIT["metre",1],
ID["EPSG",7019]]],
PRIMEM["Greenwich",0,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8901]]],
CONVERSION["unknown",
METHOD["Lambert Conic Conformal (2SP)",
ID["EPSG",9802]],
PARAMETER["Latitude of false origin",33.75,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8821]],
PARAMETER["Longitude of false origin",-79,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8822]],
PARAMETER["Latitude of 1st standard parallel",36.1666666666667,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8823]],
PARAMETER["Latitude of 2nd standard parallel",34.3333333333333,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8824]],
PARAMETER["Easting at false origin",609601.22,
LENGTHUNIT["metre",1],
ID["EPSG",8826]],
PARAMETER["Northing at false origin",0,
LENGTHUNIT["metre",1],
ID["EPSG",8827]]],
CS[Cartesian,2],
AXIS["(E)",east,
ORDER[1],
LENGTHUNIT["metre",1,
ID["EPSG",9001]]],
AXIS["(N)",north,
ORDER[2],
LENGTHUNIT["metre",1,
ID["EPSG",9001]]]]
> Sys.getenv("PROJ_LIB")
[1] "C:\\Program Files\\GRASS GIS 8.0\\share\\proj"
> 
florisvdh
commented
2 years ago I'll explore whether I can adapt rgdal to avoid resetting
PROJ_LIB. @VLucet @florisvdh @hellik : is this an extra reason to prefer terra at least until I fix rgdal (and sf) if I can?
Since I don't know about all pros and cons, I have some questions. Does it work with terra because terra doesn't overwrite PROJ_LIB if it is set already (by rgrass)? Does that aspect differ from sf and stars (which don't use rgdal either, just like terra)? IIUC this is a Windows specific aspect (in Linux, it seems that PROJ_LIB remains at ""). Further, is the aim to support sp + rgdal in rgrass only during a transition stage, or to continue their support on the long term?
rsbivand
commented
2 years ago Windows and macOS (all CRAN binaries that use bundled share/proj and share/gdal files). rgdal was updated for PROJ 6 first, to assign the path in the PROJ default context. sf and terra followed. rgdal had limited the use of the PROJ default context to PROJ 6, because at that time it was not known how the PROJ CDN would work. sf and terra did not restrict the adaptation to PROJ 6. rgdal's r1169 removes this limitation, so from this revision rgdal, sf and terra all avoid resetting the PROJ_LIB environment variable. I still need to cross-check on macOS. And no, the target remains terra for raster and vector data transfer for rgrass; changing rgdal is just fixing an omission.
florisvdh
commented
2 years ago Thank you for the clarifications. I think it's good to support terra since it handles both vector and raster data; I'm not familiar enough with stars to make a good comparison from the rgrass viewpoint. IIRC sf is supported already with use_sf().
rsbivand
commented
2 years ago It looks as though macOS is OK - most likely the GRASS binary bundle does not use PROJ_LIB as it is not set when running R inside GRASS, so no overwriting even with the present rgdal version. So the problem reported in this issue perhaps only concerned Windows, both GRASS standalone and OSGeo4W (doubt because GRASS may have been running under Rosetta, R in arm64).
rsbivand
commented
2 years ago @olenko : may I submit the fixed version of rgdal to CRAN?
olenko
commented
2 years ago @rsbivand Yes, it works as before. Thank you very much!
rsbivand
commented
2 years ago rgdal_1.5-30.tar.gz on CRAN!
Hello Roger I updated GRASS (now WinGRASS-8.0.1-1) and R (4.1.3) and all packages.
Most of the commands that I used in my lectures last year work fine, except readRAST and writeRAST.
For example I run in the GRASS console: g.region vect=precip_30ynormals res=1000 -ap Rgui
Then I run in R:
rgdal: version: 1.5-29, (SVN revision 1165M) Geospatial Data Abstraction Library extensions to R successfully loaded Loaded GDAL runtime: GDAL 3.2.1, released 2020/12/29 Path to GDAL shared files: C:/Users/aolenko/Documents/R/win-library/4.1/rgdal/gdal GDAL binary built with GEOS: TRUE Loaded PROJ runtime: Rel. 7.2.1, January 1st, 2021, [PJ_VERSION: 721] Path to PROJ shared files: C:/Users/aolenko/Documents/R/win-library/4.1/rgdal/proj PROJ CDN enabled: FALSE Linking to sp version:1.4-6 To mute warnings of possible GDAL/OSR exportToProj4() degradation, use options("rgdal_show_exportToProj4_warnings"="none") before loading sp or rgdal. Overwritten PROJ_LIB was C:/Users/aolenko/Documents/R/win-library/4.1/rgdal/proj
Similar XML error for readRAST.
I tried to change environmental variables as advised here: https://github.com/OSGeo/grass/blob/main/mswindows/env.bat but it did not solve the problem.
Could you please help to fix the issue. Thank you in advance, Andriy