Finding the real world (in^2) area of black smears on graph paper.

I have black smears on graph paper that I would like to know the area using real world in^2 units. Each square on the paper is 1/4in. x 1/4in. I ultimately would like to be able to batch process up to 6 of these images at a time where the output I would like is area of the black smears.

Any help would be appreciated.

[William Rose]

@Christopher,

Please share your inital attempt to solve this problem. Have you found the area in units of pixels? That would be a good start.

Answer 1

[Image Analyst]

See my Image Segmentation Tutorial in my File Exchange:

Image Analyst's File Exchange

It's a generic, general purpose demo of how to threshold an image to find blobs, and then measure things about the blobs, and extract certain blobs based on their areas or diameters.

But first let's solve the image capture problem. You have bad lighting uniformity and lens shading. It would be much better if you just scanned in the images on a flatbed scanner.

Now you can threshold and call regionprops.

grayImage = rgbImage(:, :, 3); % Take blue channel because the blue grid will be brightest in that channel.
mask = grayImage < 50; % Whatever value works best.
props = regionprops(mask, 'Area');
allAreasInPixels = [props.Area]

If thresholding still gives you the blue lines, see if you can use the Color Thresholding app on the Apps tab of the tool ribbon to get just the black colored dark blobs and not the blue lines. If you still have blue lines you can use imopen on the binary image mask to get rid of them.

Finally to spatially calibrate the areas you can multiply by the number of mm per pixel. You can measure this once and it will be constant for all the images you take with the scanner.

allAreasInMm = allAreasInPixels * mmPerPixel ^ 2

See my attached spatial calibration demo if you need more help on how to do spatial calibration.

To process a sequence of files, use a for loop. Code samples are in the FAQ:

https://matlab.fandom.com/wiki/FAQ#How_can_I_process_a_sequence_of_files?

Answer 2

[Image Analyst]

Here is a crude first cut at it. It's not very robust because you haven't used the suggestions in my other answer to improve the images. And I have tested only on one image and not put it into a loop over all images. Nor have I done the conversion into square millimeters. But it's a start. You might want to call histogram(allAreasInPixels) to see the distribution of blob areas.

% Demo by Image Analyst
% Initialization steps:
clc;    % Clear the command window.
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;
fontSize = 16;
markerSize = 20;

%--------------------------------------------------------------------------------------------------------
%    READ IN TEST IMAGE
folder = [];
baseFileName = 'CRS-1H 90F 10min.jpg';
fullFileName = fullfile(folder, baseFileName);
% Check if file exists.
if ~exist(fullFileName, 'file')
    % The file doesn't exist -- didn't find it there in that folder.
    % Check the entire search path (other folders) for the file by stripping off the folder.
    fullFileNameOnSearchPath = baseFileName; % No path this time.
    if ~exist(fullFileNameOnSearchPath, 'file')
        % Still didn't find it.  Alert user.
        errorMessage = sprintf('Error: %s does not exist in the search path folders.', fullFileName);
        uiwait(warndlg(errorMessage));
        return;
    end
end
rgbImage = imread(fullFileName);

%--------------------------------------------------------------------------------------------------------
% Display the image.
subplot(2, 2, 1);
imshow(rgbImage, []);
impixelinfo;
axis('on', 'image');
title('Original RGB Image', 'FontSize', fontSize, 'Interpreter', 'None');
% Update the dimensions of the image.
% numberOfColorChannels should be = 1 for a gray scale image, and 3 for an RGB color image.
[rows, columns, numberOfColorChannels] = size(rgbImage)
% Maximize window.
g = gcf;
g.WindowState = 'maximized';
g.Name = 'Demo by Image Analyst';
g.NumberTitle = 'off';
drawnow;

%--------------------------------------------------------------------------------------------
% MASK THE GRIDLINES IN THE IMAGE.  SET THEM TO WHITE.
[gridLineMask, maskedRGBImage] = createMaskOfGridLines(rgbImage);
[r,g,blueChannel] = imsplit(rgbImage);
r(gridLineMask) = 255;
g(gridLineMask) = 255;
blueChannel(gridLineMask) = 255;
% Mask the blue channel.
% Interactively and visually set a threshold on a gray scale image.
% https://www.mathworks.com/matlabcentral/fileexchange/29372-thresholding-an-image?s_tid=srchtitle
lowThreshold = 0;
highThreshold = 83;
% [lowThreshold, highThreshold] = threshold(lowThreshold, highThreshold, blueChannel);
mask = blueChannel >= lowThreshold & blueChannel <= highThreshold;
% mask = imfill(mask, 'holes');
% mask = imclearborder(mask);

%--------------------------------------------------------------------------------------------
% Assume paper is in a fixed position and field of view is fixed.
% Blacken the image outside the paper.
paperMask = false(rows, columns);
paperMask(237:2051, 464:2873) = true;
mask = mask & paperMask;
% Show the mask image.
subplot(2, 2, 2);
imshow(mask)
impixelinfo;
axis('on', 'image');
caption = sprintf('Final Mask Image');
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
drawnow;

%--------------------------------------------------------------------------------------------
% Show the masked and cropped image.
% Label each blob with 8-connectivity, so we can make measurements of it
[labeledImage, numberOfBlobs] = bwlabel(mask, 8);
% Apply a variety of pseudo-colors to the regions.
coloredLabelsImage = label2rgb (labeledImage, 'hsv', 'k', 'shuffle');
% Display the pseudo-colored image.
subplot(2, 2, 3:4);
imshow(coloredLabelsImage);
impixelinfo;
axis('on', 'image');
caption = sprintf('Final Image');
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
drawnow;

%--------------------------------------------------------------------------------------------
props = regionprops(mask, 'Area');
allAreasInPixels = [props.Area]


%===============================================================================================
function [BW,maskedRGBImage] = createMaskOfGridLines(RGB)
%createMask  Threshold RGB image using auto-generated code from colorThresholder app.
%  [BW,MASKEDRGBIMAGE] = createMask(RGB) thresholds image RGB using
%  auto-generated code from the colorThresholder app. The colorspace and
%  range for each channel of the colorspace were set within the app. The
%  segmentation mask is returned in BW, and a composite of the mask and
%  original RGB images is returned in maskedRGBImage.
% Auto-generated by colorThresholder app on 29-Oct-2024
%------------------------------------------------------
% Convert RGB image to chosen color space
I = rgb2hsv(RGB);
% Define thresholds for channel 1 based on histogram settings
channel1Min = 0.428;
channel1Max = 0.796;
% Define thresholds for channel 2 based on histogram settings
channel2Min = 0.106;
channel2Max = 1.000;
% Define thresholds for channel 3 based on histogram settings
channel3Min = 0.000;
channel3Max = 1.000;
% Create mask based on chosen histogram thresholds
sliderBW = (I(:,:,1) >= channel1Min ) & (I(:,:,1) <= channel1Max) & ...
    (I(:,:,2) >= channel2Min ) & (I(:,:,2) <= channel2Max) & ...
    (I(:,:,3) >= channel3Min ) & (I(:,:,3) <= channel3Max);
BW = sliderBW;
% Initialize output masked image based on input image.
maskedRGBImage = RGB;
% Set background pixels where BW is false to zero.
maskedRGBImage(repmat(~BW,[1 1 3])) = 0;
end