diff --git a/ShoulderCase/@Muscle/Muscle.m b/ShoulderCase/@Muscle/Muscle.m index 2585515..633c279 100644 --- a/ShoulderCase/@Muscle/Muscle.m +++ b/ShoulderCase/@Muscle/Muscle.m @@ -1,103 +1,120 @@ classdef Muscle < handle % The Muscle class is linked to a segmented muscle file (a mask) % and to the slice it has been segmented out of. % % Then, this class can measured values linked to the PCSA and the % muscle's degeneration. properties container = nan; name = ""; segmentationName = ""; sliceName = ""; PCSA = nan; atrophy = nan; fat = nan; osteochondroma = nan; degeneration = nan; centroid = nan; forceApplicationPoint = nan; forceVector = nan; end methods function obj = Muscle(musclesContainer,muscleName) obj.container = musclesContainer; obj.name = muscleName; [~, ~] = mkdir(obj.dataPath); [~, ~] = mkdir(obj.dataMaskPath); end function output = dataPath(obj) output = fullfile(obj.container.dataPath, obj.name); end function output = dataMaskPath(obj) output = fullfile(obj.dataPath, "mask"); end function setSliceName(obj, value) obj.sliceName = value; end function setSegmentationName(obj, value) obj.segmentationName = value; end function output = loadMask(obj, maskSuffix) output = logical(imread(fullfile(obj.dataMaskPath,... obj.segmentationName+"_"+maskSuffix+".png"))); end function saveMask(obj, mask, maskSuffix) imwrite(mask, fullfile(obj.dataMaskPath,... obj.segmentationName+"_"+maskSuffix+".png")); end - function output = loadPixelCoordinates(obj) - output = load(fullfile(obj.container.dataSlicesPath,... - obj.sliceName + "_PixelCoordinates.mat")).imagesPixelCoordinates; - end - + function output = loadSlice(obj, sliceSuffix) + filePattern = fullfile(obj.container.dataSlicesPath,... + obj.sliceName + "_" + sliceSuffix + "*"); + matchingFile = dir(filePattern); + assert(not(isempty(matchingFile)),... + "File matching %s not found", filePattern); + assert(length(matchingFile) <= 1,... + "Multiple matching files found. Include file extension in given suffix"); + if endsWith(matchingFile.name, ".mat") + loadedFile = load(fullfile(matchingFile.folder, matchingFile.name)); + % loaded Matlab archives are structures whose fields contain the saved + % variables. When only one variable is saved we want to return its value + % directly and avoid having to retrieve the original variable name. + loadedVariables = string(fields(loadedFile)); + if isequal(length(loadedVariables), 1) + loadedFile = loadedFile.(loadedVariables); + end + output = loadedFile; + elseif endsWith(matchingFile.name, ".png") + output = imread(fullfile(matchingFile.folder, matchingFile.name)); + end + end function output = measureFirst(obj) % Call methods that can be run after morphology() methods has been run % by all ShoulderCase objects. obj.setSliceName("rotatorCuffMatthieu"); obj.setSegmentationName("autoMatthieu"); success = Logger.timeLogExecution(... obj.name+" measure and save contour: ",... @(obj) obj.measureAndSaveContour(), obj); success = success | Logger.timeLogExecution(... obj.name+" measure and save centroid: ",... @(obj) obj.measureAndSaveCentroid(), obj); success = success | Logger.timeLogExecution(... obj.name+" PCSA: ",... @(obj) obj.measurePCSA(), obj); success = success | Logger.timeLogExecution(... obj.name+" atrophy: ",... @(obj) obj.measureAtrophy(), obj); success = success | Logger.timeLogExecution(... obj.name+" fat infiltration: ",... @(obj) obj.measureFatInfiltration(), obj); success = success | Logger.timeLogExecution(... obj.name+" osteochondroma: ",... @(obj) obj.measureOsteochondroma(), obj); success = success | Logger.timeLogExecution(... obj.name+" degeneration: ",... @(obj) obj.measureDegeneration(), obj); success = success | Logger.timeLogExecution(... obj.name+" humeral head contact point: ",... @(obj) obj.measureHumerusContactPoint(), obj); success = success | Logger.timeLogExecution(... obj.name+" force: ",... @(obj) obj.measureForceVector(), obj); output = success; end end end diff --git a/ShoulderCase/@Muscle/getMaskCoordinates.m b/ShoulderCase/@Muscle/getMaskCoordinates.m index ddac808..bced46a 100644 --- a/ShoulderCase/@Muscle/getMaskCoordinates.m +++ b/ShoulderCase/@Muscle/getMaskCoordinates.m @@ -1,10 +1,10 @@ function output = getMaskCoordinates(obj, maskSuffix) mask = obj.loadMask(maskSuffix); maskIndices = find(mask); - slicePixelCoordinates = obj.loadPixelCoordinates(); + slicePixelCoordinates = obj.loadSlice("PixelCoordinates"); maskCoordinates = [... slicePixelCoordinates.x(maskIndices),... slicePixelCoordinates.y(maskIndices),... slicePixelCoordinates.z(maskIndices)]; output = maskCoordinates; end \ No newline at end of file diff --git a/ShoulderCase/@Muscle/measureAtrophy.m b/ShoulderCase/@Muscle/measureAtrophy.m index 1a0afa7..9ceddd0 100644 --- a/ShoulderCase/@Muscle/measureAtrophy.m +++ b/ShoulderCase/@Muscle/measureAtrophy.m @@ -1,4 +1,5 @@ function measureAtrophy(obj) - measurer = MuscleMeasurer(obj); + measurer = MuscleMeasurer(obj.loadMask("Segmentation"),... + obj.loadSlice("ForMeasurements"), obj.loadSlice("PixelSpacings")); obj.atrophy = measurer.getRatioAtrophy(); end \ No newline at end of file diff --git a/ShoulderCase/@Muscle/measureDegeneration.m b/ShoulderCase/@Muscle/measureDegeneration.m index 7db9589..ad5e3d5 100644 --- a/ShoulderCase/@Muscle/measureDegeneration.m +++ b/ShoulderCase/@Muscle/measureDegeneration.m @@ -1,4 +1,5 @@ function measureDegeneration(obj) - measurer = MuscleMeasurer(obj); + measurer = MuscleMeasurer(obj.loadMask("Segmentation"),... + obj.loadSlice("ForMeasurements"), obj.loadSlice("PixelSpacings")); obj.degeneration = measurer.getRatioDegeneration(); end \ No newline at end of file diff --git a/ShoulderCase/@Muscle/measureFatInfiltration.m b/ShoulderCase/@Muscle/measureFatInfiltration.m index a5467e5..e9dd014 100644 --- a/ShoulderCase/@Muscle/measureFatInfiltration.m +++ b/ShoulderCase/@Muscle/measureFatInfiltration.m @@ -1,4 +1,5 @@ function measureFatInfiltration(obj) - measurer = MuscleMeasurer(obj); + measurer = MuscleMeasurer(obj.loadMask("Segmentation"),... + obj.loadSlice("ForMeasurements"), obj.loadSlice("PixelSpacings")); obj.fat = measurer.getRatioFat(); end \ No newline at end of file diff --git a/ShoulderCase/@Muscle/measureOsteochondroma.m b/ShoulderCase/@Muscle/measureOsteochondroma.m index e05e4a3..78d3075 100644 --- a/ShoulderCase/@Muscle/measureOsteochondroma.m +++ b/ShoulderCase/@Muscle/measureOsteochondroma.m @@ -1,4 +1,5 @@ function measureOsteochondroma(obj) - measurer = MuscleMeasurer(obj); + measurer = MuscleMeasurer(obj.loadMask("Segmentation"),... + obj.loadSlice("ForMeasurements"), obj.loadSlice("PixelSpacings")); obj.osteochondroma = measurer.getRatioOsteochondroma(); end \ No newline at end of file diff --git a/ShoulderCase/@Muscle/measurePCSA.m b/ShoulderCase/@Muscle/measurePCSA.m index 8df5381..ee11c18 100644 --- a/ShoulderCase/@Muscle/measurePCSA.m +++ b/ShoulderCase/@Muscle/measurePCSA.m @@ -1,4 +1,5 @@ function measurePCSA(obj) - measurer = MuscleMeasurer(obj); + measurer = MuscleMeasurer(obj.loadMask("Segmentation"),... + obj.loadSlice("ForMeasurements"), obj.loadSlice("PixelSpacings")); obj.PCSA = measurer.getPCSA(); end \ No newline at end of file diff --git a/ShoulderCase/@MuscleMeasurer/MuscleMeasurer.m b/ShoulderCase/@MuscleMeasurer/MuscleMeasurer.m index b4d1f54..e5b6896 100644 --- a/ShoulderCase/@MuscleMeasurer/MuscleMeasurer.m +++ b/ShoulderCase/@MuscleMeasurer/MuscleMeasurer.m @@ -1,119 +1,106 @@ classdef MuscleMeasurer < handle % Used to perfom PCSA and degeneration measurements of a Muscle object. % Slices and mask images are expected to be found at some specific places % and have specific names. % % This class is the results of extracting measurements methods from the % project of Nathan Donini. properties - muscle - segmentedArea - segmentedImage + segmentationMask + segmentationDensity + slicePixelSpacings rangeHU = double([-1000 1000]); muscleThreshold = 0; fatThreshold = 30; osteochondromaThreshold = 166; end methods - function obj = MuscleMeasurer(muscle) - obj.muscle = muscle; - obj.segmentedArea = obj.getSegmentedArea; - obj.segmentedImage = obj.getSegmentedImage; + function obj = MuscleMeasurer(segmentationMask, sliceForMeasurement, slicePixelSpacings) + obj.segmentationMask = segmentationMask; + obj.segmentationDensity =... + obj.getNormalizedImage(sliceForMeasurement).*segmentationMask; + obj.slicePixelSpacings = slicePixelSpacings; obj.normalizeThresholds; end - function output = getSegmentedArea(obj) - output = logical(imread(fullfile(obj.muscle.dataMaskPath,obj.muscle.segmentationName + "_Segmentation.png"))); - end - - function output = getSegmentedImage(obj) - load(fullfile(obj.muscle.container.dataSlicesPath,obj.muscle.sliceName + "_ForMeasurements.mat")); - normalizedImage = obj.getNormalizedImage(imageForMeasurements); - - segmentedImage = normalizedImage.*obj.segmentedArea; - - output = segmentedImage; - end - function output = getNormalizedImage(obj,image) image(image < obj.rangeHU(1)) = obj.rangeHU(1); image(image > obj.rangeHU(2)) = obj.rangeHU(2); output = (double(image)-obj.rangeHU(1))/(obj.rangeHU(2)-obj.rangeHU(1)); end function normalizeThresholds(obj) obj.muscleThreshold = ((obj.muscleThreshold - 1)-obj.rangeHU(1))/(obj.rangeHU(2)-obj.rangeHU(1)); obj.fatThreshold = ((obj.fatThreshold - 1)-obj.rangeHU(1))/(obj.rangeHU(2)-obj.rangeHU(1)); obj.osteochondromaThreshold = ((obj.osteochondromaThreshold - 1)-obj.rangeHU(1))/(obj.rangeHU(2)-obj.rangeHU(1)); end function output = getPCSA(obj) - load(fullfile(obj.muscle.container.dataSlicesPath,obj.muscle.sliceName + "_PixelSpacings.mat")); - pixelSurface = (imagesPixelSpacings(1) * imagesPixelSpacings(2)) / 100; % cm^2 - output = pixelSurface * bwarea(obj.segmentedArea); + pixelSurface = (obj.slicePixelSpacings(1) * obj.slicePixelSpacings(2)) / 100; % cm^2 + output = pixelSurface * bwarea(obj.segmentationMask); end function output = getRatioAtrophy(obj) - output = obj.getRatioAreas(obj.getAreaAtrophy,obj.segmentedArea); + output = obj.getRatioAreas(obj.getAreaAtrophy,obj.segmentationMask); end function output = getRatioFat(obj) - output = obj.getRatioAreas(obj.getAreaFat,obj.segmentedArea); + output = obj.getRatioAreas(obj.getAreaFat,obj.segmentationMask); end function output = getRatioOsteochondroma(obj) - output = obj.getRatioAreas(obj.getAreaOsteochondroma,obj.segmentedArea); + output = obj.getRatioAreas(obj.getAreaOsteochondroma,obj.segmentationMask); end function output = getRatioDegeneration(obj) output = obj.getRatioAtrophy + obj.getRatioFat + obj.getRatioOsteochondroma; end function output = getRatioAreas(obj,partialArea,totalArea) output = bwarea(partialArea)/bwarea(totalArea); end function output = getAreaMuscle(obj) % The area "Muscle" is the area of the segmented image that is not atrophied. % Looking for a better name. - areaMuscle = imbinarize(obj.segmentedImage,obj.muscleThreshold); + areaMuscle = imbinarize(obj.segmentationDensity,obj.muscleThreshold); areaMuscle = imfill(areaMuscle,'holes'); % Keep biggest island only islands = bwconncomp(areaMuscle); islandsSizes = cellfun(@numel,islands.PixelIdxList); [~,biggestIslandIndex] = max(islandsSizes); areaMuscle = false(size(areaMuscle)); areaMuscle(islands.PixelIdxList{biggestIslandIndex}) = true; output = areaMuscle; end function output = getAreaAtrophy(obj) - output = obj.segmentedArea & not(obj.getAreaMuscle); + output = obj.segmentationMask & not(obj.getAreaMuscle); end function output = getAreaFat(obj) - muscleImage = obj.segmentedImage.*obj.getAreaMuscle; + muscleImage = obj.segmentationDensity.*obj.getAreaMuscle; areaFat = obj.getAreaMuscle & not(imbinarize(muscleImage,obj.fatThreshold)); output = areaFat; end function output = getAreaOsteochondroma(obj) - muscleImage = obj.segmentedImage.*obj.getAreaMuscle; + muscleImage = obj.segmentationDensity.*obj.getAreaMuscle; areaOsteochondroma = imbinarize(muscleImage,obj.osteochondromaThreshold); areaOsteochondroma = imfill(areaOsteochondroma,'holes'); output = areaOsteochondroma; end end end