[R-sig-Geo] from x,y plot to a raster with Raster package?

[Pep SD]

Hi Robert,

Thank you for your prompt response. That was actually really helpful.

The issue comes when your x and y axis on the plot (or the data) are
actually not round values and therefore we cannot do the  rounding
('s<-round(s)'). In that case we need to reassign each value, but I do not
know how to do it efficiently

Consider this case (from your example)

# make a three column data.frame
d <- data.frame( x=c(23.2, 33.7966 , 44.3932 ), y=c(33.2 , 34.46426 ,
35.72852) , z=c(0.1,0.2,0.3))

# your two rasters, representing x and y
library(raster)
x <- y <- raster ( nrows=10, ncols=10, crs=NA)
values(x) <- runif(n=100, min=22, max=45)
values(y) <- runif(n=100, min=31, max=37)



#before attempting to perform subs with the d data.frame as a base for
substitution, we need the x and y raster to match EXACTLY the values in the
axis.

#applying it on a cell by cell basis, it is easy...we try the nearest
neighbor value.
lapply (values (x), FUN = function (i) which.min(abs(i - d$x)))

#is there a way to do it with calc? I am just saying because lapplying over
a 14M cells might be slow
v1 <- calc(x = x,fun = function (x) which.min(abs(x - d$x))) #obviously
didnt work :(



THANKS !

Pep




2015-07-19 16:35 GMT-04:00 Robert J. Hijmans <r.hijmans at gmail.com>:

> Josep,
>
> I take it you have three variables, and want to create a raster for
> the third variable, based on the other two, that are also represented
> as rasters. If so, you can
>
> # make a three column data.frame
> d <- data.frame( x=c(23,25,27), y=c(33,34,35) , z=c(0.1,0.2,0.3))
>
> # your two rasters, representing x and y
> library(raster)
> x <- y <- raster ( nrows=10, ncols=10, crs=NA)
> values(x) <- runif(n=100, min=20, max=28)
> values(y) <- runif(n=100, min=31, max=36)
>
> # combine them
> s <- stack(x, y)
>
> # round them to match to x and y in 'd'
> s <- round(s)
>
> # compute z, based on x and y using 'subs' (equivalent to base::merge)
> z <- subs(s, d, by=1:2)
>
> # test
> x <- stack(s, z)
> # should have the same values as 'd'
> unique(na.omit(values(x)))
>
> Robert
>
> On Sun, Jul 19, 2015 at 10:10 AM, Josep M Serra diaz
> <pep.bioalerts at gmail.com> wrote:
> > Hi everyone,
> >
> > I am wondering if I can speed up the process of transforming from  a
> 'plot'
> > (or matrix) to a raster.
> >
> > The idea is that I have a matrix with x,y,z values. and I have a raster
> > with spatial data for x and for y, and I want to get a raster of z.
> >
> > I would know how to do it on a cell-by-cell basis, but I was hoping there
> > was a more efficient way, since I have more than 14Milion Cells.
> >
> > I left an simple example below to see the not very efficient way in
> which I
> > was dealing with it, hoping that someone will let me know that I have a
> > better chance trying with a function in the raster package that can do
> that
> > very easily. So far, I have tried with calc, but it is saying that 'I
> > cannot use this function'.
> >
> >
> >
> >
> > ###### SIMPLE EXAMPLE
> >
> > ##dimensions of the matrix
> > mat <- list( x=c(23,25,27), y=c(33,34,35) , z=c(0.1,0.2,0.3))
> > mat2 <- matrix (data = mat$z,nrow=3, ncol=3,dimnames = mat[c(1,2)])
> >
> > #input raster x axis
> > r1 <- raster ( nrows=10, ncols=10, crs=NA)
> > r1[] <- runif(n = 100,min=20,max=28)
> >
> > #input raster y axis
> > r2 <- raster ( nrows=10, ncols=10, crs=NA)
> > r2[] <- runif(n = 100,min=31,max=36)
> >
> >
> > #AIM: I want a raster that has in the cell values the corresponding value
> > of the z dimension in mat.
> > #    That is r1 contains the values of the xaxis in mat, r2 contains
> inputs
> > values of the y axis in mat
> > #    I want the raster with the values of z in mat
> >
> >
> > #on a cell by cell basis it is 'easy'
> >
> > cell.values <- lapply (1:ncell(r1),FUN = function (i){
> >   pos.x.axis <- if (r1[i] < range(mat$x)[1] | r1[i] > range(mat$x)[2])
> {NA}
> > else{
> >     which(abs(r1[i]-mat$x) == min(abs(r1[i]-mat$x))) #it gives me the
> > closest postition in the y axis
> >     }
> >
> >   pos.y.axis <- if (r2[i] < range(mat$y)[1] | r2[i] > range(mat$y)[2])
> {NA}
> > else{
> >     which(abs(r2[i]-mat$x) == min(abs(r2[i]-mat$x))) #it gives me the
> > closest postition i nthe y axis
> >   }
> >
> >   if (is.na (pos.x.axis) ==T | is.na(pos.y.axis)==T) {out.value <- NA}
> else{
> >     out.value <- mat2[pos.x.axis, pos.y.axis]
> >   }
> >
> >   return (out.value)
> >
> > })
> >
> > raster.out <- r1 ; raster.out [] <- NA
> > raster.out[]<- unlist(cell.values)
> > plot (raster.out)
> >
> >
> > any hint on how to do it for 14Milion cells?
> >
> >         [[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]]