[R-sig-Geo] [R] Need help using GRASS within R - error when running the script for a second time

Loïc Valéry |v@|ery @end|ng |rom out|ook@|r
Tue Jul 7 21:03:59 CEST 2020


Dear Roger,

Sorry for cross-posting.
As you invite me, I'm going to post my case on https://github.com/rsbivand/rgrass7/issues.

Just a short comment to yours: if the reason you give is the right one, I don't understand why the script works the first time but not the second time. In your hypothesis, I guess it shouldn't work at all.

Best regards,
Loïc

De : Roger Bivand <Roger.Bivand using nhh.no>
Envoyé : mardi 7 juillet 2020 20:20
À : Loïc Valéry <lvalery using outlook.fr>
Cc : Micha Silver <tsvibar using gmail.com>; r-sig-geo using r-project.org <r-sig-geo using r-project.org>; grass-user <grass-user using lists.osgeo.org>
Objet : Re: [R-sig-Geo] [R] Need help using GRASS within R - error when running the script for a second time 
 
Please do not cross-post between lists, because few are likely to 
subscribe to both. Choose rather to post only to the list likely to 
include people using R and GRASS:

http://lists.osgeo.org/mailman/listinfo/grass-stats

Further, please post plain text only - I missed your earlier post to 
grass-user because it was scrubbed as HTML (or maybe Micha's reply).

Once resolved, you might report back on grass-user and R-sig-geo on the 
findings.

Briefly, it looks as though you are caught in the transition in sp/rgdal 
and sf from Proj4 to WKT2-2019 CRS, and similar changes in Grass. Grass 
7.8 and R packages may be installed with old PROJ/GDAL or new PROJ/GDAL. 
So far, nobody knows how to use g.proj and siblings sensibly in the R 
context. Most likely, using the wkt= parameter when the object in R has a 
WKT2 representation is going to be much more secure than passing a proj4= 
string, but wkt= takes a file name or standard input.

I would invite discussion in an issue on:

https://github.com/rsbivand/rgrass7/issues

with use cases (platform-independent).

Roger

On Tue, 7 Jul 2020, Loïc Valéry wrote:

> Dear all,
>
> As advised by Micha (FYI his reply is at the end of this email) , I 
> contact you to submit an inextricable problem when using rgrass7 : the 
> script (cf. just below) works the first time but, when I rerun it a 
> second time, R and windows return error messages. For this to work, I 
> must close and open R or use the rs.restartR() command from RStudio, 
> which is not very convenient because it stops the script.
>
> Many thanks in advance for your help.
> Loïc

> @Micha :  thank you very much for your reply. But I performed the tests 
> you asked me to do and, unfortunately, none of them solve the problem. 
> Also, I tried to run "execGRASS("g.proj", flags = "c", proj4 = p4str)" a 
> second time WITHOUT reinitializing GRASS but R and windows still return 
> the same error messages.
>
> ERRORS WHEN I RUN THE SCRIPT FOR A SECOND TIME :
>
>      - FROM WINDOWS : g.region.exe : the procedure entry point GEOSMakeValid_r could not be located in the dynamic link library C:\Program Files\GRASS GIS 7.8\extrabin\gdal300.dll
>
>      - FROM R :
> # Link to GRASS GIS software v.7.8.3
>> initGRASS(gisBase ="C:/Program Files/GRASS GIS 7.8",
> +           home="temp/GRASS_tmp", use_g.dirseps.exe=F,
> +           gisDbase="temp/GRASS_tmp", mapset="PERMANENT",
> +           remove_GISRC=T, override=T)
> Error in if (is.na(projstr)) uprojargs <- projstr else uprojargs <- paste(unique(unlist(strsplit(projstr,  :
>  l'argument est de longueur nulle
> De plus : There were 50 or more warnings (use warnings() to see the first 50)
>>
>> # Specifying the projection reference for the GRASS working environment
>> p4str<-sp::proj4string(seg_poly)
> Warning message:
> In sp::proj4string(seg_poly) :
>  CRS object has comment, which is lost in output
>>
>> execGRASS("g.proj", flags = "c", proj4 = p4str)
>> rgrass7::use_sp()
>>
>> # Converting the 'sp' object into a GRASS-readable file format (here, 'vec1.shp')
>> writeVECT(seg_poly,"vec1",v.in.ogr_flags=c("o", "overwrite"),
> +           driver="ESRI Shapefile")
> Error in writeVECT(seg_poly, "vec1", v.in.ogr_flags = c("o", "overwrite"),  :
>  driver %in% candDrivers is not TRUE
> De plus : Warning message:
> In system(syscmd, intern = intern, ignore.stderr = ignore.stderr,  :
>  l'exécution de la commande 'v.in.ogr.exe -f' renvoie un statut 313
>>
>>
>> # Run the 'v.generalize' function
>> execGRASS("v.generalize",flag=c("overwrite"),
> +           parameters=list(input="vec1",
> +                           output="GRASS_smooth_seg_poly",
> +                           error="GRASS_smooth_seg_poly_error",
> +                           method="distance_weighting",
> +                           threshold=1))
>>
>>
>> # Converting the shapefile 'GRASS_smooth_seg_poly.shp' into a R-readable object
>> # (here, a SpatialPolygonDataFrame named 'smooth_seg_poly')
>> smooth_seg_poly<-readVECT("GRASS_smooth_seg_poly", with_prj=T,
> +                           driver="ESRI Shapefile")
> Error in .read_vect_non_plugin(vname = vname, layer = layer, type = type,  :
>  driver %in% candDrivers is not TRUE
> De plus : Warning message:
> In system(syscmd, intern = intern, ignore.stderr = ignore.stderr,  :
>  l'exécution de la commande 'v.in.ogr.exe -f' renvoie un statut 313
>
>
>
> REPLY FROM MICHA :
> De : Micha Silver <tsvibar using gmail.com>
> Envoyé : mardi 7 juillet 2020 15:55
> À : Loïc Valéry <lvalery using outlook.fr>
> Cc : grass-user <grass-user using lists.osgeo.org>
> Objet : Re: [R] Need help using GRASS within R - error when running the script for a second time
>  
> Hello Loïc
>
> Let's keep the discussion on the list, so others can help and benefit from responses.
>
> BTW, you might want to post to the r-sig-geo list also. If you come across questions specific to R and spatial functions, that's the best place to ask.
>
> Below,inline, some guesses as to what is happening...
>
> On 07/07/2020 0:03, Loïc Valéry wrote:
> initGRASS(gisBase ="C:/Program Files/GRASS GIS 7.8",
>
> +           home="temp/GRASS_tmp", use_g.dirseps.exe=F,
> +           gisDbase="temp/GRASS_tmp", mapset="PERMANENT",
> +           remove_GISRC=T, override=T)
> Error in if (is.na(projstr)) uprojargs <- projstr else uprojargs <- paste(unique(unlist(strsplit(projstr,  :
>  l'argument est de longueur nulle
> De plus : There were 50 or more warnings (use warnings() to see the first 50)
>
> I'm not sure, but I don't think you want to rerun initGRASS() a second time. In your script, you do initGRASS only once, then if you need to perform several analyses, do them in the same (running) GRASS session.
>
>
> # Specifying the projection reference for the GRASS working environment
> p4str<-sp::proj4string(seg_poly)
>
> Warning message:
> In sp::proj4string(seg_poly) :
>  CRS object has comment, which is lost in output
>
> execGRASS("g.proj", flags = "c", proj4 = p4str)
> rgrass7::use_sp()
>
> # Converting the 'sp' object into a GRASS-readable file format (here, 'vec1.shp')
> writeVECT(seg_poly,"vec1",v.in.ogr_flags=c("o", "overwrite"),
>
> +           driver="ESRI Shapefile")
>
> The driver name is "ESRI_Shapefile" (note the underscore). Perhaps that is the problem here?
>
> Error in writeVECT(seg_poly, "vec1", v.in.ogr_flags = c("o", "overwrite"),  :
>  driver %in% candDrivers is not TRUE
> De plus : Warning message:
> In system(syscmd, intern = intern, ignore.stderr = ignore.stderr,  :
>  l'exécution de la commande 'v.in.ogr.exe -f' renvoie un statut 313
>
> # Run the 'v.generalize' function
> execGRASS("v.generalize",flag=c("overwrite"),
>
> +           parameters=list(input="vec1",
> +                           output="GRASS_smooth_seg_poly",
> +                           error="GRASS_smooth_seg_poly_error",
> +                           method="distance_weighting",
> +                           threshold=1))
>
> # Converting the shapefile 'GRASS_smooth_seg_poly.shp' into a R-readable object
> # (here, a SpatialPolygonDataFrame named 'smooth_seg_poly')
> smooth_seg_poly<-readVECT("GRASS_smooth_seg_poly", with_prj=T,
>
> +                           driver="ESRI Shapefile")
>
> Same here, the driver name for shapefile is wrong.
>
> I would also add that both R and GRASS have switched to the better Geopackage format for exchanging spatial data. You might want to try, the driver name is "GPKG"
>
> See details about Geopackage here:
> https://www.gis-blog.com/geopackage-vs-shapefile/
> and in French here:
> https://www.sigterritoires.fr/index.php/le-format-geopackage-et-qgis-3/
>
> Best, Micha
>
>
>
>
>
>
> PREVIOUS EMAILS:
>
> De : Micha Silver <tsvibar using gmail.com>
> Envoyé : mardi 9 juin 2020 23:08
> À : Loïc Valéry <lvalery using outlook.fr>
> Objet : Re: [R] Need help using GRASS within R - problem with CRS using the 'v.generalize' command
>  
>
> On 09/06/2020 16:22, Loïc Valéry wrote:
> Dear Micha,
>
> A thousand thanks for your help. It works now!
> Whoa, you'd better save some thanks for the next question. ;-)
>
> For your information and following your remark, I also updated rgdal.
> Thanks again for taking the time to answer me.
> Kind regards
> Loïc
>
> De : Micha Silver <tsvibar using gmail.com>
> Envoyé : mardi 9 juin 2020 08:44
> À : Loïc Valéry <lvalery using outlook.fr>; r-help using R-project.org <r-help using R-project.org>
> Objet : Re: [R] Need help using GRASS within R - problem with CRS using the 'v.generalize' command
>  
> Hello
>
>
> On 08/06/2020 19:52, Loïc Valéry wrote:
>
> Dear all,
>
> First of all, this is my first message on the list. Therefore, please be indulgent if my message is not perfectly formatted as it should be.
>
> I am currently encountering a difficulty with GRASS 7.8 within R when using the 'v.generalize' command to smooth the contour of polygons after a segmentation step.
>
> I tried two different ways to "call" GRASS:
>
>             1 - using the RQGIS3 package
>
> setwd("D:/test")
> library(sp)
> library(rgdal)
>
> rgdal: version: 1.4-8, (SVN revision 845)
> Geospatial Data Abstraction Library extensions to R successfully loaded
> Loaded GDAL runtime: GDAL 2.2.3, released 2017/11/20
> Path to GDAL shared files: C:/Users/toto/Documents/R/win-library/3.6/rgdal/gdal
> GDAL binary built with GEOS: TRUE
> Loaded PROJ.4 runtime: Rel. 4.9.3, 15 August 2016, [PJ_VERSION: 493]
> Path to PROJ.4 shared files: C:/Users/toto/Documents/R/win-library/3.6/rgdal/proj
> Linking to sp version: 1.4-1
>
> Please note that you have older versions of GDAL (2.2.3 here) and PROJ.4
> (4.9 here). These are currently being replaced by GDAL 3.0 and PROJ 6.3.
> You might (should) want to follow the discussion the the r-sig-geo maillist:
>
>
> https://stat.ethz.ch/pipermail/r-sig-geo/2020-June/028165.html
>
>
> ....... (skipped all the discussion regarding RQGIS3 as I think it's not
> relevant)
>
> setwd("D:/test")
> library(rgrass7)
>
> # characteristics of the SpatialPolygonsDataFrame 'seg_poly' : CRS does exist
> seg_poly
>
> class       : SpatialPolygonsDataFrame
> features    : 31
> extent      : 477371.3, 477397.6, 5631995, 5632020  (xmin, xmax, ymin, ymax)
> crs         : +proj=utm +zone=32 +datum=WGS84 +units=m +no_defs
> variables   : 1
> names       : Seg_ID
> min values  :      1
> max values  :     31
> Warning message:
> In proj4string(x) : CRS object has comment, which is lost in output
>
> # initialization of GRASS 7.8 from R
> initGRASS(gisBase ="C:/Program Files/GRASS GIS 7.8", home="temp/GRASS",gisDbase="temp/GRASS", use_g.dirseps.exe=F,remove_GISRC=T, override=T)
>
> gisdbase    temp/GRASS
> location    file19685026c56
> mapset      file196829fa7141
> rows        1
> columns     1
> north       1
> south       0
> west        0
> east        1
> nsres       1
> ewres       1
> projection  NA
>
> I suspect that the problem comes from 'projection NA' when initializing GRASS (cf. just above)
>
> What you need to do here is setup the CRS of your new location.
> Typically, you would run initGRASS and point to a *previously created* 
> LOCATION, with CRS already defined. In this case, since you are creating
> a new location, you must define it's coordinate system. (GRASS is very
> "picky" about that).
>
> Here's a GIS Stackexchange post that explains:
>
> https://gis.stackexchange.com/questions/183032/create-a-new-grass-database-in-r-with-crs-projection
>
>
> You can derive the full proj4 string from your sp object with:
>
>
> p4str = sp::proj4string(seg_poly)
>
>
> Then use that to set the project parameters for the new LOCATION, with
>
> execGRASS("g.proj", flags = "c", proj4 = p4str) Now you should be able
> to continue with...
>
>
>
> execGRASS("v.generalize",flag=c("overwrite"),parameters=list(input="vec1",
>
> +                                                              output="GRASS_smooth_seg_poly",
> +                                                              error="GRASS_smooth_seg_poly_error",
> +                                                              method="distance_weighting",
> +                                                              threshold=1))
>
>
>
> As a novice in the field, I would be very grateful for your help.
> I remain at your entire disposal for any further information you may need to help me in finding a solution to this problem.
>
> Yours sincerely,
> Loïc
> ______________________________________________
> R-help using r-project.org mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
>
>

-- 
Roger Bivand
Department of Economics, Norwegian School of Economics,
Helleveien 30, N-5045 Bergen, Norway.
voice: +47 55 95 93 55; e-mail: Roger.Bivand using nhh.no
https://orcid.org/0000-0003-2392-6140
https://scholar.google.no/citations?user=AWeghB0AAAAJ&hl=en


More information about the R-sig-Geo mailing list