Image Processing using ginput for Grain Size Analysis

The idea is the use a coin as a measurement reference in order to determine the size of the ash.

So the size of the coin isn't always the smallest or largest part of the image and I'm having trouble specifying this in my code. I know there is a way where I can use the function "ginput" so I can have the user iput/choose the reference point (the coin). Assuming this has to do with my props function (after line 85). This is important to me cause I want to also keep the measurement value (as you'll see in lines 102 and 103) so I can use that to measure other pieces of the image (which I'm still figuring out how to do). Open to suggestions for anyway to measure the ash as well.

Will post the code below and appreciate help as always. Also the image I'm using is above since I felt this was an easier one to use.

close all;
clear all;

%% Read in Image
RGB = imread('20210717_Tamarack_1524.jpg');
% Display Image
figure(1);
clf
imshow(RGB)

%% Convert Truecolor (RGB) Image into Grayscale Image
I = im2gray(RGB);
% Use Median Filter
I = medfilt2(I,[10 10]);
figure(2);
clf
imshow(I)

%% Calculate the Gradient Image and Apply a Threshold
% Use edge and Sobel operator to calculate the threshold value. Then  Tune the threshold value and use edge again to obtain a binary mask that contains the segmented cell
[~,threshold] = edge(I,'sobel');
fudgeFactor = 1;
BWs = edge(I,'sobel',threshold * fudgeFactor);
% Display resulting binary gradient mask
figure(3);
clf
imshow(BWs)

%% Dilate the Image.
% Create two perpendicular linear structuring elements by using strel function.
se90 = strel('line',3,90);
se0 = strel('line',3,0);
% Dilate the binary gradient mask using the vertical structuring element followed by the horizontal structuring element. The imdilate function dilates the image.
BWsdil = imdilate(BWs,ones(15,15));
imshow(BWsdil)
% Fill Interior Gaps
% Fill remaining holes using the imfill function
BWdfill = imfill(BWsdil,'holes');
figure(4);
clf
imshow(BWdfill)

%% Remove Connected Objects
% Remove any objects that are connected to the border of the image using
% imclearborder function. We use 4 as connectivity to remove 2D diagonal
% connections
BWnobord = imclearborder(BWdfill,4);
figure(5);
clf
imshow(BWnobord)

%% Smooth the Object
% We create the diamond structuring element using the strel function in
% order to make the object look natual/smooth
seD = strel('diamond',1);
BWfinal = imerode(BWnobord,seD);
BWfinal = imerode(BWfinal,seD);
figure(6);
clf
imshow(BWfinal)

%% Visualize the Segmentation
% Labelloverlay function allows us to display the mask over the original
% image
figure(7);
clf
imshow(labeloverlay(I,BWfinal))
%% bwareafilt function allows us to filter image, retaining (10) or n objects with
% the largest area
bwsizefilt = bwareafilt(BWfinal, 5);
% bwlabel returns the label matrix L that contains labels for the
% objects found in bwsizefilt
L = bwlabel(bwsizefilt);
figure(8);
clf
% Retains plot
hold on
% Creates pseudocolor plot using the values in matrix L
pcolor(L);
set(gca,'ydir','reverse')
% Sets color shading properties 
shading flat
% Measure properties of image regions using regionprops
% Calculate centroids and area for connected components in the image.
props = regionprops(bwsizefilt, 'Centroid', 'Area', 'EquivDiameter');
for i = 1:length(props)
    plot(props(i).Centroid(1),props(i).Centroid(2), '*k')
end

%% Find out the size of all the blobs so we know what size to use when filtering with bwareaopen.
% props = regionprops(mask, 'Area', 'EquivDiameter', 'Centroid');
allAreas = sort([props.Area]);
% Display the mask image.
% subplot(2, 1, 2);

% Coin radius, can be changed based on reference

% uiwait(msgbox('Choose Reference Point'));

% A = ginput(1);

pennyRadiusInPixels = props(3).EquivDiameter / 2;
pennyRadiusInMm = 19.05 / 2;

% Identify Points
%[x,y] = ginput(1)
%centers = [x,y];
%viscircles(centers,pennyRadiusInMm, 'Color', 'b')

% Overlay the circle over the image.

viscircles([props(3).Centroid(1), props(3).Centroid(2)], pennyRadiusInPixels, 'Color', 'b')
return
% Compute the spatial calibration factor.
mmPerPixel = pennyRadiusInMm / pennyRadiusInPixels
subplot(2, 1, 2);
caption = sprintf('The spatial calibration is %f mm per pixel', mmPerPixel)
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
% Plot the borders of all the blobs in the overlay above the original grayscale image 
% using the coordinates returned by bwboundaries().
% bwboundaries() returns a cell array, where each cell contains the row/column coordinates for an object in the image.
% Here is where we actually get the boundaries for each blob.
boundaries = bwboundaries(mask);
% boundaries is a cell array - one cell for each blob.
% In each cell is an N-by-2 list of coordinates in a (row, column) format.  Note: NOT (x,y).
% Column 1 is rows, or y.    Column 2 is columns, or x.
numberOfBoundaries = size(boundaries, 1); % Count the boundaries so we can use it in our for loop
% Here is where we actually plot the boundaries of each blob in the overlay.
subplot(2, 1, 1);
hold on; % Don't let boundaries blow away the displayed image.
for k = 1 : numberOfBoundaries
    thisBoundary = boundaries{k}; % Get boundary for this specific blob.
    x = thisBoundary(:,2); % Column 2 is the columns, which is x.
    y = thisBoundary(:,1); % Column 1 is the rows, which is y.
    plot(x, y, 'r-', 'LineWidth', 2); % Plot boundary in red.
end
hold off;
caption = sprintf('Original image with %d Outlines, from bwboundaries()', numberOfBoundaries); 
title(caption, 'FontSize', fontSize); 
axis('on', 'image'); % Make sure image is not artificially stretched because of screen's aspect ratio.
message = sprintf('Done!\nThe spatial calibration is %f mm per pixel', mmPerPixel)
uiwait(helpdlg(message));

Answer 1

[Image Analyst]

Looks like you made a common beginner mistake: doing an edge detection. Just because you can see edges in an image does not mean edge detection is the right thing to do. Beginners nearly always want to do an edge detection for some reason. Perhaps other inexperienced people told them to do that. It's usually not what is required.

First of all, I'd not use a coin. I'd get bright fluorescent spot stickers at the hobby or office supply store that are a completely different color than the soil. Or you could even use a Post-it note. Then find that color and do your spatial calibration from that.

Then segment your ash in some way that makes sense to differentiate it from the background. I don't know what that might be since I don't know what is ash and what is background in your images.

[Christopher Atkinson]

Oh I see. I think someone suggested edge detection as a first try.

Here is the main issue though (my apolies for not prefacing with this). I cannout use anything but the coin cause this image (and similar ones) were taken like a year ago. So unfortunately, bright stickers aren't an option. I guess that is why edge detection was my first method to try although all ash pieces aren't uniform, we could use a general shape.

Would love to spatial calibrate it (since my background is usually a road, concrete, etc.) but, since I don't have control over the background or reference (coin here) look like it'll be very complicated/tough for me to solve.

[Image Analyst]

Have you tried imfindcircles? It usually takes a bit of trial and error to get the input parameters correct.

[Christopher Atkinson]

I have but I feel like that's what I'm doing on line 124. Maybe it's just the beginner in me but, I believe all I need is a way to user input (to chose) what's saved on line 104 in order for line 124 to circle over the correct part of the image (the coin in this instance). Which is the reason why I was looking into ginput function. May be the long way around but just a thought I had.

[Image Analyst]

Well if you want you can use bwselect but if you're going the manual-assisted route, why not use drawcircle?

[Christopher Atkinson]

I’m going to give both a shot and aee which one I like more.

Also on a side, is there example code or discussion for finding the location of this image on google maps using Matlab? Like I want to cross reference the location of ash or background (if that makes sense). Figure it is complicated (as I would love to have the x and y coordinates, but not sure what function to use.

[Image Analyst]

@Christopher Atkinson see this:

https://www.mathworks.com/matlabcentral/fileexchange/42000-run_exiftool?s_tid=srchtitle

If you're lucky the camera put the GPS location into the metadata (EXIF or IPTC).

[Christopher Atkinson]

@Image Analyst I assume it is into metadata but is abit confusing as I ran the .exe file but it just gives random number when I drag this image to it.

[Image Analyst]

Attach the image with the paperclip icon.

[Christopher Atkinson]

Here is the image

[Walter Roberson]

That file has no useful EXIF information. Just the X Resolution and Y Resolution being the nominal "1" which tells you nothing useful.