comparing one raster to a stack and condition

Previous Topic Next Topic
 
classic Classic list List threaded Threaded
3 messages Options
Reply | Threaded
Open this post in threaded view
|

comparing one raster to a stack and condition

karsten
Hi All,
 
I am trying to compare one precipitation raster to a stack of precipitation
raster and would like to create a result raster with a count of how often
the raster value is below those of the stack.
So far I have the following code:
 
----------------------------------------------
 
library(zoo)
library(raster)
# create raster stack for Januaries
alltiffs = list.files(getwd(), pattern="*\\.tif$", full.names=TRUE)

#filter the ones with 01 in it for january
januarygrids =  alltiffs[grepl("*.01.*", alltiffs)]
 
# Create raster stack of grids
r <- stack(januarygrids, quick=TRUE)
 
# set current january layer to compare with
currentmonth <- "es_af.2017.01.tif"
currentmonthraster <- raster(currentmonth)
 
# function to count how oftetn current ratser is below values in stack,
input r and currentmonthraster
belowCurrentRaster <- function(x, y) {
  sum(x > y)
}
 
-------------------------------------------------------
 
Now I thought I could use zApply or overlay to get the count from the
belowCurrentRaster function and write those counts into a result raster.
But I could not figure out how to make this work.
Any ideas appreciated.
 
Cheers
Karsten Vennemann
Terra GIS LTD


        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Reply | Threaded
Open this post in threaded view
|

Re: comparing one raster to a stack and condition

Barry Rowlingson
Here's a way - first let's make some sample data in a stack:

 maker = function(d){raster(matrix(runif(16),4,4))}
 rains = stack(lapply(1:10, maker))

so `rains` is a stack of 10 4x4 rasters with random numbers in. Now the
raster we want to test:

 r1 = maker()

Okay, all set up. We have a raster and a stack, then:

 r1 < rains

is a stack of 1s and 0s where r1 is less than the cell in each layer of
rains, and then we can do:

 sum(r1 < rains)

to total those up.

plot(sum(r1<rains))

should map how many times r1 is less than the values in rains.

Barry



On Wed, Dec 13, 2017 at 9:13 PM, karsten <[hidden email]> wrote:

> Hi All,
>
> I am trying to compare one precipitation raster to a stack of precipitation
> raster and would like to create a result raster with a count of how often
> the raster value is below those of the stack.
> So far I have the following code:
>
> ----------------------------------------------
>
> library(zoo)
> library(raster)
> # create raster stack for Januaries
> alltiffs = list.files(getwd(), pattern="*\\.tif$", full.names=TRUE)
>
> #filter the ones with 01 in it for january
> januarygrids =  alltiffs[grepl("*.01.*", alltiffs)]
>
> # Create raster stack of grids
> r <- stack(januarygrids, quick=TRUE)
>
> # set current january layer to compare with
> currentmonth <- "es_af.2017.01.tif"
> currentmonthraster <- raster(currentmonth)
>
> # function to count how oftetn current ratser is below values in stack,
> input r and currentmonthraster
> belowCurrentRaster <- function(x, y) {
>   sum(x > y)
> }
>
> -------------------------------------------------------
>
> Now I thought I could use zApply or overlay to get the count from the
> belowCurrentRaster function and write those counts into a result raster.
> But I could not figure out how to make this work.
> Any ideas appreciated.
>
> Cheers
> Karsten Vennemann
> Terra GIS LTD
>
>
>         [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Reply | Threaded
Open this post in threaded view
|

Re: comparing one raster to a stack and condition

karsten
Thanks Barry,
 
this is ingenious ! It works really well. For the record below is my code
using this and writing the result to disk.
Karsten
 
 
library(raster)
# Use pattern arg to return a wildcard to get a list of all tifs in dir
alltiffs = list.files(getwd(), pattern="*\\.tif$", full.names=TRUE)
januarygrids =  alltiffs[grepl("*01.tif*", alltiffs)]

# Create raster stack of january rain grids
r <- stack(januarygrids, quick=TRUE)
 
# set current january layer to compare with
currentmonth <- "es_af.2017.01.tif"
currentmonthraster <- raster(currentmonth)
 
# create count for each cell where january rain is below the rain value in
the raster stack layers
resultbelow01 <- sum(r<currentmonthraster)

# write resulttif
writeRaster(resultbelow01, filename='below01.tif', overwrite=TRUE)

  _____  

From: [hidden email] [mailto:[hidden email]] On Behalf Of
Barry Rowlingson
Sent: Donnerstag, 14. Dezember 2017 09:29
To: karsten
Cc: r-sig-geo
Subject: Re: [R-sig-Geo] comparing one raster to a stack and condition


Here's a way - first let's make some sample data in a stack:

 maker = function(d){raster(matrix(runif(16),4,4))}
 rains = stack(lapply(1:10, maker))


so `rains` is a stack of 10 4x4 rasters with random numbers in. Now the
raster we want to test:


 r1 = maker()


Okay, all set up. We have a raster and a stack, then:

 r1 < rains

is a stack of 1s and 0s where r1 is less than the cell in each layer of
rains, and then we can do:

 sum(r1 < rains)

to total those up.

plot(sum(r1<rains))


should map how many times r1 is less than the values in rains.

Barry



On Wed, Dec 13, 2017 at 9:13 PM, karsten <[hidden email]> wrote:


Hi All,

I am trying to compare one precipitation raster to a stack of precipitation
raster and would like to create a result raster with a count of how often
the raster value is below those of the stack.
So far I have the following code:

----------------------------------------------

library(zoo)
library(raster)
# create raster stack for Januaries
alltiffs = list.files(getwd(), pattern="*\\.tif$", full.names=TRUE)

#filter the ones with 01 in it for january
januarygrids =  alltiffs[grepl("*.01.*", alltiffs)]

# Create raster stack of grids
r <- stack(januarygrids, quick=TRUE)

# set current january layer to compare with
currentmonth <- "es_af.2017.01.tif"
currentmonthraster <- raster(currentmonth)

# function to count how oftetn current ratser is below values in stack,
input r and currentmonthraster
belowCurrentRaster <- function(x, y) {
  sum(x > y)
}

-------------------------------------------------------

Now I thought I could use zApply or overlay to get the count from the
belowCurrentRaster function and write those counts into a result raster.
But I could not figure out how to make this work.
Any ideas appreciated.

Cheers
Karsten Vennemann
Terra GIS LTD


        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/
<https://stat.ethz.ch/mailman/listinfo/r-sig-geo> listinfo/r-sig-geo




        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo