On Thu, 11 Jun 2009, Bjarke Christensen wrote:


Hi,

I have a shapefile/SpatialPolygonsDataFrame of coastlines. Each polygon
corresponds to an island, and there are no holes. I can plot this so that
the islands are shaded by using
plot(islands, col="gray")

What I now want, is to plot the same information so that the ocean is blue
and the islands are transparent. Something like:
image(myDEM)
plot(ocean, col="lightblue", add=T)

which I would hope would allow the DEM to be visible on the islands, but
not in the ocean.

As you observe, the approach you are taking does not work, as the R graphics devices work by over-painting in layers. To see the image, it has to be painted after the enclosing rectangle. Holes are painted by re-painting the hole in a chosen background colour, which by default is "transparent", so you just get lots of blue.

Why not do something like

o <- overlay(myDEM, islands)

(untried) to get just the raster cells within the island polygons as a SpatialPixelsDataFrame object and image() that, possibly setting the background to a suitable colour (NAs will get transparent by default). This works with the graphics system, rather than trying to work round it - it doesn't "remember" that there is anything on the canvas that should be protected from overpainting, so it is safer just to paint what needs painting. I can also imagine painting first with reduced opacity (or intensity) in a different palette, then overpainting just the islands with full intensity in the target palette, which might be more visually pleasing than just flat blue sea.

Hope this helps,

Roger


The approach I have been taking does not seem to work, so I would like to
as for your suggestions. My approach is to make a new map, with one, big
rectangular polygon containing all the islands as holes.
I do this by making another shapefile, containing one polygon which
consists of the (slightly expanded) bounding box of the first shapefile. I
then import both shapefiles into GRASS, run

      v.overlay ainput=bounding...@permanent binput=isla...@permanent
output=ocean operator=not

and export the resulting vector shape to a shapefile. Plotting ocean in
GRASS

     d.vect map=oc...@permanent type=area fcolor=indigo

yields the image shown on
http://www.flickr.com/photos/39340...@n07/3616682800 (which looks exactly
like what I was looking for).

However after importing it into R with rgdal and plotting
ocean <- readOGR("myfolder", "myfilename")
plot(ocean, col="lightblue")

I get an image where the entire background is blue - not just the ocean,
but the islands as well.

Looking at the SpatialPolygonsDataFrame, I see that the polygons are not
marked as holes. But changing this does not seem to make any difference:

boundary <- which.max(sapply(oc...@polygons, function(x)
x...@polygons[[1]]@area))
for (i in ((1:length(oc...@polygons))[-boundary]))
oc...@polygons[[i]]@polygons[[...@hole <- T

It appears as if the non-hole polygon is drawn fully filled, irrespective
of any holes, and then the 'holes' are drawn on top. How might i get around
that?

Thanks in advance for any suggestions,
Bjarke Christensen

_______________________________________________
R-sig-Geo mailing list
R-sig-Geo@stat.math.ethz.ch
https://stat.ethz.ch/mailman/listinfo/r-sig-geo


--
Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Helleveien 30, N-5045 Bergen,
Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: roger.biv...@nhh.no

_______________________________________________
R-sig-Geo mailing list
R-sig-Geo@stat.math.ethz.ch
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Reply via email to