Diffusion ShoulderDB_ImportSegmentMeasureAnalyze (glenoidCap) Browse

R8218/measureSCase.m
653c5e9c2788glenoidCap

measureSCase.m
View Options

			function [status,message] = measureSCase(varargin)
			% MEASURESCASE Update the entire SCase DB with anatomical measurements.

			% The function can be run from (lbovenus) server with the following shell command
			% cd /home/shoulder/methods/matlab/database; matlab -nodisplay -nosplash -r "measureSCase;quit"

			% Progress can be checked through log file with following shell command
			% cd /home/shoulder/methods/matlab/database;tail -f log/measureSCase.log

			% It take ~30 min for 700 cases when run from server.
			% After run from server, the output file
			% /shoulder/data/Excel/xlsFromMatlab/anatomy.csv
			% must be open from Excel and saved
			% as xls at the same location.

			% The script measureSCase will replace 3 previous functions
			% 1) rebuildDatabaseScapulaMeasurements.m --> list of all SCases --> execute scapula_measure.m
			% 2) scapula_measure.m --> execute scapula_calculation and save results in Excel and SQL
			% 3) scapula_calculation.m --> perform anatomical analysis --> results varaiable

			% The current version produces correctly the entire csv when run
			% without arguments only (for the entire database)

			% Input
			% Either empty, 'N', 'P', or a SCase ('P315')

			% Output
			% File /shoulder/data/Excel/xlsFromMatlab/anatomy.csv
			% File /shoulder/data/Excel/xlsFromMatlab/anatomy.xls (only windows)
			% File /home/shoulder/data/matlab/SCaseDB.mat
			% File {SCaseId}.mat in each {SCaseId}/{CT}/matlab directory
			% logs in /log/measureSCase.log

			% Example: measureSCase

			% Author: Alexandre Terrier, EPFL-LBO
			% Creation date: 2018-07-01
			% Revision date: 2019-02-05

			% TO DO:
			% Read first last anatomy.csv and/or SCaseDB.mat
			% Fill with patient and clinical data, from excel or MySQL
			% use argument to redo or update
			% Save csv within the main loop (not sure it's a good idea)
			% Add more message in log
			% Add more test to assess validity of data
			% Update MySQL
			% Rename all objects with a capital
			% Make it faster by not loading scapula surface if already in SCaseDB.mat
			% Add argument 'load' to force reload files even if they are alerady in DB
			% Add argument 'measure' to force re-measuring event if measurement is already in DB
			% Add "scapula_addmeasure.m" --> should be coded in ShoulderCase class


			tic; % Start stopwatch timer

			%% Open log file
			logFileID = openLogFile('measureSCase.log');

			%% Set the data directory from the configuration file config.txt
			dataDir = openConfigFile('config.txt', logFileID);



			% configFile = 'config.txt';
			% if ~exist(configFile, 'file')
			% error([configFile, ' required']);
			% end
			% fileID = fopen(configFile,'r');
			% dataDir = fscanf(fileID, '%s'); % Read config file without spaces
			% fclose(fileID);
			% k = strfind(dataDir, 'dataDir='); % Find definition of dataDir
			% k = k + 8; % Remove 'dataDir='
			% dataDir = dataDir(k:end); % Assume that config file ends with dataDir content
			% % Check if dataDir exists
			% if ~exist(dataDir, 'dir')
			% fprintf(logFileID, ['Data directory not found, check ' configFile]);
			% error(['Data directory not found, check ' configFile]);
			% end


			%% Location of the XLS ShoulderDatabase
			xlsDir = [dataDir '/Excel/xlsFromMatlab'];
			matlabDir = [dataDir '/matlab'];




			%% Get list of all SCase
			% Get list of SCases from varargin

			% Delault value for inputArg
			calcGlenoidDensity = 0; % Glenoid density will not be calculated

			fprintf(logFileID, '\n\nList of SCase');
			if nargin == 0 % measureSCase called without argument
			SCaseList = listSCase;
			elseif nargin == 1
			inputArg = varargin{1};
			if strcmp(inputArg, 'glenoidDensity')
			calcGlenoidDensity = 1; % Glenoid density will be calculated
			SCaseList = listSCase; % Measure all cases
			else
			SCaseList = listSCase(inputArg);
			end
			else
			error('inputArg can be empty, N, P, or a SCaseID')
			end

			%% Load current xls database (for clinical data)
			fprintf(logFileID, '\nLoad xls database');
			addpath('XLS_MySQL_DB');
			filename = [dataDir '/Excel/ShoulderDataBase.xlsx'];
			excelRaw = rawFromExcel(filename); % cell array

			excelSCaseID = excelRaw(2:end, 1);
			excelDiagnosis = excelRaw(2:end, 4);
			excelPatientHeight = excelRaw(2:end, 23);
			excelGlenoidWalch = excelRaw(2:end, 55);

			% Lines below adapted from XLS_MySQL_DB/MainExcel2XLS
			% Excel.diagnosisList = diagnosisListFromExcel(Excel.Raw);
			% Excel.treatmentList = treatmentListFromExcel(Excel.Raw);
			% ExcelData.Patient = patientFromExcel(excelRaw); % Structure with patient data
			% Excel.shoulder = shoulderFromExcel(Excel.patient, Excel.Raw);
			% [Excel.SCase, Excel.diagnosis, Excel.treatment, Excel.outcome, Excel.study] = SCaseFromExcel(...
			% Excel.shoulder, Excel.patient, Excel.diagnosisList, Excel.treatmentList, Excel.Raw);
			% fprintf(logFileID, ': OK');

			%% Add path to ShoulderCase class
			addpath('ShoulderCase');

			%% Instance of a ShoulderCase object
			if (exist('SCase','var') > 0)
			clear SCase; % Check for delete method
			end
			SCase = ShoulderCase; % Instanciate a ShoulderCase object
			SCase.dataPath = dataDir; % Set dataDir for SCase

			%% Start loop over SCases in SCaseList
			nSCase = length(SCaseList); % Number of SCases
			for iSCaseID = 1:nSCase
			SCaseID = SCaseList(iSCaseID).id; %
			SCase(iSCaseID).id = SCaseID;
			SCase(iSCaseID).dataPath = dataDir; % Set dataDir for SCase
			percentProgress = num2str(iSCaseID/nSCase*100, '%3.1f');
			fprintf(logFileID, ['\n\nSCaseID: ' SCaseID ' (' percentProgress '%%)']);

			% There are 3 parts within this SCase loop:
			% 1) Load & analyses manual data from amira
			% 2) Load & analyses auto data from matlab (Statistical Shape Model)
			% 3) Load clinical data from Excel database, and set them to SCase
			% 4) Save SCase in file SCaseID/matlab/SCase.mat

			% 1) Load & analyses manual data from amira
			fprintf(logFileID, '\n Segmentation manual ');
			output = SCase(iSCaseID).path; % Set data path of this SCase
			if output % Continue if amira dir exists in SCase
			fprintf(logFileID, '\n Load scapula surface and landmarks');
			output = SCase(iSCaseID).shoulder.scapula.load;
			if output
			fprintf(logFileID, ': OK');
			fprintf(logFileID, '\n Set coord. syst.');
			output = SCase(iSCaseID).shoulder.scapula.coordSysSet;
			if output
			fprintf(logFileID, ': OK');
			fprintf(logFileID, '\n Load glenoid surface');
			output = SCase(iSCaseID).shoulder.scapula.glenoid.load;
			if output
			fprintf(logFileID, ': OK');
			fprintf(logFileID, '\n Measure glenoid anatomy');
			output = SCase(iSCaseID).shoulder.scapula.glenoid.morphology;
			if output
			fprintf(logFileID, ': OK');
			if calcGlenoidDensity
			fprintf(logFileID, '\n Measure glenoid density');
			output = SCase(iSCaseID).shoulder.scapula.glenoid.calcDensity;
			if output
			fprintf(logFileID, ': OK');
			else
			fprintf(logFileID, ': Could no read dicom');
			end
			end
			fprintf(logFileID, '\n Load humerus data');
			output = SCase(iSCaseID).shoulder.humerus.load;
			if output
			fprintf(logFileID, ': OK');
			fprintf(logFileID, '\n Measure humerus subluxation');
			output = SCase(iSCaseID).shoulder.humerus.subluxation;
			if output
			fprintf(logFileID, ': OK');
			fprintf(logFileID, '\n Measure acromion anatomy');
			% Should be run after SCase.shoulder.humerus (for AI)
			output = SCase(iSCaseID).shoulder.scapula.acromion.morphology;
			if output
			fprintf(logFileID, ': OK');
			else
			fprintf(logFileID, ': Error');
			end % Acromion anatmy
			else
			fprintf(logFileID, ': Error');
			end % Subluxation humerus
			else
			fprintf(logFileID, ': Error');
			end % Load humerus
			else
			fprintf(logFileID, ': Error');
			end % Anatomy glenoid
			else
			fprintf(logFileID, ': Error');
			end % Load glenoid
			else
			fprintf(logFileID, '\n No amira folder');
			end % Set coord. syst.
			end % Load scapula
			end % Amira path exist

			% 2) Load & analyses auto data from matlab (Statistical Shape Model)
			fprintf(logFileID, '\n Segmentation auto ');
			% Load scapula surface and landmarks
			fprintf(logFileID, '\n Load scapula surface and landmarks');
			output = SCase(iSCaseID).shoulder.scapulaAuto.loadAuto; % Load auto data from matlab directory
			if output
			fprintf(logFileID, ': OK');
			fprintf(logFileID, '\n Set coord. syst.');
			output = SCase(iSCaseID).shoulder.scapulaAuto.coordSysSet;
			if output
			fprintf(logFileID, ': OK');
			fprintf(logFileID, '\n Load glenoid surface');
			output = SCase(iSCaseID).shoulder.scapulaAuto.glenoid.loadAuto;
			if output
			fprintf(logFileID, ': OK');
			fprintf(logFileID, '\n Measure glenoid anatomy');
			output = SCase(iSCaseID).shoulder.scapulaAuto.glenoid.morphology;
			if output
			fprintf(logFileID, ': OK');
			else
			fprintf(logFileID, ': Error');
			end
			else
			fprintf(logFileID, ': Error');
			end
			else
			fprintf(logFileID, ': Error');
			end
			else
			fprintf(logFileID, ': Error');
			end

			% 3) Set clinical data from Excel database to SCase
			fprintf(logFileID, '\n Set clinical data from Excel');

			% Get idx of excelRaw for SCaseID
			% Use cell2struct when dots are removed from excel headers
			idx = find(strcmp(excelRaw, SCaseID)); % Index of SCaseID in excelRaw
			SCase(iSCaseID).diagnosis = excelRaw{idx,4};
			SCase(iSCaseID).treatment = excelRaw{idx,9};
			SCase(iSCaseID).patient.gender = excelRaw{idx,19};
			SCase(iSCaseID).patient.age = excelRaw{idx,21};
			SCase(iSCaseID).patient.height = excelRaw{idx,23};
			SCase(iSCaseID).patient.weight = excelRaw{idx,24};
			SCase(iSCaseID).patient.BMI = excelRaw{idx,25};
			SCase(iSCaseID).shoulder.scapula.glenoid.walch = excelRaw{idx,55};


			% Patient data --> TO DO
			% find patient id (in SQL) corresponding to iSCaseID
			%
			%
			% idxSCase = find(contains(ExcelData, SCaseID));
			% SCase(iSCaseID).patient.id = ExcelData(idxSCase).patient.patient_id; % SQL patient id
			% SCase(iSCaseID).patient.idMed = ExcelData(idxSCase).patient.IPP; % IPP to be replace by coded (anonymized) IPP from SecuTrial
			% SCase(iSCaseID).patient.gender = ExcelData(idxSCase).patient.gender;
			% SCase(iSCaseID).patient.age = []; % Age (years) at treatment, or preop CT (we might also add birthdate wi day et at 15 of the month)
			% SCase(iSCaseID).patient.ethnicity = ExcelData(idxSCase).patient.IPP;
			% SCase(iSCaseID).patient.weight = ExcelData(idxSCase).patient.weight;
			% SCase(iSCaseID).patient.height = ExcelData(idxSCase).patient.height;
			% SCase(iSCaseID).patient.BMI = ExcelData(idxSCase).pateint.BMI;
			% SCase(iSCaseID).patient.comment = ExcelData(idxSCase).patient.comment;

			% 4) Save SCase in file SCaseID/matlab/SCase.mat
			fprintf(logFileID, '\n Save SCase');
			output = SCase(iSCaseID).saveMatlab;
			if output
			fprintf(logFileID, ': OK');
			else
			fprintf(logFileID, ': Error');
			end
			end % End of loop on SCaseList

			%% Save the entire SCase array as a matlab file
			fprintf(logFileID, '\n\n Save SCase database');
			filename = 'SCaseDB';
			filename = [matlabDir '/' filename '.mat'];
			try
			SCaseDB = SCase;
			save(filename, 'SCaseDB');
			fprintf(logFileID, ': OK');
			catch
			fprintf(logFileID, ': Error');
			end

			%% Write csv file
			% This might be a function (SCase.csvSave()). The input would be a filename and a structure
			% data

			fprintf(logFileID, '\n Save csv database');

			% Replace header and data by structure. Currently not working
			%{
			txtFilename = [xlsDir, '/anatomy.txt']; % Name of the txt file

			DataStruc = struc(...
			'SCase_id', SCase.id, ...
			'shoulder_side', SCase.shoulder.side,...
			'glenoid_radius', SCase.shoulder.scapula.glenoid.radius,...
			'glenoid_sphereRMSE', SCase.shoulder.scapula.glenoid.sphereRMSE,...
			'glenoid_depth', SCase.shoulder.scapula.glenoid.depth,...
			'glenoid_width', SCase.shoulder.scapula.glenoid.width,...
			'glenoid_height', SCase.shoulder.scapula.glenoid.height,...
			'glenoid_centerPA', SCase.shoulder.scapula.glenoid.centerLocal(1),...
			'glenoid_centerIS', SCase.shoulder.scapula.glenoid.centerLocal(2),...
			'glenoid_centerML', SCase.shoulder.scapula.glenoid.centerLocal(3),...
			'glenoid_versionAmpl', SCase.shoulder.scapula.glenoid.versionAmpl,...
			'glenoid_versionOrient', SCase.shoulder.scapula.glenoid.versionOrient,...
			'glenoid_version', SCase.shoulder.scapula.glenoid.version,...
			'glenoid_inclination', SCase.shoulder.scapula.glenoid.inclination,...
			'humerus_jointRadius', SCase.shoulder.humerus.jointRadius,...
			'humerus_headRadius', SCase.shoulder.humerus.radius,...
			'humerus_GHSAmpl', SCase.shoulder.humerus.GHSAmpl,...
			'humerus_GHSOrient', SCase.shoulder.humerus.GHSOrient,...
			'humerus_SHSAmpl', SCase.shoulder.humerus.SHSAmpl,...
			'humerus_SHSOrient', SCase.shoulder.humerus.SHSOrient,...
			'humerus_SHSAngle', SCase.shoulder.humerus.SHSAngle,...
			'humerus_SHSPA', SCase.shoulder.humerus.SHSPA,...
			'humerus_SHSIS', SCase.shoulder.humerus.SHSIS,...
			'acromion_AI', SCase.shoulder.scapula.acromion.AI,...
			'acromion_CSA', SCase.shoulder.scapula.acromion.CSA,...
			'acromion_PS', SCase.shoulder.scapula.acromion.PS...
			);
			DataTable = strct2table(Data);
			writetable(DataTable,filename);
			%}

			% Header of the csv file
			dataHeader = [...
			'SCase_id,' ...
			'shoulder_side,' ...
			'glenoid_radius,' ...
			'glenoid_sphereRMSE,' ...
			'glenoid_depth,' ...
			'glenoid_width,' ...
			'glenoid_height,' ...
			'glenoid_centerPA,' ...
			'glenoid_centerIS,' ...
			'glenoid_centerML,' ...
			'glenoid_versionAmpl,' ...
			'glenoid_versionOrient,' ...
			'glenoid_version,' ...
			'glenoid_inclination,' ...
			'humerus_jointRadius,' ...
			'humerus_headRadius,' ...
			'humerus_GHSAmpl,' ...
			'humerus_GHSOrient,' ...
			'humerus_SHSAmpl,' ...
			'humerus_SHSOrient,' ...
			'humerus_SHSAngle,' ...
			'humerus_SHSPA,' ...
			'humerus_SHSIS,' ...
			'acromion_AI,' ...
			'acromion_CSA,' ...
			'acromion_PSA,'...
			'acromion_AAA'...
			];

			csvFilename = [xlsDir, '/anatomy.csv']; % Name of the csv file
			csvFileId = fopen(csvFilename, 'w');
			fprintf(csvFileId,dataHeader); % Write header
			for iSCaseID = 1:length(SCaseList) % Loop on SCaseList to write data in csv
			disp(SCase(iSCaseID).id); % For debug
			fprintf(csvFileId,['\n'...
			SCase(iSCaseID).id ','... % SCase_id
			SCase(iSCaseID).shoulder.side ','... % shoulder_side
			num2str(SCase(iSCaseID).shoulder.scapula.glenoid.radius) ','... % glenoid_radius
			num2str(SCase(iSCaseID).shoulder.scapula.glenoid.sphereRMSE) ','... % glenoid_sphereRMSE
			num2str(SCase(iSCaseID).shoulder.scapula.glenoid.depth) ','... % glenoid_depth
			num2str(SCase(iSCaseID).shoulder.scapula.glenoid.width) ','... % glenoid_width
			num2str(SCase(iSCaseID).shoulder.scapula.glenoid.height) ','... % glenoid_height
			num2str(SCase(iSCaseID).shoulder.scapula.glenoid.centerLocal.x) ','... % glenoid_centerPA
			num2str(SCase(iSCaseID).shoulder.scapula.glenoid.centerLocal.y) ','... % glenoid_centerIS
			num2str(SCase(iSCaseID).shoulder.scapula.glenoid.centerLocal.z) ','... % glenoid_centerML
			num2str(SCase(iSCaseID).shoulder.scapula.glenoid.versionAmpl) ','... % glenoid_versionAmpl
			num2str(SCase(iSCaseID).shoulder.scapula.glenoid.versionOrient) ','... % glenoid_versionOrient
			num2str(SCase(iSCaseID).shoulder.scapula.glenoid.version) ','... % glenoid_version
			num2str(SCase(iSCaseID).shoulder.scapula.glenoid.inclination) ','... % glenoid_inclination
			num2str(SCase(iSCaseID).shoulder.humerus.jointRadius) ','... % humerus_jointRadius
			num2str(SCase(iSCaseID).shoulder.humerus.radius) ','... % humerus_headRadius
			num2str(SCase(iSCaseID).shoulder.humerus.GHSAmpl) ','... % humerus_GHSAmpl
			num2str(SCase(iSCaseID).shoulder.humerus.GHSOrient) ','... % humerus_GHSOrient
			num2str(SCase(iSCaseID).shoulder.humerus.SHSAmpl) ','... % humerus_SHSAmpl
			num2str(SCase(iSCaseID).shoulder.humerus.SHSOrient) ','... % humerus_SHSOrient
			num2str(SCase(iSCaseID).shoulder.humerus.SHSAngle) ','... % humerus_SHSAgle
			num2str(SCase(iSCaseID).shoulder.humerus.SHSPA) ','... % humerus_SHSPA
			num2str(SCase(iSCaseID).shoulder.humerus.SHSIS) ','... % humerus_SHSIS
			num2str(SCase(iSCaseID).shoulder.scapula.acromion.AI) ','... % acromion_AI
			num2str(SCase(iSCaseID).shoulder.scapula.acromion.CSA) ','... % acromion_CSA
			num2str(SCase(iSCaseID).shoulder.scapula.acromion.PSA) ','... % acromion_PSA
			num2str(SCase(iSCaseID).shoulder.scapula.acromion.AAA)... % acromion_AAA
			]);
			end
			fclose(csvFileId);

			fprintf(logFileID, [': ' csvFilename]);

			%% Write xls file (Windows only)
			[~,message] = csv2xlsSCase(csvFilename);
			fprintf(logFileID, message);

			% If run from non-windows system, only csv will be produced. The xls can be
			% produced by opening the csv from Excel and save as xls.
			% Could be a function with 1 input (cvsFilenames



			%{
			xlsFilename = strrep(csvFilename, '.csv', '.xls'); % Replace .csv by .xls
			if ispc % Check if run from windows!
			% Read cvs file as table and write the table as xls
			cvsTable = readtable(csvFilename); % Get table from csv (only way to read non-numeric data)
			cvsCell = table2cell(cvsTable); % Tranform table cell array
			dataHeader = cvsTable.Properties.VariableNames; % Get header
			cvsCell = [dataHeader; cvsCell]; % Add header

			sheet = 'anatomy'; % Sheet name of the xls file

			[status,message] = xlswrite(xlsFilename, cvsCell, sheet); % Save xls
			if status
			fprintf(logFileId, ['\nSaved ' xlsFilename]);
			else
			fprintf(logFileId, ['\nCould not save ' xlsFilename]);
			fprintf(logFileId, ['\n' message.identifier]);
			fprintf(logFileId, ['\n' message.message]);
			end
			else
			fprintf(logFileId, ['\n' xlsFilename ' not saved. Open and save as xls from Excel']);
			end
			%}


			%% Close log file
			fprintf(logFileID, '\n\nSCase measured: %s', num2str(length(SCaseList))); % Number of SCases measured
			fprintf(logFileID, '\nElapsed time (s): %s', num2str(toc)); % stopwatch timer
			fprintf(logFileID, '\n');
			fclose(logFileID); % Close log file

			%% Output of the SCase if required by argument
			% SCase.output % inputArg = abaqus/density/references/

			% update mySQL
			% SCase.sqlSave


			status = 1;
			message = 'OK';

			end

R8218/measureSCase.m653c5e9c2788glenoidCap

measureSCase.mView Options

R8218/measureSCase.m
653c5e9c2788glenoidCap

measureSCase.m
View Options