[R-sig-Geo] Polygons VS MultiLineString

[Michael Sumner]

<snip>
> > I use "rasterToContour" to derive isolines for raster. After that I
> > save them as GeoJSON (or ESRI Shapefile). However, it saves contours
> > as "MultiLineString"
> >
> > Is it possible to save "contours" as "Polygons"?
> >


If your data don't build clean lines I would be inclined to classify the
raster data and then turn that into polygons.
But, that means the boundaries are not as smooth.

Using example data:

library(raster)
r <- raster(volcano)
cl <- rasterToContour(r, levels = c(110, 130))

## see how cl$level == 130 will polygonize, but not 110
library(rgeos)
p1 <- gPolygonize(subset(cl, level == 130))
p2 <- gPolygonize(subset(cl, level == 110))
plot(p1)  ## ok
is.null(p2) ## TRUE

And that is even before you decide if you want non-intersecting polygon
donuts etc.

So, forgetting contouring

## cut into our breaks
vals <- c(110, 130)
rc <- cut(r, c(vals, cellStats(r, max)))
plot(rc)
plot(cl, add = TRUE)

## another way to polygonize
rp <- rasterToPolygons(rc, dissolve = TRUE)
plot(rp, col = c("grey", "lightblue"))
rp$layer <- vals

## check all is well
plot(r)
plot(cl, add = TRUE)
contour(r, lev = vals, add = TRUE)
text(rp[2,], lab = rp$layer[2])

plot(rp, add = TRUE, lty = 2)

So there we end up with pretty rough boundaries, since the polygonizing
from cells is following pixel boundaries exactly. Maybe we can
smooth/disaggregate the raster, but it gets pretty expensive.

Another trick might be to buffer the raster with values outside the range,
and do contour individually
for each level.

So,

i <- 1
r1 <- r
r1[is.na(r1) | r1 < vals[i]] <- vals[i] - 10
r2 <- extend(r1, extent(r) + res(r) * 2, value = vals[i] - 10)

plot(r2)
plot(gPolygonize(cl1), add = TRUE, col = grey(0.5, alpha = 0.6))


Is that better, maybe - hard to generalize I imagine, and still you need to
clip your larger boundary with internal ones, and merge them into one data
set.

Thanks for the question though, it's given me some ideas to fix some things
I need to do. (Manifold, and I imagine other, GIS does this seamlessly so
it would be nice to have contour and contourLines upgraded to be able to
trace around the edges.)

Cheers, Mike.


On Wed, 4 Mar 2015 at 19:26 Frede Aakmann Tøgersen <frtog at vestas.com> wrote:

> Well, since the workhorse behind the rasterToContour() function is
> contourLines() from the lattice package don't expect the contour lines for
> a given level to form a single well defined polygon (the definition of
> which of course needs to be specified).
>
> Try this:
>
> f <- system.file("external/test.grd", package="raster")
> r <- raster(f)
> x <- rasterToContour(r)
> class(x)
> str(x)
>
> ## the levels
> levels(x at data$level)
>
> plot(r)
>
> plot(subset(x, level == 400), add=TRUE)
>
> plot(subset(x, level == 200), add=TRUE)
>
> Yours sincerely / Med venlig hilsen
>
>
> Frede Aakmann Tøgersen
> Specialist, M.Sc., Ph.D.
> Plant Performance & Modeling
>
> Technology & Service Solutions
> T +45 9730 5135
> M +45 2547 6050
> frtog at vestas.com
> http://www.vestas.com
>
> Company reg. name: Vestas Wind Systems A/S
> This e-mail is subject to our e-mail disclaimer statement.
> Please refer to www.vestas.com/legal/notice
> If you have received this e-mail in error please contact the sender.
>
> > -----Original Message-----
> > From: R-sig-Geo [mailto:r-sig-geo-bounces at r-project.org] On Behalf Of
> > Michael Sumner
> > Sent: 4. marts 2015 09:11
> > To: Antonio Rodriges; r-sig-geo at r-project.org
> > Subject: Re: [R-sig-Geo] Polygons VS MultiLineString
> >
> > Try gPolygonize from rgeos, but it *really* depends on whether your lines
> > close cleanly.
> >
> > On Wed, 4 Mar 2015 18:38 Antonio Rodriges <antonio.rrz at gmail.com>
> > wrote:
> >
> > > Hello,
> > >
> > > I use "rasterToContour" to derive isolines for raster. After that I
> > > save them as GeoJSON (or ESRI Shapefile). However, it saves contours
> > > as "MultiLineString"
> > >
> > > Is it possible to save "contours" as "Polygons"?
> > >
> > > Code:
> > >
> > > to_c - is a raster
> > > contours <- rasterToContour(to_c, levels = c(18, 33))
> > > writeOGR(contours, "file.js", layer="", driver="GeoJSON")
> > >
> > > _______________________________________________
> > > R-sig-Geo mailing list
> > > R-sig-Geo at r-project.org
> > > https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> > >
> >
> >       [[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
>

	[[alternative HTML version deleted]]