Blockproc not working as I expected

3 次查看（过去 30 天）
Jason
Jason 2019-1-4
编辑： Image Analyst 2019-1-9
Hello, I am having problems with blockproc.
I have an black image with white spots that I want to basic threshold locally rather than globally - that is to perform a basic threshold on tiles.
im127.png
IM=double(getimage(handles.axes4)); [rows,cols]=size(IM); %Get image from axes component
IM=IM-min(IM(:)); %Background subtract
% Number of tiles or chunks you want it divide up into.
xtiles =4; ytiles= 4;
%Numbers pixels per tile to maintain integer number of tiles
nx=floor(cols/xtiles); ny=floor(rows/ytiles);
Call my "LocallythresholdImage" function to define a level for binarising the image where
level = average + 2* standard deviations of the image
centroids=LocallythresholdImage(IM,nx,ny,f) %f=2 i.e. number of standard deviatins above the average
function centroids=LocallythresholdImage(IM,m,n,f)
%Uses the thresh function below and then creates a binary image using local
%threshold of size mxn. I
IM=double(IM); %make sure its a double
%Use blockproc uses my "thresh" function
myfun = @(block_struct) thresh3(block_struct.data,f); % myfun = @(block_struct) thresh2(block_struct.data,f);
%I2 = blockproc(Orig,[200 200],myfun,'BorderSize', [3 3], 'TrimBorder', true);
centroids = blockproc(IM,[m n],myfun,'BorderSize', [3 3], 'TrimBorder', true);
function [centroids]=thresh3(IM,f)
%n=str2num(get(handles.editThreshn,'String'));
%Thresholds image based on mean and S.D.
image=double(IM)-min(IM(:));
threshold=mean(image(:)) + f*std(image(:));
lessThanThreshold = image < threshold;
imageThresholded = image;
imageThresholded(lessThanThreshold) = 0;
regmax = imregionalmax(imageThresholded);
s = regionprops(regmax,double(image), 'WeightedCentroid');
format bank
centroids = cat(1, s.WeightedCentroid);
for each tile, the centroids are calculated correctly, but in my top level function "centroid=LocallythresholdImage(IM,nx,ny,f) ", centroids is empty
  1 个评论
Jason
Jason 2019-1-8
Maybe if I re-phrase my problem, someone may kindly comment.
When using blockproc, can the function fun that is applied to each block return an array of coordinates.
The aim is to plot the coordinates from each block on the original image.

请先登录，再进行评论。

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

采纳的回答

Guillaume
Guillaume 2019-1-8
Change your 'TrimBorder', true to 'TrimBorder', false. As it is you're removing a border of size 3 on each side of your output of thresh3.
Note that this thresh3 function is not well written. image should not be used as a variable name (it's a matlab function), format bank has no place in there (it only affects command line display) and the cat(1, ...) is completely wrong. I'm not sure what you're trying to do with that. There shouldn't be more than one centroid, but if there were, concatenating them all would result in an error in blockproc since the blocks would return arrays of different sizes.
  6 个评论
显示 4更早的评论
Guillaume
Guillaume 2019-1-9
When you specify 'BorderSize', matlab trims the border from the input and output unless you specify 'TrimBorder' as false. Since you were removing 6 rows from an output with only one row, you were always going to get an empty output.
Yes, you cannot use blockproc to return arrays of different sizes. Considering the discussion in IA's answer, what I'd recommend is that instead your tresh3 function return the image you were passing to regionprops. You can then do the regionprops afterward on all the blocks at once. If you forget about the local background subtraction:
function filteredblock = thresh3(IM, f)
threshold = mean(IM(:)) + f*std(IM(:));
IM(IM < threshold) = 0;
filteredblock = imregionalmax(IM);
end
function centroids = LocallythresholdImage(IM,m,n,f)
filteredimage = blockproc(double(IM), [m n], @(block_struct) thresh3(block_struct.data, f), 'BorderSize', [3 3], 'TrimBorder', false);
props = regionprops(filteredimage, IM, 'WeightedCentroid');
centroids = vertcat(props.WeightedCentroid);
end
If you do need that local background subtraction, then return a 3D output (first page: the regional max, 2nd page, the local background subtraction image:
function filteredblock = thresh3(IM, f)
IM = IM - min(IM(:));
threshold = mean(IM(:)) + f*std(IM(:));
thresholdimage = IM;
thresholdimage(IM < threshold) = 0;
filteredblock = cat(3, imregionalmax(thresholdimage), IM);
end
function centroids = LocallythresholdImage(IM, m, n, f)
filteredimage = blockproc(double(IM), [m n], @(bs) thresh3(bs.data, f), 'BorderSize', [3 3], 'TrimBorder', false);
props = regionpros(filteredimage(:, :, 1), filteredimage(:, :, 2), 'WeightedCentroid');
centroids = vertcat(props.WeightedCentroid);
end
Jason
Jason 2019-1-9
Thankyou, thats a really nice ay to get round the different size array issue - just return the binary image itself and then perform the threshold on that!
Jason

请先登录，再进行评论。

更多回答（1 个）

Image Analyst
Image Analyst 2019-1-9
"The aim is to plot the coordinates from each block on the original image."
Then don't complicate it. You can use a global threshold.
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 = 18;
folder = pwd; % Current folder.
baseFileName = 'index.png';
% % Have user browse for a file, from a specified "starting folder."
% % For convenience in browsing, set a starting folder from which to browse.
% startingFolder = pwd; % or 'C:\wherever';
% if ~exist(startingFolder, 'dir')
% % If that folder doesn't exist, just start in the current folder.
% startingFolder = pwd;
% end
% % Get the name of the file that the user wants to use.
% defaultFileName = fullfile(startingFolder, '01.png');
% [baseFileName, folder] = uigetfile(defaultFileName, 'Select a file');
% if baseFileName == 0
% % User clicked the Cancel button.
% return;
% end
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);
% Get 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)
if numberOfColorChannels > 1
% It's not really gray scale like we expected - it's color.
% Use weighted sum of ALL channels to create a gray scale image.
% grayImage = rgb2gray(rgbImage);
% ALTERNATE METHOD: Convert it to gray scale by taking only the green channel,
% which in a typical snapshot will be the least noisy channel.
grayImage = rgbImage(:, :, 3); % Take blue channel.
else
grayImage = rgbImage; % It's already gray scale.
end
% Now it's gray scale with range of 0 to 255.
% Display the image.
subplot(1, 2, 1);
imshow(grayImage, []);
title('Original Image', 'FontSize', fontSize, 'Interpreter', 'None');
axis('on', 'image');
hp = impixelinfo();
%------------------------------------------------------------------------------
% Set up figure properties:
% Enlarge figure to full screen.
set(gcf, 'Units', 'Normalized', 'OuterPosition', [0, 0.04, 1, 0.96]);
% Get rid of tool bar and pulldown menus that are along top of figure.
% set(gcf, 'Toolbar', 'none', 'Menu', 'none');
% Give a name to the title bar.
set(gcf, 'Name', 'Demo by ImageAnalyst', 'NumberTitle', 'Off')
drawnow;
% Threshold.
mask = grayImage == 255;
% Take the largest blobs only
mask = bwareaopen(mask, 100);
% Display the mask.
subplot(1, 2, 2);
imshow(mask);
hp = impixelinfo();
axis('on', 'image'); % Make sure image is not artificially stretched because of screen's aspect ratio.
title('Mask', 'FontSize', fontSize);
% Find the centroids
props = regionprops(mask, 'Centroid', 'Area');
allCentroids = vertcat(props.Centroid)
hold on;
plot(allCentroids(:, 1), allCentroids(:, 2), 'r+', 'MarkerSize', 30, 'LineWidth', 2);
0000 Screenshot.png
You could do a little better if you worked on the image without the number annotation on it so that you didn't get the 6 and 16.
  2 个评论
Jason
Jason 2019-1-9
编辑：Jason 2019-1-9
Walter / Image Analyst, thank you very much for you comments.
1: I love the bwareafilt function that Walter suggested, I had actually created my own area filtering function but using the built in function is more elegant. (its shown in the image below)
2: You both mention keep it simple and use Global thresholding. Well, I am building a GUI that allows for global thresholding & local thresholding (as well as other methods such as normalised cross correlation). The image I uploaded was a simple one where global thresholding works. However, I also need to analyse more complicated images as below which have spots much more closer together and some images have a non-uniform background, so I absolutely need a localised thresholding option.
(Sorry, there doesn't seem to be an option to shrink these images)
As you can see global thresholding is failing. Notice the bwareafilt together with imshowpair is a very nice way to show the thresholded objects that have been excluded via area classification (green blobs are excluded as too large)
I.A, in your example you nicely show bwareaopen, but often I need to go the other way and exclude the larger objects, I don't see a bwareaclose function to remove large objects from a binary image.
Also to add, the grid lines and number annotations are added after the processing and so don't affect the thresholding / binarisation.
So it still comes downto my problem. I want to break up my image into smaller tiles and do the thresholding and then reconstruct the centroids from each tile onto the main image. Im taking that this just isn't possible with blockproc which is a shame.
Thanks
Image Analyst
Image Analyst 2019-1-9
编辑：Image Analyst 2019-1-9
Jason, my code worked with your original example because you had well separated blobs that were all saturated (255), so global thresholding was sufficient. If your blobs are not saturated, then it wouldn't work.
There is no bwareaclose(), unless you wrote one. You can use bwareafilt() to explude larger blobs:
smallBlobs = bwareafilt(binaryImage, [1, maxAllowableBobSize]);
I'm attaching a program that cleans up a document by doing a local Otsu thresholding. It may be just what you want (at least the thresholding part of the program).
It looks like you could benefit by a higher resolution camera since some of those blobs are pretty blocky and only 3 pixels across.
There is an option to shrink your images. When you browse to it, turn the black twisty that says "display" downwards and select 25% or 50%, or type in something like 33%.

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 3-D Volumetric Image Processing 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by