Crack detection on concrete

Hi, i have problem in measuring the length and width the crack image. i have follow the code but it only measure straight line (endpoint). I want the measurement by following the crack pattern. Please help me on this.

      %load image
      path=imgetfile();
      originalImage1=imread(path);
      imshow(originalImage1);

      %image to grayscale
      originalImage = rgb2gray(originalImage1);
      figure
      imshow(originalImage);

      %treshold 
      thresholdValue = 100;
      binaryImage = originalImage > thresholdValue;
      binaryImage = imfill(binaryImage, 'holes');
      hold on;
      % maxYValue = ylim;
      % line([thresholdValue, thresholdValue], maxYValue, 'Color', 'r');
      % annotationText = sprintf('Thresholded at %d gray levels', thresholdValue);
      % text(double(thresholdValue + 5), double(0.5 * maxYValue(2)), annotationText, 'FontSize', 10, 'Color', [0 .5 0]);
      % text(double(thresholdValue - 70), double(0.94 * maxYValue(2)), 'Background', 'FontSize', 10, 'Color', [0 0 .5]);
      % text(double(thresholdValue + 50), double(0.94 * maxYValue(2)), 'Foreground', 'FontSize', 10, 'Color', [0 0 .5]);
      figure
      imshow(binaryImage);

      % bw = bwareaopen(binaryImage,100,8);
      bw= bwareaopen(binaryImage, round(0.1*numel(binaryImage)));
      figure
      imshow(bw)

      g=bwmorph(bw,'clean');
      figure
      imshow(g);

      % %labeled image (on hold)
      labeledImage = bwlabel(g,4); 
      figure
      imshow(labeledImage, []);

      %colored label image
      coloredLabels = label2rgb (labeledImage, 'hsv', 'k', 'shuffle');
      figure
      imshow(coloredLabels);
      axis image;

      %blobmeasurement
      blobMeasurements = regionprops(labeledImage, originalImage, 'all');
      numberOfBlobs = size(blobMeasurements, 1);
      figure
      imshow(originalImage);

      axis image; % Make sure image is not artificially stretched because of screen's aspect ratio.
      hold on;
      boundaries = bwboundaries(binaryImage);
      numberOfBoundaries = size(boundaries, 1);
      for k = 1 : numberOfBoundaries
        thisBoundary = boundaries{k};
        plot(thisBoundary(:,2), thisBoundary(:,1), 'g', 'LineWidth', 2);
      end
      hold off;

      % Measure the area
      measurements = regionprops(labeledImage, 'Area');

      boundaries = bwboundaries(binaryImage);
      numberOfBoundaries = size(boundaries, 1);
      for blobIndex = 1 : numberOfBoundaries
        thisBoundary = boundaries{blobIndex};
        x = thisBoundary(:, 2); % x = columns.
        y = thisBoundary(:, 1); % y = rows.
        
        % Find which two bounary points are farthest from each other.
        maxDistance = -inf;
        for k = 1 : length(x)
          distances = sqrt( (x(k) - x) .^ 2 + (y(k) - y) .^ 2 );
          [thisMaxDistance, indexOfMaxDistance] = max(distances);
          if thisMaxDistance > maxDistance
            maxDistance = thisMaxDistance;
            index1 = k;
            index2 = indexOfMaxDistance;
          end
          end

      % Find the midpoint of the line.
        xMidPoint = mean([x(index1), x(index2)]);
        yMidPoint = mean([y(index1), y(index2)]);
        longSlope = (y(index1) - y(index2)) / (x(index1) - x(index2))
        perpendicularSlope = -1/longSlope
        % Use point slope formula (y-ym) = slope * (x - xm) to get points
        y1 = perpendicularSlope * (1 - xMidPoint) + yMidPoint;
        y2 = perpendicularSlope * (columns - xMidPoint) + yMidPoint;
        
        % Get the profile perpendicular to the midpoint so we can find out when if first enters and last leaves the object.
        [cx,cy,c] = improfile(binaryImage,[1, columns], [y1, y2], 1000);
        % Get rid of NAN's that occur when the line's endpoints go above or below the image.
        c(isnan(c)) = 0;
        firstIndex = find(c, 1, 'first');
        lastIndex = find(c, 1, 'last');
        % Compute the distance of that perpendicular width.
        perpendicularWidth = sqrt( (cx(firstIndex) - cx(lastIndex)) .^ 2 + (cy(firstIndex) - cy(lastIndex)) .^ 2 );
        % Get the average perpendicular width.  This will approximately be the area divided by the longest length.
        averageWidth = measurements(blobIndex).Area / maxDistance;
        
        % Plot the boundaries, line, and midpoints over the two images.
        % Plot the boundary over the gray scale image
        subplot(2, 2, 3);
        plot(x, y, 'y-', 'LineWidth', 3);
        % For this blob, put a line between the points farthest away from each other.
        line([x(index1), x(index2)], [y(index1), y(index2)], 'Color', 'r', 'LineWidth', 3);
        plot(xMidPoint, yMidPoint, 'r*', 'MarkerSize', 15, 'LineWidth', 2);
        % Plot perpendicular line.  Make it green across the whole image but magenta inside the blob.
        line([1, columns], [y1, y2], 'Color', 'g', 'LineWidth', 3);  
        line([cx(firstIndex), cx(lastIndex)], [cy(firstIndex), cy(lastIndex)], 'Color', 'm', 'LineWidth', 3);
        
        % Plot the boundary over the binary image
        subplot(2, 2, 4);
        plot(x, y, 'y-', 'LineWidth', 3);
        % For this blob, put a line between the points farthest away from each other.
        line([x(index1), x(index2)], [y(index1), y(index2)], 'Color', 'r', 'LineWidth', 3);
        plot(xMidPoint, yMidPoint, 'r*', 'MarkerSize', 15, 'LineWidth', 2);
        % Plot perpendicular line.  Make it green across the whole image but magenta inside the blob.
        line([1, columns], [y1, y2], 'Color', 'g', 'LineWidth', 3);  
        line([cx(firstIndex), cx(lastIndex)], [cy(firstIndex), cy(lastIndex)], 'Color', 'm', 'LineWidth', 3);
        
        message = sprintf('The longest line is red.\nPerpendicular to that, at the midpoint, is green.\nMax distance for blob #%d = %.2f\nPerpendicular distance at midpoint = %.2f\nAverage perpendicular width = %.2f (approximately\nArea = %d', ...
          blobIndex, maxDistance, perpendicularWidth, averageWidth, measurements(blobIndex).Area);
        fprintf('%s\n', message);
        uiwait(helpdlg(message));
      end
      hold off;

[KAMOOSH BABA SHAIK]

% Use point slope formula (y-ym) = slope * (x - xm) to get points

y1 = perpendicularSlope * (1 - xMidPoint) + yMidPoint;

y2 = perpendicularSlope * (columns - xMidPoint) + yMidPoint;

in y2 what i can replace to get the answer

in matlab shows like this,

Unrecognized function or variable 'columns'.

Error in lengthwidth (line 90)

y2 = perpendicularSlope * (columns - xMidPoint) + yMidPoint;

Answer 1

[Walter Roberson]

Skeletonize and then bwdistgeodesic

[Fairuz Husna]

can you show how can i use bwdistgeodesic in the code above. i'm new with matlab.

Answer 2

[Preetham Manjunatha]

Here is the MATLAB code to calculate/estimate the crack area, width and length. In addition, this package assumes the crack is segmented either using morphological method or multiscale gradient-based or deep learning semantic segmentation methods. This package estimates the crack area, width and length (pixel scale can be provided to estimate these physical quantities).