data manipulation

Hello All,

I have got a slightly tricky matlab question that I would like to seek help from everyone out there.

I have got a logic with the size of 143x143x143. It was derived from a cubic bone sample with the size of 143x143x143 voxels (data points). It contains both "0s" and "1"s. "0s" represents marrow voxels and "1s" represents bone voxels.

What I am hoping to do is to conduct so called "resolution analysis". I will explain in detail.

Step 1: The first step that I would like to do is to determine the "overall" transit time of each "2x2" square box (slab) using the following formula along the Z axis from Z = 1 to Z = 142.

If in this 2x2 square box, I have got "1100" combination. The transit time is then going to be 0.5.

If the combination is "1111", the transit time is going to be 1

If the combination is "0000", the transit time is going to be 0

If the combination is "1110", the transit time is going to be 0.75 (the order of 1s and 0s do not matter)

If the combination is "0001", the transit time is going to be 0.25 (the order of 1s and 0s do not matter)

Step 2: Once I have done that, I would then like to determine the "final" transit time of this "2x2x142" chunck by adding all of the "overall transit" time that I have calculated in step 1 together.

Step3: Once I have done this, I would like to move on to the next "2x2x142" chunck to repeat the same process. It does not matter whether I have moved on to the one next to the "2x2x142" chunck that I have done or the one below. However, what is your suggestions on once I have reached to the boundry, can I ask the computer to move a level down or to the next line to repeat (step 1, 2 & 3)?

What I am hoping to achieve is a new matrix with "71x71 = 5041" numbers of transit time for the 2x2 resolution analysis.

Since 143/2 does not equal to a whole number? Would it possible just to treat the size of the cubic as "142x142x142"? And every single time when the size of the voxel cannot be divided by 2, I will just do (the size number -1) in order to be able to be divided by 2?

Once again, any suggestions would be greatly appreciated.

Thanks a lot in advance.

John

[Oleg Komarov]

Do you need step 1 just to calculate step 2 or for other purposes as well?

[John]

Hi Oleg,

Step 1 is essentially the boundary conditions for step 2.

Cheers!

John

Answer 1

[Teja Muppirala]

This is not so hard, and is a nice example of why MATLAB is so great for working with arrays.

% I'll make a random 143x143x143 matrix

M = randi([0 1],143,143,143,'uint8');

% Clip the dimensions to be even

S = size(M);
M2 = M(1:2*floor(S(1)/2) , 1:2*floor(S(2)/2), 1:2*floor(S(3)/2));

% Add everything in the z direction, and divide by 4.

M3 = sum(M2,3)/4;

%Do the averages

I = conv2(M3,[1 1; 1 1],'valid');

%Take every other row to get each 2x2 block

I = I(1:2:end,1:2:end);

%And "I" is the 71x71 matrix you were looking for!

[Teja Muppirala]

Or even simpler, use CONVN!

M = randi([0 1],143,143,143,'uint8');

S = size(M);

M2 = M(1:2*floor(S(1)/2) , 1:2*floor(S(2)/2), 1:2*floor(S(3)/2));

I = convn(M2,ones(2,2,size(M2,3))/4,'valid');

I = I(1:2:end,1:2:end);

[Teja Muppirala]

Sorry, one more redo...

M = randi([0 1],143,143,143,'uint8');

I = convn(M,ones(2,2,2*floor(size(M,3)/2)),'valid');

I = I(1:2:end,1:2:end,1)/4;