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. Excel file
% /shoulder/data/Excel/ShoulderDataBase.xlsx should also be open, updated, and saved.

% The script measureSCase was replacing 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-06-29

% 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);

%% 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