[R-sig-Geo] spgrass6, readVECT and overwriting ?
Rainer M Krug
Rainer at krugs.de
Wed Oct 28 16:16:02 CET 2015
Mathieu Rajerison <mathieu.rajerison at gmail.com> writes:
> Before was the code that works.
>
> Below is the code that generates issues :
>
> (sorry the multiple posts)
>
> for (id in 1:length(grids)) {
>
> (...)
>
> ## SMOOTH
> print(">> SMOOTH")
> smoothed = smooth(sieved, 9)
> values(smoothed)[values(smoothed) > 0] = 1
> values(smoothed)[values(smoothed)==0] = NA
>
> ## THIN
> # WRITE TO GRASS
> writeRAST(as(smoothed, "SpatialGridDataFrame"), "bdtopo", zcol="layer",
> flags=c("overwrite"))
> execGRASS("g.region", rast="bdtopo")
> execGRASS("r.thin", parameters=list(input="bdtopo", output="thin"),
> flags=c("overwrite"))
> thin = raster(readRAST("thin"))
>
> ## CLEAN
> execGRASS("r.to.vect", parameters=list(input="thin", output="thin",
> feature="line"), flags=c("overwrite"))
> execGRASS("v.clean", parameters=list(input="thin", tool="rmdangle",
> output="clean", thresh=20), flags=c("overwrite"))
>
> cleaned = readVECT("clean")
>
> ## EXPORT
> writeOGR(cleaned, "OUT/CLEANED", id, "ESRI Shapefile")
> }
OK. So what you are doing is in a loop
write a raster to grass
calculations in GRASS
read result vector layer into R
And the error occurs in
cleaned = readVECT("clean")
as temporary files during the readVECT are not cleaned up and not
overwritten - correct?
Then the rest is irrelevant.
You can reproduce it by doing
for (i in 1:10) {
x <- readVECT("vectorLayerInSampleDataset")
}
Could you try this and see if it works in a sample dataset?
If it does, the error comes from somewhere else.
Under grass 7 it works:
,----
| 03:37:37 {master} ~/ownCompiled/emacs-mac$ grass70
| Cleaning up temporary files...
| Starting GRASS GIS...
|
| __________ ___ __________ _______________
| / ____/ __ \/ | / ___/ ___/ / ____/ _/ ___/
| / / __/ /_/ / /| | \__ \\_ \ / / __ / / \__ \
| / /_/ / _, _/ ___ |___/ /__/ / / /_/ // / ___/ /
| \____/_/ |_/_/ |_/____/____/ \____/___//____/
|
| Welcome to GRASS GIS 7.0.1
| GRASS GIS homepage: http://grass.osgeo.org
| This version running through: Bash Shell (/usr/local/bin/bash)
| Help is available with the command: g.manual -i
| See the licence terms with: g.version -c
| If required, restart the GUI with: g.gui wxpython
| When ready to quit enter: exit
|
| Launching <wxpython> GUI in the background, please wait...
| GRASS 7.0.1 (nc_basic_spm_grass7):~/ownCompiled/emacs-mac > R
|
| R version 3.2.2 (2015-08-14) -- "Fire Safety"
| Copyright (C) 2015 The R Foundation for Statistical Computing
| Platform: x86_64-apple-darwin14.5.0 (64-bit)
|
| R is free software and comes with ABSOLUTELY NO WARRANTY.
| You are welcome to redistribute it under certain conditions.
| Type 'license()' or 'licence()' for distribution details.
|
| Natural language support but running in an English locale
|
| R is a collaborative project with many contributors.
| Type 'contributors()' for more information and
| 'citation()' on how to cite R or R packages in publications.
|
| Type 'demo()' for some demos, 'help()' for on-line help, or
| 'help.start()' for an HTML browser interface to help.
| Type 'q()' to quit R.
|
| > library(rgrass7)
| Loading required package: sp
| Loading required package: XML
| GRASS GIS interface loaded with GRASS version: GRASS 7.0.1 (2015)
| and location: nc_basic_spm_grass7
| > x <- readVECT("streams")
| Exporting 8554 features...
| 100%
| v.out.ogr complete. 8554 features (Line String type) written to <streams>
| (ESRI_Shapefile format).
| OGR data source with driver: ESRI Shapefile
| Source: "/Users/rainerkrug/nc_basic_spm_grass7/PERMANENT/.tmp/Rainers-MacBook-Pro.local", layer: "streams"
| with 8554 features
| It has 14 fields
| > x <- readVECT("streams")
| Exporting 8554 features...
| 100%
| v.out.ogr complete. 8554 features (Line String type) written to <streams>
| (ESRI_Shapefile format).
| OGR data source with driver: ESRI Shapefile
| Source: "/Users/rainerkrug/nc_basic_spm_grass7/PERMANENT/.tmp/Rainers-MacBook-Pro.local", layer: "streams"
| with 8554 features
| It has 14 fields
| >
`----
and also under grass-64:
,----
| > library(spgrass6)
| Loading required package: sp
| Loading required package: XML
| GRASS GIS interface loaded with GRASS version: GRASS 6.4.4 (2014)
| and location: spearfish60
| > x_readVECT("soils")
| Error: could not find function "x_readVECT"
| > x <- readVECT("soils")
| with_c: argument reversed from version 0.7-11 and in GRASS 6
| Exporting 737 areas (may take some time)...
| 100%
| WARNING: 1 features without attributes were written
| v.out.ogr complete. 737 features written to <soils> (ESRI_Shapefile).
| OGR data source with driver: ESRI Shapefile
| Source: "/Users/rainerkrug/spearfish60/user1/.tmp/Rainers-MacBook-Pro.local", layer: "soils"
| with 737 features
| It has 2 fields
| Warning message:
| In execGRASS("v.out.ogr", flags = flags, input = vname, type = type, :
| The command:
| v.out.ogr -e -c input=soils type=area layer=1 dsn=/Users/rainerkrug/spearfish60/user1/.tmp/Rainers-MacBook-Pro.local olayer=soils format=ESRI_Shapefile
| produced at least one warning during execution:
| Exporting 737 areas (may take some time)...
| 100%
| WARNING: 1 features without attributes were written
| v.out.ogr complete. 737 features written to <soils> (ESRI_Shapefile).
| > plot(x)
| > plot(x)
| > x <- readVECT("soils")
| with_c: argument reversed from version 0.7-11 and in GRASS 6
| Exporting 737 areas (may take some time)...
| 100%
| WARNING: 1 features without attributes were written
| v.out.ogr complete. 737 features written to <soils> (ESRI_Shapefile).
| OGR data source with driver: ESRI Shapefile
| Source: "/Users/rainerkrug/spearfish60/user1/.tmp/Rainers-MacBook-Pro.local", layer: "soils"
| with 737 features
| It has 2 fields
| Warning message:
| In execGRASS("v.out.ogr", flags = flags, input = vname, type = type, :
| The command:
| v.out.ogr -e -c input=soils type=area layer=1 dsn=/Users/rainerkrug/spearfish60/user1/.tmp/Rainers-MacBook-Pro.local olayer=soils format=ESRI_Shapefile
| produced at least one warning during execution:
| Exporting 737 areas (may take some time)...
| 100%
| WARNING: 1 features without attributes were written
| v.out.ogr complete. 737 features written to <soils> (ESRI_Shapefile).
| >
`----
Also: if you want to export it to a shape file, why do you import it
into R and save it afterwards? Why not write it directly from GRASS as a
shape file?
Cheers,
Rainer
>
> 2015-10-28 14:18 GMT+01:00 Mathieu Rajerison <mathieu.rajerison at gmail.com>:
>
>> Mr Krug,
>>
>>
>> It is difficult to reproduce an example from the spearfish dataset as it
>> doesn't natively contain data in it close to the one I used.
>>
>> To give you th context, my goal is to detect all blue features : rivers
>> from a topographic map, ie 3-band RGB raster.
>>
>> With RStoolbox, I perform a Spectral Angle Classification of my data using
>> sample points located on (near)blue pixels
>>
>> With GRASS, I then thin, vectorize and clean the result.
>>
>> As the process must be performed on a very large scale, I split it using
>> grids where sample RGB parts of the topo map are extracted with
>> gdal_translate. Everything is executed within a loop. At the end of the
>> loop, I export the resulted cleant shapefile in a folder. At the end, I
>> have as many shapefiles as processing grids.
>>
>> That's why I use temporary vector datasets inside GRASS that are
>> overwritten at each loop instead of storing each individual result. I
>> noticed that in my process, the vector data read in the grass database was
>> not correct : always the first one. Because in fact, it didn't manage to
>> delete the temporary file before writing the new one..
>>
>> The best I can do, is give you part of the code as as you can see the
>> context.
>>
>> Thanks for the code improvements with file.remove and unlink !
>>
>> Cheers,
>>
>> Mathieu
>>
>> for (id in 1:length(grids)) {
>>
>> (...)
>>
>> ## SMOOTH
>> print(">> SMOOTH")
>> smoothed = smooth(sieved, 9)
>> values(smoothed)[values(smoothed) > 0] = 1
>> values(smoothed)[values(smoothed)==0] = NA
>>
>> ## THIN
>> # WRITE TO GRASS
>> writeRAST(as(smoothed, "SpatialGridDataFrame"), "bdtopo", zcol="layer",
>> flags=c("overwrite"))
>> execGRASS("g.region", rast="bdtopo")
>> execGRASS("r.thin", parameters=list(input="bdtopo", output="thin"),
>> flags=c("overwrite"))
>> thin = raster(readRAST("thin"))
>>
>> ## CLEAN
>> execGRASS("r.to.vect", parameters=list(input="thin", output="thin",
>> feature="line"), flags=c("overwrite"))
>> execGRASS("v.clean", parameters=list(input="thin", tool="rmdangle",
>> output="clean", thresh=20), flags=c("overwrite"))
>>
>> # REMOVE FILES
>> G = gmeta()
>> tmpDir = file.path(G$GISDBASE, G$LOCATION_NAME, G$MAPSET, ".tmp")
>> file.remove( list.files(path = tmpDir, Pattern ="^clean.*$", full.names
>> = TRUE))
>>
>> # READ TO SPDF
>> execGRASS("v.out.ogr", parameters=list(input="clean", type="line",
>> dsn=tmpDir, olayer="clean", format="ESRI_Shapefile"))
>> cleaned = readOGR(tmpDir, "clean")
>>
>> ## EXPORT
>> writeOGR(cleaned, "OUT/CLEANED", id, "ESRI Shapefile")
>> }
>>
>> 2015-10-27 14:11 GMT+01:00 Rainer M Krug <Rainer at krugs.de>:
>>
>>> Mathieu Rajerison <mathieu.rajerison at gmail.com> writes:
>>>
>>> > Dear Mr Krug,
>>> >
>>>
>>> I'll Cc the R-sig-geo list in to keep the conversation there for future
>>> reference and for other users who might have the same problem.
>>>
>>> > Sorry in advance for my english.
>>>
>>> Don't worry abut your English - it is perfectly fine.
>>>
>>> >
>>> > My question was not about debugging, but about a potentiel generic
>>> > behaviour of spgrass6.
>>> >
>>> > I should have reformulated my question in a simpler manner : does
>>> spgrass6
>>> > allow overwriting an existing dataset when using readVECT function ?
>>>
>>> The functions readVect / RAST are using temporary files - you are right,
>>> but they should delete these upon exit. So if they do nbot do this, it
>>> is a bug.
>>>
>>> >
>>> > 1) no version information about R and GRASS (there are several GRASS 6
>>> > versions)
>>> > I've just put it at the end of this post.
>>> > 2) no session information i.e. package versions used
>>> > the latest spgrass6
>>> > 3) no information about the OS (I guess it is windows based on the
>>> > Source path above)
>>> > you're right
>>>
>>> OK - thanks.
>>>
>>> > 4) no reproducible example which could be easily done by using the GRASS
>>>
>>> As I said - if you manually have to delete the temporary files, then
>>> there is a bug in the package. Could you please provide a reproducible
>>> example, using one of the example data sets from GRASS, so that I can
>>> look at it?
>>>
>>> > 5) as far as I can see, in the newest version of spgrass6 there is no
>>> > function readVECT()
>>> > Yes, there is, as you will see on page 14 of the documentation.
>>>
>>> Sorry - overlooked it.
>>>
>>> >
>>> > Finally, I solved the problem by using the following hint, "clean" being
>>> > the name of my temporary shapefile dataset :
>>> > l = list.files("D:/GRASSDB/paca/mapset/.tmp", "^clean.*$"); l =
>>> > file.path("D:/GRASSDB/paca/mapset/.tmp", l)
>>> > lapply(l, file.remove)
>>> > execGRASS("v.out.ogr", parameters=list(input="clean", type="line",
>>> > dsn="D:/GRASSDB/paca/mapset/.tmp", olayer="clean",
>>> format="ESRI_Shapefile"))
>>>
>>> OK - reformated, the code looks as follow:
>>>
>>> --8<---------------cut here---------------start------------->8---
>>> l = list.files("D:/GRASSDB/paca/mapset/.tmp", "^clean.*$")
>>> l = file.path("D:/GRASSDB/paca/mapset/.tmp", l)
>>> lapply(l, file.remove)
>>> execGRASS( "v.out.ogr",
>>> parameters=list(input="clean",
>>> type="line",dsn="D:/GRASSDB/paca/mapset/.tmp", olayer="clean",
>>> format="ESRI_Shapefile"))
>>> --8<---------------cut here---------------end--------------->8---
>>>
>>> as a side note, you should be able to do
>>>
>>> file.remove( list.files(path = "D:/GRASSDB/paca/mapset/.tmp", Pattern =
>>> "^clean.*$", full.names = TRUE) )
>>>
>>> or
>>>
>>> unlink("D:/GRASSDB/paca/mapset/.tmp/clean.*")
>>>
>>> instead of the first three commands.
>>>
>>> Now - in your example, you don't use the readVECT() function?
>>>
>>> If you check in the GRASS help for v.out.ogr, you can specify the
>>> overwrite flag and the exported layer will be overwritten.
>>>
>>> But I have the feeling you are doing something in a to complicated way.
>>>
>>> What is it you want to do, and e=what is the problem you have or the
>>> error message you get?
>>>
>>> Cheers,
>>>
>>> Rainer
>>>
>>> >
>>> >
>>> >> R.Version()
>>> > $platform
>>> > [1] "x86_64-w64-mingw32"
>>> > $arch
>>> > [1] "x86_64"
>>> > $os
>>> > [1] "mingw32"
>>> >
>>> > $system
>>> > [1] "x86_64, mingw32"
>>> > (...)
>>> > $version.string
>>> > [1] "R version 3.1.2 (2014-10-31)"
>>> > $nickname
>>> > [1] "Pumpkin Helmet"
>>> >
>>> > I will happily look int this if you can provide the necessary info.
>>> >
>>> > Cheers,
>>> >
>>> > Mathieu
>>> >
>>> > 2015-10-24 13:34 GMT+02:00 Rainer M Krug <Rainer at krugs.de>:
>>> >
>>> >> Mathieu Rajerison <mathieu.rajerison at gmail.com> writes:
>>> >>
>>> >> > Hi,
>>> >> >
>>> >> > I'm using spgrass6 and I use readVECT function in a loop.
>>> >> >
>>> >> > Using it causes an error because it doesn't overwrite the the
>>> temporary
>>> >> > shapefile..
>>> >> >
>>> >> >> thin = readVECT(vname="zoneBDTOPOthin")
>>> >> > with_c: argument reversed from version 0.7-11 and in GRASS 6
>>> >> > ERROR :Layer <zoneBDTO> already exists in OGR data source
>>> >> > 'D:/GRASSDB/paca/mapset/.tmp'
>>> >> > OGR data source with driver: ESRI Shapefile
>>> >> > Source: "D:/GRASSDB/paca/mapset/.tmp", layer: "zoneBDTO"
>>> >> > with 770 features
>>> >> > It has 3 fields
>>> >> > Warning message:
>>> >> > running command 'v.out.ogr.exe -e -c input=zoneBDTOPOthin2 type=line
>>> >> > layer=1 dsn=D:/GRASSDB/paca/mapset/.tmp olayer=zoneBDTO
>>> >> > format=ESRI_Shapefile' had status 1
>>> >> >
>>> >> >
>>> >> > I wonder how to specify to spgrass that I want to overwrite the file
>>> in
>>> >> > readVECT, although I know using execGRASS and v.out.ogr with the
>>> >> > appropriate overwrite flags would do the job.
>>> >> >
>>> >> > Thanks in advance for the answer..
>>> >>
>>> >> One reason why probably nobody replied is that
>>> >>
>>> >> 1) no version information about R and GRASS (there are several GRASS 6
>>> >> versions)
>>> >> 2) no session information i.e. package versions used
>>> >> 3) no information about the OS (I guess it is windows based on the
>>> >> Source path above)
>>> >> 4) no reproducible example which could be easily done by using the
>>> GRASS
>>> >> 5) as far as I can see, in the newest version of spgrass6 there is no
>>> >> function readVECT()
>>> >>
>>> >> So I would suggest to provide the info and give a reproducible example
>>> >> using a GRASS sample dataset
>>> >> (see https://grass.osgeo.org/download/sample-data/ - the Spearfish
>>> >> dataset is probably the best)
>>> >>
>>> >> I will happily look int this if you can provide the necessary info.
>>> >>
>>> >> Cheers,
>>> >>
>>> >> Rainer
>>> >>
>>> >> >
>>> >> > Mathieu
>>> >> >
>>> >> > [[alternative HTML version deleted]]
>>> >> >
>>> >> > _______________________________________________
>>> >> > R-sig-Geo mailing list
>>> >> > R-sig-Geo at r-project.org
>>> >> > https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>> >>
>>> >> --
>>> >> Rainer M. Krug, PhD (Conservation Ecology, SUN), MSc (Conservation
>>> >> Biology, UCT), Dipl. Phys. (Germany)
>>> >>
>>> >> Centre of Excellence for Invasion Biology
>>> >> Stellenbosch University
>>> >> South Africa
>>> >>
>>> >> Tel : +33 - (0)9 53 10 27 44
>>> >> Cell: +33 - (0)6 85 62 59 98
>>> >> Fax : +33 - (0)9 58 10 27 44
>>> >>
>>> >> Fax (D): +49 - (0)3 21 21 25 22 44
>>> >>
>>> >> email: Rainer at krugs.de
>>> >>
>>> >> Skype: RMkrug
>>> >>
>>> >> PGP: 0x0F52F982
>>> >>
>>>
>>> --
>>> Rainer M. Krug, PhD (Conservation Ecology, SUN), MSc (Conservation
>>> Biology, UCT), Dipl. Phys. (Germany)
>>>
>>> Centre of Excellence for Invasion Biology
>>> Stellenbosch University
>>> South Africa
>>>
>>> Tel : +33 - (0)9 53 10 27 44
>>> Cell: +33 - (0)6 85 62 59 98
>>> Fax : +33 - (0)9 58 10 27 44
>>>
>>> Fax (D): +49 - (0)3 21 21 25 22 44
>>>
>>> email: Rainer at krugs.de
>>>
>>> Skype: RMkrug
>>>
>>> PGP: 0x0F52F982
>>>
>>
>>
--
Rainer M. Krug, PhD (Conservation Ecology, SUN), MSc (Conservation Biology, UCT), Dipl. Phys. (Germany)
Centre of Excellence for Invasion Biology
Stellenbosch University
South Africa
Tel : +33 - (0)9 53 10 27 44
Cell: +33 - (0)6 85 62 59 98
Fax : +33 - (0)9 58 10 27 44
Fax (D): +49 - (0)3 21 21 25 22 44
email: Rainer at krugs.de
Skype: RMkrug
PGP: 0x0F52F982
-------------- next part --------------
A non-text attachment was scrubbed...
Name: not available
Type: application/pgp-signature
Size: 454 bytes
Desc: not available
URL: <https://stat.ethz.ch/pipermail/r-sig-geo/attachments/20151028/f2602873/attachment.bin>
More information about the R-sig-Geo
mailing list