How to subtract images and get the largest blob from the subtracted image

I have few images from which I want to subtract the first image from the subsequent images. The first image is used as a reference point which would be subtracted from the other images. Please find the code below:

dirpath = uigetdir(pwd);
dirList = dir(dirpath);
dirList = {dirList.name};
find_regexp=regexp(dirList,'[0-9]F+[0-9]');
find_regex_indxs=find(~cellfun(@isempty,find_regexp));
dirList=dirList(find_regex_indxs);
dirList=natsort(dirList);

num_files = length(dirList);
imgs = cell(1,num_files);
for val = 1:num_files 
    imgs{val} = imread(dirList{val}); 
end
%%
final_segments = cell(1,num_files);
for idx = 1:num_files
    %% convert rgb to grayscale
    L = rgb2gray(imgs{idx});
%     figure, imshow(L)
    
    %% Thresholding multi-level
    thresh = multithresh(imgs{idx},2);
    seg_I = imquantize(imgs{idx},thresh);
%     figure, imshow(seg_I,[]);
    
    %% Removing the background
    seg_img = seg_I;
    seg_img(seg_I ~= 2) =  0;
%     figure, imshow(seg_img,[]);
    
    %% filtering using imclose
    struct_elem = strel('diamond',4);
    seg_img_close = imclose(seg_img,struct_elem);
%     figure, imshow(seg_img_close,[]);
    %% Get largest object

    bw_keep = bwareafilt(logical(seg_img_close),1);
%     figure, imshow(bw_keep,[]);
    %% fill holes within
    fill_bw_keep = imfill(bw_keep,'holes');
%     figure, imshow(fill_bw_keep,[]);
    
    %% Count nnz or area
    area_fill_bw_keep = nnz(fill_bw_keep);
    sprintf('area of segmented region: %d', area_fill_bw_keep)
    %% map the image back to original pixel values
    img_class = class(imgs{1});
    fill_bw_keep = cast(fill_bw_keep, img_class);
    map_img = fill_bw_keep.*imgs{idx};
%     figure, imshow(map_img,[]);
    
    %% save to cell
    final_segments{idx} = map_img;
end


%% disp final segments
for idx = 1:num_files
%     figure(1), imshow(final_segments{idx},[]);
    pause(1)
end

%% subtract 2 to n from 1st image and display
sub_final_segments = cell(1,num_files-1); % initialisation
for idx = 2:num_files
    sub_final_segments{idx-1} = final_segments{idx} - final_segments{1};
    figure(2), imshow(sub_final_segments{idx-1},[]);
    pause(1)
end

g = cell(1,num_files);
for idx = 2:num_files
    g = rgb2gray(sub_final_segments{idx-1});
%     figure, imshow(g,[]);
    pause(1)
end
%% subtraction areas
% areas_subtraction = [];
% for idx = 1:num_files-1
%     areas_subtraction = [areas_subtraction, nnz(sub_final_segments{idx})];
% end
% sprintf('Areas after subtraction: '), disp(areas_subtraction)

% Standard Deviation
subtraction_standard_deviation = [];
for idx = 1:num_files-1
    subtraction_standard_deviation = [subtraction_standard_deviation, std2(sub_final_segments{idx})];
end
sprintf('Standard Deviation after subtraction: '), disp(subtraction_standard_deviation)
    
g_deviation = [];
for idx = 1:num_files-1
    g_deviation = [g_deviation, std2(g(idx))];
end
sprintf('Standard Deviation after subtraction: '), disp(g_deviation)

%average density
subtraction_average_denstiy = [];
for idx = 1:num_files-1
    subtraction_average_denstiy = [subtraction_average_denstiy, mean2(sub_final_segments{idx})];
end
sprintf('Average Density after subtraction: '), disp(subtraction_average_denstiy)
    

g_denstiy = [];
for idx = 1:num_files-1
    g_denstiy = [g_denstiy, mean2(g(idx))];
end
sprintf('Average Density after subtraction: '), disp(g_denstiy)

However, I get the following error message:

Error in warid_codes (line 37)

bw_keep = bwareafilt(logical(seg_img_close),1);

I am using the bwareafilt function to take the largest blob from the subtracted image. But I am getting confused about the error message. I have attached few images. Any suggestions would be very much appreciated. Thank you.

[Matt J]

You have not shown us the complete error message.

[Image Analyst]

Is seg_img_close already a logical image?

[Warid Islam]

Hi @Matt,

I am sorry that I have missed the full error message. Please find it below:

Error using bwpropfilt

Expected input number 1, BW, to be two-dimensional.

Error in bwpropfilt>parse_inputs (line 121)

validateattributes(bw, {'logical'}, {'nonsparse', '2d'}, mfilename, 'BW', 1);

Error in bwpropfilt (line 58)

[bw, I, attrib, p, direction, conn] = parse_inputs(args{:});

Error in bwareafilt (line 34)

bw2 = bwpropfilt(bw, 'area', p, direction, conn);

Error in warid_codes (line 37)

bw_keep = bwareafilt(logical(seg_img_close),1);

Thank you.

[Warid Islam]

Hi @Image Analyst,

Yes, seg_img_close is already a logical image.

Answer 1

[Sid Singh]

Hi, bwareafilt expects a 2-D logical image. In your code, you are using RGB (3-D) images for preprocessing and the input to bwareafilt is a 3-D logical image. That's why you are getting the error. Check your preprocessing steps.

    %% convert rgb to grayscale
    L = rgb2gray(imgs{idx});
%     figure, imshow(L)
    
    %% Thresholding multi-level
    thresh = multithresh(imgs{idx},2);
    seg_I = imquantize(imgs{idx},thresh);
%     figure, imshow(seg_I,[]);
    

[Warid Islam]

Hi Sid,

I changed couple of lines of code above. But I am getting the same error.

   %% convert rgb to grayscale
%     L = rgb2gray(imgs{idx});
%     figure, imshow(L)
    
    %% Thresholding multi-level
    thresh = multithresh(imgs{idx},2);
    seg_I = imquantize(imgs{idx},thresh);
%     figure, imshow(seg_I,[]);
    
    %% Removing the background
    seg_img = seg_I;
    seg_img(seg_I ~= 2) =  0;
%     figure, imshow(seg_img,[]);
    
    %% filtering using imclose
    struct_elem = strel('diamond',4);
    seg_img_close = imclose(seg_img,struct_elem);
%     figure, imshow(seg_img_close,[]);
    %% Get largest object

    bw_keep = bwareafilt(logical(seg_img_close),1);
%     figure, imshow(bw_keep,[]);

Previously, I converted the image from rgb to gray but my image here is not rgb. So I commented out that part of the code. But the error still persists.