[R-sig-Geo] trying to add points to a raster

Robert J. Hijmans r.hijmans at gmail.com
Tue Apr 13 18:37:55 CEST 2010


Hi Roman,

The RasterLayer rst.p is correct, I think. The problem has only to do
with the plotting. The error/warning occurs because you have (a) a
relatively large RasterLayer and (b) hardly any cells have values.
Large RasterLayers are sampled (regularly) for plotting. In this case,
none of the sampled cells get a value.  You can stop this from
happening by adding a maxpixels argument to plot:

plot(rst.p, maxpixels=ncell(rst.p))

Not that it helps much; I was not able to see the difference because
there are so few cells with values (and they are still resampled to
match the screen resolution).

For visualization (only) you could make a coarser grid:

res(rst) = 10
rst.p <- pointsToRaster(raster = rst,  xy = r,   values = r at data$walker)
plot(rst.p)

Note that the default function 'fun' is length (ie. count the number
of points per cell). As you provide a 'values' argument, I assume you
want something else. To transfer the value (of the last point that is
found for a cell) you can use

rst.p <- pointsToRaster(raster = rst,  xy = r,   values =
r at data$walker,  fun=function(x){x} )

# the mean would be fun=mean, but as 'walker'  is a factor:
rst.p <- pointsToRaster(raster = rst,  xy = r,   values =
r at data$walker,  fun=function(x){mean(as.numeric(x)) } )


I am sorry for your head.

Best,
Robert

# This is the command I use to add points but for some reason, the
# values are not added (strictly speaking, only one value is added).
# Can someone point out what I'm doing wrong?
rst.p <- pointsToRaster(raster = rst,
       xy = r,
       values = r at data$walker)
rst.p



> rst.p <- pointsToRaster(raster = rst2,
+        xy = r,
+        values = r at data$walker)
> plot(rst.p)
> ?plot



On Tue, Apr 13, 2010 at 4:45 AM, Roman Luštrik <roman.lustrik at gmail.com> wrote:
> # Dear list,
> #
> # I'm apologizing in advance if my question is rather basic, but I've
> # been banging my head against the wall over this all day and the
> # seismology department has already contacted me to stop it.
> #
> # I'm trying to add points (SpatialPointsDataFrame object r) to a
> # raster object (rst).
>
> library(raster)
> library(sp)
>
> # object with points I'm trying to add
>
> r <- new("SpatialPointsDataFrame"
>        , data = structure(list(walker = structure(c(3L, 3L, 3L, 3L, 3L,
> 12L,
>                                        12L, 12L, 12L, 12L, 19L, 19L, 19L,
> 19L, 19L, 20L, 20L, 20L, 20L,
>                                        20L), .Label = c("1", "2", "3", "4",
> "5", "6", "7", "8", "9",
>                                        "10", "11", "12", "13", "14", "15",
> "16", "17", "18", "19", "20"
>                                ), class = "factor")), .Names = "walker",
> class = "data.frame", row.names = c(212L,
>                        216L, 232L, 294L, 279L, 1134L, 1124L, 1173L, 1195L,
> 1123L, 1877L,
>                        1819L, 1810L, 1802L, 1899L, 1918L, 1965L, 1977L,
> 1952L, 1915L
>                ))
>        , coords.nrs = 2:3
>        , coords = structure(c(-23.2432695470495, -42.8570787452707,
> -25.9987883825460,
>                        -64.2282194183847, 15.9572133404169,
> 841.345031514424, 846.560380022882,
>                        858.356702340203, 829.290173246473,
> 855.204942248044, -153.570529391965,
>                        -222.615649663872, -189.494697286148,
> -150.755239868264, -172.405589596216,
>                        700.62999593643, 714.319660569588, 769.937984828572,
> 721.180459026309,
>                        688.932274938503, 118.294179820811,
> 114.236959667182, 141.137206092922,
>                        55.7030239659073, 136.785318376395,
> 707.681240819362, 691.95292905712,
>                        595.302752910708, 590.191085932507,
> 671.318331940053, -545.390752528143,
>                        -593.352201166765, -582.711927913932,
> -577.494929150682, -579.851251631267,
>                        -169.417426401479, -199.140015976785,
> -195.301934322738, -130.762717546069,
>                        -126.873035982826), .Dim = c(20L, 2L), .Dimnames =
> list(NULL,
>                        c("x", "y")))
>        , bbox = structure(c(-222.615649663872, -593.352201166765,
> 858.356702340203,
>                        707.681240819362), .Dim = c(2L, 2L), .Dimnames =
> list(c("x",
>                                "y"), c("min", "max")))
>        , proj4string = new("CRS"
>                , projargs = NA_character_
>        )
> )
>
> # This is the raster object used as a recipient
>
> rst <- raster(nrow = 2000, ncol = 2000, xmn = -1000, xmx = 1000, ymn =
> -1000, ymx = 1000)
>
> # This is the command I use to add points but for some reason, the
> # values are not added (strictly speaking, only one value is added).
> # Can someone point out what I'm doing wrong?
> rst.p <- pointsToRaster(raster = rst,
>        xy = r,
>        values = r at data$walker)
> rst.p
>
> # plot(rst.p)
> # Error in if (del == 0 && to == 0) return(to) :
> #   missing value where TRUE/FALSE needed
> # In addition: Warning messages:
> # 1: In min(x) : no non-missing arguments to min; returning Inf
> # 2: In max(x) : no non-missing arguments to max; returning -Inf
> # 3: In min(x) : no non-missing arguments to min; returning Inf
> # 4: In max(x) : no non-missing arguments to max; returning -Inf
>
>
> # Cheers,
> # Roman
>
>
>
> --
> In God we trust, all others bring data.
>
>        [[alternative HTML version deleted]]
>
> _______________________________________________
> 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