raster icon indicating copy to clipboard operation
raster copied to clipboard

Using Variable Window Sizes within Focal Function

Open azh3 opened this issue 5 years ago • 1 comments

So I am attempting to implement a maximum filter on a CHM (Canopy Height Model) with a variable window size here, 'Tree_Heights' is a data frame containing the fitted values of a linear regression model to predict tree diameter, the window size is supposed to match the lower prediction interval for the radius of a tree at a given height.

RoundOdd <- function(x) {2*floor(x/2)+1}

fitlwr <- function(x){ RoundOdd(Tree_Heights[Tree_Heights$Tree_Height == x, "fit.lwr"]/2) }

m <- raster::focalWeight(x = CMM, d = fitlwr(), type = "circle")

CMM <- raster::focal(x = CMM, w = m, fun = max)

This returns the following error:

Error in [.data.frame(Tree_Heights, Tree_Heights$Tree_Height == x, "fit.lwr") : argument "x" is missing, with no default 6.[.data.frame(Tree_Heights, Tree_Heights$Tree_Height == x, "fit.lwr") 5.Tree_Heights[Tree_Heights$Tree_Height == x, "fit.lwr"] 4.RoundOdd(Tree_Heights[Tree_Heights$Tree_Height == x, "fit.lwr"]/2) 3.fitlwr() 2..circular.weight(x, d[1]) 1.raster::focalWeight(x = CMM, d = fitlwr(), type = "circle")

If I try instead to use the function in the argument for window size, I get this error:

Error in .local(x, ...) : is.matrix(w) is not TRUE 5. stop(simpleError(msg, call = if (p <- sys.parent(1L)) sys.call(p))) 4. stopifnot(is.matrix(w)) 3. .local(x, ...) 2. raster::focal(x = CMM, w = fitlwr, fun = max)

  1. raster::focal(x = CMM, w = fitlwr, fun = max)

Is there no way to use a variable window size inherently in this function? Are there any other options for doing this in the Raster package?

azh3 avatar Oct 30 '20 01:10 azh3

I think I found a temporary

First I converted my raster to a matrix, removed NAs, and round the values (this is unique to my use case and not necessary for the algorithm to function.

    X <- raster::as.matrix(Z)
    X <- round(X, digits = 0)
    X[is.na(X)] <- 0

This is to calculate the maximum for a variable size rectangular moving window:

Y <- X
  for (i in 1:nrow(X)){ 
     for (j in 1:ncol(X)){ 
        N <- fitlwr(X[i,j])
        Y[i,j] = max(X[max(1, i-N):min(nrow(X), i+N), max(1, j-N):min(ncol(X), j+N)]) 
    }
  }

fitlwr() is a custom function that calls a linear model that matches the value of a cell to the expected radius of the moving window.

And here is for a circular moving window:

   Y <- X for (i in 1:nrow(X)){     for (j in 1:ncol(X)){ 
          N = fitlwr(X[i,j])
          M = X[max(1, i-N):min(nrow(X), i+N), max(1, j-N):min(ncol(X), j+N)]
          W = reshape2::melt(M)
          W$d2 = sqrt((W$Var1-mean(W$Var1))^2 + (W$Var2-mean(W$Var2))^2)
          Y[i,j] = max(X[i,j], max(subset(W, d2 <= N, select = value)))}

I then write the values back to my raster, this maintains the CRS, projection, etc.

raster::values(Z) <- Y

I imagine this would be more efficient if it could be implemented within the raster package itself and think there are many potential use-cases for being able to specify a variable window size within the focal function.

azh3 avatar Nov 03 '20 07:11 azh3