Location of centers of Circular regions

I have a square region digitized to 40x40 grid. Each cell is denoted either by 1 or by 2. The 1's show circular regions for which I would like to determine the location of center. for e.g. 12 units to left and 9 units from bottom up. and so on. The region has only 4 cicular subregions and it gets a bit complicated as these try to cross the border. They get mirrored to the opposite side. I am attaching a figure so that this gets clear. Any help or hint to program this is appreciated.

[darova]

Can you use bwlabel and regionprops?

[Pappu Murthy]

I can. I have the image tool box. But not sure how I use them to find centers

Answer 1

[Image Analyst]

Please attach your original image, not some colorized screenshot with gridlines and a huge white surround. Also, you do not have any circles in there. Well, maybe one that's roughly circular in the middle, but the other ones are not circles, but only partial circles. If you want an estimate of the centroids of those, as if you had the complete blob, you'd have to get the perimeter that is in the image (and not on the border) using bwboundaries() and pass it into the code to estimate the circle from the FAQ:

https://matlab.fandom.com/wiki/FAQ#How_can_I_fit_a_circle_to_a_set_of_XY_data.3F

It's trivial for either you or I to do that once we have the actual image, not a pseudocolored representation of it. Here's a start using the screenshot image you attached. Obviously there is some stuff about cropping and color segmentation that you won't have to do if you simply start with the original, non-colorized, non-screenshot version of your image.

% Demo by Image Analyst.
clc;    % Clear the command window.
close all;  % Close all figures (except those of imtool.)
clearvars;
workspace;  % Make sure the workspace panel is showing.
format long g;
format compact;
fontSize = 16;
fprintf('Beginning to run %s.m ...\n', mfilename);

%-----------------------------------------------------------------------------------------------------------------------------------
% Read in image.  This is a horrible image.  NEVER use JPG format for image analysis.  Use PNG, TIFF, or BMP instead.
folder = [];
baseFileName = 'example.png';
fullFileName = fullfile(folder, baseFileName);
% Check if file exists.
if ~exist(fullFileName, 'file')
	% The file doesn't exist -- didn't find it there in that folder.
	% Check the entire search path (other folders) for the file by stripping off the folder.
	fullFileNameOnSearchPath = baseFileName; % No path this time.
	if ~exist(fullFileNameOnSearchPath, 'file')
		% Still didn't find it.  Alert user.
		errorMessage = sprintf('Error: %s does not exist in the search path folders.', fullFileName);
		uiwait(warndlg(errorMessage));
		return;
	end
end
rgbImage = imread(fullFileName);
[rows, columns, numberOfColorChannels] = size(rgbImage)
% Display the test image.
subplot(2, 2, 1);
imshow(rgbImage, []);
axis('on', 'image');
hp = impixelinfo(); % Set up status line to see values when you mouse over the image.
caption = sprintf('Image : "%s", Final Mask Image : %d rows by %d columns', baseFileName, rows, columns);
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
drawnow;
hp = impixelinfo(); % Set up status line to see values when you mouse over the image.

% Set up figure properties:
% Enlarge figure to full screen.
hFig1 = gcf;
hFig1.Units = 'Normalized';
hFig1.WindowState = 'maximized';
% Get rid of tool bar and pulldown menus that are along top of figure.
% set(gcf, 'Toolbar', 'none', 'Menu', 'none');
% Give a name to the title bar.
hFig1.Name = 'Demo by Image Analyst';

% Create a binary image of the green spots.
threshold = 128;
binaryImage = rgbImage(:, :, 1) <= threshold & rgbImage(:, :, 2) >= threshold & rgbImage(:, :, 3) <= threshold;
% Crop away the huge white frame that surrounds the actual ROI:
binaryImage = binaryImage(41:466, 150:576);
% Update the size because we cropped the image.
[rows, columns] = size(binaryImage)
% Display the image.
subplot(2, 2, 2);
imshow(binaryImage, []);
hp = impixelinfo(); % Set up status line to see values when you mouse over the image.
axis('on', 'image');
caption = sprintf('Initial, Cropped Mask Image : %d rows by %d columns', rows, columns);
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
drawnow;

% Get rid of grid lines and holes caused by red spots.
se = strel('disk', 9, 0);
binaryImage = imclose(binaryImage, se);
col1 = binaryImage(:, 1);
binaryImage(:, 1) = true;
% Fill holes
binaryImage = imfill(binaryImage, 'holes');

% Display the image.
subplot(2, 2, 3);
imshow(binaryImage, []);
hp = impixelinfo(); % Set up status line to see values when you mouse over the image.
axis('on', 'image');
caption = sprintf('Final Mask Image : %d rows by %d columns', rows, columns);
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
drawnow;

hold on;
boundaries = bwboundaries(binaryImage);
% Get an array that says whether the spot's center is in the image.
centerIsInImage = false(numel(boundaries), 1);
for k = 1 : numel(boundaries)
	thisBoundary = boundaries{k};
	x = thisBoundary(:, 2);
	y = thisBoundary(:, 1);
	% Get rid of values on the left or right edge of the image
	edgeIndexes = x == 1 | x == columns;
	x(edgeIndexes) = [];
	y(edgeIndexes) = []; % Delete them
	% Get rid of values on the top or bottom edge of the image
	edgeIndexes = y == 1 | y == rows;
	x(edgeIndexes) = [];
	y(edgeIndexes) = []; % Delete them
	% Plot the boundaries
	plot(x, y, 'r-', 'LineWidth', 2);
	% Fit the circle to what's left.
	[xCenter(k), yCenter(k), radius(k), a] = circlefit(x, y);
	if (xCenter(k) >= 1) && (xCenter(k) <= columns) && (yCenter(k) > 1) && (yCenter(k) <= rows)
		plot(xCenter(k), yCenter(k), 'r.', 'MarkerSize', 40);
		str = sprintf('  #%d is at (%.1f, %.1f)', k, xCenter(k), yCenter(k));
		text(xCenter(k), yCenter(k), str, 'FontWeight', 'bold', 'Color', 'r', 'FontSize', 12);
		centerIsInImage(k) = true;
	end	
end
% Set radii that are outside the image to Nan
radius(~centerIsInImage) = nan;
% Get the area(s) of the blob(s) from the final mask image.
allAreas = pi * radius .^ 2
subplot(2, 2, 4);
bar(allAreas, 1);
grid on;
title('Spot Areas', 'FontSize', fontSize, 'Interpreter', 'None');
xlabel('Circle Number', 'FontSize', fontSize, 'Interpreter', 'None');
ylabel('Area in Pixels', 'FontSize', fontSize, 'Interpreter', 'None');

msgbox('Done!');

%==========================================================================================
function [xCenter, yCenter, radius, a] = circlefit(x, y)
% circlefit(): Fits a circle through a set of points in the x - y plane.
% USAGE :
% [xCenter, yCenter, radius, a] = circlefit(X, Y)
% The output is the center point (xCenter, yCenter) and the radius of the fitted circle.
% "a" is an optional output vector describing the coefficients in the circle's equation:
%     x ^ 2 + y ^ 2 + a(1) * x + a(2) * y + a(3) = 0
% by Bucher Izhak 25 - Oct - 1991

numPoints = numel(x);
xx = x .* x;
yy = y .* y;
xy = x .* y;
A = [sum(x),  sum(y),  numPoints;
	sum(xy), sum(yy), sum(y);
	sum(xx), sum(xy), sum(x)];
B = [-sum(xx + yy) ;
	-sum(xx .* y + yy .* y);
	-sum(xx .* x + xy .* y)];
a = A \ B;
xCenter = -.5 * a(1);
yCenter = -.5 * a(2);
radius  =  sqrt((a(1) ^ 2 + a(2) ^ 2) / 4 - a(3));
end

[Pappu Murthy]

Thanks much for the elaborate answer and the accompanying matlab code. It does solve the problem. Regarding original image that you enquired, I really do not have image probably I did not clarify properly in my original question. I only have a square matrix with a bunch of 1's and 2's (40x40) and I made an image out of it which is why you see those two colors. My aim was to find the centers of the region of "2's" which I believe is programmatically possible. But your solution miraculously captured the correct centers. I am wondering whether there is another way without using the image tool box, to get to the centers. Again, your help is much appreciated.

[Image Analyst]

It was the square matrix that I wanted, not "the image you made out of it". The matrix is already an image. Displaying it with pcolor() or something merely adds a bunch of rendering/visualization stuff that makes it harder to analyze.

You can do it without the image processing toolbox. Just write versions of imclose(), imfill(), and bwboundaries() yourself. You should be able to replicate imclose with movmin() and movmax(). The others will take more work but there is probably open source code out there to do that.

Answer 2

[darova]

Here is an example

I = imread('image.png');        % your image
B = im2bw(I);                   % binarize your image
A = regionprops(B,'centroid');  % find centroid of each region

[Pappu Murthy]

This does not work for me. please check the attached picture. My image is 40x40 cells. The centers after visual inspection are marked and they are located at: (28,2), (23,18),(39,29),(12,40) ... The x coordinates counted from left to right from top.. that is the very first point is 1 and and so on therefore there are 41 points. Similarly y is counted from top to bottom. the first 1 and the last one 41. The centers are marked and shown in the picture as well. I got them by trial and error and after some visual inspection. Looks complicated to me I agree. But I am just wondering whether there is a slick way of doing this programmatically. Hope that explanation clarifies the problem a bit better. My fault and I should have explained this better up front.