Isolate and measure diameter of mouse pupil in MATLAB

I have several different images of eyes of mice where I would like to isolate the pupil and measure the diameter. How would I create a mask for that? How does one include a function in a script? I've attached an example image. Would really appreciate some guidance. Thank you!!

[Scott Satinover]

Have you tried using the imfindcircles function?

[Madiha Irshad]

I can try that but I'm not sure what range to use. Is there a way for me to figure that out?

[Scott Satinover]

You need to know how many pixels the eye is in diameter. Assuming you have the Image Processing Toolbox, you can use the Image Viewer for an initial guess:

https://www.mathworks.com/help/images/measure-distance-between-pixels-in-image-viewer-app.html

I am also assuming that the images aren't all going to have the same sized pupil. Without getting too complex, just playing with the range for this function with some manual guess and check should narrow this down enough. Because the pupil has a fair amount of contrast, this should be pretty easy to detect.

[Madiha Irshad]

Thank you so much for the suggestion. I tried that and the outputs for [centers,radii] = imfindcircles(A,radiusRange) are empty after having put in a suitable range. Do you perhaps have any other ideas?

Answer 1

[Image Analyst]

Madiha, try this:

% Demo by Image Analyst.  July 2021
clc;    % Clear the command window.
close all;  % Close all figures (except those of imtool.)
clear;  % Erase all existing variables. Or clearvars if you want.
workspace;  % Make sure the workspace panel is showing.
format long g;
format compact;
fontSize = 20;

%--------------------------------------------------------------------------------------------------------
%    READ IN IMAGE
folder = pwd;
baseFileName = 'Mouse Eye.png';
grayImage = imread(baseFileName);
% Get the dimensions of the image.
% numberOfColorChannels should be = 1 for a gray scale image, and 3 for an RGB color image.
[rows, columns, numberOfColorChannels] = size(grayImage)
if numberOfColorChannels > 1
	% It's not really gray scale like we expected - it's color.
	% Use weighted sum of ALL channels to create a gray scale image.
	grayImage = min(grayImage, [], 3);
end

%--------------------------------------------------------------------------------------------------------
% Display the image.
subplot(2, 2, 1);
imshow(grayImage, []);
axis('on', 'image');
title('Original Image', 'FontSize', fontSize, 'Interpreter', 'None');
impixelinfo;
hFig = gcf;
hFig.WindowState = 'maximized'; % May not work in earlier versions of MATLAB.

%--------------------------------------------------------------------------------------------------------
% Get the histogram so we can see where to threshold it.
subplot(2, 2, 2);
imhist(grayImage);
grid on;
xlabel('Gray Level', 'FontSize', fontSize, 'Interpreter', 'None');
ylabel('Pixel Count', 'FontSize', fontSize, 'Interpreter', 'None');
title('Histogram of Original Image', 'FontSize', fontSize, 'Interpreter', 'None');

%--------------------------------------------------------------------------------------------------------
% Threshold to get the pupil.
% Use interactive thresholding app at https://www.mathworks.com/matlabcentral/fileexchange/29372-thresholding-an-image?s_tid=srchtitle
lowThreshold = 0.0;	% Set an initial guess.
highThreshold = 48;	% Set an initial guess.  (Depends on neighborhood size used for stdfilt().)
if ~isempty(which('threshold'))
	[lowThreshold, highThreshold] = threshold(lowThreshold, highThreshold, grayImage)
else
	message = sprintf('Please download interactive threshold app at https://www.mathworks.com/matlabcentral/fileexchange/29372-thresholding-an-image?s_tid=srchtitle');
	uiwait(helpdlg(message));
end
% Put lines of threshold on the histogram.
xline(lowThreshold, 'Color', 'r', 'LineWidth', 2);
xline(highThreshold, 'Color', 'r', 'LineWidth', 2);

pupilImage = (grayImage >= lowThreshold) & (grayImage <= highThreshold);
% Display the image.
subplot(2, 2, 3);
imshow(pupilImage, []);
axis('on', 'image');
title('Binary (Thresholded) Image', 'FontSize', fontSize, 'Interpreter', 'None');
caption = sprintf('Binary Image Thresholded at %d.', highThreshold);
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
impixelinfo;
drawnow;

%--------------------------------------------------------------------------------------------------------
% Further clean up.
% Remove blobs touching edge of the image
% grayImage = imclearborder(grayImage);
% Take largest blob only.
pupilImage = bwareafilt(pupilImage, 1);
% Fill holes due to specular reflections.
pupilImage = imfill(pupilImage, 'holes');
subplot(2, 2, 4);
imshow(pupilImage, []);
axis('on', 'image');
title('Final Pupil Mask', 'FontSize', fontSize, 'Interpreter', 'None');
impixelinfo;

%--------------------------------------------------------------------------------------------------------
% Measure areas
props = regionprops(pupilImage, 'Area', 'Centroid');
pupilArea = props.Area; % A simple pixel ocount.
% Measure area an alternative way that takes into account the shape of the boundary of the pupil.
weightedArea = bwarea(pupilImage)
caption = sprintf('Final Pupil Mask.  The pupil area is %d pixels (or %.1f pixels)', pupilArea, weightedArea);
subplot(2, 2, 4);
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');

%--------------------------------------------------------------------------------------------------------
% Plot the boundaries and center over the original image.
subplot(2, 2, 1);
boundaries = bwboundaries(pupilImage);
boundaries = boundaries{1};
x = boundaries(:, 2);
y = boundaries(:, 1);
hold on;
plot(x, y, 'r-', 'LineWidth', 2); % Plot outline (boundary)
xCenter = props.Centroid(1);
yCenter = props.Centroid(2);
plot(xCenter, yCenter, 'r+', 'LineWidth', 2, 'MarkerSize', 40); % Plot crosshairs at centroid
title('Original Image with Pupil Outlined in Red', 'FontSize', fontSize, 'Interpreter', 'None');

message = sprintf('The pupil area is %d pixels (or %.1f pixels)', pupilArea, weightedArea);
fprintf('%s\n', message);
message = sprintf('Done!\nThe pupil area is %d pixels (or %.1f pixels).', pupilArea, weightedArea);
uiwait(helpdlg(message));

[Madiha Irshad]

Thank you so so much!! I'll try this and let you know how it goes. Really appreciate the help!!

[Madiha Irshad]

If i have a video of the pupil instead, How would i modify this script to get the same result? I have seperate videos of the pupil at baseline, while contracting and while expanding.

[Image Analyst]

See my attached demo for dealing with a video.

[Madiha Irshad]

Thank you!!