% Detection of parasites by means of granulometry and regional maxima.
+function [THS, c] = DetectionOfParasites(imgRGB, gt)
-%% Parameters
-clear;
+ %% Extracting of H, S and V components
+ imgHSV = rgb2hsv(imgRGB);
-imageFolder = '1305121398';
-imageNumber = '0001';
-scaleFactor = 0.6;
+ imgH = imgHSV(:, :, 1) * 255;
+ % We shift the hue to have the greatest value for the nuclei.
+ hueShiftValue = 195; % 100
+ imgH = uint8(mod(255 - imgH + hueShiftValue, 256));
-% Load the image and its ground truth.
-[imgRGB, gt] = loadImg(imageFolder, imageNumber);
+ imgS = uint8(imgHSV(:, :, 2) * 255);
+ % imgV = uint8(imgHSV(:, :, 3) * 255); % The value component isn't used.
+ % imgH = 255 - imgV; % We use the value component instead of the hue one (just for testing)
-%% Extracting of H, S and V components
+ % p. 136: median filter to smooth the noise and
+ % area closing to enhance the bright objects and make flatter, darker and cleaner the image background.
+ imgFiltered{1} = mmareaclose(medfilt2(imgH, [5, 5]), 400); % Hue.
+ imgFiltered{2} = mmareaclose(medfilt2(imgS, [5, 5]), 400); % Saturation.
-% Resample the image and its ground truth.
-imgRGBResampled = imresize(imgRGB, scaleFactor);
-gtResampled = imresize(gt, scaleFactor);
+ % Shading correction with a Top-hat transformation (p. 673).
+ % imgHFiltered = mmopenth(imgHFiltered, mmsedisk(80, '2D', 'OCTAGON'));
+ % imgSFiltered = mmopenth(imgSFiltered, mmsedisk(80, '2D', 'OCTAGON'));
-imgHSV = rgb2hsv(imgRGBResampled);
-imgH = imgHSV(:, :, 1) * 255;
-% We shift the hue to have the greatest value for the nuclei.
-hueShiftValue = 195; % 100
-imgH = uint8(mod(255 - imgH + hueShiftValue, 256));
+ %% Granulometry
+ % We use the saturation component to find the red cells mean size.
-imgS = uint8(imgHSV(:, :, 2) * 255);
+ redCellMaxSize = 42; % Radius [px].
+ funVolume = @(m) sum(sum(m));
-% imgV = uint8(imgHSV(:, :, 3) * 255); % The value component isn't used.
-% imgH = 255 - imgV; % We use the value component instead of the hue one (just for testing)
+ volImg = funVolume(imgFiltered{2});
-% p. 136: median filter to smooth the noise and
-% area closing to enhance the bright objects and make flatter, darker and cleaner the image background.
-imgFiltered{1} = mmareaclose(medfilt2(imgH, [5, 5]), 400); % Hue.
-imgFiltered{2} = mmareaclose(medfilt2(imgS, [5, 5]), 400); % Saturation.
+ sizeDistribution = zeros(redCellMaxSize, 1);
+ parfor k = 1:redCellMaxSize
+ SE = mmsedisk(k, '2D', 'OCTAGON');
+ imgOpened = mmopen(imgFiltered{2}, SE);
+ A = funVolume(imgOpened);
+ N = 1 - A / volImg;
+ sizeDistribution(k) = N;
+ end
-% Shading correction with a Top-hat transformation (p. 673).
-% imgHFiltered = mmopenth(imgHFiltered, mmsedisk(80, '2D', 'OCTAGON'));
-% imgSFiltered = mmopenth(imgSFiltered, mmsedisk(80, '2D', 'OCTAGON'));
+ patternSpectrum = zeros(redCellMaxSize - 1, 1);
+ for k = 1:redCellMaxSize - 1
+ patternSpectrum(k) = abs(sizeDistribution(k + 1) - sizeDistribution(k));
+ end
+ % The paper choose the biggest red cell size among the possible red cell sizes. (FIXME)
+ [m, c] = max(patternSpectrum);
+ c = c + 3;
+ nucleiRadius = 5; % Find a way to extract this information from the pattern spectrum histogram. (FIXME)
-%% Granulometry
-% We use the saturation component to find the red cells mean size.
+ bar(patternSpectrum);
-redCellMaxSize = 42; % Radius [px].
-funVolume = @(m) sum(sum(m));
-volImg = funVolume(imgFiltered{2});
+ %% Regional extrema
+ % The size of the SE (c) must be equal to the size of red cells (see granulometry below).
-sizeDistribution = zeros(redCellMaxSize, 1);
-parfor k = 1:redCellMaxSize
- SE = mmsedisk(k, '2D', 'OCTAGON');
- imgOpened = mmopen(imgFiltered{2}, SE);
- A = funVolume(imgOpened);
- N = 1 - A / volImg;
- sizeDistribution(k) = N;
-end
-
-patternSpectrum = zeros(redCellMaxSize - 1, 1);
-for k = 1:redCellMaxSize - 1
- patternSpectrum(k) = abs(sizeDistribution(k + 1) - sizeDistribution(k));
-end
-
-% The paper choose the biggest red cell size among the possible red cell sizes. (FIXME)
-[m, c] = max(patternSpectrum);
-c = c + 3;
-nucleiRadius = 5; % Find a way to extract this information from the pattern spectrum histogram. (FIXME)
-
-bar(patternSpectrum);
+ regionalMaxSE = mmsedisk(c, '2D', 'OCTAGON');
+ parfor i = 1:2
+ % This operation is very time consuming!
+ M{i} = mmregmax(imgFiltered{i}, regionalMaxSE);
+ end
-%% Regional extrema
-% The size of the SE (c) must be equal to the size of red cells (see granulometry below).
-
-regionalMaxSE = mmsedisk(c, '2D', 'OCTAGON');
-
-parfor i = 1:2
- % This operation is very time consuming!
- M{i} = mmregmax(imgFiltered{i}, regionalMaxSE);
-end
-
-nucleiSE = mmsedisk(nucleiRadius, '2D', 'OCTAGON');
-MHS = mmdil(M{1}, nucleiSE) & mmdil(M{2}, nucleiSE);
+ nucleiSE = mmsedisk(nucleiRadius, '2D', 'OCTAGON');
+ MHS = mmdil(M{1}, nucleiSE) & mmdil(M{2}, nucleiSE);
-volumeMHS = funVolume(~MHS);
-muH = funVolume(uint8(~MHS) .* imgFiltered{1}) / volumeMHS;
-muS = funVolume(uint8(~MHS) .* imgFiltered{2}) / volumeMHS;
+ volumeMHS = funVolume(~MHS);
+ muH = funVolume(uint8(~MHS) .* imgFiltered{1}) / volumeMHS;
+ muS = funVolume(uint8(~MHS) .* imgFiltered{2}) / volumeMHS;
-TH = im2bw(imgFiltered{1}, double(muH) / 255);
-TS = im2bw(imgFiltered{2}, double(muS) / 255);
-THS = TH & TS;
+ TH = im2bw(imgFiltered{1}, double(muH) / 255);
+ TS = im2bw(imgFiltered{2}, double(muS) / 255);
+
+ THS = TH & TS; % Parasites of all type and WBC.
-%% Output
+ %% Output
-mkdir('../output') % Just in case the 'output' directory doesn't exist.
+ mkdir('../output') % Just in case the 'output' directory doesn't exist.
-imwrite(imgRGBResampled, '../output/imgRGBResampled.png')
-imwrite(gtResampled, '../output/gtResampled.png')
+ imwrite(imgRGB, '../output/imgRGB.png')
+ imwrite(gt, '../output/gt.png')
-imwrite(imgH, '../output/imgH.png')
-imwrite(imgS, '../output/imgS.png')
+ imwrite(imgH, '../output/imgH.png')
+ imwrite(imgS, '../output/imgS.png')
-imwrite(M{1}, '../output/MH.png')
-imwrite(imgFiltered{1}, '../output/imgHFiltered.png');
+ imwrite(M{1}, '../output/MH.png')
+ imwrite(imgFiltered{1}, '../output/imgHFiltered.png');
-imwrite(M{2}, '../output/MS.png')
-imwrite(imgFiltered{2}, '../output/imgSFiltered.png');
+ imwrite(M{2}, '../output/MS.png')
+ imwrite(imgFiltered{2}, '../output/imgSFiltered.png');
-imwrite(MHS, '../output/MHS.png')
+ imwrite(MHS, '../output/MHS.png')
-imwrite(THS, '../output/THS.png')
+ imwrite(THS, '../output/THS.png')
-%% Display
+ %% Display
-figure('Position', [100 100 1600 800])
-colormap(gray);
-
-subplot(2, 4, 1);
-imagesc(imgRGBResampled);
-title(['Original: ', imageFolder, '/', imageNumber]);
-
-subplot(2, 4, 2);
-imagesc(imgFiltered{1});
-title('Hue component');
-
-subplot(2, 4, 3);
-imagesc(imgFiltered{2});
-title('Saturation component');
-
-subplot(2, 4, 4);
-imagesc(M{1});
-title('MH');
-
-subplot(2, 4, 5);
-imagesc(M{2});
-title('MS');
-
-subplot(2, 4, 6);
-imagesc(MHS);
-title('MHS');
-
-subplot(2, 4, 7);
-imagesc(THS);
-title('THS');
+% figure('Position', [100 100 1600 800])
+% colormap(gray);
+%
+% subplot(2, 4, 1);
+% imagesc(imgRGB);
+% title(['Original: ', imageFolder, '/', imageNumber]);
+%
+% subplot(2, 4, 2);
+% imagesc(imgFiltered{1});
+% title('Hue component');
+%
+% subplot(2, 4, 3);
+% imagesc(imgFiltered{2});
+% title('Saturation component');
+%
+% subplot(2, 4, 4);
+% imagesc(M{1});
+% title('MH');
+%
+% subplot(2, 4, 5);
+% imagesc(M{2});
+% title('MS');
+%
+% subplot(2, 4, 6);
+% imagesc(MHS);
+% title('MHS');
+%
+% subplot(2, 4, 7);
+% imagesc(THS);
+% title('THS');
+end
% print("-f1"\92, "-dpng", "Toto.png")