Resize 3D binary image by averaging

Hi, I want to reduce the size of a binary image stack by averaging the values in 3D. For example if the image is 1260*1260*6. I wrote the code below to process the image to average the values in every 4*4*2 cube thereby creating a new matrix that is 315*315*3 in size.

X(1:315)=4;

Y(1:315)=4;

Z(1:3)=2;

A=mat2cell(imagestack_1,[X],[Y],[Z]);

average=cellfun(@mean2,A);

The problem with this code is that its very slow. Please suggest a faster way. My actual image is super big 2000*2000*2000, so please consider the size when suggesting solutions.

Also, the object in the binary image whose size I want to reduce is cylindrical in shape (like a bunch of circles with equal size stacked together). Since I don't know how to average only the values in the circle, I am first cropping the largest square from the circle and then performing the average. If anyone can think of a way to average only the values in the circle.

I am not using imresize or interp3 because I don't any kind of interpolation, I just want a simple arithmetic average as 1 in the object represents solid and 0 represents void space, so I want to know the mean of solid and void space.

[Walter Roberson]

How much memory do you have available?

2000*2000*2000 is at least 8 gigabytes of uint8 (or logical). If you were totaling over a 4 x 4 x 2 block then that would be a maximum count of 32, would could fit into uint8, requiring about 250 megabytes. If you need to convert that into floating point instead of totals, you would have to go to at least single precision, 4 bytes per location, which would take you up to a gigabyte again.

Some of the faster calculation routines would end up converting to single or double along the way, which could be a massive memory problem for you, as it would take you up to 64 gigabytes for your image.

I can think of a vectorized mechanism that would do the additions fairly efficiently, but the indexing involved in it could take up to 2000^3 * 8 bytes * 24 copies if not carefully managed. I will wait for you to indicate your available memory before I work on algorithms that make the best trade-off between speed and memory.

[Xiaoli]

Hey, thanks, here are the memory stats from matlab

>> memory

Maximum possible array: 62253 MB (6.528e+10 bytes) *

Memory available for all arrays: 62253 MB (6.528e+10 bytes) *

Memory used by MATLAB: 811 MB (8.501e+08 bytes)

Physical Memory (RAM): 32692 MB (3.428e+10 bytes)

• Limited by System Memory (physical + swap file) available.

Answer 1

[Image Analyst]

To convolve in 3-D, which can compute a local average, use the convn() function.

[Xiaoli]

Hi, thanks. I am not familiar with convolution. I tried using the convn function but the averages are overlapping. Can you please provide a simple example like the code I wrote above without the averages overlapping.

[Image Analyst]

You can use blockproc() if you want to move in "jumps" without overlapping. Sorry I can't attach m-files until the Mathworks fixes the bug they introduced where you can no longer attach m-files.

[Walter Roberson]

blockproc() is only documented for 2D images, not for 3D volumes. Testing briefly, I see that you can pass in a volume, but only two parameters for the block size, and that it processes all of the data planes at once. You can output a 3D array.

This suggests that you could at least blockproc to get 4 x 4 x 2000 volumes, and your routine to do the averaging could then produce 1000 results

blockproc(data, [4 4], @(block) uint8( block.data(:,:,1:2:end) + block.data(:,:,2:2:end) ) )

This would probably not be the most efficient arrangement, but it should be functional.

[Xiaoli]

Ok, thank you both!

[zahra ghadery]

hi

one time you asked this question

i have a bunch (thousands) of 2D CT scan slices(grayscale) of a single rock core.The slices are from the top of the core to the bottom. I want to combine the images into one 3D image to create a digital core

did you understad the solution?

please help me if you know

[Image Analyst]

Sorry, I don't remember asking that question.

But you asked one. But it doesn't seem related to this one, so ask it in a new question after you read this link.