MATLAB Answers


Help with the bwtraceboundary function

Asked by mounim
on 9 Dec 2012


i am a totally new beginner in matlab, and i hope someone can help me with this code, what i want to do is finding the intersection point between the edge of the bone and the horizontal line, i tried to modify the code in the demo section and i got it to work for the left side of the bone, i want to do the same thing for the other side of the bone as well.

please how can i achieve that ?

( i attached my code, and the output image )

      %%%%%%%% LOADING THE IMAGE %%%%%%%%%
RGB = imread('1.jpg');
%%%%%%%% EXTRACT THE ROI %%%%%%%%%%%%%%%
% you can obtain the coordinates of the rectangular region using
% pixel information displayed by imtool
start_row = 1;
start_col = 1;
start_row2 = 1100;
start_col2 = 1;
cropRGB = RGB(start_row:115, start_col:1210, :);
cropRGB2 = RGB(start_row2:1265, start_col2:1210, :);
% Store (X,Y) offsets for later use; subtract 1 so that each offset will
% correspond to the last pixel before the region of interest
offsetX = start_col-1;
offsetY = start_row-1;
offsetX2 = start_col2-1;
offsetY2 = start_row2-1;
%%%%% THRESHOLD THE IMAGE %%%%%%%%%
I = rgb2gray(cropRGB);
II = rgb2gray(cropRGB2);
threshold = graythresh(I);
BW = im2bw(I,threshold);
%BW = ~BW;  % complement the image (objects of interest must be white)
threshold2 = graythresh(II);
BW2 = im2bw(II,threshold2);
%BW = ~BW;  % complement the image (objects of interest must be white)
%%%%%%%%%% INITIAL POINT OF BOUNDERY %%%%%%%%%%%%%%%
dim = size(BW);
dim2 = size(BW2);
% horizontal beam
col1 = 420;
row1 = min(find(BW(:,col1)));
col11 = 420;
row11 = min(find(BW2(:,col11)));
% angled beam
row2 = 12;
col2 = min(find(BW(row2,:)))
row22 = 60;
col22 = min(find(BW2(row22,:)))
%%%%%%%%% TRACING THE BOUNDERY %%%%%%%%%%%
boundary1 = bwtraceboundary(BW, [row1, col1], 'N', 8, 70);
boundary11 = bwtraceboundary(BW2, [row11, col11], 'N', 8, 70);
% set the search direction to counterclockwise, in order to trace downward.
boundary2 = bwtraceboundary(BW, [row2, col2], 'E', 8, 90,'counter');
boundary22 = bwtraceboundary(BW2, [row22, col22], 'E', 8, 90,'counter');
%imshow(RGB); hold on;
% apply offsets in order to draw in the original image
%%%%%%%%%% FIT LINES TO THE BOUNDERY %%%%%%%
ab1 = polyfit(boundary1(:,2), boundary1(:,1), 1);
ab2 = polyfit(boundary2(:,2), boundary2(:,1), 1);
ab11 = polyfit(boundary11(:,2), boundary11(:,1), 1);
ab22 = polyfit(boundary22(:,2), boundary22(:,1), 1);
vect1 = [1 ab1(1)]; % create a vector based on the line equation
vect2 = [1 ab2(1)];
dp = dot(vect1, vect2);
vect11 = [1 ab11(1)]; % create a vector based on the line equation
vect22 = [1 ab22(1)];
dp2 = dot(vect11, vect22);
% compute vector lengths
length1 = sqrt(sum(vect1.^2));
length2 = sqrt(sum(vect2.^2));
% compute vector lengths
length11 = sqrt(sum(vect11.^2));
length22 = sqrt(sum(vect22.^2));
% obtain the larger angle of intersection in degrees
angle = 180-acos(dp/(length1*length2))*180/pi
angle1 = 180-acos(dp2/(length11*length22))*180/pi
intersection = [1 ,-ab1(1); 1, -ab2(1)] \ [ab1(2); ab2(2)];
% apply offsets in order to compute the location in the original,
% i.e. not cropped, image.
intersection = intersection + [offsetY; offsetX]
intersection1 = [1 ,-ab11(1); 1, -ab22(1)] \ [ab11(2); ab22(2)];
% apply offsets in order to compute the location in the original,
% i.e. not cropped, image.
intersection1 = intersection1 + [offsetY2; offsetX2]
%%%%%%%%%%%% PLOT THE RESULTS %%%%%%%%%
inter_x = intersection(2);
inter_y = intersection(1);
inter_x1 = intersection1(2);
inter_y1 = intersection1(1);
% draw an "X" at the point of intersection
text(inter_x-100, inter_y-30, [sprintf('%1.3f',angle),'{\circ}'],...
%interString = sprintf('(%2.1f,%2.1f)', inter_x, inter_y);
text(inter_x1-100, inter_y1-30, [sprintf('%1.3f',angle1),'{\circ}'],...
%interString = sprintf('(%2.1f,%2.1f)', inter_x1, inter_y1);


2 Answers

Answer by Image Analyst
on 9 Dec 2012
Edited by Image Analyst
on 9 Dec 2012
 Accepted answer

Don't use bwtraceboundary() - it's almost never necessary. I think you're making this way way too complicated. Just scan down the left and right columns to find the line, then find the first and last point on the row where it's white. First convert the image to binary.

rightColumn = binaryImage(:, end); % Get right column
topLine = find(rightColumn, 1, 'first'); % Find row where line is
% Extract the row above
oneRow = binaryImage(topLine-1, :);
% Find first and last places - the bone edges
leftEdge1 = find(oneRow, 1, 'first');
rightEdge1 = find(oneRow, 1, 'last');
bottomLine = find(rightColumn, 1, 'last'); % Find row where line is
% Extract the row above
oneRow = binaryImage(bottomLine -3, :); % 3 (or whatever) lines up.
% Find first and last places - the bone edges
leftEdge2 = find(oneRow, 1, 'first');
rightEdge2 = find(oneRow, 1, 'last');


i am kind of confused, that's the closest i got :

for k=50:100
oneRow{k} = binaryImage(k, :);
leftEdge{k} = find(oneRow{k}, 1, 'first');
rightEdge{k} = find(oneRow{k}, 1, 'last');
hold on;
midpoints(k) = (leftEdge{k} + rightEdge{k})/2;
plot(leftEdge{k}, k, 'rx','LineWidth',2);
plot(rightEdge{k}, k, 'rx','LineWidth',2);
coeffs = polyfit(midpoints,Y,1);

i am sorry i dont think i am getting how to fit the polyfit in. i tried looking at the documentation. but i dont know.

That looks about right except that you're still using cell arrays for oneRow and leftEdge and rightEdge when you don't have to, but it doesn't hurt anything other than making it more coplicated. Just get rid of the {k} from all of them. Other than that, what's the problem? It looks like it should work once you flip Y and midpoints, since Y is your "x" or independent axis in this case.

coeffs = polyfit(Y, midpoints, 1);

Thanks again, well i just cannot draw that line, what i am doing wrong ?. Here is my code :

for k=50:100
oneRow = binaryImage(k, :);
leftEdge = find(oneRow, 1, 'first');
rightEdge = find(oneRow, 1, 'last');
hold on;
midpoints(k) = (leftEdge + rightEdge)/2;
plot(leftEdge, k, 'rx','LineWidth',2);
plot(rightEdge, k, 'rx','LineWidth',2);
coeffs = polyfit(X,midpoints, 1);
Y = coeffs(1) .* X + coeffs(2);
plot(X,Y, '-','LineWidth',2)

And the output in matlab ( i used a white rectangle, to be sure that things are clean )

Answer by mounim
on 14 Dec 2012

It is okay, i figured that out !! ... i had to remove the zeros - from the array in my for-loop. so now it is working.

thanks for your help


Join the 15-year community celebration.

Play games and win prizes!

Learn more
Discover MakerZone

MATLAB and Simulink resources for Arduino, LEGO, and Raspberry Pi

Learn more

Discover what MATLAB® can do for your career.

Opportunities for recent engineering grads.

Apply Today

MATLAB Academy

New to MATLAB?

Learn MATLAB today!