[R-sig-Geo] Needing to speed up a process involving calc() and cover() raster functions
Mathieu Rajerison
mathieu.rajerison at gmail.com
Mon Jan 4 11:51:13 CET 2016
Hi,
And thanks for all the answers.
I've tested Robert's elegant solution successfully.
I haven't yet tested your solution, Banjamin, but it seems a really good
improvement for speeding all this up and I'll soon give it a try. Thanks
for your awesome RSToolbox package
To answer your question "why setting all pixels to 1", it's just that
afterwards, I smooth the result with a mean-shift focal function. This way,
it gives all the extracted pixels the same value/weight.
Thanks again and happy new yeaR
Mathieu
2015-12-25 6:46 GMT+01:00 Robert J. Hijmans <r.hijmans at gmail.com>:
> How about this:
>
> b <- brick(system.file("external/rlogo.grd", package="raster"))
> classified <- b / 255
>
> threshs <- c(.1, .3, .5)
> x <- classified < threshs
> y <- max(x, na.rm=TRUE)
> z <- reclassify(y, cbind(0, NA))
>
>
> Robert
>
>
> On Tue, Dec 22, 2015 at 1:57 AM, Mathieu Rajerison
> <mathieu.rajerison at gmail.com> wrote:
> > Hi,
> >
> >
> > I use RSToolBox to classify a RGB raster.
> >
> > I have a resulting RasterBrick which has as many layer as end members, in
> > my case 3 for different tones of blue.
> >
> > I reclassify each band with calc to extract the pixels which have a small
> > angle mapping value. The threshold used is different depending on the
> > endmember layer.
> >
> > I finally assembly all the bands with the cover function.
> >
> > I needed to increase the memory limit assigned to R to have it worked. I
> > suspect that my code could be optimized, but I don't know in which way.
> >
> > Here is the part of my code, that I think, could be optilmized, if you
> want
> > to have a look and give some advice :
> >
> > # RECLASSIFY
> > # classified is the classified RasterStack
> > # here I change the values of each band to 1 or NA depending on the
> > spectral angle mapping value.
> > # Is calc() slower than reclassify() for this purpose as I have only one
> > threshold value ?
> >
> > threshs = c(.1,.2,.1)
> > for (i in 1:nlayers(classified )) {
> >
> > clas = classified[[i]]
> > thresh=threshs[i]
> >
> > out[[i]] = calc(clas, function(x) {x[x >= thresh] = NA;
> > x[x < thresh] = 1;
> > return(x)})
> > }
> >
> > # COVERING
> > r = out[[1]]
> > for(i in 2:length(out)) {
> > r = cover(out[[i]], r) ## I cover by iteration
> > }
> >
> > plot(r) # r is the final combined raster
> >
> > ================================================================
> > My sessionInfo() :
> >> sessionInfo()
> > R version 3.1.2 (2014-10-31)
> > Platform: x86_64-w64-mingw32/x64 (64-bit)
> >
> > locale:
> > [1] LC_COLLATE=French_France.1252 LC_CTYPE=French_France.1252
> > LC_MONETARY=French_France.1252
> > [4] LC_NUMERIC=C LC_TIME=French_France.1252
> >
> > attached base packages:
> > [1] stats graphics grDevices utils datasets methods base
> >
> > other attached packages:
> > [1] R.utils_2.1.0 R.oo_1.19.0 R.methodsS3_1.7.0
> > igraph_1.0.1 scatterplot3d_0.3-36
> > [6] gdalUtils_2.0.1.7 spdep_0.5-88 Matrix_1.2-2
> > maptools_0.8-34 spgrass6_0.8-8
> > [11] XML_3.98-1.3 rgeos_0.3-8 FNN_1.1
> > rgdal_0.9-2 RStoolbox_0.1.1
> > [16] raster_2.3-40 sp_1.0-17
> >
> > loaded via a namespace (and not attached):
> > [1] boot_1.3-13 car_2.0-25 caret_6.0-57 coda_0.18-1
> > codetools_0.2-9 colorspace_1.2-6
> > [7] deldir_0.1-9 digest_0.6.8 doParallel_1.0.8
> foreach_1.4.2
> > foreign_0.8-61 geosphere_1.4-3
> > [13] ggplot2_1.0.1 grid_3.1.2 gtable_0.1.2
> > iterators_1.0.7 lattice_0.20-29 LearnBayes_2.15
> > [19] lme4_1.1-10 magrittr_1.5 MASS_7.3-35
> > MatrixModels_0.4-1 mgcv_1.8-3 minqa_1.2.4
> > [25] munsell_0.4.2 nlme_3.1-118 nloptr_1.0.4 nnet_7.3-8
> > parallel_3.1.2 pbkrtest_0.4-2
> > [31] plyr_1.8.3 proto_0.3-10 quantreg_5.19 Rcpp_0.12.0
> > reshape2_1.4.1 scales_0.2.5
> > [37] SparseM_1.7 splines_3.1.2 stats4_3.1.2
> stringi_0.5-5
> > stringr_1.0.0 tools_3.1.2
> >
> > ============================================================
> > Best,
> >
> > Mathieu
> >
> > [[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]]
More information about the R-sig-Geo
mailing list