Slow Find Loop no vectorisation

Hi,

I have a 3d image where each colour is a unique labelled cell. I'm trying to find the centre of mass of each of these individual colours.

        for colour = 1:2^16
            disp(100*colour/(2^16));
            coords = find(image == colour);
            [y_pos,x_pos,z_pos] = ind2sub(size(image),coords);
            data(colour,:) = mean([x_pos,y_pos,z_pos]);
        end

This loop does the job but it is incredible slow across all 2^16 colours. Is there a way of vectorising this problem?

Thanks,

Craig

Answer 1

[Andrei Bobrov]

Please try it:

co = (1:2^16)';

[l0,i0] = ismember(image1(:),co);
s = size(image1);
[ii jj k] = ndgrid(1:s(1),1:s(2),1:s(3));
coor = [ii(:),jj(:),k(:)];
i1 = (1:numel(image1))';
v = accumarray(i0(l0),i1(l0),[numel(co) 1],@(x){mean(coor(x,:),1)});
out = [co(~cellfun(@isempty,v)), cell2mat(v)];

[Craig Russell]

This works! I don't understand how, but it works.

Answer 2

[Teja Muppirala]

FOR loops are not necessarily slow. This is in R2016a.

%%Make some data...

rng(0);
image = randi(2^16,[30 40 50]);

%% 1. Original method takes 7 seconds. (By the way, it should be MEAN([...],1) not just MEAN([...])

tic
numColors = 2^16;
data = zeros(numColors,3);
for colour = 1:numColors
            %disp(100*colour/(2^16));
            coords = find(image == colour);
            [y_pos,x_pos,z_pos] = ind2sub(size(image),coords);
            data(colour,:) = mean([x_pos,y_pos,z_pos],1);
end
toc    

% Elapsed time is 7.447380 seconds

%% Using a for loop more efficiently takes less than 0.1 seconds.

tic
data2 = zeros(numColors,4);
for n1 = 1:size(image,1)
    for n2 = 1:size(image,2)
        for n3 = 1:size(image,3)
            data2(image(n1,n2,n3),:) = data2(image(n1,n2,n3),:) + [1 n2 n1 n3];
        end
    end
end

data2 = bsxfun(@rdivide,data2(:,2:4),data2(:,1));
toc

% Elapsed time is 0.091517 seconds

%% Show they are the same.

isequaln(data,data2)

% ans =
%
%     1

[Andrei Bobrov]

+1

Answer 3

[Guillaume]

I don't think you can get rid of the colour loop, but you can speed up the loop code with:

[y_pos, x_pos, z_pos] = ndgrid(1:size(image, 2), 1:size(image, 1), 1:size(image, 3)); 
for colour = 1 : pow2(16)
  iscolour = image == colour; %use logical instead of find
  data(colour, :) = mean([x_pos(iscolour), y_pos(iscolour), z_pos(iscolour)]);
end

edit: I've just thought that you can replace the loop with accumarray and I see that Andrei answered with that in the meantime.

edit-edit: Here is a (more readable?) alternative to andrei's answer:

[y_pos, x_pos, z_pos] = ndgrid(1:size(image, 2), 1:size(image, 1), 1:size(image, 3)); 
datarows = [pow2(16), 1];
data = [accumarray(image(:), x_pos(:), datarows, @mean), ...
        accumarray(image(:), y_pos(:), datarows, @mean), ...
        accumarray(image(:), z_pos(:), datarows, @mean)];