function I=CalcAdvancedFeatures(Scan, mask)
% This function calculates advanced features from medical images using
% Matlab's Radiomics tools.
%
% Inputs:
%   Scan - A cell array containing paths to Scan images. This is generated
%       via the SortDicomFiles_OUH and LoadDicomData_OUH pipeline in the
%       CollectDataForModel function.
%   mask - A cell array containing a.o. Voxel Data for cooresponding maths.
%       This mask is also generated in the CollectDataForModel pathway.
% 
% Outputs:
%   I - A structure containing various properties extraced from the image,
%       including shape, intensity and teture features
%
% Please refer to Matlab's documentation about Radiomics for more
% information: https://se.mathworks.com/help/medical-imaging/ref/radiomics.html#mw_61595c22-92ef-46b1-b2f2-82febab77a71
%
% Note:
% Matlab's radiomics is very new and has only been included since R2023b.
% The radiomics is (imo) rather slow and therefore this code will increase
% the calculation time of the data collection quite a lot. Therefore the
% CalcAdvancedFeature is at the moment not used as a default.

imageFilepaths = strings(length(Scan{1,1}), 1);
% Loop through each filename in the first element of the Scan cell array and
% assign the filename to the imageFilepaths array to get a total list of
% all Dicom files making up an image
for iFilenames = 1:length(Scan{1,1})
    imageFilepaths(iFilenames) = Scan{1,1}{1,iFilenames}.Filename;
end

% Load the medical volume image using the filenames stored in
% imageFilepaths. Note that we will actually not use this medicalVolume for
% the radiomics calculation and only as a reference for registration of
% scan and image matrix that are already existend in the Scan and mask
% variable.
medVolImage = medicalVolume(imageFilepaths);

% Create the reference information object based on the size and voxel
% distances of the loaded medical volume image. Note that this is a
% medicalref3d object and not a imref3d object.
referenceInfo = medicalref3d(size(medVolImage.Voxels), medVolImage.VolumeGeometry.Position, medVolImage.VolumeGeometry.VoxelDistances);

% Now we create the actual medical Volumes that will be used in the
% radiomics calculation. The lines at the end of the function can help you
% to verify that scan and mask were registered correctly.
medVolImage = medicalVolume(Scan{1,2}, referenceInfo);
medVolStruct = medicalVolume(mask{1,2}, referenceInfo);

radiomicsObj = radiomics(medVolImage, medVolStruct);

I.shapeProperties = shapeFeatures(radiomicsObj);
I.intensityProperties = intensityFeatures(radiomicsObj);
I.textureProperties = textureFeatures(radiomicsObj);

% % Uncomment the following lines if you want to visualize a slice of the
% % image/mask combination.
% slice = 130;
% imshowpair(medVolImage.Voxels(:,:,slice), medVolStruct.Voxels(:,:,slice))

end