[R-sig-Geo] 3rd+ rasters added using "plot(r, add=TRUE)" are misregistered

Michael Sumner mdsumner at gmail.com
Fri Dec 6 02:26:47 CET 2013


Interesting and difficult!  Thanks for pointing out the useRaster
difference, I had mixed up TRUE/FALSE in my original report (offlist)
on that one. That might help me find the solution, though sorry I have
nothing for your issue. :)

Cheers, Mike.

On Fri, Dec 6, 2013 at 11:58 AM, Josh O'Brien <joshmobrien at gmail.com> wrote:
> Hi Michael,
>
> Thanks for those thoughts.
>
> The issue you raise (which has bit me before as well) appears to be
> unrelated. In particular, your resizing issue disappears if one sets
> useRaster=FALSE, whereas my problem persists whether or not rasterImage() is
> used to draw the raster. Try this to get a version of the volcano and
> contour that resize/stretch in sync:
>
>     plot(raster(v), xaxs = "i", useRaster=FALSE)
>     lines(cl)
>     ## Check that resizing works correctly. (It does, presumably b/c
>     ## image(..., useRaster=FALSE) doesn't attempt to maintain aspect ratio
> of the raster.
>
> The problem I described seems related to resetting of elements of `par` by
> the first call to `plot(r, ..., add=TRUE)`, as can be seen by doing:
>
>
>     plot(r, col="yellow", legend=FALSE)
>     p0 <- par()
>
>     plot(r, col="green", legend=FALSE, add=TRUE)
>     p1 <- par()
>
>     plot(r, col="grey", legend=FALSE, add=TRUE)
>     p2 <- par()
>
>     plot(r, col="yellow", legend=FALSE, add=TRUE)
>     p3 <- par()
>
>     ## Find elements of par that are changed by first call to plot(r, ...,
> add=TRUE)
>     names(p1)[!mapply(identical, p0,p1)]
>     # [1] "cxy" "pin" "plt"
>
>     ## Show that par does not change thereafter
>     all(identical(p1,p2), identical(p2,p3))
>     # [1] TRUE
>
> As might be expected if par is to blame, frequently resetting par to its
> 'original' value makes the plots come out right:
>
>     plot(r, col="yellow", legend=FALSE)
>     p0 <- par(no.readonly=TRUE)
>
>     plot(r, col="green", legend=FALSE, add=TRUE)
>     plot(r, col="blue", legend=FALSE, add=TRUE)
>     par(p0)
>
>     plot(r, col="red", legend=FALSE, add=TRUE)
>     plot(r, col="grey", legend=FALSE, add=TRUE)
>     par(p0)
>     plot(r, col="dodgerblue", legend=FALSE, add=TRUE)
>
> The piece that is most puzzling to me is that, although cxy, pin, & plt get
> reset after the first added layer is plotted, it isn't until _two_ more
> layers are added that the problem manifests itself
>
> Seems like an interesting puzzle here, which might be solvable by doing
> `debug(raster:::.imageplot)` and stepping through it, paying careful
> attention to how par gets set and reset during each of the calls to plot().
>
> Cheers,
>
> Josh
>
>
> On Thu, Dec 5, 2013 at 2:02 PM, Michael Sumner <mdsumner at gmail.com> wrote:
>>
>> I just wonder if this is related to an existing resize problem when
>> overplotting? I see this on the windows() graphics device, you need to
>> try these two sets of plot commands and resize the window (especially
>> stretch it in one direction) to see the problem at ## 1 and ## 2
>>
>> v <- list(x = seq(0, 1, length = nrow(volcano)), y = seq(0, 1, length
>> = ncol(volcano)), z = volcano)
>> cl <- contourLines(v, levels = 160)[[2]]
>>
>> image(v, xaxs = "i")
>> lines(cl)
>> ## (1) now resize the window, it's all good - the contour stays in the
>> ## right place
>>
>>
>> ## however, resize after adding the lines to a raster plot
>> ## and the relationship is broken
>> library(raster)
>>
>> plot(raster(v), xaxs = "i")
>>
>> lines(cl)
>>
>> ## (2)
>>
>> May be totally unrelated but I thought it worth pointing out. It's
>> hard to trace through the setting up that raster's plot does (well, I
>> find it hard to do - I've tried a few times to figure out where these
>> problems are happening).
>>
>> Cheers, Mike.
>>
>> On Wed, Dec 4, 2013 at 5:57 AM, Josh O'Brien <joshmobrien at gmail.com>
>> wrote:
>> > Hi all,
>> >
>> > Briefly, running the following code will not produce a single grey ring,
>> > as
>> > would be expected. The first three raster layers all stack up neatly,
>> > but
>> > the 4th and subsequent raster layers are all shifted/skewed to the
>> > right.
>> >
>> > It appears that the plotting engine loses track of the fact that some
>> > space
>> > in the original plot was allocated to a legend, and just begins plotting
>> > as
>> > if a wider canvas were available to it.
>> >
>> > I'm not sure whether this is ultimately an issue with raster or sp or
>> > something else, which is why I'm noting it here on R-sig-geo. (I'm also
>> > fully aware that there are other and better ways to add several raster
>> > layers to a plot, but still think the following behavior qualifies as a
>> > bug
>> > that I ought to report.)
>> >
>> >   library(maptools) ## Only needs to be installed for example data
>> >   library(raster)
>> >   library(rgeos)
>> >
>> >   ## Create an example raster
>> >   p <- shapefile(system.file("shapes/co37_d90.shp", package="maptools"))
>> >   p <- p[31,]  ## A tall narrow county polygon
>> >   pr <- gDifference(gBuffer(p, width=.01), p)
>> >   r <- rasterize(pr, raster(extent(pr), ncol=100, nrow=100))
>> >
>> >   ## These three are properly registered on one another
>> >   plot(r, col="yellow", legend=FALSE)
>> >   plot(r, col="green", legend=FALSE, add=TRUE)
>> >   plot(r, col="grey", legend=FALSE, add=TRUE)
>> >   ## All subsequent "layers" are improperly shifted/skewed to right
>> >   plot(r, col="yellow", legend=FALSE, add=TRUE)
>> >   plot(r, col="blue", legend=FALSE, add=TRUE)
>> >   plot(r, col="red", legend=FALSE, add=TRUE)
>> >   plot(r, col="grey20", legend=FALSE, add=TRUE)
>> >   ## Following the above, SpatialPolygons are also shifted/skewed
>> >   plot(p)
>> >
>> >
>> >
>> > My session info, in case it matters, is:
>> >
>> >> sessionInfo()
>> > R version 3.0.2 (2013-09-25)
>> > Platform: i386-w64-mingw32/i386 (32-bit)
>> >
>> > locale:
>> > [1] LC_COLLATE=English_United States.1252
>> > [2] LC_CTYPE=English_United States.1252
>> > [3] LC_MONETARY=English_United States.1252
>> > [4] LC_NUMERIC=C
>> > [5] LC_TIME=English_United States.1252
>> >
>> > attached base packages:
>> > [1] stats     graphics  grDevices utils     datasets  methods   base
>> >
>> > other attached packages:
>> > [1] rgdal_0.8-14  rgeos_0.3-2   raster_2.1-66 sp_1.0-14
>> >
>> > loaded via a namespace (and not attached):
>> > [1] compiler_3.0.2  grid_3.0.2      lattice_0.20-24 tools_3.0.2
>> >
>> >         [[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
>>
>>
>>
>> --
>> Michael Sumner
>> Hobart, Australia
>> e-mail: mdsumner at gmail.com
>
>



-- 
Michael Sumner
Hobart, Australia
e-mail: mdsumner at gmail.com



More information about the R-sig-Geo mailing list