optimization of scobel image processing - better approach needed

Hello,

i have pictures with bands. The bands can be tilted. There can be 1 upto 5 bands on a image. I need the width of atleast one band!

My script just work for 50% of the images and almost always i need some manual adjustments for different types of the pattern. Im looking for a generalized solution without manual adjustments.

Thank You.

My Code:

%% User variable
image_name = '1';
image_extension = '.png';

% two lines will be detected at top and bottom of the imagee. That is why
% you will see +/- 2 in the script below at number of places. if band is not
% detected fully this number can be increased.
number_of_lines = 10;  

% possible value (0.1 - 0.5) This should be increased if band is present in lighter region.
intensity_factor = 0.1; 

% This is filter to segment the line. If line is bolder then 2 in (2,20) should be
% increased to appropriate value. for these images, 2 is fine.
% if there are more feature in the band boundary then 20 in (2,20) can be
% adjusted.
line_filter = [ones(2,20);zeros(2,20);-ones(2,20)] ;


%% Line detection
im1 = imread([image_name image_extension]) ;
im2 = (im2double(rgb2gray(im1)));

imProcess = adapthisteq(im2);
fimg = imfilter((imProcess),line_filter ); 

bw = edge(fimg,'sobel');
[H,T,R] = hough(bw);
P  = houghpeaks(H,number_of_lines+2,'threshold',ceil(intensity_factor*max(H(:))));
lines = houghlines(fimg,T,R,P,'FillGap',10,'MinLength',100);

f_1 = figure(1);
clf
imshow(im1), hold on
mid_point = zeros(length(lines)-2,2) ;

for k = 1:length(lines)-2
    xy = [lines(k+2).point1; lines(k+2).point2];
    mid_point(k,:) = mean(xy) ; 
    plot(xy(:,1),xy(:,2),'LineWidth',2,'Color','green');
end

%% Bandwidth calculation
dist_points = triu(pdist2(mid_point,mid_point)) ;
dist_points(dist_points==0) = nan;
bandwidth = min(dist_points,[],'all') ;
[row,col] = find(dist_points==bandwidth);
plot(mid_point([row(1) col(1)],1),mid_point([row(1) col(1)],2),'LineWidth',2,'Color','red') ;
title(['Band Width = ', num2str(bandwidth), ' pixels'])
drawnow
%% File save
save([image_name,'.mat'],'bandwidth')
print(f_1,'-dpng',[image_name,'_line'])

Answer 1

[Image Analyst]

Here is what I have so far:

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

%===============================================================================
% Read in gray scale image.
folder = pwd;
baseFileName = '4.png';
% Get the full filename, with path prepended.
fullFileName = fullfile(folder, baseFileName);
if ~exist(fullFileName, 'file')
	% Didn't find it there.  Check the search path for it.
	fullFileName = baseFileName; % No path this time.
	if ~exist(fullFileName, 'file')
		% Still didn't find it.  Alert user.
		errorMessage = sprintf('Error: %s does not exist.', fullFileName);
		uiwait(warndlg(errorMessage));
		return;
	end
end
[rgbImage, map] = imread(fullFileName);
% Get the dimensions of the image.  numberOfColorBands should be = 1.
[rows, columns, numberOfColorBands] = size(rgbImage);
% If it's RGB instead of grayscale, convert it to gray scale.
if numberOfColorBands > 1
	grayImage = rgbImage(:, :, 1); % Take red channel.
else
	grayImage = rgbImage;
end
% Display the original image.
subplot(2, 3, 1);
imshow(grayImage, []);
axis on;
impixelinfo; % Let user mouse around and see gray level.
caption = sprintf('Original Gray Scale Image : %s', baseFileName);
title(caption, 'FontSize', fontSize, 'Interpreter', 'none');
impixelinfo;
% Enlarge figure to full screen.
g = gcf;
g.WindowState = 'maximized';
g.NumberTitle = 'off';
g.Name = 'Histogram Analysis';
drawnow;

mkdir('results');
startingImage = grayImage;
kx = 4;
ky = 24;
% for kx = 2 : 2 : 24
% 	for ky = 2 : 2 : 24
		grayImage = adapthisteq(startingImage, 'numTiles', [kx, ky]);
		% Display the image.
		subplot(2, 3, 2);
		imshow(grayImage, []);
		axis on;
		impixelinfo; % Let user mouse around and see gray level.
		caption = sprintf('/results/%d %d.png', kx, ky);
		fn = fullfile(pwd, caption);
		title(caption, 'FontSize', fontSize, 'Interpreter', 'none');
		impixelinfo;
		drawnow;
% 		imwrite(grayImage, fn);
% 	end
% end

%===============================================================================
% Display the histogram
subplot(2, 3, 3);
[pixelCounts, grayLevels] = imhist(grayImage);
bar(grayLevels, pixelCounts);
grid on;
title('Image Histogram', 'FontSize', fontSize);
xlabel('Gray Level', 'FontSize', fontSize);
ylabel('Pixel Count', 'FontSize', fontSize);

% verticalSDProfile = std(double(grayImage), [], 2);
% verticalSDProfile = sgolayfilt(verticalSDProfile, 2, 31);
verticalProfile = mean(double(grayImage), 2);
verticalProfile = sgolayfilt(verticalProfile, 3, 29);

% Display the vertical profile.
subplot(2, 3, 4);
plot(verticalProfile, '-', 'Linewidth', 2);
grid on;
impixelinfo; % Let user mouse around and see gray level.
title('Mean Profile', 'FontSize', fontSize, 'Interpreter', 'none');
impixelinfo;


% Display the vertical profile.
% subplot(2, 3, 6);
% plot(verticalSDProfile, '-', 'Linewidth', 2);
% grid on;
% impixelinfo; % Let user mouse around and see gray level.
% title('SD Profile', 'FontSize', fontSize, 'Interpreter', 'none');
% impixelinfo;

As you can see the bands are quite well found from the vertical profile. So what you need to do is to take my attached radon transform demo and rotate your images so that the bands are perfectly horizontal. Then use the code above to find the starting and stopping lines for each band.

Answer 2

[Image Analyst]

What happens if you try stdfilt(), or look at the vertical profile by summing the image horizontally?

[Image Analyst]

What are these things? What do you really want to know?

1. A line fitted along the top and bottom of each "band"?
2. The angle of that band?
3. The average width of the band (why would it vary)?
4. Do you want the outline to be not a line but "hug" the regions somewhat closely (how closely)?
5. The width of each band as a function of column in the image?
6. The width of each band perpendicular to the main axis of each band?
7. Have you tried activecontour() or Chan-Vese segmentation?
8. Have you tried SegNet (segmentation deep learning)?