imgH = uint8(mod(255 - imgH + hueShiftValue, 256));
imwrite(imgH, '../output/imgH.png')
- imgS = uint8(imgHSV(:, :, 2) * 255);
+ imgS = uint8(255 - 255 * imgHSV(:, :, 2));
imwrite(imgS, '../output/imgS.png')
- imgV = uint8(255 - imgHSV(:, :, 3) * 255);
+ imgV = uint8(255 - 255 * imgHSV(:, :, 3));
imwrite(imgV, '../output/imgV.png')
% 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.
medianFilterWindow = [5, 5];
- imgFiltered{1} = mmareaclose(medfilt2(imgH, medianFilterWindow), 800); % Hue.
- % We use an opening fot the saturation because the parasites are darker than the rest of the image.
- imgFiltered{2} = mmareaopen(medfilt2(imgS, medianFilterWindow), 800); % Saturation.
- imgFiltered{3} = mmareaclose(medfilt2(imgV, medianFilterWindow), 800); % Value. This component isn't used.
+ imgFiltered{1} = mmareaclose(medfilt2(imgH, medianFilterWindow), 200); % Hue.
+ imgFiltered{2} = mmareaclose(medfilt2(imgS, medianFilterWindow), 200); % Saturation.
+ imgFiltered{3} = mmareaclose(medfilt2(imgV, medianFilterWindow), 200); % Value. This component isn't used.
imwrite(imgFiltered{1}, '../output/imgHFiltered.png');
imwrite(imgFiltered{2}, '../output/imgSFiltered.png');
imwrite(imgFiltered{3}, '../output/imgVFiltered.png');
-
%% Granulometry
% We use the saturation component to find the red cells mean size.
saturationComponent = imgFiltered{2};
- redCellMaxSize = round(size(imgRGB, 2) / 35); % Radius [px], depending of the input size.
+ redCellMaxSize = uint32(size(imgRGB, 2) / 35); % Radius [px], depending of the input size.
funVolume = @(m) sum(sum(m));
volImg = funVolume(saturationComponent);
sizeDistribution = zeros(redCellMaxSize, 1);
parfor k = 1:redCellMaxSize
SE = mmsedisk(k, '2D', 'OCTAGON'); % 'EUCLIDEAN' is more precise.
- imgOpened = mmopen(saturationComponent, SE);
- A = funVolume(imgOpened);
+ imgClosed = mmclose(saturationComponent, SE); % Here we use a closing because the cells are darker than the background.
+ A = funVolume(imgClosed);
N = 1 - A / volImg;
sizeDistribution(k) = N;
end
- % TODO : voir diff
+ % TODO : voir diff.
patternSpectrum = zeros(redCellMaxSize - 1, 1);
for k = 1:redCellMaxSize - 1
patternSpectrum(k) = abs(sizeDistribution(k + 1) - sizeDistribution(k));
% The paper chooses the biggest red cell size among the possible red cell sizes. (FIXME)
[~, c] = max(patternSpectrum);
- nucleiRadius = c / 8; % We admit the size of a parasite is eight time smaller than a red blood cell. Find a way to extract this information from the pattern spectrum histogram. (FIXME)
+ nucleiRadius = c / 5; % We admit the size of a parasite is five time smaller than a red blood cell. Find a way to extract this information from the pattern spectrum histogram. (FIXME)
bar(patternSpectrum);
% The functions 'mmregmax' and 'mmregmin' are very time consuming!
parfor i = 1:3
- if i == 2
- M{i} = mmregmin(imgFiltered{i}, regionalExtremumSE); % Special case for the saturation component.
- else
- M{i} = mmregmax(imgFiltered{i}, regionalExtremumSE);
- end
+ M{i} = mmregmax(imgFiltered{i}, regionalExtremumSE);
end
nucleiSE = mmsedisk(nucleiRadius, '2D', 'OCTAGON');
muS = funVolume(uint8(MHS) .* imgFiltered{2}) / volumeMHS;
TH = im2bw(imgFiltered{1}, double(muH) / 255);
- TS = ~im2bw(imgFiltered{2}, double(muS) / 255);
+ TS = im2bw(imgFiltered{2}, double(muS) / 255);
THS = TH & TS;
imwrite(M{2}, '../output/MS.png')
imwrite(MHS, '../output/MHS.png')
imwrite(THS, '../output/THS.png')
-
end
% c: (TODO: should be the average red cell size)
function [schizontsSmoothed] = DetectionOfSchizonts(THS, c)
- schizontsConnectivityMin = 3; % Minimum cluster size when detecting schizonts.
+ schizontsConnectivityMin = 6; % Minimum cluster size when detecting schizonts.
redCellsSE = mmsedisk(c, '2D', 'OCTAGON');
schizonts = zeros(size(THS), 'uint8');
end
end
- %% Fill 'schizonts' by testing the connectivity of each cluster set.
+ %% Fill the variable 'schizonts' by testing the connectivity of each cluster set.
[~, C] = graphconncomp(adjacencyMatrix, 'Directed', false, 'Weak', true);
n = 1;
while 1
% Dilation + Morphological smoothing (Hole filling).
schizontsSmoothed = logical(mmclohole(mmdil(schizonts, mmsedisk(4))));
- imwrite(schizontsSmoothed, '../output/schizonts.png')
+ WriteImageGradient(schizontsSmoothed, '../output/schizonts.png', [255, 30, 30]);
end
\ No newline at end of file
% c: The typical red cell radius [px].
function [WBCSmoothed] = DetectionOfWhiteCells(THS, c)
% We use a radius near the smallest red blood cell size (80% of the radius)
- redCellsSE = mmsedisk(c * 0.8, '2D', 'OCTAGON');
+ redCellsSE = mmsedisk(c * 0.5, '2D', 'OCTAGON');
% The white cells may not be correctly marked in the case they have a very special shape (oblong).
WBCMarker = mmero(THS, redCellsSE); % Erosion to keep only fragments of white cells.
% Morphological smoothing (Hole filling).
WBCSmoothed = mmclohole(WBC);
- imwrite(WBCSmoothed, '../output/whiteCells.png')
+ WriteImageGradient(WBCSmoothed, '../output/whiteCells.png', [150, 255, 150]);
end
\ No newline at end of file
function [infectedRedCells] = IdentificationOfInfectedRedCells(redCells, parasites)
infectedRedCells = mminfrec(parasites, redCells);
- imwrite(infectedRedCells, '../output/infected red cells.png')
+
+
+ WriteImageGradient(infectedRedCells, '../output/infected red cells.png', [50, 50, 255]);
end
\ No newline at end of file
% imageFolder = '1405022890'; % 4.
% imageFolder = '1409191647'; % 9.
imageFolder = '1412151257'; % 10. (Gamma)
+% imageFolder = '13.05.2015';
%% Gamma
-% imageNumber = 'Gamma-0.8-0002';
-% imageNumber = 'Gamma-0.8-0005';
-% imageNumber = 'Gamma-0.8-0008';
+% imageName = '1412151257-Gamma-0.8-0002.png';
+% imageName = '1412151257-Gamma-0.8-0005.png';
+% imageName = '1412151257-Gamma-0.8-0008.png';
-% imageNumber = 'Gamma-1-0001';
-% imageNumber = 'Gamma-1-0004';
-imageNumber = 'Gamma-1-0007';
+% imageName = '1412151257-Gamma-1-0001.png';
+% imageName = '1412151257-Gamma-1-0004.png';
+imageName = '1412151257-Gamma-1-0007.png';
-% imageNumber = 'Gamma-1.2-0003';
-% imageNumber = 'Gamma-1.2-0006';
-% imageNumber = 'Gamma-1.2-0009';
+% imageName = 'Gamma-1.2-0003';
+% imageName = 'Gamma-1.2-0006';
+% imageName = 'Gamma-1.2-0009';
+
+%% 13 mai 2015
+
+% imageName = '1412151257-100x-Teinte30-Saturation3-0013.tif';
+% imageName = '1412151257-100x-Teinte30-Saturation3-0014.tif'; % Good result.
+% imageName = '1412151257-100x-Teinte0-Saturation0-0015.tif'; % Bad result.
+
+segmentationMethod = SegmentationMethod.WatershedByMorphologicalGradient;
%% Main
-fprintf('Folder: %s, Number: %s\n', imageFolder, imageNumber);
+fprintf('Folder: %s, Image name: %s\n', imageFolder, imageName);
-scaleFactor = 1;
+scaleFactor = 0.6;
% Load the image and its ground truth.
-[imgRGB, gt] = loadImg(imageFolder, imageNumber);
+imgRGB = loadImg(imageFolder, imageName);
% Resample the image and its ground truth.
imgRGBResampled = imresize(imgRGB, scaleFactor);
-gtResampled = imresize(gt, scaleFactor);
imwrite(imgRGBResampled, '../output/imgRGB.png')
-imwrite(gtResampled, '../output/gt.png')
disp('1) Detection of parasites ...')
[THS, c, nucleiRadius] = DetectionOfParasites(imgRGBResampled);
disp('2) Detection of white cells ...')
[WBC] = DetectionOfWhiteCells(THS, c);
+THSWithoutWBC = THS & ~WBC; % We remove the white cells from THS.
disp('3) Detection of schizonts ...')
-[schizonts] = DetectionOfSchizonts(THS, c);
-[imgRGBRedCells] = RemoveWBCAndSchizonts(imgRGBResampled, WBC, schizonts);
+[schizonts] = DetectionOfSchizonts(THSWithoutWBC, c);
+THS(schizonts) = 0; % We remove the white cells from THS.
+THSWithoutWBCandSchizonts = THSWithoutWBC & ~schizonts; % We remove the shizonts from THS.
disp('4) Segmentation of red cells ...')
-[redCells] = SegmentationOfRedCells(imgRGBRedCells, c);
+[redCells] = SegmentationOfRedCells(imgRGBResampled, c, WBC, schizonts, segmentationMethod);
disp('5) Finding infected red cells ...')
-parasites = THS & ~WBC & ~schizonts;
-[infectedRedCells] = IdentificationOfInfectedRedCells(redCells, parasites);
+[infectedRedCells] = IdentificationOfInfectedRedCells(redCells, THSWithoutWBCandSchizonts);
+disp('6) Counting red cells and infected red cells ...')
nbRedCells = mmstats(mmlabel(redCells), 'max');
nbInfectedRedCells = mmstats(mmlabel(infectedRedCells), 'max');
infectionPercentage = 100 * nbInfectedRedCells / nbRedCells;
-fprintf('Percentage of infected cell: %.1f%%\n', infectionPercentage)
+fprintf('Percentage of infected cell: %.1f%% (%d/%d)\n', infectionPercentage, nbInfectedRedCells, nbRedCells)
disp('Finished')
toc
\ No newline at end of file
-function [redCells] = SegmentationOfRedCells(imgRGB, c)
- cMin = uint32(c * 0.8); % 80%.
+function [segmentedCells] = SegmentationOfRedCells(imgRGB, c, WBC, schizonts, segmentationMethod)
+ cMin = uint32(c * 0.75); % 75%.
cMinArea = pi * cMin ^ 2;
- gray = mmareaopen(rgb2gray(imgRGB), 2 * cMinArea); % The gray image is better than the green component alone.
+ gray = imgRGB(:,:,2);
+ gray(WBC) = 255;
+ gray(schizonts) = 255;
imwrite(gray, '../output/red cells segmentation - gray.png')
- % Thresholding using the Otsu's method .
- threshold = ~im2bw(gray, graythresh(gray));
+ grayFlat = mmareaopen(gray, 2 * cMinArea, mmsebox(1)); % The gray image is better than the green component alone.
+ imwrite(grayFlat, '../output/red cells segmentation - gray flat.png')
+
+ % Thresholding using the Otsu's method.
+ % TODO : montrer l'histogramme
+ threshold = ~im2bw(grayFlat, graythresh(grayFlat));
% We remove all objects with a smaller area than a red blood cell.
threshold = mmareaopen(threshold, cMinArea);
% We set the background to 0 by using the previous threshold.
- redCellsFiltered = gray;
+ redCellsFiltered = grayFlat;
redCellsFiltered(~threshold) = 0;
imwrite(redCellsFiltered, '../output/red cells segmentation - input.png')
%% Segmentation
- opened = mmopen(redCellsFiltered, mmsedisk(cMin, '2D', 'EUCLIDEAN', 'NON-FLAT', 0));
+
+ % To smooth the borders and to remove little elements.
+ opened = mmopen(redCellsFiltered, mmsedisk(cMin, '2D', 'EUCLIDEAN', 'NON-FLAT'));
imwrite(opened, '../output/red cells segmentation - opened.png')
- openedGradient = mmgradm(opened);
- imwrite(openedGradient, '../output/red cells segmentation - opened - gradient.png')
-
- openedGradientThreshold = im2bw(openedGradient, 0.02);
- imwrite(openedGradientThreshold, '../output/red cells segmentation - opened - gradient - threshold.png')
-
- holesClosed = mmareaclose(openedGradientThreshold, 10 * cMinArea);
- imwrite(holesClosed, '../output/red cells segmentation - holes closed.png')
-
- % Watershed algorithm from here: http://mmorph.com/mxmorph/html/mmdemos/mmdcells.html.
- d = holesClosed;
- e1 = mmdist(d, mmsebox,'EUCLIDEAN');
- e3 = mmregmax(e1);
- e = mmdil(e3, mmsedisk(5, '2D', 'EUCLIDEAN'));
- f = mmneg(e1);
- g = mmcwatershed(f, e, mmsebox);
- h = mmintersec(d, mmneg(g));
- i = mmedgeoff(h);
- j = mmareaopen(i, cMinArea / 3);
- imwrite(j, '../output/red cells segmentation - individual.png')
-
- redCells = j;
+ openedThreshold = im2bw(opened, 0.02);
+ imwrite(openedThreshold, '../output/red cells segmentation - opened - threshold.png')
+
+ switch segmentationMethod
+ case SegmentationMethod.WatershedByDistanceTransform
+ %% Simplified version.
+
+ % Remove cells which are touching the edges.
+ segmentedCells = mmedgeoff(WatershedsByDistanceTransform(openedThreshold));
+
+ case SegmentationMethod.WatershedByMorphologicalGradient
+
+ %% Don't know the meaning. (end of the right column of page 139).
+% openedGradient = mmgradm(opened);
+% imwrite(openedGradient, '../output/red cells segmentation - opened - gradient.png')
+
+% openedGradientThreshold = im2bw(openedGradient, 0.02);
+% imwrite(openedGradientThreshold, '../output/red cells segmentation - opened - gradient - threshold.png')
+
+% holesClosed = mmareaclose(openedGradientThreshold, 10 * cMinArea);
+% imwrite(holesClosed, '../output/red cells segmentation - holes closed.png')
+
+ % Watershed by distance transform.
+ watershed = mmedgeoff(WatershedsByDistanceTransform(openedThreshold));
+ imwrite(watershed, '../output/red cells segmentation - watershed.png')
+
+ watershedLabeled = mmlabel(watershed);
+ axisProperties = regionprops(watershedLabeled, 'MajorAxisLength', 'MinorAxisLength');
+
+ % The area opening remove the little artifacts.
+ watershedSingles = mmareaopen(watershed, cMinArea / 5);
+
+ for k = 1:length(axisProperties)
+ % We remove each composite cell from 'watershedSingle'.
+ if axisProperties(k).MajorAxisLength / axisProperties(k).MinorAxisLength >= 1.3
+ watershedSingles(watershedLabeled == k) = 0;
+ end
+ end
+ imwrite(watershedSingles, '../output/red cells segmentation - watershed singles.png')
+
+ % We remove the single cells from the initial image to obtain the composite cells.
+ composites = mmedgeoff(and(openedThreshold, ~watershedSingles));
+ imwrite(composites, '../output/red cells segmentation - composites.png')
+
+ compositesOpened = grayFlat;
+ compositesOpened(~composites) = 0;
+ compositesOpened = mmopen(compositesOpened, mmsedisk(cMin, '2D', 'EUCLIDEAN', 'FLAT'));
+
+ imwrite(compositesOpened, '../output/red cells segmentation - composites opened.png')
+
+ compositesOpenedGradient = mmgradm(compositesOpened);
+ imwrite(compositesOpenedGradient, '../output/red cells segmentaion - composites opened gradient.png');
+
+ compositesOpenedGradientThreshold = im2bw(compositesOpenedGradient, 0.00001);
+ imwrite(compositesOpenedGradientThreshold, '../output/red cells segmentation - composites opened gradient threshold.png')
+
+ compositesOpenedGradientThresholdHolesFilled = ~mmareaopen(~compositesOpenedGradientThreshold, cMinArea / 2);
+ imwrite(compositesOpenedGradientThresholdHolesFilled, '../output/red cells segmentation - composites opened gradient threshold holes filled.png')
+
+ compositesThinning = mmthin(compositesOpenedGradientThresholdHolesFilled);
+ compositesThinningFlooded = mmareaclose(compositesThinning, 5 * cMinArea);
+ compositesThinningFlooded(compositesThinning) = 0;
+ compositesThinningFlooded = mmero(compositesThinningFlooded); % We erode by a 3x3 cross to avoid the composites to touch the singles when fusionning later.
+ imwrite(compositesThinningFlooded, '../output/red cells segmentation - composites thinning.png')
+
+ % Compute the ratio of the major axis over the minor axis and keep
+ % cells with a ration greater than or equal to 1.3.
+ segmentedCells = or(compositesThinningFlooded, watershedSingles);
+ end
+
+ WriteImageGradient(segmentedCells, '../output/red cells segmentation - individual.png', [255, 255, 255]);
end
\ No newline at end of file
%% Data
% RGB = imread ('/Users/michel.kocher/michel/Data/Movable/imgs_corrMK/Gamma/1412151257-Gamma-0.8-0002.png');
-RGB = imread('../imgs_corrMK/Gamma/1412151257-Gamma-0.8-0002.png');
+RGB = imread('../imgs_corrMK/1412151257/1412151257-Gamma-0.8-0002.png');
%% HSV
HSV = rgb2hsv (RGB);
linkaxes (ax)
%% Print 2 file
-% print ('-f1', '-dpng', 'THue.png')
\ No newline at end of file
+print ('-f1', '-dpng', 'THue.png')
\ No newline at end of file
-% Load the image and the ground-truth from the given foldername and suffix.
-% If the ground-truth isn't found then it is put to white.
-function [img, gt] = loadImg(name, num)
-
+% Load the image from the given foldername and filename.
+function img = loadImg(dirname, filename)
imagesFolder = '../imgs_corrMK';
- groundTruthFolder = '../gt';
- filename = [name, '-', num, '.png'];
- filepath_img = [imagesFolder, '/', name, '/', filename];
- filepath_gt = [groundTruthFolder, '/', name, '/', filename];
-
- img = imread(filepath_img);
-
- try
- gt = imread(filepath_gt);
- catch exception
- fprintf('Can''t load the ground-truth : %s\n', exception.message)
- gt = ones(size(img));
- end
+ filepath_img = [imagesFolder, '/', dirname, '/', filename];
+ img = imread(filepath_img);
end
projects / malaria.git / commitdiff