Skip to content

How to upscale (reduce the spatial resolution) an image using a Gaussian filter in R?

4 messages · Nikolaos Tziokas, Michael Sumner, Greg Snow

#
Hi,

I want to resample a raster from 15m to 460m (pixel size) using a Gaussian
filter.

*The goal*

I am having a coarse image which I want to downscale. I also have a fine
resolution band to assist the downscaling. The downscaling method I am
using is called geographically weighted area-to-point regression Kriging
(GWATPRK). The method consists of two steps:

   1. GWR
   2. area-to-point Kriging on the GWR's residuals

In order to perform GWR using raster data, those needs to have the same
pixel size. This means that, my fine resolution image needs to be upscaled
to match the spatial resolution of the coarse band. This upscaling of the
fine band needs to be done using a Gaussian kernel with ?=0.5(i.e., the
PSF).

How can I upscale (reduce the spatial resolution) a satellite image using a
Gaussian kernel (i.e., point spread function)?

For reference, I am following the paper The effect of point spread function
on downscaling continua
<https://www.sciencedirect.com/science/article/pii/S092427162030229X> where
the authors at p.253 in Eq (9) mention:

*the coarse image produced by upscaling the corresponding fine band k using
a PSF.*

I googled how I can achieve that but unfortunately I couldn't find any
solution. So to do this, how can I use this Gaussian filter to change the
resolution (pixel size) of my image?

Here is the image
<https://drive.google.com/drive/folders/18_1Kshb8WbT04gwOw4d_xhfQenULDXdB?usp=sharing>
I
am trying to convolve.


Many thanks
#
Hi,  I would construct the kernel and push it through terra::focal.

It should be easy enough to find the right kernel matrix but (below, as a
kind of strange alt -  to give you something to compare with) -  you might
invoke gauss resampling with overviews:

Make a dummy copy with GDAL

gdal_translate pan15.tif dummy.tif -outsize 200% 200%

then add overviews with gauss resampling, then your 2nd level will have the
same resolution as your input tif

gdaladdo dummy.tif -r gauss

then extract that level (you can't do this with gdal_translate until GDAL
3.6 but it can be arranged in R, and if you gdal_translate to the exact
outsize of the 2nd level you probably will get it exactly):

You otherwise can't control the parameters  with the app, but you possibly
could with a computed band in VRT.

Cheers, Mike

On Wed, Sep 7, 2022 at 2:46 AM Nikolaos Tziokas <nikos.tziokas at gmail.com>
wrote:

  
    
#
The `image_resize` function in the `magick` package may be a simple
option.  A gaussian filter is one of the options.  This uses the
ImageMagick program behind the scenes, but is a quick and easy
interface within R for images (raster data).

On Tue, Sep 6, 2022 at 10:46 AM Nikolaos Tziokas
<nikos.tziokas at gmail.com> wrote:

  
    
17 days later
#
I managed to solve the issue using the OpenImageR package and the function
down_sample_image.

library(raster)
library(OpenImageR)

r = raster("path/tirs.tif")

m = as.matrix(r)

psf = down_sample_image(m,
                        factor = 4.6, # zoom factor
                        gaussian_blur = T,
                        gauss_sigma = 0.5 * 100) # sigma * pixel size

e <- extent(r)

m2r <- raster(psf)
extent(m2r) <- e
raster::crs(m2r) <- "EPSG:7767"
res(m2r)

writeRaster(m2r, "path/tirs460.tif")

The function seems to do exactly what the creator of the area-to-point
regression Kriging, Qunming Wang, did with in his Matlab code (
https://github.com/qunmingwang/Code-for-S2-fusion/blob/master/PSF_template.m).
That is, he multiplied the width of the PSF with scale factor.


???? ???, 7 ??? 2022, 15:46 ? ??????? Greg Snow <538280 at gmail.com> ??????: