Deal all, I'm merging several rasters representing patches of forest, and the pixel values represent the area of the patches. The patches' area should be summed along the edges of the merged rasters. I provide a reproducible example. Any help is welcome. Hugo library(raster) r1<-raster(nrows=10, ncols=10, xmn=0, xmx=1, ymn=0, ymx=1) r2<-raster(nrows=10, ncols=10, xmn=0, xmx=1, ymn=1, ymx=2) r1[]<-c(rep(24,24), rep(0,76)) r2[]<-c(rep(0,86), rep(14,14)) par(mfrow=c(2,2)); plot(r1,main="r1"); plot(r2,main="r2") r3<-merge(r1,r2) ; plot(r3, main="merged rasters") # do something r3[r3>0]<-38 ; plot(r3, main="desired result")
How to sum area of adjacent patches?
3 messages · Ben Tupper, Hugo Costa
Hi, I think you may need to use raster::clump() and table() to label your regions. Check out an old message about clump... https://stat.ethz.ch/pipermail/r-sig-geo/2017-January/025346.html It's sort of cumbersome, but I think the following works... library(raster) r1<-raster(nrows=10, ncols=10, xmn=0, xmx=1, ymn=0, ymx=1) r2<-raster(nrows=10, ncols=10, xmn=0, xmx=1, ymn=1, ymx=2) r1[]<-c(rep(24,24), rep(0,76)) r2[]<-c(rep(0,86), rep(14,14)) par(mfrow=c(3,2)); plot(r1,main="r1"); plot(r2,main="r2") r3<-merge(r1,r2) ; plot(r3, main="merged rasters") r3 <- r3 r4[r4>0]<-38 ; plot(r4, main="desired result") r5 <- clump(r3) ; plot(r5, main = 'clumped') r5[is.na(r5)] <- 0 r5t <- table(getValues(r5)) for (i in names(r5t)[-1]){ id <- as.numeric(i) print(id) r5[r5 == id] <- r5t[i] } plot(r5, main = 'relabeled with area') Cheers, Ben P.S. Years ago I used a language called IDL. It's histogram function has an output keyword called reverse_indices (more here https://www.harrisgeospatial.com/docs/HISTOGRAM.html#H_835179117_677188 and here http://www.idlcoyote.com/tips/histogram_tutorial.html ) - super useful! Your situation just begs for that same solution in R.
On Nov 23, 2018, at 5:13 PM, Hugo Costa <hugoagcosta at gmail.com> wrote: Deal all, I'm merging several rasters representing patches of forest, and the pixel values represent the area of the patches. The patches' area should be summed along the edges of the merged rasters. I provide a reproducible example. Any help is welcome. Hugo library(raster) r1<-raster(nrows=10, ncols=10, xmn=0, xmx=1, ymn=0, ymx=1) r2<-raster(nrows=10, ncols=10, xmn=0, xmx=1, ymn=1, ymx=2) r1[]<-c(rep(24,24), rep(0,76)) r2[]<-c(rep(0,86), rep(14,14)) par(mfrow=c(2,2)); plot(r1,main="r1"); plot(r2,main="r2") r3<-merge(r1,r2) ; plot(r3, main="merged rasters") # do something r3[r3>0]<-38 ; plot(r3, main="desired result") [[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
Ben Tupper Bigelow Laboratory for Ocean Sciences 60 Bigelow Drive, P.O. Box 380 East Boothbay, Maine 04544 http://www.bigelow.org Ecological Forecasting: https://eco.bigelow.org/
Hi Ben, thank you for your help. Clump is a possibility, however I was looking for something different. In my real data, 9 rasters have to be merged (1 central raster and the 8 surrounding rasters), and running clump after merging 9 rasters is slow and memory consuming. And this has to be repeated thousands of times (hopefully in parallel with the rasters in memory). I didn't explain these details to make the message short. I'll continue searching for more alternatives, and all ideas are welcome. Thanks Hugo Ben Tupper <btupper at bigelow.org> escreveu no dia sexta, 23/11/2018 ?(s) 23:58:
Hi, I think you may need to use raster::clump() and table() to label your regions. Check out an old message about clump... https://stat.ethz.ch/pipermail/r-sig-geo/2017-January/025346.html It's sort of cumbersome, but I think the following works... library(raster) r1<-raster(nrows=10, ncols=10, xmn=0, xmx=1, ymn=0, ymx=1) r2<-raster(nrows=10, ncols=10, xmn=0, xmx=1, ymn=1, ymx=2) r1[]<-c(rep(24,24), rep(0,76)) r2[]<-c(rep(0,86), rep(14,14)) par(mfrow=c(3,2)); plot(r1,main="r1"); plot(r2,main="r2") r3<-merge(r1,r2) ; plot(r3, main="merged rasters") r3 <- r3 r4[r4>0]<-38 ; plot(r4, main="desired result") r5 <- clump(r3) ; plot(r5, main = 'clumped') r5[is.na(r5)] <- 0 r5t <- table(getValues(r5)) for (i in names(r5t)[-1]){ id <- as.numeric(i) print(id) r5[r5 == id] <- r5t[i] } plot(r5, main = 'relabeled with area') Cheers, Ben P.S. Years ago I used a language called IDL. It's histogram function has an output keyword called reverse_indices (more here https://www.harrisgeospatial.com/docs/HISTOGRAM.html#H_835179117_677188 and here http://www.idlcoyote.com/tips/histogram_tutorial.html ) - super useful! Your situation just begs for that same solution in R.
On Nov 23, 2018, at 5:13 PM, Hugo Costa <hugoagcosta at gmail.com> wrote:
Deal all,
I'm merging several rasters representing patches of forest, and the pixel
values represent the area of the patches. The patches' area should be
summed along the edges of the merged rasters. I provide a reproducible
example. Any help is welcome.
Hugo
library(raster)
r1<-raster(nrows=10, ncols=10, xmn=0, xmx=1, ymn=0, ymx=1)
r2<-raster(nrows=10, ncols=10, xmn=0, xmx=1, ymn=1, ymx=2)
r1[]<-c(rep(24,24), rep(0,76))
r2[]<-c(rep(0,86), rep(14,14))
par(mfrow=c(2,2)); plot(r1,main="r1"); plot(r2,main="r2")
r3<-merge(r1,r2) ; plot(r3, main="merged rasters")
# do something
r3[r3>0]<-38 ; plot(r3, main="desired result")
[[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
Ben Tupper Bigelow Laboratory for Ocean Sciences 60 Bigelow Drive, P.O. Box 380 East Boothbay, Maine 04544 http://www.bigelow.org Ecological Forecasting: https://eco.bigelow.org/