[R-sig-Geo] spgrass6, readVECT and overwriting ?

Mathieu Rajerison mathieu.rajerison at gmail.com
Wed Oct 28 14:20:47 CET 2015


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")
}

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
>>
>
>

	[[alternative HTML version deleted]]



More information about the R-sig-Geo mailing list