[R-sig-Geo] Converting GRID to SHAPE in R (was: Re: [Fwd: Re: rgdal: problem with grid data])

Agustin Lobo alobolistas at gmail.com
Wed Mar 24 10:35:18 CET 2010


Thanks for all the inputs I've got. I've
tried the following solution, but the steps
with as.SpatialPolygons.GridTopology() and
unionSpatialPolygons() take too long a time, had to abort the
second one after 1h, so
the procedure does not seem to be a really working
solution (note than I'm using the landcover map of a tiny country as example)

#Data from http://www.diva-gis.org/data/cov/ISR_cov.zip

isrlc <-
readGDAL("/media/Transcend/MASTER_ICTA2007_2008/GEODATA_miniprojs/Kefah/ISR_cov/isr_cov.gri")

 > class(isrlc)
[1] "SpatialGridDataFrame"
attr(,"package")
[1] "sp"

isrlcSP <- as.SpatialPolygons.GridTopology(getGridTopology(isrlc))

 > class(isrlcSP)
[1] "SpatialPolygons"
attr(,"package")
[1] "sp"

isrlcSPU <- unionSpatialPolygons(isrlcSP, isrlc at data[,1])
Aborted after 1h

A second solution is saving as gtif, then gdal_polygonize:

gdal_translate isr_cov.gri isr_cov.tif
gdal_polygonize.py -f 'ESRI Shapefile' isr_cov.tif isr_cov

which takes longer for writing the commands than the actual execution time!
Then we still have to attach the legend, which can be easily done in R. So I 
think the best is writing a function using system() for the above 2 commands, then
reading the grid and the shape, attach the legend to the SPDF and save again as
shape.

So, problem solved. I'll post the function to the list as soon as I make it.

Agus


Michael Sumner wrote:
> Hi, I'm not sure on the license status of gpclib in maptools - not for
> commercial use now - but I'm sure Roger will comment on that on the
> list. This seems to work:
> 
> library(maptools)
> 
> gdata <- image2Grid(list(x = 1:nrow(volcano), y = 1:ncol(volcano), z = volcano))
> 
> SpP_grd <- as.SpatialPolygons.GridTopology(getGridTopology(gdata))
> 
> 
> unionpoly <-   unionSpatialPolygons(SpP_grd, gdata$z)
> 
> 
> trdata <- data.frame(z = sort(unique(gdata$z)), row.names =
> sapply(slot(unionpoly, "polygons"), function(i) slot(i, "ID")))
> 
> 
> SpPDF <- SpatialPolygonsDataFrame(unionpoly, trdata)
> 
> 
> I checked the result in GIS and it seems to work perfectly.
> 
> Cheers, Mike.
> 
> On Tue, Mar 23, 2010 at 6:01 PM, Agustin Lobo <alobolistas at gmail.com> wrote:
>>
>> -------- Original Message --------
>> Subject: Re: [R-sig-Geo] rgdal: problem with grid data
>> Date: Tue, 23 Mar 2010 07:53:38 +0100
>> From: Agustin Lobo <aloboaleu at gmail.com>
>> Reply-To: Agustin.Lobo at ija.csic.es
>> To: Michael Sumner <mdsumner at gmail.com>
>> CC: sig-geo <r-sig-geo at stat.math.ethz.ch>
>> References: <4BA86074.8010709 at gmail.com>
>> <522664f81003222349g39cc0c4ya756547d80c1d2f6 at mail.gmail.com>
>>
>> ok, I see, I have to convert to Spatial Polygons data frame first, right?
>> Ideally, I should be able of clumping those cells with the same value,
>> so this is a bit more complicated than I expected. But should be possible.
>> I'll try and let the list know.
>> Thanks
>>
>> Agus
>>
>> Michael Sumner wrote:
>>> What are you expecting to get - a point for every raster cell?  If so,
>>> use gridded(isrlc) <- FALSE to convert to a SPDF first - then the
>>> write will work.
>>>
>>> More complicated conversions are possible, such as with
>>> as.SpatialPolygons.GridTopology.
>>>
>>> Cheers, Mike.
>>>
>>>
>>>
>>>
>>>
>>> On Tue, Mar 23, 2010 at 5:32 PM, Agustin Lobo <alobolistas at gmail.com>
>>> wrote:
>>>> Hi!
>>>>
>>>>  I read a grid file into R with
>>>> isrlc <-
>>>>
>>>> readGDAL("/media/Transcend/MASTER_ICTA2007_2008/GEODATA_miniprojs/Kefah/ISR_cov/isr_cov.gri")
>>>>
>>>> While isrlc seems correct,  I'm getting an empty shapefile when I attempt
>>>> to
>>>> write:
>>>>
>>>>
>>>> writeOGR(isrlc,dsn="/media/Transcend/MASTER_ICTA2007_2008/GEODATA_miniprojs/Kefah",
>>>> layer="ISR_cov2",driver="ESRI Shapefile")
>>>>
>>>> Data from
>>>> http://www.diva-gis.org/data/cov/ISR_cov.zip
>>>>
>>>> Thanks
>>>>
>>>> Agus
>>>>
>>>> _______________________________________________
>>>> R-sig-Geo mailing list
>>>> R-sig-Geo at stat.math.ethz.ch
>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>>
>> --
>> Dr. Agustin Lobo
>> Institut de Ciencies de la Terra "Jaume Almera" (CSIC)
>> LLuis Sole Sabaris s/n
>> 08028 Barcelona
>> Spain
>> Tel. 34 934095410
>> Fax. 34 934110012
>> email: Agustin.Lobo at ija.csic.es
>> http://www.ija.csic.es/gt/obster
>>
>> _______________________________________________
>> R-sig-Geo mailing list
>> R-sig-Geo at stat.math.ethz.ch
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
>



More information about the R-sig-Geo mailing list