MATLAB: Gradient for fibers orientation – Local orientation map

colormapfibershistogram of orientationsimage processingImage Processing Toolboxorientation

Hello,
From an original I have which is representating muring skin fibers, I have done a the max of openings by a rotating linear segments. What I am now trying to do, is to get a colormap gradient of the fiber on my image depending on their orientation (with a colorbar map on the right). Also I am trying to get an orientation histogram.
If anyone could help 🙂 thanks !!
close all; clear all; clc
workDir='';
addpath(workDir)
%Open the image containing the collagen fibers of murine skin.
im = imread('fibres.png');
%im = imshow(im);
subplot(131); imshow(im);colormap gray; axis image;
% Compute the Local Orientation Map per pixel using the max of openings by a
% rotating linear segment.
im_supopen = zeros(size(im),'like',im);
for theta=1:180
SE = strel('line',20,theta);
im_supopen = max(im_supopen, imopen(im,SE));
end;
subplot(132); imshow(im_supopen); colormap gray; axis image;
title('sup opening by rotating segment')

Best Answer

That's a really cool technique. I don't think I've heard of that before but it gives a nicely denoised output image. I've improved it to give the histogram of fiber orientations by keeping track of what angle the max signal occurred at.
% Demo to get a histogram of fiber angles.
clc;    % Clear the command window.
clear all;
close all;
workspace;  % Make sure the workspace panel is showing.
format long g;
format compact;
fontSize = 15;
folder = pwd;
fullFileName = fullfile(folder, 'fibres.png')
% Open the image containing the collagen fibers of murine skin.
grayImage = imread(fullFileName);
% Convert to gray scale
if ndims(grayImage) > 1
	grayImage = rgb2gray(grayImage);
end
subplot(2, 3, 1); 
imshow(grayImage);
colormap gray; 
axis('on', 'image');
caption = sprintf('Original Image');
title(caption, 'fontSize', fontSize)
drawnow;
g = gcf;
g.Units = 'normalized';
g.Position = [0.1, 0.2, 0.8, 0.7]
% Compute the Local Orientation Map per pixel 
% using the max of openings by a rotating line segment
% as a structuring element for a morphological opening.
maxOfOpenings = zeros(size(grayImage), 'like', grayImage);
% Create an image where we're going to keep track of what angle was max for each pixel.
angleOfMax = zeros(size(grayImage));
lengthOfLine = 30; % 30 pixel long line that is the structuring element.
for theta = 1 : 0.25 : 180 % Degrees.
	fprintf('Theta = %.2f.\n', theta);
	% Get a structuring element of a line at this angle.
	SE = strel('line', lengthOfLine, theta);
	thisImage = imopen(grayImage, SE);
	maxOfOpenings = max(maxOfOpenings, thisImage);
	% Find out where this image is the max.
	thisIsMax = (thisImage == maxOfOpenings);
	% That is a "map" of where this angle was the max.
	% Store the angle into angleOfMax at those locations.
	angleOfMax(thisIsMax) = theta;
	if rem(theta, 10) == 1
		% Show the opening at this one particular angle.
		subplot(2, 3, 2);
		imshow(thisImage, []);
		axis('on', 'image');
		caption = sprintf('Opened with line at %d degrees', theta);
		title(caption, 'fontSize', fontSize)
		% Show the output image which is the max over all angles so far.
		subplot(2, 3, 3);
		imshow(maxOfOpenings, []); 
		axis('on', 'image');
		caption = sprintf('Output After theta = %d', theta);
		title(caption, 'fontSize', fontSize)
		drawnow;
		% Show the output image which keeps track of the max angle so far at each pixel.
		subplot(2, 3, 4);
		imshow(angleOfMax, []); 
		axis('on', 'image');
		caption = sprintf('Angle Of Max');
		title(caption, 'fontSize', fontSize)
		impixelinfo; % Let user mouse around and see the angle that was max at every pixel.
		drawnow;
	end
end
axis('on', 'image');
% Get a distribution of angles so we can see what angles 
% aligned with the structuring element the most.
subplot(2, 3, 5:6);
histogram(angleOfMax);
grid on;
title('Distribution of Angle of Max', 'fontSize', fontSize);
xlabel('Angle of Max', 'fontSize', fontSize);
ylabel('Count', 'fontSize', fontSize);
Related Question