function [logTransformRoiPCAData]=CalculateFurtherParameterValuesForInference(StudyOverview, averageDirection, Model, mapFilepath, removeRoiName)

% Calculate further parameters for our PCA Outlier model.
% This includes: Euler Angles, PCA Ratio and log transforms.

index = 1;
nStudy = length(StudyOverview);
for iStudy =  1:nStudy
    % Creates a new struct that complies with the structure
    % needed for follow-up methods.
    % The initial information (i.e. StudyOverview) stems from the 
    % BasicFeature calculation.
    % C.f.: CollectDataForModel.m
    if isempty(StudyOverview(iStudy).ROIInformation)
        warning('Study in row %i of StudyOverview is empty. Skipping...', iStudy);
        continue
    end
    ROIs = fieldnames(StudyOverview(iStudy).ROIInformation); % Name of all the ROIs contained in RTStruct
    nROI = length(ROIs); % Number of ROIs contained in RTStruct
    for iRoi = 1:nROI
        validFieldName = matlab.lang.makeValidName(ROIs{iRoi});
        roiInfo = StudyOverview(iStudy).ROIInformation.(validFieldName);
        if isstruct(roiInfo)
            try
                roiInfo.ROIName = validFieldName;
            catch ME
                fprintf('iStudy %d caused issues for ROI Nr. %d \n', iStudy, iRoi);
            end
            roiInfo.Study = StudyOverview(iStudy).StudyInstanceUID;
            roiInfo.Series = StudyOverview(iStudy).ScanInformation.SeriesInstanceUID;
            roiInfo.StdOverviewRowNumber = iStudy;
            roiInfo.MRN = StudyOverview(iStudy).PatientInformation.PatientID;
            roiInfo.PtNumber = iStudy;
    
            patientInfoStruct = StudyOverview(iStudy).PatientInformation;
            fieldnamesPatient = fieldnames(patientInfoStruct);
        
            for iPatientFields = 1:numel(fieldnamesPatient)
                roiInfo.(fieldnamesPatient{iPatientFields}) =  patientInfoStruct.(fieldnamesPatient{iPatientFields});
            end
            
            scanInfoStruct = StudyOverview(iStudy).ScanInformation;
            fieldnamesScan = fieldnames(scanInfoStruct);
            for iScanFields = 1:numel(fieldnamesScan)
                roiInfo.(fieldnamesScan{iScanFields}) =  scanInfoStruct.(fieldnamesScan{iScanFields});
            end

            fileInfoStruct = StudyOverview(iStudy).FileInformation;
            for iFiles = 1:size(fileInfoStruct,1)
                roiInfo.(fileInfoStruct{iFiles,1}) =  fileInfoStruct{iFiles,2};
            end
            try
                roiData(index) = roiInfo;
            catch ME
                roiData
                disp(roiInfo)
                rethrow(ME);
            end
        end
        index = index + 1;
    end    
end

% The ROIs in roiPCA struct contains the following parameters in the following order
% Volume
% xcm
% ycm
% zcm
% lengthPCA1
% lengthPCA2
% lengthPCA3
% transformedCoord
% CoeffPCA
% Area
% Compactness
% CMToSurface
% MaxLength
if nargin > 3
  roiData = correctRoiNames(roiData, mapFilepath);
  if removeRoiName
    roiData = removeRoiNames(roiData, mapFilepath);
  end
end

%%
% In some instances CalcBasicFeatures returns a ROI with Volume 0.
% This either happens because the ROI is defined in the RTStruct and is
% indeed empty (i.e. no delineation) or because the ROI could not be found
% in an RTStruct despite being requested by the calcBasicFeatures function.
% We therefore filter out all these cases.
keepIndex = ~cellfun(@isempty, {roiData.Volume});
roiData = roiData(keepIndex);

% The following code takes the CoeffPCA and calculates the Euler Angles 
% This code is essential identical to Carsten's CalcPCAAxis code
% with the only change being to rename the temp variables ;-)
% This is not needed during inference.
% averageDirection = CalcPCAAxis(roiData);

roiData = AddPCADirection(roiData, averageDirection);
roiData = AddPCARatio(roiData);

% Temp2=RestructureRoiData(roiData); done with applyLogTransform
% Temp2=RemoveModelOutliers(Temp2);
% Temp2=TranslationCorrection(Temp2);
% FinalModelData=ImputePositionMatrix(Temp2);

roiData = CalcIslandPoissonProbability(roiData, Model);
roiData = RelativePos(roiData, Model);

% Define the fieldnames of roiData you want to plot
% TODO: add patient Parameters and scan Parameters that we want to keep
% NOTE: We remove xcm, ycm and zcm here because we will rely on relative
% positions instead.
fieldsToKeep = {'PatientID', 'StudyInstanceUID', 'SeriesInstanceUID', 'StdOverviewRowNumber', 'StudyDate', 'ROIName', ...
    'Volume', ...
    'lengthPCA1', 'lengthPCA2', 'lengthPCA3', ...
    'Area', 'Compactness', 'CMToSurface', 'MaxLength', ...
    'PCAalpha', 'PCAbeta', 'PCAgamma', 'RatioPCA21', 'RatioPCA31', 'Solidity', 'IslandsLikelihood', 'IslandsLambda', 'Islands', 'Scan', 'RTStruct'};

roiData = FlattenRelativePosition(roiData, true);

availableFieldNames = fieldnames(roiData);
fieldsToKeep = [fieldsToKeep, availableFieldNames(startsWith(availableFieldNames, "relPos_"))']; % add all the relPos Parameters


% Create a new struct to hold the copied fields
roiPCAData = struct();

% Loop through the fields to keep and copy them to the new struct
for i = 1: length(roiData)
    for j =  1:numel(fieldsToKeep)
        fieldName = fieldsToKeep{j};
        if isfield(roiData(i), fieldName)
            roiPCAData(i).(fieldName) = roiData(i).(fieldName);
        end
    end
end

% Create a map to hold unique PatientIDs and their corresponding PtNumber
ptMap = containers.Map('KeyType', 'char', 'ValueType', 'double');

% Initialize PtNumber
ptNumber = 1;

% Loop through each struct in the array
for i = 1:length(roiPCAData)
    % Check if the PatientID is already in the map
    if ~isKey(ptMap, roiPCAData(i).PatientID)
        % If not, add it to the map with the current PtNumber
        ptMap(roiPCAData(i).PatientID) = ptNumber;
        % Increment PtNumber for the next unique PatientID
        ptNumber = ptNumber + 1;
    end
    % Add the PtNumber field to the current struct
    roiPCAData(i).PtNumber = ptMap(roiPCAData(i).PatientID);
end

% Create a map to hold unique StudyIDs and their corresponding StudyNumber
stdMap = containers.Map('KeyType', 'double', 'ValueType', 'double');

% Initialize PtNumber
stdNumber = 1;

% Loop through each struct in the array and see if a RTStruct is unique
for i = 1:length(roiPCAData)
    % Check if the RTStruct is already in the map.
    % We can't unfortunately not just look at the StudyInstanceUID or
    % SeriesInstanceUID as there might be weird datasets where we actually
    % have two RTStructs that are NOT identical but part of the same
    % SeriesInstance. So that's why we create a hash over some parameters
    % and the SeriesInstanceUID
    uniqueIdentifierGroup = [roiPCAData(i).StudyInstanceUID, roiPCAData(i).SeriesInstanceUID, roiPCAData(i).StdOverviewRowNumber];
    uniqueIdentifierHash = keyHash(uniqueIdentifierGroup);
    if ~isKey(stdMap, uniqueIdentifierHash)
        % If not existend, add it to the map with the current StdNumber
        stdMap(uniqueIdentifierHash) = stdNumber;
        % Increment StdNumber for the next unique RTStruct
        stdNumber = stdNumber + 1;
    end
    % Add the StdNumber field to the current struct
    roiPCAData(i).StdNumber = stdMap(uniqueIdentifierHash);
end

% Restructure Roi Data and log transform
logTransformRoiPCAData = applyLogTransformation(roiPCAData);

% Convert the struct to a table
% roiPCAData = struct2table(roiPCAData);
logTransformRoiPCAData = struct2table(logTransformRoiPCAData);
end

function [roiData] = correctRoiNames(roiData, mapFilepath)
% CORRECTROINAMES
% Based on a mapping datastructure obtained via a excel spreadsheet
% roiNames in the provided data are corrected.
% (e.g. Prostata -> Prostate).
% The excel spreadsheet needs to be structures as follows:
% Each row is ROI specific with the first column containing the "desired"
% ROI name, while the other columns contain the "existent but not desired"
% ROI names. The "existent but not desired" ROI names are changed into the
% "desired" ROI name.
% NOTE: As opposed to the inference ROI name correction via the
% Damerau-Levenshtein distance, this method is not as robust to tiny errors
% in the mapping. This function ONLY corrects what you tell it to correct.
% 
% E.g.:
% Prostate_AI	Prostata	Prostate	Prosttate	Prostata_AI
% FemoralHead_AI	FemoralHead	Femoral_Head		
% Sigmoid_AI	Sigmo	Sigmi	Sigmoaid	
%
% This will ONLY correct Prostata, Prostate, Prosttate, Prostata_AI to
% Prostate_AI. Same goes for FemoralHead_AI and Sigmoid_AI.
% Thus you have to be aware of all the deviations in the namegiving of your
% data yourself!


  mapROITable = readcell(mapFilepath);
  wrongRoiNames = mapROITable(1:end,2:end);
  wrongRoiNames(cellfun(@(x) any(ismissing(x)), wrongRoiNames)) = {''};
  
  for iRow = 1:length(roiData)
    [Lia,Locb] = ismember(roiData(iRow).ROIName, wrongRoiNames);
    if Lia
      rowIndex = mod(Locb,length(wrongRoiNames(:,1)));
      roiData(iRow).ROIName = mapROITable{rowIndex,1};
    end
  end
end

function [keepStruct] = removeRoiNames(roiData, mapFilepath)
% REMOVEROINAMES
% This removes any ROI Name that is not part of the excel table of
% mapFilepath
mapROITable = readcell(mapFilepath);
roiNamesToKeep = mapROITable(:,1);

for iRow = 1:length(roiData)
  [Lia,Locb] = ismember(roiData(iRow).ROIName, roiNamesToKeep);
  if Lia && iRow == 1
    keepStruct = roiData(iRow);
  end
  if Lia && iRow > 1
    keepStruct(end+1) = roiData(iRow);
  end
end
C
end