Raster moving window analysis

I am quite new to Matlab so I apologize for being naive and am looking for guidance on what direction I should go in/where to look. I am hoping to use some sort of a moving window analysis on a raster (.tif) that is an output of a model of elevation data. With that in mind I would like to calculate local summary statistics (e.g. standard deviation) of the raster within a moving window of a user defined size (e.g. 20 m) and shape (e.g rectangle, square, circle) and then move that window at some increment in the x and y on the raster that is also user defined (e.g. 10 m). All of these values would then be output into a new raster that I could use for statistical analysis.

My thought was to overlay some sort of a grid onto the raster image that is the size of the desired increment so that each xy intersect/node is a point where the window stops to calculate the summary value. This will ultimately reduce the resolution of the original raster because the output raster resolution will be equal to the user defined x y increment (e.g. 10mm). I've looked into Blockproc, but that isn't quite what I'm looking for since it functions on the pixel scale and only works as a rectangle. It also appears that movstd is similar in this way. Any help would be greatly appreciated!

Answer 1

[Nate]

Hi Brian,

Your plan overall seems sound.

Blockproc should work well for this procedure, or at least be able to handle the square and rectangle case. In order to do this however, you will need to convert the distance scale to pixels when you use the function.

I am curious about the use case for circle, however. Circles do not tessellate so I would expect that you would lose some information at random. If you still want to pursue this, it would be possible to write a function that iterates over the TIFF data, taking blocks of pixels in a circular shape each time.

One way to get this data is by using the "read" function of the Tiff class:

Tiff read documentation

Let me know if you have any questions.

[Brian Danhoff]

Hi Nate,

I appreciate your help with this. I think you are right in that the circle would potentially lose some information (along the edges, correct?). I played with Blockproc today and visually, it appears to perform the calculation I am looking for. However, I realized that it is not actually producing the output I am looking for. It wasn't until I used the example image "moon.tif" from the Blockproc help page (https://www.mathworks.com/help/images/ref/blockproc.html) that I realized the there are no fewer pixels in my output.

I believe I was misunderstanding what Blockproc does. My initial thought was that it applies some function to all pixels within a window and then outputs one value, but in reality it appears that blockproc applies the calculated function value to all the pixels within the window and keeps all those pixels. For example, the moon.tif image is a 537x358 pixel image. After using the example provided on the help page (standard deviation within a 32x32 window) the output is still a 537x358 pixel image.

Would this require some type of resampling at an interval equal to my block spacing to get the output I am looking for? Or is there some other tool that might be more appropriate for this.

Thank you for taking the time to help with this!

[Nate]

I was assuming that the circles would not overlap, so that the spaces in between them are never actually inspected. This image from Wolfram may help clarify:

Circle Packing Image

If your goal is for the output image to be of smaller resolution, you can simply use the "imresize" function during or after the "blockproc" call. Calling "imresize" afterwards may be simpler.