MATLAB/Octave: cut a lot of circles from a image

Question

I have an matrix (image) and information about interesting part within circles (center corrdinates and radii given). I want to cut for all the circles the parts of the matrix in order to do some more calculations for each circle. Or at least I want to have a bitmask with all the circle.

I use Octave (but could also use MATLAB but it would be difficult because of licence iusses) and have the following script with some hints from stackoverflow. I have information of 20 circles and it takes about 0.7 s on my Core i5 using Octave:

% image
dim_x = 1000;
dim_y = 1000;
A=rand(dim_x,dim_y);

% center positions and ...
c = [222 111; 878 112; 81 718; 89 112; 222 111; 878 112; 81 718; 89 112; 222 111; 878 112; 81 718; 89 112; 222 111; 878 112; 81 718; 89 112; 222 111; 878 112; 81 718; 89 112];
%...  radii of the circles
r = [10 33 55 2 22 10 33 55 2 22 10 33 55 2 22 10 33 55 2 22];

tic;
for i=1:size(c,1)
    % create a bitmask ...
    mask = bsxfun(@plus, ((1:dim_y) - c(i,1)).^2, (transpose(1:dim_x) - c(i,2)).^2) < r(i)^2;
    % ... cut the circles out of the image
    B=A.*mask;
end;
toc;

Do you know a more performant solution since I want to have about 600 circles.

Thanks in advance

Answer 1

Try

mask = bsxfun(@lt, ((1:dim_y) - c(i,1)).^2,  r(i)^2 - ((1:dim_x).' - c(i,2)).^2);

According to my MATLAB profiler, this is about 4 times faster than your version. Also, the B = A.*mask line takes about the same amount of time as the original mask = ... line. Not sure there's much you can do about that.

Answer 2

You may want to look into Matlab's strel (not sure about Octave availability, in Matlab it is part of the image processing toolbox).

radius = 10;
center = [320 240];
nn = 0; 
se = strel('disk', radius, nn);
px = se.getneighbors;
px = px + repmat(center, [length(px) 1]);

The nn parameter effects performance. Making it 4, 6, or 8 will improve performance at the cost of your mask not being precisely a circle.

You can also squeeze some performance out of it by rewriting the repmat bit using a bsxfun.

Answer 3

There's several things you can do to make your code more efficient, though some of them depend on exactly what you want in the end, and on what assumptions you can make about the circles (e.g. can they overlap? Are there many similar radii?).

Below is a solution that assumes that there is very little repetition among radii, and center coordinates are always integer pixel values.

%# image
dim_x = 1000;
dim_y = 1000;
A=rand(dim_x,dim_y);

%# center positions and ...
c = [222 111; 878 112; 81 718; 89 112; 222 111; 878 112; 81 718; 89 112; 222 111; 878 112; 81 718; 89 112; 222 111; 878 112; 81 718; 89 112; 222 111; 878 112; 81 718; 89 112];
%#...  radii of the circles
r = [10 33 55 2 22 10 33 55 2 22 10 33 55 2 22 10 33 55 2 22];

%# find the largest circle...
rMax = max(r);
%#... and create a distance array
distFromCenterSquared = bsxfun(@plus,(-rMax:rMax).^2,transpose(-rMax:rMax).^2);

%# now we can loop over the radii to create the logical mask for all circles
mask = false(dim_x,dim_y); %# initialize inside the loop if you want one circle at a time
for i=1:length(r)

    %# create logical mini-circle mask
    miniMask = distFromCenterSquared(rMax-r(i)+1:end-(rMax-r(i)),rMax-r(i)+1:end-(rMax-r(i)))...
       < r(i)^2;

    %# add to the mask. The ranges need to be fixed, obviously, if 
    %# circles can be only partially inside the image
    %# also, the "or" is only necessary if you're adding to
    %# a mask, instead of recreating it each iteration
    mask(c(i,1)-r(i):c(i,1)+r(i),c(i,2)-r(i):c(i,2)+r(i)) = ...
       mask(c(i,1)-r(i):c(i,1)+r(i),c(i,2)-r(i):c(i,2)+r(i)) | ...
       miniMask;

end

By the way: If you have non-overlapping circles, you can use bwlabel after the loop (or use find and sub2ind to write i into the individual circles), so that you can process all circles in one go using accumarray.

Answer 4

I am going to suggest using POLY2MASK function from the MATLAB Image Processing Toolbox (also available in the Image package for Octave). Check the "algorithm" section to see how it handles the discrete pixels.

Here is an example to test performance:

%# image
dim_x = 1000;
dim_y = 1000;
A = rand(dim_x,dim_y);

%# center positions and radii of the circles
c = [222 111; 878 112; 81 718; 89 112; 222 111; 878 112; 81 718; 89 112; 222 111; 878 112; 81 718; 89 112; 222 111; 878 112; 81 718; 89 112; 222 111; 878 112; 81 718; 89 112];
r = [10 33 55 2 22 10 33 55 2 22 10 33 55 2 22 10 33 55 2 22];

%# lets make them 600 circles
c = repmat(c,30,1);
r = repmat(r,1,30);

%# zero-centered unit circle
t = linspace(0,2*pi,50);
ct = cos(t);
st = sin(t);

%# compute binary mask for each circle
tic
for i=1:numel(r)
    %# scale and shift scale circle, and use to get mask
    BW = poly2mask(r(i).*ct + c(i,1), r(i).*st + c(i,2), dim_x, dim_y);

    %# use the mask ...
end
toc

On my laptop, this finishes in:

Elapsed time is 4.864494 seconds.