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@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@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@r-project.org
> > https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
>
>
> --
> Michael Sumner
> Hobart, Australia
> e-mail: mdsumner@gmail.com
>

	[[alternative HTML version deleted]]

