Concrete Crack Width Measurement

Hi,

I have many images of cracked concrete samples. I have to calculate the average widths of cracks. My supervisor suggested me to use ImageJ for this issue. However I have never used imageJ before. I will be really appreciate if anyone help me and offer me a method and so on. Thanks in advance.

Also I attached the image of cracked concrete sample.

[Walter Roberson]

https://www.mathworks.com/matlabcentral/answers/?term=tag%3A%22crack%22

Answer 1

[Image Analyst]

See my Image Segmentation Tutorial. My File Exchange

What I might do is to use bwskel() to get the length of the crack. Then use regionprops() to get the area. Then divide the two to get the mean width. Also search the forum for cracks - it's been seen here several times before.

[Ilgin Sandalci]

Thanks for your help.

Answer 2

[Image Analyst]

Ilgin, your last post does not exactly scream "I know I can do this!" So I'm giving you a start with this little snippet:

% Code to detect dark cracks in a concrete image and compute the average width of the cracks.
% By Image Analyst
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 = 22;

%--------------------------------------------------------------------------------------------------------
%    READ IN IMAGE
folder = pwd;
baseFileName = 'concrete.jpeg';
% 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');

% Do a morphological opening filter
se = strel('disk', 9, 0);
filteredImage = imopen(grayImage, se);
subplot(2, 2, 3);
imshow(filteredImage);
impixelinfo;
title('Opened Image', 'FontSize', fontSize, 'Interpreter', 'None');

subplot(2, 2, 2);
imhist(filteredImage);
grid on;
title('Histogram of Morphologically Opened 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
% 60 seems good for this image.
lowThreshold = 0;
highThreshold = 60;
% [lowThreshold, highThreshold] = threshold(lowThreshold, highThreshold, filteredImage)
mask = filteredImage > lowThreshold & filteredImage < highThreshold;
mask = imfill(mask, 'holes');
minAcceptableArea = 10000;
mask = bwareafilt(mask, [minAcceptableArea, inf]);
subplot(2, 2, 4);
imshow(mask, []);
impixelinfo;
caption = sprintf('Binary Image of Cracks larger than %d', minAcceptableArea);
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
hFig = gcf;
hFig.WindowState = 'maximized'; % May not work in earlier versions of MATLAB.
drawnow;

% Measure the areas
props = regionprops(mask, 'Area');
allAreas = [props.Area]
out = bwferet(mask)
averageLengths = sort(allAreas ./ out.MaxDiameter, 'descend')
message = sprintf('The average width = %.2f pixels.', averageLengths(1));
fprintf('%s\n', message);
caption = sprintf('Binary Image of Cracks larger than %d\nAverage width of largest = %.2f pixels', minAcceptableArea, averageLengths(1));
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
uiwait(helpdlg(message));

[Ilgin Sandalci]

Thank you so much. Actually I have not enough information about Matlab and image analysis. It is a really new topic for me. Your tutorial and snippet will be good and effective beginning for me. Again thanks for your kindness.

[yasmin ismail]

what do you mean by:

minAcceptableArea=10000

and if change what is different?

[Image Analyst]

@yasmin ismail when you threshold you will get a lot of small noise blobs, like tiny black grains or pebbles. These might be so small that you do not want to consider them as cracks. So what is the size of the smallest blob that you consider to be a legitimate crack? One pixel? Probably not? 10 pixels - again too small. Whatever it is, you must specify it so that it does not get counted. You can get the areas of all the blobs before doing bwareafilt by doing this:

% Measure the areas
props = regionprops(mask, 'Area');
allAreas = sort([props.Area], 'descend')

This will give you an idea of what particle sizes you have in your image after thresholding.

[yasmin ismail]

you mean these functions it will give me the whole area of the noise (blobs) without the areas of cracks?

[yasmin ismail]

when i tried these two functions before using bwareafilt and click on run I got like this :

allAreas =

Columns 1 through 7

0 0 0 0 0 0 1

Columns 8 through 14

5 10 9 13 14 43 130

Columns 15 through 21

345 647 623 729 795 922 1010

Columns 22 through 28

1023 1140 1069 1215 1413 2129 5782

Columns 29 through 35

15098 25281 32779 34692 34480 27301 19644

Columns 36 through 42

14152 10536 8292 6419 5051 4379 3350

Columns 43 through 49

2921 2899 2995 2539 2206 1778 1151

Columns 50 through 56

1132 986 805 714 761 533 346

Columns 57 through 63

460 299 349 245 285 208 158

Columns 64 through 70

179 156 154 151 109 99 83

Columns 71 through 77

92 79 79 115 69 86 87

Columns 78 through 84

81 71 76 67 50 68 58

Columns 85 through 91

62 97 35 63 57 74 47

Columns 92 through 98

68 58 61 48 45 58 44

Columns 99 through 105

38 64 59 46 41 48 64

Columns 106 through 112

41 65 45 23 59 33 54

Columns 113 through 119

45 43 48 57 34 66 64

Columns 120 through 126

59 70 50 25 62 71 61

Columns 127 through 133

50 47 53 45 102 62 38

Columns 134 through 140

38 91 55 98 47 97 54

Columns 141 through 147

65 79 51 51 101 60 80

Columns 148 through 154

68 48 70 67 75 74 56

Columns 155 through 161

91 54 61 104 118 53 106

Columns 162 through 168

19 43 64 68 53 32 46

Columns 169 through 175

26 31 35 22 39 17 44

Columns 176 through 182

16 16 15 30 20 24 23

Columns 183 through 189

36 25 18 14 10 16 15

Columns 190 through 196

14 11 28 14 8 10 7

Columns 197 through 203

19 8 17 2 4 16 14

Columns 204 through 210

4 17 7 5 4 5 4

Columns 211 through 217

7 6 2 5 3 2 5

Columns 218 through 220

3 2 1

>>

[Image Analyst]

Those aren't sorted like I had in my code. Try this:

% Measure the areas
props = regionprops(mask, '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?

[yasmin ismail]

clc

clear all

fontSize = 10;

I=imread('crack2.jpg');

II=rgb2gray(I);

I2 = II(1:640,350:end);

subplot(2,2,1);imshow(I2);title('Original Image');

subplot(2,2,2);histogram(I2);title('histogram Gray image');

III=I2>90;

mask = imfill(III, 'holes');

Bw=imcomplement(mask);

se=strel('sphere',2);

filteredImage = imopen(Bw,se);

subplot(2, 2, 3);imshow(filteredImage);title('Opened Image');

minAcceptableArea =10000;

mask = bwareafilt(filteredImage,[minAcceptableArea]);

subplot(2, 2, 4);imshow(mask);

props = regionprops(mask, 'Area');

allAreas = [props.Area];

out = bwferet(mask);

averageLengths = sort(allAreas / out.MaxDiameter, 'descend');

message = sprintf('The average width = %.2f pixels.', averageLengths(1));

fprintf('%s\n', message);

caption = sprintf('Binary Image of Cracks larger than %d\nAverage width of largest = %.2f pixels', minAcceptableArea, averageLengths(1));

title(caption, 'FontSize', fontSize, 'Interpreter', 'None');

uiwait(helpdlg(message));

I attached the image that I used and;

I found the crack area between 10000 and 15000 wich is the less frequency in histogram as I understand, when i put minAcceptablearea=10000 the width not apear, only at 1 but unfortunatley some parts of cracks lost and samething happened with other image I dont know why?

[Image Analyst]

This works fine for this particular image:

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 = 14;

%----------------------------------------------------------------------------
% Read in image and convert to gray scale.
rgbImage = imread('crack2.jpg');
grayImage = rgb2gray(rgbImage);

%----------------------------------------------------------------------------
% Crop image.
I2 = grayImage(1:640, 350:end);
subplot(2,2,1);
imshow(I2);
title('Original Image');

%----------------------------------------------------------------------------
% Take histogram.
subplot(2,2,2);
histogram(I2);
grid on;
title('histogram Gray image');

%----------------------------------------------------------------------------
% Binarize image.
threshold = 90;
xline(threshold, 'Color', 'r', 'LineWidth',2)
mask=I2 < threshold;
mask = imfill(mask, 'holes');

%----------------------------------------------------------------------------
% Filter image.
se=strel('sphere',2);
filteredImage = imopen(mask, se);
subplot(2, 2, 3);
imshow(filteredImage);
title('Opened Image');
minAcceptableArea = 10000;
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)
% 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]

first of all thanks alot

second I have few questions because i have litle confusing as following:

1- Is mask=I2 < threshold;

mask = imfill(mask, 'holes');

same to mask = imfill(III, 'holes');

Bw=imcomplement(mask);

I didnt understand mask=I2 < threshold; why this function? is it same to imcomplement function?

2- the second point is variable mask, its repeated many times starting from binarize

i usually when start new funtion I used different variable name, can i use the same variable name?what I mean as following:

mask=I2 < threshold;

mask = imfill(mask, 'holes'); here 2 masks I get confusing

mask = bwareafilt(filteredImage, [minAcceptableArea, inf]); related to which mask?

[Image Analyst]

We want to find dark blobs. So rather than finding blobs brighter than 90 and inverting it, we can just get there directly by thresholding to find blobs darker than 90.

imfill() fills dark holes in a white mask so that it is solid.

I don't see any real need to keep the intermediate versions of mask. Why would you want to save those? I mean you could if you want, like if you're going to need them later on but if it's just an intermediate step in a pipeline of operations then normally there is no need to keep all the intermediate masks. It just takes up memory for no reason. If you want to see what that operation did, you can call imshow() right after you get the new mask.

[yasmin ismail]

So, forthe second part if I want to display the mask at the following :

(1)

threshold = 90;

xline(threshold, 'Color', 'r', 'LineWidth',2)

mask=I2 < threshold;

mask = imfill(mask, 'holes')

imshow(mask)

and then later

(2)

minAcceptableArea = 10000;

mask = bwareafilt(filteredImage, [minAcceptableArea, inf]);

subplot(2, 2, 4);

imshow(mask);

so, i have to use another variable in (2) not to use same variable (mask) if i want to display both , otherwise its fine to use same variable name ,right?

[yasmin ismail]

another question what is (inf) which used in the following

mask = bwareafilt(filteredImage,[minAcceptableArea,inf]);

[Image Analyst]

@yasmin ismail did you look in the help? The second value is the largest area you want to allow. If it's inf, then there is no upper limit on the size of the blobs to extract. If that's the case you can also use bwareaopen:

mask = bwareaopen(filteredImage, minAcceptableArea);

bwareaopen just gives you an image with blobs minAcceptableArea pixels or larger in area.

To display mask at different stages you can simply display it as soon as it's generated. If mask is changed, then it will NOT automatically change/update any prior version of it that you displayed. Of course if you want to show both at the same time you'll need to have them in different suplots, otherwise the later one will overwrite the earlier one and you'll just see the latest on.

[yasmin ismail]

and can you explain the following

props = regionprops(mask, 'Area'); this one to return area measurment

allAreas = [props.Area]; this one not clear

out = bwferet(mask) this also not clear

[Image Analyst]

@yasmin ismail I explained the bwareafilt in my last comment.

bwferet is explained in the help but it is basically the caliper width of blobs.

props is a structure and putting brackets around one of the fields is a shortcut way to extract and concatenate all the structure fields, from each structure in the structure array, into a 1-D list of all the numbers.

regionprops is the main function in the Image Processing Toolbox to make various measurements on the segmented (binary/logical) mask image.

Have you gone through my Image Segmentation Tutorial?

Image Analyst's File Exchange

Answer 3

[yasmin ismail]

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('7002-9.jpg');
I = rgb2gray(rgbImage);
subplot(2,2,1);
imshow(I);
title('7001-21');
% 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('sphere',1);
filteredImage = imopen(mask, se);
subplot(2, 2, 3);
imshow(filteredImage);
title('Opened 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));

I got the answer of crack width as shown in attachment. My question is, as you see after filteration, there are discontinuous lines of cracks. How can I connect between them to be continuous?

[Image Analyst]

@yasmin ismail if you look at the places where there are gaps in your crack, you can see that the crack is very faint or nonexistent there. A little bit of gray breaks the crack into two cracks. You can try adjusting the threshold and see if you can raise it enough to consider the gray parts in the gaps as cracks.

[yasmin ismail]

I noticed something that when I used

se=strel('sphere',1);

or

se = strel("line",1,100);

I got different width in each, I dont know why and how can I fix it to choose the correct one?

[Image Analyst]

Of course, because you're using different shapes to eat away at your crack. I don't even know why you're using imopen on your mask anyway. That will just make the gaps worse (more open). You should be using something like imclose

se = strel('disk', 5, 0);
mask = imopen(mask, se);

This will close up the gaps. Change the 5 to various even or odd numbers and see how they affect the mask. Use the value that produces the final crack shape you like best.

[yasmin ismail]

you write in code (imopen) you mean

mask= imclose(mask,se);

right?

[Image Analyst]

Sorry, yes. Use imclose. (Copy and paste error).

[yasmin ismail]

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

[Image Analyst]

@yasmin ismail, so we don't keep bothering @Ilgin Sandalci with emails on your problem, please start your own discussion thread and I'll help you there when I get more time.

[yasmin ismail]

I am very sorry @llgin Sandalci I didnt know that my questions are sending to him

and I have already post my questions in new discussion.

Hope @Image Analyst have time to answer them

How to measure crack width - MATLAB Answers - MATLAB Central (mathworks.com)

Sorry again