How to get the length of n curves present in an image in Matlab?

Question

Currently I have been working on obtaining the length of a curve, with the following code I have managed to get the length of a curve present in an image.

test image one curve

Then I paste the code that I used to get the length of the curve of a simple image. What I did is the following:

1. I got the columns and rows of the image
2. I got the columns in x and the rows in y
3. I obtained the coefficients of the curve, based on the formula of the parable
4. Build the equation

5. Implement the arc length formula to obtain the length of the curve

   grayImage = imread(fullFileName);
   [rows, columns, numberOfColorBands] = size(grayImage);
   if numberOfColorBands > 1
     grayImage = grayImage(:, :, 2); % Take green channel.
   end
   subplot(2, 2, 1);
   imshow(grayImage, []);
   
   % Get the rows (y) and columns (x).
   [rows, columns] = find(binaryImage);
   
   coefficients = polyfit(columns, rows, 2); % Gets coefficients of the    formula.
   
   % Fit a curve to 500 points in the range that x has.
   fittedX = linspace(min(columns), max(columns), 500);
   
   % Now get the y values.
   fittedY = polyval(coefficients, fittedX);
   
   % Plot the fitting:
   
   subplot(2,2,3:4);
   plot(fittedX, fittedY, 'b-', 'linewidth', 4);
   grid on;
   xlabel('X', 'FontSize', fontSize);
   ylabel('Y', 'FontSize', fontSize);
   % Overlay the original points in red.
   hold on;
   plot(columns, rows, 'r+', 'LineWidth', 2, 'MarkerSize', 10)
   
   formula = poly2sym([coefficients(1),coefficients(2),coefficients(3)]);
   % formulaD = vpa(formula)
   df=diff(formula);
   df = df^2;
   
   f= (sqrt(1+df));
   
   i = int(f,min(columns),max(columns));
   j = double(i);
   disp(j);
   

Now I have the image 2 which has n curves, I do not know how I can do to get the length of each curve

test image n curves

Answer 1

I suggest you to look at Hough Transformation: https://uk.mathworks.com/help/images/hough-transform.html

You will need Image Processing Toolbox. Otherwise, you have to develop your own logic.

https://en.wikipedia.org/wiki/Hough_transform

Update 1

I had a two-hour thinking about your problem and I'm only able to extract the first curve. The problem is to locate the starting points of the curves. Anyway, here is the code I come up with and hopefully will give you some ideas for further development.

clc;clear;close all;

grayImage = imread('2.png');
[rows, columns, numberOfColorBands] = size(grayImage);
if numberOfColorBands > 1
    grayImage = grayImage(:, :, 2); % Take green channel.
end

% find edge.
bw = edge(grayImage,'canny');
imshow(bw);

[x, y] = find(bw == 1);
P = [x,y];

% For each point, find a point that is of distance 1 or sqrt(2) to it, i.e.
% find its connectivity.
cP = cell(1,length(x));
for i = 1:length(x)
    px = x(i);
    py = y(i);
    dx = x - px*ones(size(x));
    dy = y - py*ones(size(y));
    distances = (dx.^2 + dy.^2).^0.5;
    cP{i} = [x(distances == 1), y(distances == 1);
        x(distances == sqrt(2)), y(distances == sqrt(2))];
end

% pick the first point and a second point that is connected to it.
fP = P(1,:);
Q(1,:) = fP;
Q(2,:) = cP{1}(1,:);
m = 2;
while true
    % take the previous point from point set Q, when current point is
    % Q(m,1)
    pP = Q(m-1,:);
    % find the index of the current point in point set P.
    i = find(P(:,1) == Q(m,1) & P(:,2) == Q(m,2));
    % Find the distances from the previous points to all points connected
    % to the current point.
    dx = cP{i}(:,1) - pP(1)*ones(length(cP{i}),1);
    dy = cP{i}(:,2) - pP(2)*ones(length(cP{i}),1);
    distances = (dx.^2 + dy.^2).^0.5;
    % Take the farthest point as the next point.
    m = m+1;
    p_cache = cP{i}(find(distances==max(distances),1),:);
    % Calculate the distance of this point to the first point.
    distance = ((p_cache(1) - fP(1))^2 + (p_cache(2) - fP(2))^2).^0.5;
    if distance == 0 || distance == 1
        break;
    else
        Q(m,:) = p_cache;
    end
end
% By now we should have built the ordered point set Q for the first curve.
% However, there is a significant weakness and this weakness prevents us to
% build the second curve.

---

Update 2

Some more work since the last update. I'm able to separate each curve now. The only problem I can see here is to have a good curve fitting. I would suggest B-spline or Bezier curves than polynomial fit. I think I will stop here and leave you to figure out the rest. Hope this helps.

Note that the following script uses Image Processing Toolbox to find the edges of the curves.

clc;clear;close all;

grayImage = imread('2.png');
[rows, columns, numberOfColorBands] = size(grayImage);
if numberOfColorBands > 1
    grayImage = grayImage(:, :, 2); % Take green channel.
end

% find edge.
bw = edge(grayImage,'canny');
imshow(bw);

[x, y] = find(bw == 1);
P = [x,y];

% For each point, find a point that is of distance 1 or sqrt(2) to it, i.e.
% find its connectivity.
cP =[0,0]; % add a place holder
for i = 1:length(x)
    px = x(i);
    py = y(i);
    dx = x - px*ones(size(x));
    dy = y - py*ones(size(y));
    distances = (dx.^2 + dy.^2).^0.5;
    c = [find(distances == 1); find(distances == sqrt(2))];
    cP(end+1:end+length(c),:) = [ones(length(c),1)*i, c];
end
cP (1,:) = [];% remove the place holder

% remove duplicates
cP = unique(sort(cP,2),'rows');

% seperating curves
Q{1} = cP(1,:);
for i = 2:length(cP)
    cp = cP(i,:);
    % search for points in cp in Q.
    for j = 1:length(Q)
        check = ismember(cp,Q{j});
        if ~any(check) && j == length(Q) % if neither has been saved in Q
            Q{end+1} = cp;
            break;
        elseif sum(check) == 2 % if both points cp has been saved in Q
            break;
        elseif sum(check) == 1 % if only one of the points exists in Q, add the one missing.
            Q{j} = [Q{j}, cp(~check)];
            break;
        end
    end
    % review sets in Q, merge the ones having common points
    for j = 1:length(Q)-1
        q = Q{j};
        for m = j+1:length(Q)
            check = ismember(q,Q{m});
            if sum(check)>=1 % if there are common points
                Q{m} = [Q{m}, q(~check)]; % merge
                Q{j} = []; % delete the merged set
                break;
            end
        end
    end
    Q = Q(~cellfun('isempty',Q)); % remove empty cells;
end
% each cell in Q represents a curve. Note that points are not ordered.

figure;hold on;axis equal;grid on;
for i = 1:length(Q)
    x_ = x(Q{i});
    y_ = y(Q{i});

    coefficients = polyfit(y_, x_, 3); % Gets coefficients of the    formula.

    % Fit a curve to 500 points in the range that x has.
    fittedX = linspace(min(y_), max(y_), 500);

    % Now get the y values.
    fittedY = polyval(coefficients, fittedX);
    plot(fittedX, fittedY, 'b-', 'linewidth', 4);
    % Overlay the original points in red.
    plot(y_, x_, 'r.', 'LineWidth', 2, 'MarkerSize', 1)

    formula = poly2sym([coefficients(1),coefficients(2),coefficients(3)]);
    % formulaD = vpa(formula)
    df=diff(formula);
    lengthOfCurve(i) = double(int((sqrt(1+df^2)),min(y_),max(y_)));
end

Result: