Detection of Particle Orientation

How can I detect the orientation of the bright region with respect to the center of the particle using matlab?? I have done a rough analysis using paint and simple trignometry and it turns out to be around 45 degree (considering the center of the bright region and the center of the particle).

Answer 1

[Image Analyst]

Use regionprops():

% Code to detect angle of a blob.
% By Image Analyst
clc;    % Clear the command window.
fprintf('Beginning to run %s.m ...\n', mfilename);
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 = 15;

%--------------------------------------------------------------------------------------------------------
%    READ IN IMAGE
folder = pwd;
baseFileName = 'ZoomIn_100f.jpg';
% Get the full filename, with path prepended.
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
grayImage = imread(fullFileName);

% 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 = rgb2gray(grayImage);
	% ALTERNATE METHOD: Convert it to gray scale by taking only the green channel,
	% which in a typical snapshot will be the least noisy channel.
	grayImage = grayImage(:, :, 1); % Take red channel.
end
% Crop it to the right lines.
% grayImage = grayImage(1750:2250, 800:end);
subplot(2, 2, 1);
imshow(grayImage, []);
impixelinfo;
title('Red Channel Image', 'FontSize', fontSize, 'Interpreter', 'None');
hFig = gcf;
hFig.WindowState = 'maximized'; % May not work in earlier versions of MATLAB.
drawnow;

subplot(2, 2, 2);
imhist(grayImage);
grid on;
title('Histogram of Original Image', 'FontSize', fontSize, 'Interpreter', 'None');

%--------------------------------------------------------------------------------------------------------
%    SEGMENTATION OF IMAGE
% Get a binary image by interactively thresholding using the function at
% https://www.mathworks.com/matlabcentral/fileexchange/29372-thresholding-an-image
% 91 seems good for this image.
lowThreshold = 91;
highThreshold = 255;
% [lowThreshold, highThreshold] = threshold(lowThreshold, highThreshold, grayImage)
mask = grayImage > lowThreshold & grayImage < highThreshold;
mask = imclearborder(mask);
mask = imfill(mask, 'holes');
% Take largest blobs.
mask = bwareafilt(mask, 1);
subplot(2, 2, 3);
imshow(mask, []);
impixelinfo;
caption = sprintf('Mask of White SemiCircles');
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');

% Measure centroid and orientation.
props = regionprops(mask, 'Centroid', 'Orientation')
xCenter = props.Centroid(1);
yCenter = props.Centroid(2);
hold on;
plot(xCenter, yCenter, 'r+', 'MarkerSize', 40, 'LineWidth', 2);

% Display the masked image.
subplot(2, 2, 4);
overlayImage = imoverlay(grayImage, mask, 'r');
imshow(overlayImage, []);
impixelinfo;
caption = sprintf('Masked Image.  Angle = %.2f degrees.', props.Orientation);
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
hFig = gcf;
hFig.WindowState = 'maximized'; % May not work in earlier versions of MATLAB.
drawnow;

% Make a line at that angle at the centroid.
slope = tan(props.Orientation)
x = 1 : columns;
y = slope * (x - xCenter) + yCenter;
hold on;
plot(x, y, 'r-', 'LineWidth', 2);
fprintf('Done running %s.m ...\n', mfilename);

[Qaisar Latif]

Hi Image Analyst. Thanks a lot for your help and the code. If I may ask, was the blob angle calculated considering the center of the particle itself?

I actually used the method shown below to calculate the angle of blob with respect to the center of the particle, i.e x = ?, y = ? and then getting theta.

[Image Analyst]

No, I didn't see that you asked to do that. So what I did was to find the bright blob inside the dark circle and fit that shape to an ellipse (via regionprops) and then get the angle of the major axis of the fitted ellipse to the horizontal axis. If you want to do what you said, you'd first have to find the centroid of the black blob, then the centroid of the white blob inside it, then compute the angle. Pretty easy once you see what I've done already, right? I think you can manage that, right?

[Qaisar Latif]

Hi Image Analyst. Thanks for your prompt response. To get the centroid for dark circle do I have to start from the image segmentation (that's what I understood using your demo_script.m) as you did to get the white blob.

[Qaisar Latif]

Hi Image Analyst, I did manage to find the centroids and orientation of the blob and particle. Is it possible to overlay the point on top of each other and draw a line between two center points and then getting the angle?

I did the following changes to your code;

%--------------------------------------------------------------------------------------------------------
%    SEGMENTATION OF IMAGE
% Get a binary image by interactively thresholding using the function at
% https://www.mathworks.com/matlabcentral/fileexchange/29372-thresholding-an-image
% 91 seems good for this image.

lowThreshold = 91;
highThreshold = 140;

% For Blob segmentation

% [lowThreshold, highThreshold] = threshold(lowThreshold, highThreshold, grayImage)
mask = grayImage > lowThreshold & grayImage < highThreshold;

mask = imclearborder(mask);
mask = imfill(mask, 'holes');

% Take largest blobs.
mask = bwareafilt(mask, 1);
subplot(2, 2, 1);
imshow(mask, []);
impixelinfo;
caption = sprintf('Mask of White SemiCircles');
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');

% Measure centroid and orientation.

props = regionprops(mask, 'Centroid', 'Orientation')
xCenter = props.Centroid(1);
yCenter = props.Centroid(2);
hold on;
plot(xCenter, yCenter, 'r+', 'MarkerSize', 40, 'LineWidth', 2);

% For Particle segmentation

mask = grayImage < lowThreshold & grayImage < highThreshold;

mask = imclearborder(mask);
mask = imfill(mask, 'holes');

mask = bwareafilt(mask, 1);
subplot(2, 2, 2);
imshow(mask, []);
impixelinfo;
caption = sprintf('Mask of particle');
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');

% Measure centroid and orientation.
props1 = regionprops(mask, 'Centroid', 'Orientation')
xCenter1 = props1.Centroid(1);
yCenter1 = props1.Centroid(2);
hold on;
plot(xCenter1, yCenter1, 'r+', 'MarkerSize', 40, 'LineWidth', 2);


% Display the masked image.

subplot(2, 2, 3);

overlayImage = imoverlay(grayImage, mask, 'r');
imshow(overlayImage, []);
impixelinfo;
caption = sprintf('Masked Image.  Angle = %.2f degrees.', props.Orientation);
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
hFig = gcf;
hFig.WindowState = 'maximized'; % May not work in earlier versions of MATLAB.
drawnow;

% Make a line at that angle at the centroid.
slope = tan(props.Orientation)
x = 1 : columns;
y = slope * (x - xCenter) + yCenter;
hold on;
plot(x, y, 'r-', 'LineWidth', 2);
fprintf('Done running %s.m ...\n', mfilename);


% Display the masked image.

subplot(2, 2, 4);
overlayImage = imoverlay(grayImage, mask, 'r');
imshow(overlayImage, []);
impixelinfo;
caption = sprintf('Masked Image.  Angle = %.2f degrees.', props1.Orientation);
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
hFig = gcf;
hFig.WindowState = 'maximized'; % May not work in earlier versions of MATLAB.
drawnow;

% Make a line at that angle at the centroid.
slope = tan(props1.Orientation)
x = 1 : columns;
y = slope * (x - xCenter1) + yCenter1;
hold on;
plot(x, y, 'r-', 'LineWidth', 2);
fprintf('Done running %s.m ...\n', mfilename);

and this the output:

I don't understand why in case of blob the whole particle is filled up with red, though the tangent line does pass through the blob centroid.

[Image Analyst]

Qaisar: try this:

% Code to detect angle of a blob.
% By Image Analyst
clc;    % Clear the command window.
fprintf('Beginning to run %s.m ...\n', mfilename);
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 = 15;

%--------------------------------------------------------------------------------------------------------
%    READ IN IMAGE
folder = pwd;
baseFileName = 'ZoomIn_100f.jpg';
% Get the full filename, with path prepended.
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
grayImage = imread(fullFileName);

% 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 = rgb2gray(grayImage);
	% ALTERNATE METHOD: Convert it to gray scale by taking only the green channel,
	% which in a typical snapshot will be the least noisy channel.
	grayImage = grayImage(:, :, 1); % Take red channel.
end
% Crop it to the right lines.
% grayImage = grayImage(1750:2250, 800:end);
subplot(2, 3, 1);
imshow(grayImage, []);
impixelinfo;
title('Red Channel Image', 'FontSize', fontSize, 'Interpreter', 'None');
hFig = gcf;
hFig.WindowState = 'maximized'; % May not work in earlier versions of MATLAB.
drawnow;

hHist = subplot(2, 3, 2:3);
imhist(grayImage);
grid on;
title('Histogram of Original Image', 'FontSize', fontSize, 'Interpreter', 'None');

%--------------------------------------------------------------------------------------------------------
%    SEGMENTATION OF IMAGE : Large Dark Circle
% Get a binary image by interactively thresholding using the function at
% https://www.mathworks.com/matlabcentral/fileexchange/29372-thresholding-an-image
% 91 seems good for this image.
lowThreshold = 0;
highThreshold = 108;
% [lowThreshold, highThreshold] = threshold(lowThreshold, highThreshold, grayImage)
% Mark threshold on the histogram.
xline(hHist, highThreshold, 'Color', 'b', 'LineWidth', 3);
mask = grayImage > lowThreshold & grayImage < highThreshold;
mask = imclearborder(mask);
mask = imfill(mask, 'holes');
% Take largest blobs.
mask = bwareafilt(mask, 1);
subplot(2, 3, 4);
imshow(mask, []);
impixelinfo;
caption = sprintf('Mask of Dark Circle');
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');

% Measure centroid and orientation.
props1 = regionprops(mask, 'Centroid', 'Orientation')
xCenter1 = props1.Centroid(1);
yCenter1 = props1.Centroid(2);
hold on;
plot(xCenter1, yCenter1, 'b+', 'MarkerSize', 40, 'LineWidth', 2);

%--------------------------------------------------------------------------------------------------------
%    SEGMENTATION OF IMAGE : Small oval bright blob within larger black circle.
% Get a binary image by interactively thresholding using the function at
% https://www.mathworks.com/matlabcentral/fileexchange/29372-thresholding-an-image
% 91 seems good for this image.
lowThreshold = 91;
highThreshold = 255;
% [lowThreshold, highThreshold] = threshold(lowThreshold, highThreshold, grayImage)
% Mark threshold on the histogram.
xline(hHist, lowThreshold, 'Color', 'r', 'LineWidth', 3);
mask = grayImage > lowThreshold & grayImage < highThreshold;
mask = imclearborder(mask);
mask = imfill(mask, 'holes');
% Take largest blobs.
mask = bwareafilt(mask, 1);
subplot(2, 3, 5);
imshow(mask, []);
impixelinfo;
caption = sprintf('Mask of White SemiCircles');
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');

% Measure centroid and orientation.
props2 = regionprops(mask, 'Centroid', 'Orientation')
xCenter2 = props2.Centroid(1);
yCenter2 = props2.Centroid(2);
hold on;
plot(xCenter2, yCenter2, 'r+', 'MarkerSize', 40, 'LineWidth', 2);

% Display the masked image.
subplot(2, 3, 6);
% overlayImage = imoverlay(grayImage, mask, 'r');
imshow(grayImage, []);
axis('on', 'image');
impixelinfo;

% Make a line at that angle at the centroid.
deltay = yCenter2 - yCenter1
deltax = xCenter2 - xCenter1
angle = -atan2d(deltay, deltax)
slope = (yCenter2 - yCenter1) / (xCenter2 - xCenter1)
x = 1 : columns;
y = slope * (x - xCenter2) + yCenter2;
hold on;
plot(xCenter1, yCenter1, 'b.', 'MarkerSize', 30);
plot(xCenter2, yCenter2, 'r.', 'MarkerSize', 30);
plot(x, y, 'r-', 'LineWidth', 2);

caption = sprintf('Gray Image.  Angle = %.2f degrees.', angle);
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
hFig = gcf;
hFig.WindowState = 'maximized'; % May not work in earlier versions of MATLAB.
drawnow;

fprintf('Done running %s.m ...\n', mfilename);

[Qaisar Latif]

Hi Image Analyst. Thanks a lot. This is exactly what I was looking for..

[Image Analyst]

You're welcome. Could you then "Accept this answer"? Thanks in advance.

[Venkata Manikantah Sai Ganesh Tanuku]

Hi Image analyst,

   I want to track the same orientation orientation in a video, where one record theta with respect to the deviation from intital position (X0,Y0).

First instance define a unit vector at (X0,Y0) and as particle move to the position (X1 Y1) record the deviation of orientation angle.

[Image Analyst]

@Tanuku please start a new discussion thread after you read this:

TUTORIAL: How to ask a question (on Answers) and get a fast answer

Attach your video.