[R] expand gridded matrix to higher resolution
Anthoni, Peter (IMK)
peter.anthoni at kit.edu
Thu Jul 6 08:42:18 CEST 2017
Hi Jeff,
thanks, the raster package disaggregate will do the trick as well.
library(raster)
rmm <- raster(ncols=5, nrows=3)
rmm[] <- matrix(1:15,nrow=3,byrow = T)
xrmm <- disaggregate(rmm, fact=c(3, 3))
> > as.matrix(rmm)
> [,1] [,2] [,3] [,4] [,5]
> [1,] 1 2 3 4 5
> [2,] 6 7 8 9 10
> [3,] 11 12 13 14 15
> > as.matrix(xrmm)
> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12] [,13] [,14] [,15]
> [1,] 1 1 1 2 2 2 3 3 3 4 4 4 5 5 5
> [2,] 1 1 1 2 2 2 3 3 3 4 4 4 5 5 5
> [3,] 1 1 1 2 2 2 3 3 3 4 4 4 5 5 5
> [4,] 6 6 6 7 7 7 8 8 8 9 9 9 10 10 10
> [5,] 6 6 6 7 7 7 8 8 8 9 9 9 10 10 10
> [6,] 6 6 6 7 7 7 8 8 8 9 9 9 10 10 10
> [7,] 11 11 11 12 12 12 13 13 13 14 14 14 15 15 15
> [8,] 11 11 11 12 12 12 13 13 13 14 14 14 15 15 15
> [9,] 11 11 11 12 12 12 13 13 13 14 14 14 15 15 15
the disaggregate as a bit faster than the tapply.
mmb=matrix(1:259200,nrow=720,ncol=360)
rmmb <- raster(ncols=360, nrows=720)
rmmb[] <- mmb[]
system.time(for(i in 1:10) {xmm=matrix(NA,nrow=nrow(mmb)*3,ncol=ncol(mmb)*3)
for(icol in 1:ncol(mmb)) {
for(irow in 1:nrow(mmb)) {
xicol=(icol-1)*3 +c(1:3)
xirow=(irow-1)*3 +c(1:3)
xmm[xirow,xicol]=mmb[irow,icol]
}
}
})
system.time(for(i in 1:10) {apply(t(apply(mmb,1,rep,each=3)),2,rep,each=3)}) #ca. 10x faster
system.time(for(i in 1:10) {xrmmb <- disaggregate(rmmb, fact=c(3, 3))})
> > system.time(for(i in 1:10) {xmm=matrix(NA,nrow=nrow(mmb)*3,ncol=ncol(mmb)*3)
> + for(icol in 1:ncol(mmb)) {
> + for(irow in 1:nrow(mmb)) {
> + xicol=(icol-1)*3 +c(1:3)
> + xirow=(irow-1)*3 +c(1:3)
> + xmm[xirow,xicol]=mmb[irow,icol]
> + }
> + }
> + })
> user system elapsed
> 8.297 0.048 8.364
> > system.time(for(i in 1:10) {apply(t(apply(mmb,1,rep,each=3)),2,rep,each=3)}) #ca. 10x faster
> user system elapsed
> 0.785 0.093 0.881
> > system.time(for(i in 1:10) {xrmmb <- disaggregate(rmmb, fact=c(3, 3))})
> user system elapsed
> 0.583 0.147 0.731
cheers
Peter
> On 5. Jul 2017, at 16:57, Jeff Newmiller <jdnewmil at dcn.davis.ca.us> wrote:
>
> You probably ought to be using the raster package. See the CRAN Spatial Task View.
> --
> Sent from my phone. Please excuse my brevity.
>
> On July 5, 2017 12:20:28 AM PDT, "Anthoni, Peter (IMK)" <peter.anthoni at kit.edu> wrote:
>> Hi all,
>> (if me email goes out as html, than my email client don't do as told,
>> and I apologies already.)
>>
>> We need to downscale climate data and therefore first need to expand
>> the climate from 0.5deg to the higher resolution 10min, before we can
>> add high resolution deviations. We basically need to have the original
>> data at each gridcell replicated into 3x3 gridcells.
>> A simple for loop can do this, but I could need a faster procedure.
>> Anybody know a faster way? Is there package than can do what we need
>> already?
>> I tried matrix with rep, but I am missing some magic there, since it
>> doesn't do what we need.
>> replicate might be promising, but then still need to rearrange the
>> output into the column and row format we need.
>>
>> A simple example:
>> mm=matrix(1:15,nrow=3,byrow = T)
>> xmm=matrix(NA,nrow=nrow(mm)*3,ncol=ncol(mm)*3)
>> for(icol in 1:ncol(mm)) {
>> for(irow in 1:nrow(mm)) {
>> xicol=(icol-1)*3 +c(1:3)
>> xirow=(irow-1)*3 +c(1:3)
>> xmm[xirow,xicol]=mm[irow,icol]
>> }
>> }
>> mm
>>>> mm
>>> [,1] [,2] [,3] [,4] [,5]
>>> [1,] 1 2 3 4 5
>>> [2,] 6 7 8 9 10
>>> [3,] 11 12 13 14 15
>>>
>> xmm
>>>> xmm
>>> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12]
>> [,13] [,14] [,15]
>>> [1,] 1 1 1 2 2 2 3 3 3 4 4 4
>> 5 5 5
>>> [2,] 1 1 1 2 2 2 3 3 3 4 4 4
>> 5 5 5
>>> [3,] 1 1 1 2 2 2 3 3 3 4 4 4
>> 5 5 5
>>> [4,] 6 6 6 7 7 7 8 8 8 9 9 9
>> 10 10 10
>>> [5,] 6 6 6 7 7 7 8 8 8 9 9 9
>> 10 10 10
>>> [6,] 6 6 6 7 7 7 8 8 8 9 9 9
>> 10 10 10
>>> [7,] 11 11 11 12 12 12 13 13 13 14 14 14
>> 15 15 15
>>> [8,] 11 11 11 12 12 12 13 13 13 14 14 14
>> 15 15 15
>>> [9,] 11 11 11 12 12 12 13 13 13 14 14 14
>> 15 15 15
>>
>> I tried various rep with matrix, but don't get the right result.
>> xmm2=matrix(rep(rep(mm,each=3),times=3),nrow=nrow(mm)*3,ncol=ncol(mm)*3,byrow
>> = F)
>>> identical(xmm,xmm2)
>> [1] FALSE
>>
>> rr=replicate(3,rep(t(mm),each=3))
>> rr
>>>> rr
>>> [,1] [,2] [,3]
>>> [1,] 1 1 1
>>> [2,] 1 1 1
>>> [3,] 1 1 1
>>> [4,] 2 2 2
>>> [5,] 2 2 2
>>> [6,] 2 2 2
>>> [7,] 3 3 3
>>> ...
>> identical(xmm,matrix(rr,ncol=15,nrow=9,byrow=T))
>>>> identical(xmm,matrix(rr,ncol=15,nrow=9,byrow=T))
>>> [1] FALSE
>>
>> Many thanks for any advice.
>>
>> cheers
>> Peter
>>
>> ______________________________________________
>> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide
>> http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
More information about the R-help
mailing list