function createRotatorCuffSlicesMatthieu(obj,sliceName)
  % This function and its results might be reffered as 
  % Matthieu algorithm and Matthieu slices in some studies.
  %
  % This method has been written to simplify, improve, and make reusable 
  % the Nathan's slicing algorithm written in createRotatorCuffSlicesNathan(). 
  % This means:
  %     - taking advantage of the MATLAB obliqueslice() method wich wasn't available
  %       when Nathan's algorithm has been written.
  %     - encapsulating the slicing data and methods in a DicomVolumeSlicer object.
  %     - being able to orient a slice according to dicom points coordinates.
  %     - being able to crop the slice according to dicom points coordinates
  %       and cropping dimensions given in the dicom set units. 
  %
  % Images and pixel spacings are saved at:
  % shoulder/dataDev/*/*/*/SCase-IPP/CT-SCase-IPP/matlab/shoulder/muscles/oblique_slice/

  assert(not(obj.shoulder.scapula.coordSys.isEmpty),'Scapula coordinate system not measured.');
  [imageForSegmentation,imageForMeasurements,imagesPixelSpacings,imagesPixelCoordinates] = getSlices(obj);
  saveImages(obj,imageForSegmentation,imageForMeasurements,sliceName);
  saveImagesPixelSpacings(obj,imagesPixelSpacings,sliceName);
  saveImagesPixelCoordinates(obj,imagesPixelCoordinates,sliceName);
end



function [slice8Bit,slice16Bit,slicePixelSpacings,slicePixelCoordinates] = getSlices(obj)
  SCase = obj.shoulder.SCase;
  scapula = obj.shoulder.scapula;
  path = SCase.dataDicomPath;
  try
    rotatorCuffSlicer = DicomVolumeSlicer(SCase.getSmoothDicomPath);
  catch
    rotatorCuffSlicer = DicomVolumeSlicer(SCase.dataDicomPath);
  end

  % The 3D volume must be normalised for the slices not to be distorted
  rotatorCuffSlicer.setMinimalResolutionDividedBy(2);
  rotatorCuffSlicer.normaliseVolume;
  rotatorCuffSlicer.slice(scapula.coordSys.origin,scapula.coordSys.ML);

  % The following lines set the slice to display a 300x300 mm area around the
  % scapula with the "Y" shape of the scapula oriented up
  rotatorCuffSlicer.orientSliceUpwardVector(scapula.coordSys.origin,scapula.coordSys.IS);
  rotatorCuffSlicer.crop(scapula.coordSys.origin-100*scapula.coordSys.IS,300,300);

  % Rotator cuff segmentation code requires 1024x1024 pixels images
  if size(rotatorCuffSlicer.sliced,1) > size(rotatorCuffSlicer.sliced,2)
    rotatorCuffSlicer.resize([1024 nan]);
  else
    rotatorCuffSlicer.resize([nan 1024]);
  end
  rotatorCuffSlicer.addEmptyBackgroundToSlice(1024,1024);


  rawSlice = rotatorCuffSlicer.sliced;

  rotatorCuffSlicer.rescaleSliceToUint8;
  slice8Bit = rotatorCuffSlicer.sliced;
  rotatorCuffSlicer.sliced = rawSlice;

  rotatorCuffSlicer.rescaleSliceToInt16;
  slice16Bit = rotatorCuffSlicer.sliced;

  % Extract the pixels coordinates in space
  slicePixelIndices = arrayfun(@(x,y,z) [max(1,round(x)),max(1,round(y)),max(1,round(z))],...
    rotatorCuffSlicer.slicedX,...
    rotatorCuffSlicer.slicedY,...
    rotatorCuffSlicer.slicedZ,...
    "UniformOutput", false);
  % Pixels of points outside the volume are arbitrarily set to be [1 1 1]. The final
  % coordinates will be set to nan. 
  slicePixelIndices(not(rotatorCuffSlicer.pointsInVolume)) = {[1 1 1]};
  slicePixelCoordinates = cellfun(@(pointIndices) rotatorCuffSlicer.getPointCoordinates(pointIndices),...
    slicePixelIndices, "UniformOutput", false);
  slicePixelCoordinates(not(rotatorCuffSlicer.pointsInVolume)) = {nan};

  % Flip the slices for right and left shoulders to give similar slices
  if obj.shoulder.side=='L'
    slice8Bit=rot90(flipud(slice8Bit),2);
    slice16Bit=rot90(flipud(slice16Bit),2);
    slicePixelCoordinates=rot90(flipud(slicePixelCoordinates),2);
  end

  slicePixelSpacings = rotatorCuffSlicer.slicedPixelSpacings;
end



function saveImages(obj,imageForSegmentation,imageForMeasurements,sliceName)
  imwrite(imageForSegmentation,...
    fullfile(obj.slicesDataPath,sliceName + "_ForSegmentation.png"));
  save(fullfile(obj.slicesDataPath,sliceName + "_ForMeasurements.mat"),...
    'imageForMeasurements');
end



function saveImagesPixelSpacings(obj,imagesPixelSpacings,sliceName)
  save(fullfile(obj.slicesDataPath,sliceName + "_PixelSpacings.mat"),...
    'imagesPixelSpacings');
end


function saveImagesPixelCoordinates(obj,imagesPixelCoordinates,sliceName)
  save(fullfile(obj.slicesDataPath,sliceName + "_PixelCoordinates.mat"),...
    'imagesPixelCoordinates');
end