Random non-overlapping circles in MATLAB
For a recent project I had to generate a stimulus set consisting of random non-overlapping circles in MATLAB. My solution was to keep a list of the outside of a gradually increasing outer ring of circles, and add new ones, expanding the ring. This requirement is similar to circle-packing, except I needed to have gaps between the circles.
If you use especially large gaps, or the ratio between your maximum and minimum circle sizes is great, I detected the occasional overlap. But for the parameters I needed, I tested over 10000 tests of 400 circles, with no overlaps.
Here is how it works:
% Parameters for the circles setup.min_circle_radius = 15; setup.max_circle_radius = 20; setup.min_gap = 10; setup.max_gap = 20; setup.n_circles = 400;
First I set up some parameters - minimum and maximum circle radiuses and gaps, plus the number of circles I want to generate.
function [ circles ] = CirclePacking( setup ) % Start with a circle in the middle first_circle = generate_random_radius_circle( setup ); first_circle.position = [ 0, 0 ]; % Add second circle at random angle second_circle = generate_random_radius_circle( setup ); gap = randi([setup.min_gap, setup.max_gap]); total_distance = first_circle.radius + gap + second_circle.radius; [x, y] = pol2cart(2 * pi * rand, total_distance); second_circle.position = [x, y];
Next I begin with two circles, one in the centre, and another at a random angle from the first.
function [circle] = generate_random_radius_circle( setup ) circle.radius = randi([setup.min_circle_radius, setup.max_circle_radius]); end
I represent circles as structs, first defined solely by radius, but then by position. All circles are first generated using this function which makes a random radius.
Back to the main function:
first_circle.before = 2; first_circle.after = 2; second_circle.before = 1; second_circle.after = 1; % Setup doubly linked list first_circle.id = 1; circles(1) = first_circle; % have to set first circle manually circles = save_circle(second_circle, circles);
I set up a doubly linked list, that is, each circle keeps track of the circle before and after it. With just two, before and after are the same.
function [circles, circle] = save_circle(circle, circles) if (~isfield(circle, 'id')) circle.id = size(circles,2) + 1; end circles(circle.id) = circle; end
'save_circle' adds or saves a circle to the main circles array. I couldn't work out after a google search how to initialise an empty array of structs, so I didn't waste time doing so.
Back to main function:
current_position = 1; for i = 1:setup.n_circles [circles, current_position] = add_new_circle( setup, circles, current_position ); end
To finish the main function, I just call 'add_new_circle', maintaining a variable 'current_position' which points to a particular circle in the circles array.
This is the end of the main function, most of what happens is going on inside 'add_new_circle', but before I go into it, I will show a few more important helper functions:
function circles_overlap = circles_overlap(circle1, circle2, min_gap) dist = edist(circle1.position, circle2.position); circles_overlap = dist < circle1.radius + circle2.radius + min_gap; end function edist = edist(a, b) % euclidean distance edist = sqrt((a(1)-b(1))^2+(a(2)-b(2))^2); end
The above functions are self explanatory, just checking whether circles overlap, and a function for euclidean distance.
function next_circle = next_circle(circle, circles) next_circle = circles(circle.after); end function previous_circle = previous_circle(circle, circles) previous_circle = circles(circle.before); end
The above functions allow me to get the next and previous circles - this is not necessary really, but I find it's easier on the eyes using these functions.
function new_circle = poisiton_new_circle(new_circle, circle1, circle2, gap1, gap2) distance1 = new_circle.radius + gap1 + circle1.radius; distance2 = new_circle.radius + gap2 + circle2.radius; distance_between = edist(circle1.position, circle2.position); % Calculate angle to new circle angle = get_law_cosines_angle_a(distance1, distance2, distance_between); % Calculate angle of first two circles vector_first_to_second = circle2.position - circle1.position; angle_first_to_second = atan2(vector_first_to_second(2), vector_first_to_second(1)); angle_of_new_from_first = mod(angle_first_to_second + angle, 2 * pi); [x, y] = pol2cart(angle_of_new_from_first, distance1); new_circle.position = [x, y] + circle1.position; end function angle = get_law_cosines_angle_a(ab,bc,ca) num = (ab * ab) + (ca * ca) - (bc * bc); den = 2 * ab * ca; angle = acos(num / den); end
The 'position_new_circle' function puts a circle next to circle1 and circle2, taking into account the gaps passed in too. It always results in the anticlockwise order circle1, new_circle, circle2. It is using the law of cosines to do this - I won't explain the maths here.
function [circles, first_circle_index] = add_new_circle( setup, circles, first_circle_index ) % circles must contain 2 or more circles gap_with_first = randi([setup.min_gap, setup.max_gap]); gap_with_second = randi([setup.min_gap, setup.max_gap]); first_circle = circles(first_circle_index); second_circle = next_circle(first_circle, circles); new_circle = generate_random_radius_circle(setup); % Check backwards for overlaps for i = 1:8; new_circle = poisiton_new_circle(new_circle, first_circle, second_circle, gap_with_first, gap_with_second); before_circle = previous_circle(first_circle, circles); number_to_check = min(size(circles,2) ,8); overlaps = false; for j=1:number_to_check % Check if circle overlaps if (before_circle.id == second_circle.id) break; end if (circles_overlap(new_circle, before_circle, setup.min_gap)) first_circle = before_circle; first_circle_index = first_circle.id; overlaps = true; break; end before_circle = previous_circle(before_circle, circles); end if (~overlaps) break; end end % Update linked list circles = update_list(circles, first_circle, second_circle, new_circle); end
This is the main part of the algorithm that works out where to put a new circle. It starts by setting the gaps, and taking the key circle, and the circle after that. 'first_circle_index' indexes the key circles, all circles are generated around this circle, until overlaps. Often, a new circle will overlap one of the previous circles in the chain. When this occurs I make the overlapping circle the new key circle, and generate circles around that.
I found that checking 8 circles back was more than enough (on the idea that you can fit roughly 6 circles around another circle - if they are similar size) - but this can be modified to check all circles, at a performance hit.
So after positioning a circle between the key circle and the next one, I check circles before the key (making sure I don't check the same circle). If the circles overlap, I set the key circle to the one that overlaps.
If there are no overlaps, then I update the linked list with that circle.
function circles = update_list(circles, circle1, circle2, new_circle) new_circle.before = circle1.id; new_circle.after = circle2.id; [circles, new_circle] = save_circle(new_circle, circles); circle1.after = new_circle.id; circles = save_circle(circle1, circles); circle2.before = new_circle.id; circles = save_circle(circle2, circles); end
Note this could be massively optimised up - but remembering that 'premature optimization is the root of all evil' - I have not yet needed to so.