cluster analysis from edge detection

2020 5 3
2 个回答
回答已采纳
更新时间：2020 5 14
3 次查看（30 天）

0 个投票

Hello,
I am looking to perform a cluster analysis (K-means) for an image. I've converted the original image to a binary image and performed an edge detection, which should help with the clustering, but am stuck on how to implement a cluster analysis from the binary image. Any suggestions are greatly appreciated! The code for the binary edge detection is below and the image is attached. The end goal here is to identify patterns of wood orientation.
clear
RGB=imread('DJI_0022.jpg'); %inputs image
I=rgb2gray(RGB); %convers to grayscale
figure
imshow(I)
BW1 = edge(I,'Canny',0.6);
imshow(BW1)

请先登录，再回答此问题。

 采纳的回答

Image Analyst
Image Analyst 2020-5-3

0 个投票

Assuming you have a bunch of blobs that are the tree edges, I'd call regionprops() and ask for 'Orientation'. This will give you the average angle of every blob. You can then histogram the angles if you want.
props = regionprops(edgeImage, 'Orientation');
allAngles = [props.Orientation];
histogram(allAngles);
grid on;
xlabel('Angle', 'FontSize', 20);
ylabel('Count', 'FontSize', 20);

16 个评论
显示 14更早的评论 隐藏 14更早的评论

Anna Marshall
Anna Marshall 2020-5-3
Thank you! I have calculated the orientations of the wood pieces, but actually am filtering out the trees because I only want to focus on the wood orientations for clustering. With that being said, is it possible to do a similar thing as what you mentioned above, but only look at the wood orientations based on the edge detection output (filtered out to remove trees). I'm curious as to the interpretation of the histogram. Is that grouping similar angles within the image? Thanks again!
Image Analyst
Image Analyst 2020-5-3
You can filter out the trees by doing a color segmentation. Find green and mask it away.
I do a histogram because that's what makes sense. You can choose the bin width if you want. Why use kmeans()? That makes no sense at all to me. Why cluster the angles into some predetermined number of angular groups, where the angles might change from image to image? I can't imagine any use case where that would be desireable and better than a histogram, which will accurately give a distribution of angles.
Anna Marshall
Anna Marshall 2020-5-3
Thank you! I'm following most of what you suggested, but don't have a clear understanding of what the histogram in the code you suggested above is providing? My research goals are trying to cluster wood orientations in the image by orientation so I was thinking that some sort of cluster analysis would achieve this, although without ever having done one. I don't want to cluster the angles into some predetermined number of groups, but would rather like to have the analysis identify where those grouping patterns are.
Image Analyst
Image Analyst 2020-5-3
You can get the angles from the histogram of the orientations. You can look at various peaks in the histogram. Look at it this way: let's say the orientations were completely random from 0 to 360 so that the histogram was flat or almost flat. What would you say k should be? 2? 3? 10? Whatever you pick is arbitrary, and the groupings would be aribtrary and changing from one image to the next. Sure if you had trees only at 30, 60, and 90 degrees, you might say k=3, but what is the likelihood you'll have 3 groupings? Chances are there will be trees at in between angles, so then where would kmeans put the centroids? Maybe at 25, 64, 87 one time, and for another image at 29, 51, 89 degrees. Why would you want that? I still maintain you would want to look at how many trees there are at every angle or in every 10 degree range, depending on how many bins you told it to use. You could use findpeaks() on the histogram counts if you want to find local groupings. But at least with a histogram you have control over the angles get counted. With kmeans you don't have that control. All you can control is how many groups there are. Even if there are no groupings, if you tell kmeans to find 5 angle groups, it will be forced to find 5 whether or not those are valid groups.
Anna Marshall
Anna Marshall 2020-5-3
That makes a lot of sense! Thank you for the thorough explanation! I understand what the histogram is showing now and agree with the approach you suggested, but it doesn't quite reach my objective. I can use the histogram to see the number of pieces that fall into each angle catagory, but in thinking through this more, I think what I am looking for is something that is slightly more visual (i.e. an overlay on the original image with an outline of the angle groupings) so that I can see how those groupings appear within the image. Does this seem like a possible task in matlab?
Image Analyst
Image Analyst 2020-5-3
编辑：Image Analyst 2020-5-3
Yes, you can use quiver to show little arrows over each tree, if that meets your needs. I haven't used quiver but I think you can color each arrow according to its angle if you want.
This may help you: http://blogs.mathworks.com/steve/2010/07/30/visualizing-regionprops-ellipse-measurements/
Also you could draw a short line at the blob's angle over the centroid into a binary image. Then us imoverlay() to overlay the lines on the original image.
Anna Marshall
Anna Marshall 2020-5-3
Very neat! Quiver sounds like what I might be looking for. I'll look into it. Thanks again!
Image Analyst
Image Analyst 2020-5-3
Anna, see this code where I drew a line over the gray scale image where the line color corresponds to the angle:
% Initialization steps. Brute force cleanup of everything currently existing to start with a clean slate.
% Just for the demo. If you put this into your own function you'll want to get rid of the close and clear commands.
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;
rgbImage = imread('DJI_0022.jpg'); %inputs image
hFig = figure;
subplot(2, 2, 1);
imshow(rgbImage)
axis('on', 'image');
impixelinfo;
title('Original RGB Image', 'FontSize', 20);
grayImage = rgb2gray(rgbImage); %convers to grayscale
subplot(2, 2, 2);
imshow(grayImage)
axis('on', 'image');
impixelinfo;
title('Gray Scale Image', 'FontSize', 20);
edgeImage = edge(grayImage, 'Canny', 0.6);
subplot(2, 2, 3);
imshow(edgeImage)
axis('on', 'image');
title('Canny Edge Image', 'FontSize', 20);
hFig.WindowState = 'maximized';
drawnow;
% Get coordinates and plot a line over them.
props = regionprops(edgeImage, 'PixelList', 'Orientation', 'Centroid');
lineLength = 100; % Whatever you want.
subplot(2, 2, 2); % Switch to the gray scale image.
hold on;
lineColorMap = jet(181); % One color for every 1 degrees.
for k = 1 : length(props)
% thisList = props(k).PixelList;
% fprintf('\nFor blob #%d of %d:\n ', k, length(props))
% for k2 = 1 : size(thisList, 1)
% thisx = thisList(k2, 1);
% thisy = thisList(k2, 2);
% fprintf('(%d, %d), ', thisx, thisy);
% end
% Get the centroid
x = props(k).Centroid(1);
y = props(k).Centroid(2);
% Get the angle.
angle = props(k).Orientation;
% Get the line endpoints.
x1 = x + lineLength * cosd(angle);
x2 = x - lineLength * cosd(angle);
y1 = y + lineLength * sind(angle);
y2 = y - lineLength * sind(angle);
thisColor = lineColorMap(round(angle + 90) + 1, :);
fprintf('Drawing line from (%.1f to %.1f) to (%.1f, %.1f)\n', x1, y1, x2, y2);
line([x1,x2], [y1,y2], 'Color', thisColor);
if mod(k, 100) == 0
drawnow;
end
end
% Plot the histogram of angles.
subplot(2, 2, 4)
allAngles = [props.Orientation];
histObject = histogram(allAngles);
grid on;
xlabel('Angle', 'FontSize', 20);
ylabel('Count', 'FontSize', 20);
title('Histogram of angles', 'FontSize', 20);
Anna Marshall
Anna Marshall 2020-5-3
This is too cool!! Wow! I'm excited to play around with this more, but is exactly what I was trying to figure out a way to do. Thank you!
Image Analyst
Image Analyst 2020-5-3
You're welcome. Thanks for accepting the answer. I actually thought it was pretty cool myself. However you might want to zoom in and make sure that the Canny lines actually line up well with the logs, because to me they look like a bunch of little dots. They might -- it's just worth zooming in and checking.
Anna Marshall
Anna Marshall 2020-5-3
It's very cool! One thing I did notice when zooming in is that the lines don't actually match the orientations of the wood pieces. This might be either an issue with my canny edge detection or orientations from the histogram. Are the colored lines based on the canny detection or directly from region.prop?
Image Analyst
Image Analyst 2020-5-3
Doesn't look like they align that well. Try using imoverlay() to put the Canny image over the gray scale image. If it still doesn't align well, then maybe something's wrong. You might try a COSFIRE filter, Hessian filter, or a Frangi filter which look for ridges. When you think about it, a ridgeline-finding filter should be better than an edge filter. I think they're are some in the File Exchange.
Anna Marshall
Anna Marshall 2020-5-3
Thanks for the suggestion! I'll give that a try.
Anna Marshall
Anna Marshall 2020-5-4
Interestingly, it does look like the Canny filter lines up well with the wood orientations in the image overlay. Looking back at the code, I'm wondering if the code to find the centroid angles before plotting those lines is what is creating a discrepancy in the colored lines?
Image Analyst
Image Analyst 2020-5-4
You might try getting the PixelList, like I already showed you in your other post, and then put those x,y values into polyfit to get the angle. Something like
% Get coordinates and fit to a line.
props = regionprops(edgeImage, 'PixelList', 'Centroid');
for k = 1 : length(props)
thisList = props(k).PixelList;
fprintf('\nGetting angle for blob #%d of %d.\n ', k, length(props))
thisx = thisList(k2, 1);
thisy = thisList(k2, 2);
coefficients = polyfit(thisx, thisy, 1); % Fit to a line
angles(k) = atand(coefficients(1));
end

请先登录，再进行评论。

更多回答（1 个）

Anna Marshall
Anna Marshall 2020-5-4

0 个投票

I'll give it a try! Thanks!

6 个评论
显示 4更早的评论 隐藏 4更早的评论

Anna Marshall
Anna Marshall 2020-5-4
One more question for you! I modified the code slightly to include k2=, but am no longer able to get the orientation lines to show up.
% Get coordinates and fit to a line.
props = regionprops(edgeImage, 'PixelList', 'Centroid');
for k = 1 : length(props)
thisList = props(k).PixelList;
fprintf('\nGetting angle for blob #%d of %d.\n ', k, length(props));
k2 = 1 : size(thisList, 1);
thisx = thisList(k2, 1);
thisy = thisList(k2, 2);
coefficients = polyfit(thisx, thisy, 1); % Fit to a line
angles(k) = atand(coefficients(1));
end
Image Analyst
Image Analyst 2020-5-4
I think it should be
% Get coordinates and fit to a line.
props = regionprops(edgeImage, 'PixelList', 'Centroid');
for k = 1 : length(props)
thisList = props(k).PixelList;
fprintf('\nGetting angle for blob #%d of %d.\n ', k, length(props));
thisx = thisList(:, 1);
thisy = thisList(:, 2);
coefficients = polyfit(thisx, thisy, 1); % Fit to a line
angles(k) = atand(coefficients(1));
end
Anna Marshall
Anna Marshall 2020-5-4
Getting rid of the k2 makes more sense! Thanks. The funny thing is, I keep running into an error in the props=regionprops line now that 'orientation' is removed. I'm using the following as code for this section, but the error is in props = regionprops(edgeImage, 'PixelList', 'Centroid ')
% Get coordinates and plot a line over them.
imshow(grayImage)
axis('on', 'image');
title('Line Grouping', 'FontSize', 20);
hFig.WindowState = 'maximized';
drawnow;
props = regionprops(edgeImage, 'PixelList', 'Centroid ');
lineLength = 20; % Whatever you want.
hold on;
lineColorMap = jet(181); % One color for every 1 degrees.
for k = 1 : length(props)
thisList = props(k).PixelList;
fprintf('\nFor blob #%d of %d:\n ', k, length(props))
thisx = thisList(:, 1);
thisy = thisList(:, 2);
coefficients = polyfit(thisx, thisy, 1); % Fit to a line
angles(k) = atand(coefficients(1 ));
thisColor = lineColorMap(round(angle + 90) + 1, :);
fprintf('Drawing line from (%.1f to %.1f) to (%.1f, %.1f)\n', x1, y1, x2, y2);
line([x1,x2], [y1,y2], 'Color', thisColor);
if mod(k, 100) == 0
drawnow;
end
end
Image Analyst
Image Analyst 2020-5-4
Anna, this code works but it looks like you're still missing quite a few edges with that filter.
% Initialization steps. Brute force cleanup of everything currently existing to start with a clean slate.
% Just for the demo. If you put this into your own function you'll want to get rid of the close and clear commands.
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;
rgbImage = imread('DJI_0022.jpg'); %inputs image
hFig = figure;
subplot(2, 2, 1);
imshow(rgbImage)
axis('on', 'image');
impixelinfo;
title('Original RGB Image', 'FontSize', 20);
grayImage = rgb2gray(rgbImage); %convers to grayscale
subplot(2, 2, 2);
imshow(grayImage)
axis('on', 'image');
impixelinfo;
title('Gray Scale Image', 'FontSize', 20);
edgeImage = edge(grayImage, 'Canny', 0.6);
subplot(2, 2, 3);
imshow(edgeImage)
axis('on', 'image');
title('Canny Edge Image', 'FontSize', 20);
hFig.WindowState = 'maximized';
drawnow;
% Get coordinates and plot a line over them.
imshow(grayImage)
axis('on', 'image');
title('Line Grouping', 'FontSize', 20);
hFig.WindowState = 'maximized';
drawnow;
props = regionprops(edgeImage, 'PixelList', 'Centroid');
lineLength = 20; % Whatever you want.
hold on;
lineColorMap = jet(181); % One color for every 1 degrees.
angles = zeros(1, length(props));
for k = 1 : length(props)
thisList = props(k).PixelList;
fprintf('\nFor blob #%d of %d:\n ', k, length(props))
thisx = thisList(:, 1);
thisy = thisList(:, 2);
coefficients = polyfit(thisx, thisy, 1); % Fit to a line
angles(k) = atand(coefficients(1 ));
% Get the centroid
x = props(k).Centroid(1);
y = props(k).Centroid(2);
% Get the line endpoints.
x1 = x + lineLength * cosd(angles(k));
x2 = x - lineLength * cosd(angles(k));
y1 = y + lineLength * sind(angles(k));
y2 = y - lineLength * sind(angles(k));
thisColor = lineColorMap(round(angles(k) + 90) + 1, :);
fprintf('Drawing line from (%.1f to %.1f) to (%.1f, %.1f)\n', x1, y1, x2, y2); line([x1,x2], [y1,y2], 'Color', thisColor);
if mod(k, 100) == 0
drawnow;
end
end
% Plot the histogram of angles.
subplot(2, 2, 4)
histObject = histogram(angles);
grid on;
xlabel('Angle', 'FontSize', 20);
ylabel('Count', 'FontSize', 20);
title('Histogram of angles', 'FontSize', 20);
Anna Marshall
Anna Marshall 2020-5-4
Thank you!! That works great. I'm going to play around with different filtersm particularly the ridgeline-finding filters and see if there is one that might work a bit better than the canny filter. Thanks again for all the suggestions!!
Anna Marshall
Anna Marshall 2020-5-14
@Image Analyst- I have another question to throw your way! I've been playing around with this code and one thing that I'm looking to try and do is add borders that group together "clusters" of the same colors aka same angles. For example, something like the attached sketch. I've tried the boundary function, but can't get it to work quite right. Do you have any ideas on what might work?

请先登录，再进行评论。

类别

在 帮助中心 和 File Exchange 中查找有关 Object Analysis 的更多信息

产品

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by