Reduce noise in a raster (DEM) with R

Question

I have several DEMs obtained from NASA-SRTM (https://www2.jpl.nasa.gov/srtm/).

I try reduce the noise infor these applying the function noise.filter ofwvtool package and receive the next message:

> library(raster)
> DEM <- raster("S04W075.hgt")
> DEMnoise <- noise.filter(DEM, n= 5, method = "gaussian")
Warning message:
In noise.filter(DEM, n = 5, method = "gaussian") :
  data must be grayscale image
> DEMnoise
[1] "data must be grayscale image"

I have to process several rasters, so I need the solution in R.

Answer 1

This error is resulting from an error check within the function.

if (length(dim(x)) > 2) {
  warning("data must be grayscale image")
}

Given this condition length(dim(x)), a single band rasterLayer class object will always return 3 representing row, column and nlayers dimensions.

Whereas this function, in theory, operates on a matrix object it is not quite as simple as it would seem. There is attribute encoding that is necessary for the function to run. So, it is not as simple as just coercing a raster to a matrix. You can see this attribution in the example data where the class is matrix but there are additional attributes included. The function at hand requires these attributes to run and I do not see a function in the package that calculates/adds these attributes for a matrix object.

library(wvtools)
data(camphora)
  class(camphora)
  attributes(camphora)

I dropped these attributes in the camphora object and found that you receive an error. You can brute force the function by manually adding the require attributes. You could also just use a Guassian kernel for smoothing within the raster focal function. Something that I notice in comparing results is that the noise.filter function forces an 8-bit output.

library(raster)
library(spatialEco)
library(wvtool)

This coerces the camphora matrix to a raster object and then to a matrix, adding the require attributes for the noise.filter function.

data(camphora) 
r <- raster(wvtool::crop(camphora,200,200))
 rmat <- as.matrix(r)
   attributes(rmat) <- list( dim=c(nrow(rmat),ncol(rmat)),  
                             bits.per.sample=8, 
                             samples.per.pixel=1)

We then run the noise.filter and focal with gaussian.kernel functions. Since the sigma of the Gaussian kernel is not an argument in the noise.filter function it is difficult to directly compare but, the out rescaling to 8-bit is a little problematic. Not sure if this is something controllable by changing $bits.per.sample but, since this is not indicated in the functions help it is a moot point.

( r.filt <- raster(noise.filter(rmat,3,"median")) )  
( r.smooth <- focal(r, w=gaussian.kernel(sigma=2, n=3), 
                    fun=median) )

And, plot results.

par(mfrow=c(2,2))
 plot(r, main="original")
 plot(r.filt, main="wvtool") 
 plot(r.smooth, main="gauss") 

Answer 2

You can try the function 'blur' in the package 'spatstat' to apply a kernel based blur (Gaussian or otherwise) to your image. You can perform multiple iterations or adjust the parameters of the function until you achieve the desired result.

library(spatstat)
DEM_blurred <- blur(DEM, sigma=0.5, kernel="gaussian")
#adjust sigma value to change the dimensions of your smoothing kernel

https://rdrr.io/cran/spatstat/man/blur.html#heading-1

Also I recommend reading the answer by @whuber on this post:What raster smoothing/generalization tools are available?