Homec4science
Diffusion ShoulderDB_ImportSegmentMeasureAnalyze (aterrier-plotUpdate) Browse

R8218/muscleDegeneration.m
30aa777b9d8baterrier-plotUpdate

Last Change
Raw File

muscleDegeneration.m
View Options

function muscle = muscleDegeneration(caseID,muscleName,HUThresholdMuscle,HUThresholdOsteochondroma,muscleDirectory)
%%MUSCLEDEGENERATION computes degeneration parameters from muscle contours
%
% This function returns muscle atrophy, fat infiltration, osteochondroma
% and degeneration values for a given muscle
%
% USE: muscleDegeneration(caseID,muscleName,HUThresholdMuscle,HUThresholdOsteochondroma,output)
%
% INPUT caseID: patient number
% muscleName: name of the rotator cuff muscle being measured (SS,IS, SC or TM)
% HUThresholdMuscle: HU threshold used for muscle segmentation (muscle-fat limit)
% HUThresholdMuscle:HU threshold used for osteochondroma segmentation (muscle-osteochondroma limit)
% muscleDirectory: Path to 'muscles' directory where contours and output images are stored
%
% OUTPUT muscle: structure containing fields 'atrophy', 'infiltration',
% 'osteochondroma' and 'degeneration'
%
% REMARKS An image of the segmentation of atrophied muscle, fat infiltration
% and osteochondroma is saved in muscles directory, in the
% 'muscleName' folder.
%
% created with MATLAB ver.: 9.2.0.556344 (R2017a) on Windows 7
% Author: EPFL-LBO-VMC
% Date: 07-Jun-2017
%
inputFile = ['ssContours' caseID '.mat'];
contouredMuscle = load([muscleDirectory '\' muscleName '\' inputFile]);
% --Atrophied muscle segmentation--
% MATLAB does not accept to threshold an image with a negative value, as is
% the case when using HU scale. CT scanners usually have a CT range from
% -1000 HU to 1000 HU or -2000 HU to 2000 HU for modern scanners. In this
% case, we're not using HU values lower than -1000 and higher than +1000 so
% all pixel intensity values of the input image< -1000 HU are set to -1000
% HU, and values > 1000 HU are set to 1000 HU. Then, the input image is
% mapped to the [0 1] range -1000->0 1000->1, as are the thresholds for
% muscle and osteochondroma.
% Normalize image and thresholds
contouredMuscle.inputImage(contouredMuscle.inputImage < -1000) = -1000;
contouredMuscle.inputImage(contouredMuscle.inputImage > 1000) = 1000;
rangeHU = double([-1000 1000]);
normalizedThresholdMuscle = (HUThresholdMuscle-rangeHU(1))/(rangeHU(2)-rangeHU(1));
normalizedThresholdOsteochondroma = (HUThresholdOsteochondroma-rangeHU(1))/(rangeHU(2)-rangeHU(1));
normalizedInputImage = (double(contouredMuscle.inputImage)-rangeHU(1))/(rangeHU(2)-rangeHU(1));
% Mask image using muscle fossa contours
maskedImage = normalizedInputImage;
maskedImage(~contouredMuscle.binaryImage1) = 0;
% Segment using muscle threshold
myMuscle = imbinarize(maskedImage,normalizedThresholdMuscle);
myMuscle = imfill(myMuscle, 'holes'); % Fill holes
% --Atrophied muscle without fat infiltration segmentation--
% Mask image using atrophied muscle segmentation
maskedImage = normalizedInputImage;
maskedImage(~myMuscle) = 0;
% Segment using muscle threshold
myMuscleNoFat = imbinarize(maskedImage,normalizedThresholdMuscle);
% --Osteochondroma segmentation--
%Segment using osteochondroma threshold
myOsteo = imbinarize(maskedImage,normalizedThresholdOsteochondroma);
myOsteo = imfill(myOsteo, 'holes'); % Fill holes
% --Measure atrophy and degeneration so that Stot = Sa + Sm + Si + So--
Stot = bwarea(contouredMuscle.binaryImage1); % Stot: total area of the muscle fossa (manually delimited contours)
Sa = bwarea(contouredMuscle.binaryImage1) - bwarea(myMuscle); % Sa: area of the atrophy (Stot - Sm)
So = bwarea(myOsteo); % So: area of osteochondroma
Sm = bwarea(myMuscleNoFat) - bwarea(myOsteo); % Sm: area of the degenerated muscle, without atrophy, fat infiltration and osteochondroma
Si = bwarea(myMuscle) - Sm - So; % Si: area of fat infiltration
muscle.atrophy = Sa/Stot; % ratio of muscle atrophy (area atrophy over total area of muscle fossa)
muscle.infiltration = Si/Stot; % ratio of fat infiltration in the muscle (area fat infiltration over total area of muscle fossa)
muscle.osteochondroma = So/Stot; % ratio of osteochondroma in the muscle (area osteochondroma over total area of muscle fossa)
muscle.degeneration = (Stot-Sm)/Stot; % the ratio of degeneration of the muscle (sum of atrophy, fat infiltration and osteochondroma over total area of muscle fossa)
% --Plot CT slice with all contours--
fontSize = 14;
visualizationWindow = [-100 200];
imshow(contouredMuscle.inputImage,visualizationWindow);
title(sprintf('Muscle threshold: %.0fHU, Muscle atrophy: %.2f%%, Muscle infiltration: %.2f%%, Osteochondroma: %.2f%%', ...
HUThresholdMuscle,muscle.atrophy*100,muscle.infiltration*100,muscle.osteochondroma*100), 'FontSize', fontSize);
set(gcf, 'Position', get(0,'Screensize'), 'Visible', 'on'); % Maximize figure.
hold on
% Plot muscle fossa contours.
structBoundaries = bwboundaries(contouredMuscle.binaryImage1);
xy=structBoundaries{1}; % Get n by 2 array of x,y coordinates.
x = xy(:, 2); % Columns.
y = xy(:, 1); % Rows.
plot(x, y, 'color','green');
% Plot atrophied muscle contours
structBoundaries = bwboundaries(myMuscle);
xy=structBoundaries{1}; % Get n by 2 array of x,y coordinates.
x = xy(:, 2); % Columns.
y = xy(:, 1); % Rows.
plot(x, y, 'color','blue');
% Plot fat infiltration contours
structBoundaries = bwboundaries(myMuscleNoFat);
for i = 1:size(structBoundaries,1)
xy = structBoundaries{i}; % Get n by 2 array of x,y coordinates.
x = xy(:, 2); % Columns.
y = xy(:, 1); % Rows.
hold on; % Keep the image.
plot(x, y,'color','red');
end
% Plot osteochondroma contours
structBoundaries = bwboundaries(myOsteo);
for i = 1:size(structBoundaries,1)
xy = structBoundaries{i}; % Get n by 2 array of x,y coordinates.
x = xy(:, 2); % Columns.
y = xy(:, 1); % Rows.
hold on; % Keep the image.
plot(x, y,'color','magenta');
end
% --Save images--
outputFile = sprintf('ssDegeneration%s.%d',char(caseID),HUThresholdMuscle);
saveas(gcf,sprintf('%s%s.png',[muscleDirectory '\' muscleName '\'],outputFile))
close all
end
c4science · Help