How to measure crack width

I noticed that when I reduced the minimum acceptable area the crack width decreased although the appearance of blobs in this situation, while when i increased minimum acceptable area and the blobs removed and the width increased

for example:

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 = 10;
% Read in image and convert to gray scale.
rgbImage = imread('7004-46.jpg');
I = rgb2gray(rgbImage);
subplot(2,2,1);
imshow(I);
title('7004-46');
% Take histogram.
subplot(2,2,2);
histogram(I);
grid on;  %to make grid line
title('histogram Gray image');
% Binarize image.
threshold = 145;
xline(threshold, 'Color', 'r', 'LineWidth',2)
mask=I < threshold;
mask = imfill(mask, 'holes');
% Filter image.
se = strel('line',1, 0);
filteredImage = imclose(mask, se);
subplot(2, 2, 3);
imshow(filteredImage);
title('closed Image');
%%% Measure the areas to determin minAcceptableArea show the histogram of
%%% area
% props = regionprops(filteredImage, 'Area');
% allAreas = sort([props.Area], 'descend');
% histogram(allAreas)
%%Look at the histogram.  What area do you think is the minimum size to be a valid crack?
minAcceptableArea =20;
mask = bwareafilt(filteredImage, [minAcceptableArea, inf]);
subplot(2, 2, 4);
imshow(mask);
% 
% % Measure size of crack.
props = regionprops(mask, 'Area');
allAreas = [props.Area];
out=bwferet(mask)
% Measure the areas to know the area not to be considered
% props = regionprops(mask, 'Area');
% allAreas = sort([props.Area], 'descend')
% Get width = area / maximum length
averageWidths = sort(allAreas ./ out.MaxDiameter, 'descend');
message = sprintf('The average width = %.2f pixels.', averageWidths(1));
fprintf('%s\n', message);
caption = sprintf('Binary Image of Cracks larger than %d\nAverage width of largest = %.2f pixels', minAcceptableArea, averageWidths(1));
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
uiwait(helpdlg(message));

[yasmin ismail]

another thing when i added the whole areas and divided by the maximum diamter I got different answer for example if minimum accpetable area is 100 the output is :

3×3 table

MaxDiameter MaxAngle MaxCoordinates

________________ _________________ ______________

136.297468795279 -108.833752912659 {2×2 double}

45.6946386351835 113.198590513648 {2×2 double}

69.3541635375988 -95.7927964950321 {2×2 double}

allAreas =

523 253 230

The average width = 11.45 pixels.

if I take (523+253+230)/136.29= 7.46 not 11.45

[yasmin ismail]

@Image Analyst

Answer 1

[Image Analyst]

That doesn't look like the right algorithm if you want to find crack widths. Did you see my post where I found the mean width of a blob (image of lightning)? Basically you want to

1. Threshold/mask your image
2. Get the Euclidean Distance Transform (EDT) of your mask
3. Get the skeleton of your mask
4. Multiply your EDT image by your skeleton image to get the radius of the blob at each point along it's skeleton;
5. Double the radii to get the widths.
6. Display the histogram to see the distribution, or take the mean.

See attached demo.

[yasmin ismail]

@Image Analyst In the demo you didnt mentioned the value of the threshold .

I used your demo with another image but I didnt get the answer and I dont know where the values must be changed to let this demo be matched with other images

for example the attached image

[Image Analyst]

Try this, but first download my interactive thresholding app:

https://www.mathworks.com/matlabcentral/fileexchange/29372-thresholding-an-image?s_tid=srchtitle

% Program to compute the mean width of a blob in an image.
clearvars;
close all;
clc;
fontSize = 15;
% Read in original image, with white lightning on black background.
baseFileName = '7004-46.jpg';
fullFileName = fullfile(pwd, baseFileName);
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(rgbImage);
    % 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(:, :, 2); % Take green channel.
else
    grayImage = grayImage; % It's already gray scale.
end
% Now it's gray scale with range of 0 to 255.
subplot(2, 3, 1);
imshow(grayImage, [])
impixelinfo; % Let user mouse around and see values in the status line at the lower right.
title('Original Image', 'FontSize', fontSize);
% Binarize the image.
% mask = imbinarize(grayImage);
lowThreshold = 0;
highThreshold = 160;

% Interactively and visually set a threshold on a gray scale image.
% https://www.mathworks.com/matlabcentral/fileexchange/29372-thresholding-an-image?s_tid=srchtitle
[lowThreshold, highThreshold] = threshold(lowThreshold, highThreshold, grayImage);
mask = grayImage >= lowThreshold & grayImage <= highThreshold;
% Fill holes.
mask = imfill(mask, 'holes');
% Take largest blob only.
mask = bwareafilt(mask, 1);
subplot(2, 3, 2);
imshow(mask)
impixelinfo; % Let user mouse around and see values in the status line at the lower right.
title('Mask', 'FontSize', fontSize)
% Compute the skeleton
skelImage = bwskel(mask);
subplot(2, 3, 3);
imshow(skelImage)
impixelinfo; % Let user mouse around and see values in the status line at the lower right.
title('Thinned', 'FontSize', fontSize)
% Enlarge figure to full screen.
g = gcf;
g.WindowState = 'maximized';
% Compute the Euclidean distance image.
edtImage = bwdist(~mask);
subplot(2, 3, 4);
imshow(edtImage, [])
title('Distance Transform Image', 'FontSize', fontSize);
impixelinfo; % Let user mouse around and see values in the status line at the lower right.
% Multiply them to get an image where the only pixels in the image
% are along the skeleton and their value is the radius.
% Multiply radius image by 2 to get diameter image.
diameterImage = 2 * edtImage .* single(skelImage);
subplot(2, 3, 5);
imshow(diameterImage, [])
title('Diameter Image', 'FontSize', fontSize);
impixelinfo; % Let user mouse around and see values in the status line at the lower right.
% Get the widths.  These will be where the image is not zero.
widths = diameterImage(diameterImage > 0);
% Show histogram of widths.
subplot(2, 3, 6);
histogram(widths);
grid on;
xlabel('Width in Pixels', 'FontSize', fontSize);
ylabel('Count', 'FontSize', fontSize);
% Compute the mean width
meanWidth = mean(widths)
% Put a line on the histogram at the mean width
xline(meanWidth, 'LineWidth', 2, 'Color', 'r');
caption = sprintf('Histogram of Widths.  Mean Width = %.1f Pixels', meanWidth);
title(caption, 'FontSize', fontSize);
message = sprintf('Mean Width = %.1f Pixels', meanWidth);
msgbox(message);

[yasmin ismail]

What I understood that the demo of threshold will automatically show the min and mx threshold from threshold interface right?

what about option of otsu ,shall i use it or no need?

[Image Analyst]

Yes, I used my utility to determine the threshold. You might be able to use an automatic way like graythresh or otsuthresh, or my attached triangle threshold.

[yasmin ismail]

when i used traingle threshold I got warning dialog:

Undefined function'triangle-thrshold' for input arrgument of type 'double'

[Image Analyst]

You can also use the triangle threshold in the interactive thresholding app I attached. It's a button in the lower right of the window.

I'm pretty sure triangle threshold works for double inputs (the counts array from the histogram are always double I believe). You'd need to attach your script that calls the triangle threshold function so I can see what's wrong.

[yasmin ismail]

1-I uploaded the image that i got warning function of triangle thresholding

2-The other thing I noticed , that if i used automatic threshold in the result image I attached for the width, you will find that part of crack not considered at the left side in (diameter image) not whole cracked appear?

3- also in width histogram the average width is 3.8 although its zero frequency in histogram, so why it is not 4 which is at middle and highest frequency?

[yasmin ismail]

I have also a question, how can I get the maximum width instead of mean?

@Image Analyst

[Image Analyst]

maxWidths = max(widths)

[yasmin ismail]

I am confusing when choosing the threshold for example the image i attached:

when I run automatically the interface show min and max threshold 79 and 160, when I click ok the result shown as i attached and you will find that part of crack wont apear as I marked in the attached image and maybe this part include max width, so can you clarify this point and how to fix it to let all crack appear.

Thanks alot for your help

I attached the original image and the result of the run

@Image Analyst

[yasmin ismail]

@Image Analyst

I am confusing when choosing the threshold for example the image i attached:

when I run automatically the interface show min and max threshold 79 and 160, when I click ok the result shown as i attached and you will find that part of crack wont apear as I marked in the attached image and maybe this part include max width, so can you clarify this point and how to fix it to let all crack appear.

Thanks alot for your help

I attached the original image and the result of the run

[Image Analyst]

You need to take out the call to bwareafilt() where it extracts only the largest skeleton and use bwareaopen

% Program to compute the mean width of a blob in an image.
clearvars;
close all;
clc;
fontSize = 15;
% Read in original image, with white lightning on black background.
baseFileName = '7049-198.jpg';
fullFileName = fullfile(pwd, baseFileName);
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(rgbImage);
    % 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(:, :, 2); % Take green channel.
else
    grayImage = grayImage; % It's already gray scale.
end
% Now it's gray scale with range of 0 to 255.
subplot(2, 3, 1);
imshow(grayImage, [])
impixelinfo; % Let user mouse around and see values in the status line at the lower right.
title('Original Image', 'FontSize', fontSize);
% Binarize the image.
% mask = imbinarize(grayImage);
lowThreshold = 0;
highThreshold = 160;
% Interactively and visually set a threshold on a gray scale image.
% https://www.mathworks.com/matlabcentral/fileexchange/29372-thresholding-an-image?s_tid=srchtitle
[lowThreshold, highThreshold] = threshold(lowThreshold, highThreshold, grayImage);
mask = grayImage >= lowThreshold & grayImage <= highThreshold;
% Fill holes.
mask = imfill(mask, 'holes');
subplot(2, 3, 2);
imshow(mask)
impixelinfo; % Let user mouse around and see values in the status line at the lower right.
title('Mask', 'FontSize', fontSize)

% Take largest blob only.
% mask = bwareafilt(mask, 1);
% Compute the skeleton
skelImage = bwskel(mask);
% Find the areas of all the skeletons.
props = regionprops(skelImage, 'Area');
allAreas = sort([props.Area])
% Extract only skeletons longer than 60 pixels.
skelImage = bwareaopen(skelImage, 60);
subplot(2, 3, 3);
imshow(skelImage)
impixelinfo; % Let user mouse around and see values in the status line at the lower right.
title('Thinned', 'FontSize', fontSize)
% Enlarge figure to full screen.
g = gcf;
g.WindowState = 'maximized';

% Compute the Euclidean distance image.
edtImage = bwdist(~mask);
subplot(2, 3, 4);
imshow(edtImage, [])
title('Distance Transform Image', 'FontSize', fontSize);
impixelinfo; % Let user mouse around and see values in the status line at the lower right.

% Multiply them to get an image where the only pixels in the image
% are along the skeleton and their value is the radius.
% Multiply radius image by 2 to get diameter image.
diameterImage = 2 * edtImage .* single(skelImage);
subplot(2, 3, 5);
imshow(diameterImage, [])
title('Diameter Image', 'FontSize', fontSize);
impixelinfo; % Let user mouse around and see values in the status line at the lower right.

% Get the widths.  These will be where the image is not zero.
widths = diameterImage(diameterImage > 0);
% Show histogram of widths.
subplot(2, 3, 6);
histogram(widths);
grid on;
xlabel('Width in Pixels', 'FontSize', fontSize);
ylabel('Count', 'FontSize', fontSize);
% Compute the mean width
meanWidth = mean(widths)
% Put a line on the histogram at the mean width
xline(meanWidth, 'LineWidth', 2, 'Color', 'r');
caption = sprintf('Histogram of Widths.  Mean Width = %.1f Pixels', meanWidth);
title(caption, 'FontSize', fontSize);
message = sprintf('Mean Width = %.1f Pixels', meanWidth);
msgbox(message);

[yasmin ismail]

But how did you choose 60 ?

% Extract only skeletons longer than 60 pixels.

skelImage = bwareaopen(skelImage, 60);

@Image Analyst

[Image Analyst]

I just picked it. There is lots of noise in the image. There are cracks of all different lengths from just a few pixels to thousands of pixels. What is YOUR definition of how long a legitimate crack should be? Do you want cracks just 4 pixels long? A hundred long? It's very complicated. There is no one single overall, grand "right" answer for the minimum size of valid cracks, so you just have to define one. Do you care about analyzing a crack if it's only 2 pixels long? I think not. So what is your minimum length that you care about?

[yasmin ismail]

Thanks alot @Image Analyst

[yasmin ismail]

@Image Analyst

I have question about determining threshold:

In the inteface of threshold there are three options to determine it ;

1)otsu

2)triangle

3)manually by detriming the lower and upper threshold

My question is what you call option No.3, is it greythresh or what you call it exactly?

[Image Analyst]

I call it "manual" because the thresholds are what you set by moving the sliders (scrollbars). The threshold is not selected by an automatic method such as Otsu or triangle or graythresh, it's selected by you manually.

[yasmin ismail]

@Image Analyst

what do you mean by

impixelinfo; % Let user mouse around and see values in the status line at the lower right.

which values you mean?

[Image Analyst]

The (x,y) location and RGB color values of the pixel underneath your cursor.

[yasmin ismail]

yes I got it.

But why you asked for it twice, for diameter image and distance transform image, wht is the important to see thier values in both?

[Image Analyst]

I never know where or when I might need it so I usually just call impixelinfo as a matter of habit after I call imshow(). Better to have it if you need it than to not have it. It doesn't hurt to have it, so why not? If you're sure you'll never need it, you can delete that line of code.

[yasmin ismail]

@Image Analyst thnx lot for your clarification

[yasmin ismail]

@Image Analyst How can I measure the following:

1- length of crack and

2-its angle at different points

[Image Analyst]

1. Sum up the skeleton image (will count all cracks in the image if there is more than one). This will give you a pixel count.
2. Not sure what you mean but atand or atan2d will give you angles, given a delta y and delta x.

[yasmin ismail]

@Image Analyst another question:

To find width,if I used (bwperim) instead of (bwskel) will be more accurate?

[yasmin ismail]

@Image Analyst

I am trying to find max width of this image which includes many cracks but I have difficulities with determing thresholding with same above technique

thus, can you help me to fix it please?

[Image Analyst]

Try adapthisteq to flatten the image so that you can use a global threshold. Or use imbinarize with the locally adaptive option. Once you have the binary image, you can take the max of the Euclidean Distance Transform to find the radius of the max crack or blob.

[yasmin ismail]

@Image Analyst

I did both option

for first option using adapthisteq I used following code and the result is attached image called (adapthistic)

I = imread('real2.jpg');

II=rgb2gray(I);

J = adapthisteq(II);

level = graythresh(J);

imshowpair(II,J,'montage')

for second option the result also attached called (imbinarize with local adaptive) , the following code:

T=imread('real2.jpg');

T2=rgb2gray(T);

BW2=imbinarize(T2,'adaptive');

imshowpair(T2,BW2,'montage')

are that what do you mean?

[Image Analyst]

Yeah, basically, though you may need to change some of the parameters, like the threshold level or, in the case of adapthisteq, the tile size.

[yasmin ismail]

@Image Analyst

I tried to find mean and max width fo attached image i tried many options but i couldnt fix it, can you help me please

% Program to compute the mean width of a blob in an image.

clearvars;

close all;

clc;

fontSize = 15;

% Read in original image, with white lightning on black background.

baseFileName = '44.jpg';

fullFileName = fullfile(pwd, baseFileName);

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(rgbImage);

% 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(:, :, 2); % Take green channel.

else

grayImage = grayImage; % It's already gray scale.

end

% Now it's gray scale with range of 0 to 255.

subplot(2, 3, 1);

imshow(grayImage, [])

impixelinfo; % Let user mouse around and see values in the status line at the lower right.

title('Original Image', 'FontSize', fontSize);

% Binarize the image.

% mask = imbinarize(grayImage);

lowThreshold = 0;

highThreshold = 160;

% Interactively and visually set a threshold on a gray scale image.

% https://www.mathworks.com/matlabcentral/fileexchange/29372-thresholding-an-image?s_tid=srchtitle

[lowThreshold, highThreshold] = threshold(lowThreshold, highThreshold, grayImage);

mask = grayImage >= lowThreshold & grayImage <= highThreshold;

% Fill holes.

mask = imfill(mask, 'holes');

subplot(2, 3, 2);

imshow(mask)

impixelinfo; % Let user mouse around and see values in the status line at the lower right.

title('Mask', 'FontSize', fontSize)

% Take largest blob only.

% mask = bwareafilt(mask, 1);

% Compute the skeleton

skelImage = bwskel(mask);

% Find the areas of all the skeletons.

props = regionprops(skelImage, 'Area');

allAreas = sort([props.Area])

% Extract only skeletons longer than 60 pixels.

skelImage = bwareaopen(skelImage, 60);

subplot(2, 3, 3);

imshow(skelImage)

impixelinfo; % Let user mouse around and see values in the status line at the lower right.

title('Thinned', 'FontSize', fontSize)

% Enlarge figure to full screen.

g = gcf;

g.WindowState = 'maximized';

% Compute the Euclidean distance image.

edtImage = bwdist(~mask);

subplot(2, 3, 4);

imshow(edtImage, [])

title('Distance Transform Image', 'FontSize', fontSize);

impixelinfo; % Let user mouse around and see values in the status line at the lower right.

% Multiply them to get an image where the only pixels in the image

% are along the skeleton and their value is the radius.

% Multiply radius image by 2 to get diameter image.

diameterImage = 2 * edtImage .* single(skelImage);

subplot(2, 3, 5);

imshow(diameterImage, [])

title('Diameter Image', 'FontSize', fontSize);

impixelinfo; % Let user mouse around and see values in the status line at the lower right.

% Get the widths. These will be where the image is not zero.

widths = diameterImage(diameterImage > 0);

% Show histogram of widths.

subplot(2, 3, 6);

histogram(widths);

grid on;

xlabel('Width in Pixels', 'FontSize', fontSize);

ylabel('Count', 'FontSize', fontSize);

% Compute the mean width

meanWidth = mean(widths)

% Put a line on the histogram at the mean width

xline(meanWidth, 'LineWidth', 2, 'Color', 'r');

caption = sprintf('Histogram of Widths. Mean Width = %.1f Pixels', meanWidth);

title(caption, 'FontSize', fontSize);

message = sprintf('Mean Width = %.1f Pixels', meanWidth);

msgbox(message);

[Image Analyst]

That is a pretty variable image. Could you try a Hessian, Frangi, or COSFIRE-B filter on it? There are ones in the File Exchange.

[yasmin ismail]

@Image Analyst so, what you mean is I have to use filteration function first then saved image after filtrarion then I can applied above codes to get width,right?

[Image Analyst]

Yes, try various spatial filters so that you can get as much contrast between the crack and the background, and have the background be fairly flat (same gray level range). Then you can apply a fixed global threshold to discriminate crack from background.

[yasmin ismail]

@Image Analyst

I used filters as you advised me before and I could determine the threshold as shown in attached image but still the diameter image and distance transform the lines of crack is too small not clear enough, although I reduced (skeletons longer than 1 pixels)

but still not clear , can you help to fix it please

% Program to compute the mean width of a blob in an image.

clearvars;

close all;

clc;

fontSize = 15;

% Read in original image, with white lightning on black background.

baseFileName = 'A2.jpg';

fullFileName = fullfile(pwd, baseFileName);

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(rgbImage);

% 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(:, :, 2); % Take green channel.

else

grayImage = grayImage; % It's already gray scale.

end

% Now it's gray scale with range of 0 to 255.

subplot(2, 3, 1);

imshow(grayImage, [])

impixelinfo; % Let user mouse around and see values in the status line at the lower right.

title('Original Image', 'FontSize', fontSize);

% Binarize the image.

% mask = imbinarize(grayImage);

lowThreshold = 0;

highThreshold = 40;

% Interactively and visually set a threshold on a gray scale image.

% https://www.mathworks.com/matlabcentral/fileexchange/29372-thresholding-an-image?s_tid=srchtitle

[lowThreshold, highThreshold] = threshold(lowThreshold, highThreshold, grayImage);

mask = grayImage >= lowThreshold & grayImage <= highThreshold;

% Fill holes.

mask = imfill(mask, 'holes');

subplot(2, 3, 2);

imshow(mask)

impixelinfo; % Let user mouse around and see values in the status line at the lower right.

title('Mask', 'FontSize', fontSize)

% Take largest blob only.

% mask = bwareafilt(mask, 1);

% Compute the skeleton

skelImage = bwskel(mask);

% Find the areas of all the skeletons.

props = regionprops(skelImage, 'Area');

allAreas = sort([props.Area])

% Extract only skeletons longer than 60 pixels.

skelImage = bwareaopen(skelImage, 1);

subplot(2, 3, 3);

imshow(skelImage)

impixelinfo; % Let user mouse around and see values in the status line at the lower right.

title('Thinned', 'FontSize', fontSize)

% Enlarge figure to full screen.

g = gcf;

g.WindowState = 'maximized';

% Compute the Euclidean distance image.

edtImage = bwdist(~mask);

subplot(2, 3, 4);

imshow(edtImage, [])

title('Distance Transform Image', 'FontSize', fontSize);

impixelinfo; % Let user mouse around and see values in the status line at the lower right.

% Multiply them to get an image where the only pixels in the image

% are along the skeleton and their value is the radius.

% Multiply radius image by 2 to get diameter image.

diameterImage = 2 * edtImage .* single(skelImage);

subplot(2, 3, 5);

imshow(diameterImage, [])

title('Diameter Image', 'FontSize', fontSize);

impixelinfo; % Let user mouse around and see values in the status line at the lower right.

% Get the widths. These will be where the image is not zero.

widths = diameterImage(diameterImage > 0);

% Show histogram of widths.

subplot(2, 3, 6);

histogram(widths);

grid on;

xlabel('Width in Pixels', 'FontSize', fontSize);

ylabel('Count', 'FontSize', fontSize);

% Compute the mean width

meanWidth = mean(widths)

% Put a line on the histogram at the mean width

xline(meanWidth, 'LineWidth', 2, 'Color', 'r');

caption = sprintf('Histogram of Widths. Mean Width = %.1f Pixels', meanWidth);

title(caption, 'FontSize', fontSize);

message = sprintf('Mean Width = %.1f Pixels', meanWidth);

msgbox(message);

% Compute the max width

maxWidth = max(widths)

% Put a line on the histogram at the max width

xline(maxWidth, 'LineWidth', 2, 'Color', 'r');

caption = sprintf('Histogram of Widths. Max Width = %.1f Pixels', maxWidth);

title(caption, 'FontSize', fontSize);

message = sprintf('Max Width = %.1f Pixels', maxWidth);

msgbox(message);

[yasmin ismail]

@Image Analyst

can you help me to know how to let image 1 has same gray level as 2, I tried using imadjust but not working?

I attached the 2 images 7010-17 and CC1-3

[Image Analyst]

Let's say that we want the mean of all images to be 150 and the cracks are small enough that they don't affect the mean of the entire image much. So, before you apply the threshold, just scale each image by dividing by it's mean and multiplying by 150.

image1 = 150 * double(image1) / mean2(image1);
image2 = 150 * double(image2) / mean2(image2);

[yasmin ismail]

@Image Analyst first of all thanks

I got it but I couldnit display the image with new mean, how can I save it and enter it to demo of thresholding to find threshold

I=imread('CC1-3.jpg');

meanIntensity = mean(I(:))

image1 = 186.1056 * double(meanIntensity) / mean2(meanIntensity)

% Program to compute the mean width of a blob in an image.

clearvars;

close all;

clc;

fontSize = 15;

% Read in original image, with white lightning on black background.

baseFileName = '7009-30.jpg';

fullFileName = fullfile(pwd, baseFileName);

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(rgbImage);

% 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(:, :, 2); % Take green channel.

else

grayImage = grayImage; % It's already gray scale.

end

% Now it's gray scale with range of 0 to 255.

subplot(2, 3, 1);

imshow(grayImage, [])

impixelinfo; % Let user mouse around and see values in the status line at the lower right.

title('Original Image', 'FontSize', fontSize);

% Binarize the image.

% mask = imbinarize(grayImage);

lowThreshold = 0;

highThreshold = 160;

% Interactively and visually set a threshold on a gray scale image.

% https://www.mathworks.com/matlabcentral/fileexchange/29372-thresholding-an-image?s_tid=srchtitle

[lowThreshold, highThreshold] = threshold(lowThreshold, highThreshold, grayImage);

mask = grayImage >= lowThreshold & grayImage <= highThreshold;

[Image Analyst]

I don't think you actually thought about what you were doing. First of all meanIntensity is a scalar, so when you do mean2(meanIntensity) you get the same scalar. Let's say meanIntensity was 186.1056. Well, mean(186.1056) is simply just 186.1056, right? Secondly, you then call clearvars after that so it clears everything you just computed!

You need to do something like I told you:

% Program to compute the mean width of a blob in an image.
clearvars;
close all;
clc;
fontSize = 15;

%--------------------------------------------------------------------------------------------
% Compute mean intensity that we will set all images to.
refImage = imread('CC1-3.jpg');
meanRefIntensity = mean(refImage, 'all')

%--------------------------------------------------------------------------------------------
% Read in test image, with white lightning on black background.
baseFileName = '7009-30.jpg';
fullFileName = fullfile(pwd, baseFileName);
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(rgbImage);
    % 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(:, :, 2); % Take green channel.
end
% Now it's gray scale with range of 0 to 255.

%--------------------------------------------------------------------------------------------
% Scale this intensity so that its mean intensity is the same as the
% reference image
grayImage = meanRefIntensity * double(grayimage) / mean2(grayImage);
subplot(2, 3, 1);
imshow(grayImage, [])
impixelinfo; % Let user mouse around and see values in the status line at the lower right.
title('Original Image', 'FontSize', fontSize);

%--------------------------------------------------------------------------------------------
% Binarize the image by thresholding.
lowThreshold = 0;
highThreshold = 160;
% Interactively and visually set a threshold on a gray scale image.
% https://www.mathworks.com/matlabcentral/fileexchange/29372-thresholding-an-image?s_tid=srchtitle
[lowThreshold, highThreshold] = threshold(lowThreshold, highThreshold, grayImage);
mask = grayImage >= lowThreshold & grayImage <= highThreshold;

[yasmin ismail]

@Image Analyst

unfortunately it is not what I want, let me exaplain

If you see the gray level in image 2 , the concave down of grey level is between 140 to 200

while image (Gray CC1-3) is between 40 to 160 which is affecting when i tried to detect and find maximum crack width, the result is unaccurate for image ( CC1-3) you wil find the crack is not detected well and discontinuity as shown in attachment (Result CC1-3) the image of (Diameter Image) not detected the crack well, while the attached No.4 the crack is clear in (Diameter Image). Thus, I though to reduced the range of gray image and let it between 40 to 160,because when I teseted the image with same range I got clear and accurate detectection otherwise not. So, how to fix it?

The following code I used to detect and find max crack width for image CC1-3

% Program to compute the mean width of a blob in an image.

clearvars;

close all;

clc;

fontSize = 15;

% Read in original image, with white lightning on black background.

baseFileName = 'CC1-3.jpg';

fullFileName = fullfile(pwd, baseFileName);

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(rgbImage);

% 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(:, :, 3); % Take Blue channel.

else

grayImage = grayImage; % It's already gray scale.

end

% Now it's gray scale with range of 0 to 255.

subplot(2, 3, 1);

imshow(grayImage, [])

impixelinfo; % Let user mouse around and see values in the status line at the lower right.

title('Original Image', 'FontSize', fontSize);

% Binarize the image.

% mask = imbinarize(grayImage);

lowThreshold = 0;

highThreshold = 40;

% Interactively and visually set a threshold on a gray scale image.

% https://www.mathworks.com/matlabcentral/fileexchange/29372-thresholding-an-image?s_tid=srchtitle

[lowThreshold, highThreshold] = threshold(lowThreshold, highThreshold, grayImage);

mask = grayImage >= lowThreshold & grayImage <= highThreshold;

% Fill holes.

mask = imfill(mask, 'holes');

subplot(2, 3, 2);

imshow(mask)

impixelinfo; % Let user mouse around and see values in the status line at the lower right.

title('Mask', 'FontSize', fontSize)

% Take largest blob only.

% mask = bwareafilt(mask, 1);

% Compute the skeleton

skelImage = bwskel(mask);

% Find the areas of all the skeletons.

props = regionprops(skelImage, 'Area');

allAreas = sort([props.Area])

% Extract only skeletons longer than 60 pixels.

skelImage = bwareaopen(skelImage, 10);

subplot(2, 3, 3);

imshow(skelImage)

impixelinfo; % Let user mouse around and see values in the status line at the lower right.

title('Thinned', 'FontSize', fontSize)

% Enlarge figure to full screen.

g = gcf;

g.WindowState = 'maximized';

% Compute the Euclidean distance image.

edtImage = bwdist(~mask);

subplot(2, 3, 4);

imshow(edtImage, [])

title('Distance Transform Image', 'FontSize', fontSize);

impixelinfo; % Let user mouse around and see values in the status line at the lower right.

% Multiply them to get an image where the only pixels in the image

% are along the skeleton and their value is the radius.

% Multiply radius image by 2 to get diameter image.

diameterImage = 2 * edtImage .* single(skelImage);

subplot(2, 3, 5);

imshow(diameterImage, [])

title('Diameter Image', 'FontSize', fontSize);

impixelinfo; % Let user mouse around and see values in the status line at the lower right.

% Get the widths. These will be where the image is not zero.

widths = diameterImage(diameterImage > 0);

% Show histogram of widths.

subplot(2, 3, 6);

histogram(widths);

grid on;

xlabel('Width in Pixels', 'FontSize', fontSize);

ylabel('Count', 'FontSize', fontSize);

% Compute the mean width

meanWidth = mean(widths)

% Put a line on the histogram at the mean width

xline(meanWidth, 'LineWidth', 2, 'Color', 'r');

caption = sprintf('Histogram of Widths. Mean Width = %.1f Pixels', meanWidth);

title(caption, 'FontSize', fontSize);

message = sprintf('Mean Width = %.1f Pixels', meanWidth);

msgbox(message);

% Compute the max width

maxWidth = max(widths)

% Put a line on the histogram at the max width

xline(maxWidth, 'LineWidth', 2, 'Color', 'r');

caption = sprintf('Histogram of Widths. Max Width = %.1f Pixels', maxWidth);

title(caption, 'FontSize', fontSize);

message = sprintf('Max Width = %.1f Pixels', maxWidth);

msgbox(message);

[Image Analyst]

OK, it looks like you just refuse to standardize the gray levels to use a fixed threshold. So in that case you're going to be stuck with trying to use an automatic threshold. From the looks of your histogram you might try a triangle threshold on the dark/left side. There is a button for that in my utility, or you can call it directly with the attached function. See if that selects a good threshold. But there may be no global threshold that is perfect. If your crack is not solid then changing the threshold to make it more solid may introduce other noise blobs, which you may be able to get rid of with morphological filtering methods based on area or shape. Or if you have no background noise then your crack may be a series of separated smaller cracks. You have to decide if you REALLY need just one giant crack or if you're still able to make some sort of decision on the quality of the material based on analyzing the smaller cracks. Realize that just because you think you want one single big long crack, this may not actually be necessary to analyze the quality or suitability of the material for use in some building/construction application. For example, maybe the area fraction of all the cracks, regardless of how many or what size they are, is a good indicator of material quality, and maybe that's even better than if you had to go through special image processing gymnastics to get just one big long crack.

[yasmin ismail]

@Image Analyst If I would like to apply anytype of filters at the shown place in the code to be input for next stage to find the threshold how can I make it (how can I type it)

% Program to compute the mean width of a blob in an image.

clearvars;

close all;

clc;

fontSize = 15;

% Read in original image, with white lightning on black background.

baseFileName = 'bb.jpg';

fullFileName = fullfile(pwd, baseFileName);

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(rgbImage);

% 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(:, :, 3); % Take Blue channel.

else

grayImage = grayImage; % It's already gray scale.

end

% Now it's gray scale with range of 0 to 255.

subplot(2, 3, 1);

imshow(grayImage, [])

impixelinfo; % Let user mouse around and see values in the status line at the lower right.

title('Original Image', 'FontSize', fontSize);

% apply filter for example gaussian or any thing else

fi = (imgray);

figure, imshow(fi)

subplot(2, 3, 2);

imshow(fi)

impixelinfo; % Let user mouse around and see values in the status line at the lower right.

title('filtered Image', 'FontSize', fontSize);

% Binarize the image.

% mask = imbinarize(grayImage);

lowThreshold = 0;

highThreshold = 230;

% Interactively and visually set a threshold on a gray scale image.

% https://www.mathworks.com/matlabcentral/fileexchange/29372-thresholding-an-image?s_tid=srchtitle

[lowThreshold, highThreshold] = threshold(lowThreshold, highThreshold, fi);

mask = fi >= lowThreshold & fi <= highThreshold

% Fill holes.

mask = imfill(mask, 'holes');

subplot(2, 3, 2);

imshow(mask)

impixelinfo; % Let user mouse around and see values in the status line at the lower right.

title('Mask', 'FontSize', fontSize)

[Image Analyst]

I have no idea what filters you want to use, but just go ahead and use them.

"how can I make it" <== it totally depends on what kind of filter you want. If it's a filter that is built-in to MATLAB, then just read the documentation for it.

[yasmin ismail]

@Image Analyst assume I want to use gaussian filter and the output be the input for next step for binarize to find threshold by the interface , how to do it?

[Image Analyst]

Isn't "how to do it" given in the help for imgaussfilt?

help imgaussfilt

 IMGAUSSFILT 2-D Gaussian filtering of images
    B = imgaussfilt(A) filters image A with a 2-D Gaussian smoothing kernel
    with standard deviation of 0.5. A can have any number of dimensions.
 
    B = imgaussfilt(A,SIGMA) filters image A with a 2-D Gaussian smoothing
    kernel with standard deviation specified by SIGMA. SIGMA can be a
    scalar or a 2-element vector with positive values. If sigma is a
    scalar, a square Gaussian kernel is used.
 
    B = imgaussfilt(___,Name,Value,...) filters image A with a 2-D
    Gaussian smoothing kernel with Name-Value pairs used to control aspects
    of the filtering. 
 
    Parameters include:
 
    'FilterSize'    -   Scalar or 2-element vector, of positive, odd 
                        integers that specifies the size of the Gaussian
                        filter. If a scalar Q is specified, then a square
                        Gaussian filter of size [Q Q] is used.
 
                        Default value is 2*ceil(2*SIGMA)+1.
 
    'Padding'       -   String or character vector or numeric scalar that
                        specifies padding to be used on image before
                        filtering. 
 
                        If a scalar (X) is specified, input image values
                        outside the bounds of the image are implicitly
                        assumed to have the value X. 
 
                        If a string is specified, it can be 'replicate',
                        'circular' or 'symmetric'. These options are
                        analogous to the padding options provided by
                        imfilter. 
 
                        'replicate'
                        Input image values outside the bounds of the image
                        are assumed equal to the nearest image border
                        value.
 
                        'circular'
                        Input image values outside the bounds of the image
                        are computed by implicitly assuming the input image
                        is periodic.
 
                        'symmetric'
                        Input image values outside the bounds of the image
                        are computed by mirror-reflecting the array across
                        the array border. 
 
                        Default value is 'replicate'.
 
    'FilterDomain'  -   String or character vector that specifies domain in
                        which filtering is performed. It can be 'spatial',
                        'frequency' or 'auto'. For 'spatial', convolution
                        is performed in the spatial domain, for 'frequency',
                        convolution is performed in the frequency domain
                        and for 'auto', convolution may be performed in
                        spatial or frequency domain based on internal
                        heuristics.
 
                        Default value is 'auto'.
 
 
    Class Support
    -------------
    The input image A must be a real, non-sparse matrix of any dimension of
    the following classes: uint8, int8, uint16, int16, uint32, int32,
    single or double.
 
 
    Notes
    -----
    1. If the image A contains Infs or NaNs, the behavior of imgaussfilt
       for frequency domain filtering is undefined. This can happen when
       the 'FilterDomain' parameter is set to 'frequency' or when it is set
       to 'auto' and frequency domain filtering is used internally. To
       restrict the propagation of Infs and NaNs in the output in a manner
       similar to imfilter, consider setting the 'FilterDomain' parameter
       to 'spatial'.
 
    2. When 'FilterDomain' parameter is set to 'auto', an internal
       heuristic is used to determine whether spatial or frequency domain
       filtering is faster. This heuristic is machine dependent and may
       vary for different configurations. For optimal performance, consider
       comparing 'spatial' and 'frequency' to determine the best filtering
       domain for your image and kernel size.
 
    3. When the 'Padding' parameter is not specified, 'replicate' padding 
       is used as the default, which is different from the default used in
       imfilter.
 
 
    Example
    ---------
    % Smooth an image with Gaussian filters of increasing standard deviations
        I = imread('cameraman.tif');
 
        subplot(2,2,1), imshow(I), title('Original Image');
 
        Iblur = imgaussfilt(I, 2);
        subplot(2,2,2), imshow(Iblur)
        title('Gaussian filtered image, \sigma = 2')
 
        Iblur = imgaussfilt(I, 4);
        subplot(2,2,3), imshow(Iblur)
        title('Gaussian filtered image, \sigma = 4')
 
        Iblur = imgaussfilt(I, 6);
        subplot(2,2,4), imshow(Iblur)
        title('Gaussian filtered image, \sigma = 6')
 
    See also imgaussfilt3, imfilter, fspecial.

    Documentation for imgaussfilt
       doc imgaussfilt

    Other uses of imgaussfilt

       gpuArray/imgaussfilt

Open the documentation in your own version of MATLAB (with the doc imgaussfilt command) for more complete documentation.

Of course, this function will blur the image so I hope that is what you want.

[yasmin ismail]

@Image AnalystI didnt mean how to apply gaussian , I know how to apply it

but I want it to be part of the above code , I want the output of gaussian automatically enter the next step as shown above , I put fi as symbol for the output of gaussian to enter it to get threshold

I tried it, but it doesnt work I have to save the output of gaussian filter separetly and then put it name by function imread(' imge.jpg)

[Image Analyst]

@yasmin ismail If you want the variable to be used in your other function, just pass the (poorly named) variable into the other function, for example

fi = imgaussfilt(...............)
outputImage = SomeOtherFunction(fi); % Pass fi to the other function.

If you want to have your other function read it in from a disk file, then you need to call imwrite to write fi to a disk file, and then tell the other function what the name of the file is. By the way, NEVER use jpg format for images you're doing image analysis on because it's lossy compression and you won't get the same image back that you saved after reading in a JPG image.

[yasmin ismail]

@Image Analyst what the format i have to use (png or tif) or what? and shall I have to use the one of those formats from first step or when I used filter i have to save image as tif or png then continue to binarization and other steps with one of those format?

[Image Analyst]

For saving to disk, use PNG format. Of course once it's read back into an array in MATLAB, the format it had on disk doesn't matter, but you want to make sure you use a lossless format such as PNG or TIF to make sure you get back exactly the same data you saved out..

[yasmin ismail]

@Image Analyst so when I take photo by mobile , I have to save it in labtop as PNG, then start use it in matlab, and also the result output from process by matlb save it as png

is it correct?

The thing that I noticed when I send captured photo from my mobile to labtop automatically formated as jpg ,so its entered matlab as it jpg, so its already jpg. The only thing I can do if my model has several process for example start with filteration, the output could be saved as png then input it to next step for binarization.

Thus , what do you think?

[Image Analyst]

@yasmin ismail youi don't have to use PNG, though it's best. If your image capture system takes JPG and you can't improve your image capture system then you're stuck with JPG. Hopefully the degradation is not too bad. But it's not worth converting the files into a PNG file format - that won't improve the quality one tiny bit. Just leave them as JPG format images.

And I don't see any need to save intermediate images to disk unless you want to for some reason. Just keep processing, using the images in memory, until you have the final result.

Also for spatial calibration you need to make sure your camera is the same distance from the scene. If not you'd have to have some calibration object in the scene and you'd have to detect that and calibrate to that embedded object.

Answer 2

[yasmin ismail]

@Image Analyst can you explain to me why the mean width = 1.7mm as shown in the result image No.1

althought there is no blue column on 1.7 mm and the mean width is 0.9 mm althougth also no blue column on 0.9 mm as shown in the histogram .

The code used:

% Program to compute the mean width of a blob in an image.
clearvars;
close all;
clc;
fontSize = 15;
% Read in original image, with white lightning on black background.
baseFileName = '7002-9.jpg';
fullFileName = fullfile(pwd, baseFileName);
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(rgbImage);
    % 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(:, :, 3); % Take Blue channel.
else
    grayImage = grayImage; % It's already gray scale.
end
% Now it's gray scale with range of 0 to 255.
subplot(2, 3, 1);
imshow(grayImage, [])
impixelinfo; % Let user mouse around and see values in the status line at the lower right.
title('Original Image', 'FontSize', fontSize);
% Binarize the image.
% mask = imbinarize(grayImage);
lowThreshold = 0;
highThreshold = 150;
% Interactively and visually set a threshold on a gray scale image.
% https://www.mathworks.com/matlabcentral/fileexchange/29372-thresholding-an-image?s_tid=srchtitle
[lowThreshold, highThreshold] = threshold(lowThreshold, highThreshold, grayImage);
mask = grayImage >= lowThreshold & grayImage <= highThreshold;
% Bridges unconnected pixels (morphological operation).
mask = bwmorph(mask, 'bridge');
% Fill holes.
mask = imfill(mask, 'holes');
subplot(2, 3, 2);
imshow(mask)
impixelinfo; % Let user mouse around and see values in the status line at the lower right.
title('Mask', 'FontSize', fontSize)
% Take largest blob only.
% mask = bwareafilt(mask, 1);
% Compute the skeleton
skelImage = bwskel(mask);
% Find the areas of all the skeletons.
props = regionprops(skelImage, 'Area');
allAreas = sort([props.Area])
% Extract only skeletons longer than 60 pixels.
skelImage = bwareaopen(skelImage, 10);
subplot(2, 3, 3);
imshow(skelImage)
impixelinfo; % Let user mouse around and see values in the status line at the lower right.
title('Thinned', 'FontSize', fontSize)
% Enlarge figure to full screen.
g = gcf;
g.WindowState = 'maximized';
% Compute the Euclidean distance image.
edtImage = bwdist(~mask);
subplot(2, 3, 4);
imshow(edtImage, [])
title('Distance Transform Image', 'FontSize', fontSize);
impixelinfo; % Let user mouse around and see values in the status line at the lower right.
% Multiply them to get an image where the only pixels in the image
% are along the skeleton and their value is the radius.
% Multiply radius image by 2 to get diameter image.
diameterImage = 2 * edtImage .* single(skelImage);
subplot(2, 3, 5);
imshow(diameterImage, [])
title('Diameter Image', 'FontSize', fontSize);
impixelinfo; % Let user mouse around and see values in the status line at the lower right.
% Get the widths.  These will be where the image is not zero.
widths = diameterImage(diameterImage > 0);
% Show histogram of widths.
subplot(2, 3, 6);
histogram(widths);
grid on;
xlabel('Width in mm', 'FontSize', fontSize);
ylabel('Count', 'FontSize', fontSize);
% Compute the mean width in pixels
meanWidthPixls = mean(widths)
%%Calaculate Mean crack width in mm
calibration_width=25.4;   
calibration_pixels=96;
Mean_Width=meanWidthPixls;
crack_Mean_width=(Mean_Width *calibration_width)/calibration_pixels
% Put a line on the histogram at the mean width
xline(crack_Mean_width, 'LineWidth', 2, 'Color', 'r');
caption = sprintf('Histogram of Widths.  Mean Width = %.1f mm', crack_Mean_width);
title(caption, 'FontSize', fontSize);
message = sprintf('Mean Width = %.1f mm', crack_Mean_width);
msgbox(message);
% Compute the max width in pixels
maxWidthPixels = max(widths)
%%Calaculate Max crack width in mm
calibration_width=25.4;   
calibration_pixels=96;
Max_Width=maxWidthPixels;
crack_width=(Max_Width *calibration_width)/calibration_pixels

% Put a line on the histogram at the max width
xline(crack_width, 'LineWidth', 2, 'Color', 'r');
caption = sprintf('Histogram of Widths.  Max Width = %.1f mm', crack_width);
title(caption, 'FontSize', fontSize);
message = sprintf(' Max crack width = %.1f mm', crack_width);
msgbox(message);

[Image Analyst]

@yasmin ismail you plotted the histogram of widths in pixels. But you then calibrated to millimeters with this:

%%Calaculate Mean crack width in mm
calibration_width=25.4;   
calibration_pixels=96;
Mean_Width=meanWidthPixls;
crack_Mean_width=(Mean_Width *calibration_width)/calibration_pixels

then you used xline to put a line up at the mm width but the x axis was in pixel widths, so the red vertical lines won't be in the right place because the vertical red lines are in calibrated units while the histogram was in pixel units.

[yasmin ismail]

@Image Analyst so, how can I fix it? I want histogram in millimeters not pixels

can you help me to modify it, please?

[yasmin ismail]

@Image Analyst so, how can I fix it? I want histogram in millimeters not pixels

can you help me to modify it, please?

[Image Analyst]

Try this:

% Program to compute the mean width of a crack in an image.
% Initialization steps.
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 original image, with white lightning on black background.
baseFileName = '7002-9.jpg';
fullFileName = fullfile(pwd, baseFileName);
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(rgbImage);
    % 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(:, :, 3); % Take Blue channel.
else
    grayImage = grayImage; % It's already gray scale.
end
% Now it's gray scale with range of 0 to 255.
subplot(2, 3, 1);
imshow(grayImage, [])
impixelinfo; % Let user mouse around and see values in the status line at the lower right.
title('Original Image', 'FontSize', fontSize);
% Binarize the image.
% mask = imbinarize(grayImage);
lowThreshold = 0;
highThreshold = 150;
% Interactively and visually set a threshold on a gray scale image.
% https://www.mathworks.com/matlabcentral/fileexchange/29372-thresholding-an-image?s_tid=srchtitle
[lowThreshold, highThreshold] = threshold(lowThreshold, highThreshold, grayImage);
mask = grayImage >= lowThreshold & grayImage <= highThreshold;
% Bridges unconnected pixels (morphological operation).
mask = bwmorph(mask, 'bridge');
% Fill holes.
mask = imfill(mask, 'holes');
subplot(2, 3, 2);
imshow(mask)
impixelinfo; % Let user mouse around and see values in the status line at the lower right.
title('Mask', 'FontSize', fontSize)
% Take largest blob only.
% mask = bwareafilt(mask, 1);
% Compute the skeleton
skelImage = bwskel(mask);
% Find the areas of all the skeletons.
props = regionprops(skelImage, 'Area');
allAreas = sort([props.Area])
% Extract only skeletons longer than 60 pixels.
skelImage = bwareaopen(skelImage, 10);
subplot(2, 3, 3);
imshow(skelImage)
impixelinfo; % Let user mouse around and see values in the status line at the lower right.
title('Thinned', 'FontSize', fontSize)
% Enlarge figure to full screen.
g = gcf;
g.WindowState = 'maximized';
% Compute the Euclidean distance image.
edtImage = bwdist(~mask);
subplot(2, 3, 4);
imshow(edtImage, [])
title('Distance Transform Image', 'FontSize', fontSize);
impixelinfo; % Let user mouse around and see values in the status line at the lower right.
% Multiply them to get an image where the only pixels in the image
% are along the skeleton and their value is the radius.
% Multiply radius image by 2 to get diameter image.
diameterImage = 2 * edtImage .* single(skelImage);
subplot(2, 3, 5);
imshow(diameterImage, [])
title('Diameter Image', 'FontSize', fontSize);
impixelinfo; % Let user mouse around and see values in the status line at the lower right.

% Get the widths.  These will be where the image is not zero.
widthsInPixels = diameterImage(diameterImage > 0);
% Calculate crack widths in mm
calibration_width=25.4;   
calibration_pixels=96;
widthsInMms =(widthsInPixels * calibration_width) / calibration_pixels;

% Show histogram of widths.
subplot(2, 3, 6);
histogram(widthsInMms);
grid on;
xlabel('Width in mm', 'FontSize', fontSize);
ylabel('Count', 'FontSize', fontSize);
% Compute the mean width in pixels
meanWidthInPixels = mean(widthsInPixels)
% Compute the mean width in mm
meanWidthInMms = mean(widthsInMms)

% Put a line on the histogram at the mean width
xline(meanWidthInMms, 'LineWidth', 2, 'Color', 'r');
caption = sprintf('Histogram of Widths.  Mean Width = %.1f mm', meanWidthInMms);
title(caption, 'FontSize', fontSize);

% Compute the max width in pixels
maxWidthPixels = max(widthsInPixels)
% Compute the max width in mm
maxWidthMm = max(widthsInMms)

% Put a line on the histogram at the max width
xline(maxWidthMm, 'LineWidth', 2, 'Color', 'r');
caption = sprintf('Histogram of Widths.\nMean Width = %.2f mm,  Max Width = %.2f mm', meanWidthInMms, maxWidthMm);
title(caption, 'FontSize', fontSize);
message = sprintf('Mean Crack Width = %.2f mm\nMax Crack Width = %.2f mm', meanWidthInMms, maxWidthMm);
uiwait(msgbox(message));

[yasmin ismail]

@Image Analyst yes its working thanks alot, but can you explain why mean on 0.91 mm although no blue column too?

[yasmin ismail]

@Image Analyst also can I get width size for each part in the crack image in the same above code not only the max and mean?

[Image Analyst]

The widths of each crack are given in the variable widthsInMms.

Why do you think there should be a crack with exactly the mean crack width? For a simple example, let's say your data set is two numbers 10 and 20. The mean value is 15 but the mean value is not one of the values in the dataset.

data = [10, 20];
meanOfData = mean(data)

meanOfData = 15

histogram(data);
xline(meanOfData, 'LineWidth', 3, 'Color', 'r');
grid on;

[yasmin ismail]

@Image Analyst for the mean width I got it thanks

but my second question I mean if i divided crack into 10 parts how to get width in each part like in attached image 2 , but I want to apply it in image 7002-9

[Image Analyst]

Sorry for the delay - you've probably solved it by now. I've just been busy with other stuff and couldn't spend time on Answers recently. So with your skeleton image, you need to label it with bwlabel, then use ismember to loop over all blobs and extract one at a time.

[labeledImage, numberOfCracks] = bwlabel(skelImage);
allCrackWidths10 = zeros(numberOfCracks, 10); % 10 widths for every crack.
for k = 1 : numberOfCracks
    oneSingleCrack = ismember(labeledImage, k);
end

If you want the histogram of that, that is easy. If you insist on traversing down the skeleton, that is not so easy since there is no way to give you all those coordinates in order going along the curve. So you can use bwmorph to identify endpoints.

endPointImage = bwmorph(oneSingleCrack, 'endpoints');

Then use regionprops to get all the pixel coordinates.

props = regionprops(oneSingleCrack, 'PixelList');

Then go though all the pixel coordinates to see which one is closest to the last one. This is the hardest part. Unfortunately I don't have time to do it for you so I'm leaving that part to you.

Then use those coordinates and the edtImage to get the width at every coordinate in order. This loop should be inside the k loop over all cracks.

% x and y are the column and row coordinates of the single crack.
widthsOfSingleCrack = zeros(1, numel(x));
for k2 = 1 : numel(x);
    row = y(k2);
    col = x(k2);
    widthsOfSingleCrack(k2) = edtImage(row, col);
end

Then, after the loop, in another loop, use linspace and interp1 to get 10 distances between the starting point and the ending point. This loop is inside the k loop and after the k2 loop.

for k3 = 1 : numberOfCracks
    % Get 10 indexes from first element to last element.
    xCrack = linspace(1, numel(widthsOfSingleCrack), 10);
    % Get the 10 widths at those xCrack indexes.
    allCrackWidths10(k3, :) = interp1(1:numel(widthsOfSingleCrack), widthsOfSingleCrack, xCrack);
end

Give it a try. The code is untested so do the part that I didn't do and debug it. I'll check back in a day or two.

[yasmin ismail]

@Image Analyst first of all i would thanks and appreciated your helps

I am civil engineer and I donnt have a good experince in this topic, thus I tried to test the code above as followig but I didnt get the answer as shwon below:

% Program to compute the mean width of a crack in an image.

% Initialization steps.

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 original image, with white lightning on black background.

baseFileName = '7006-134.jpg';

fullFileName = fullfile(pwd, baseFileName);

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(rgbImage);

% 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(:, :, 3); % Take Blue channel.

else

grayImage = grayImage; % It's already gray scale.

end

% Now it's gray scale with range of 0 to 255.

subplot(2, 3, 1);

imshow(grayImage, [])

impixelinfo; % Let user mouse around and see values in the status line at the lower right.

title('Original Image', 'FontSize', fontSize);

% Binarize the image.

% mask = imbinarize(grayImage);

lowThreshold = 0;

highThreshold = 150;

% Interactively and visually set a threshold on a gray scale image.

% https://www.mathworks.com/matlabcentral/fileexchange/29372-thresholding-an-image?s_tid=srchtitle

[lowThreshold, highThreshold] = threshold(lowThreshold, highThreshold, grayImage);

mask = grayImage >= lowThreshold & grayImage <= highThreshold;

% Bridges unconnected pixels (morphological operation).

mask = bwmorph(mask, 'bridge');

% Fill holes.

mask = imfill(mask, 'holes');

subplot(2, 3, 2);

imshow(mask)

impixelinfo; % Let user mouse around and see values in the status line at the lower right.

title('Mask', 'FontSize', fontSize)

% Take largest blob only.

% mask = bwareafilt(mask, 1);

% Compute the skeleton

skelImage = bwskel(mask);

% Find the areas of all the skeletons.

props = regionprops(skelImage, 'Area');

allAreas = sort([props.Area])

% Extract only skeletons longer than 60 pixels.

skelImage = bwareaopen(skelImage, 10);

subplot(2, 3, 3);

imshow(skelImage)

impixelinfo; % Let user mouse around and see values in the status line at the lower right.

title('Thinned', 'FontSize', fontSize)

[labeledImage, numberOfCracks] = bwlabel(skelImage);

allCrackWidths10 = zeros(numberOfCracks, 10); % 10 widths for every crack.

for k = 1 : numberOfCracks

oneSingleCrack = ismember(labeledImage, k);

end

endPointImage = bwmorph(oneSingleCrack, 'endpoints');

props = regionprops(oneSingleCrack, 'PixelList');

% x and y are the column and row coordinates of the single crack.

widthsOfSingleCrack = zeros(1, numel(x));

for k2 = 1 : numel(x);

row = y(k2);

col = x(k2);

widthsOfSingleCrack(k2) = edtImage(row, col);

end

for k3 = 1 : numberOfCracks

% Get 10 indexes from first element to last element.

xCrack = linspace(1, numel(widthsOfSingleCrack), 10);

% Get the 10 widths at those xCrack indexes.

allCrackWidths10(k3, :) = interp1(1:numel(widthsOfSingleCrack), widthsOfSingleCrack, xCrack);

end

So, could you help me to fix it, please?

[Image Analyst]

I told you the k2 and k3 loops should be INSIDE the k loop, not outside it.

[yasmin ismail]

@Image Analyst I dont have well experience in matlab can you show me how to do it pleases?

[yasmin ismail]

@Image Analyst

[Image Analyst]

How about just a list of widths for each crack? I don't have the time to sort the widths in 10 different locations along the midline. If you need it you'll have to figure out that part for yourself. But I don't see why you'd need that anyway since you have the width at every center pixel - just click on it or nover the mouse to see what the value is. Or just write all the widths out to a text file or something.