ENH: Support 5D files in NRRD IO

Several modern imaging devices and processing algorithms now creates 5D data sets, including:
- time sequence of 3D displacement fields (result of deformable registration of 4D images)
- time sequence of 3D velocity fields (acquired by 4D MRI flow images)
- time sequence of 3D color (RGB or RGBA) images (result of colorized volumes provided by AI segmentation)

NRRD IO can now read and write all these data types, and in general any 5D arrays (2 or 3 spatial dimensions + scalar/vector pixel type + optional list axis).

1. Writing

Previously the NRRD writer handled one component axis, and all the other dimensions were space.
This commit allows specifying one of the axes as list (by setting "NRRD_kinds[i]=list" metadata key).

Example:

File header before this commit:

  dimension: 5
  space dimension: 4
  sizes: 4 320 320 1 5
  space directions: none (1,0,0,0) (0,1,0,0) (0,0,1,0) (0,0,0,1)
  kinds: RGBA-color domain domain domain domain
  encoding: gzip
  space origin: (0,0,0,1)

File header after this commit:

  dimension: 5
  space: left-posterior-superior
  sizes: 4 320 320 1 5
  space directions: none (-1,0,0) (0,-1,0) (0,0,1) none
  kinds: vector domain domain domain list
  encoding: gzip
  space origin: (0,0,0)

2. Reading

NRRD IO now allows having more than one non-spatial axis.
The non-spatial axes of the NRRD file can be used in the ITK image as pixel components or as additional image dimensions.
The default strategy (UseAnyRangeAxisAsPixel) implements the same behavior as before (first non-spatial NRRD axis always becomes pixel component).
The new UseNonListRangeAxisAsPixel option allows "list" dimension (e.g., time points) to always mapped to additional image dimension.
This new option simplifies application implementation, as time sequences are loaded the same way for scalar and vector volumes.

3. Add support for axis unit in NRRD IO; Save label and unit for range axis

NRRD IO did not support the "units" metadata for axes, now it is added both for reading and writing.
Per-axis metadata (kinds, labels, units) was only saved for spatial axes, now it is added for all axes.
For pixel components, "pixel" is used instead of axis ("NRRD_labels[pixel]" and "NRRD_units[pixel]" fields).

Added many tests for reading/writing 3D, 4D, 5D images in NRRD format.

Co-authored-by: Andras Lasso <[email protected]>

Author: 2 people

Modules/IO/NRRD/include/itkNrrdImageIO.h

@@ -27,6 +27,29 @@ struct NrrdEncoding_t;
 
 namespace itk
 {
+/** \class NrrdImageIOEnums
+ * \brief Contains all enum classes used by NrrdImageIO class.
+ * \ingroup ITKIONRRD
+ */
+class NrrdImageIOEnums
+{
+public:
+  /** \class AxesReorder
+   * \ingroup ITKIONRRD
+   * Enum type that specifies how to reorder axes in NRRD file when reading it into an ITK image.  */
+  enum class AxesReorder : uint8_t
+  {
+    UseAnyRangeAxisAsPixel = 0,
+    UseNonListRangeAxisAsPixel = 1,
+    UseScalarPixel = 2
+  };
+};
+// Define how to print enumeration
+extern ITKIONRRD_EXPORT std::ostream &
+                        operator<<(std::ostream & out, const NrrdImageIOEnums::AxesReorder value);
+
+using AxesReorderEnum = NrrdImageIOEnums::AxesReorder;
+
 /**
  * \class NrrdImageIO
  *
@@ -109,6 +132,73 @@ class ITKIONRRD_EXPORT NrrdImageIO : public ImageIOBase
   void
   Write(const void * buffer) override;
 
+  /** Set/Get the strategy for axis reordering.
+   * ITK image can be an Image object that can store scalar, RGB, CovariantVector, etc. pixels in an arbitrarily high
+   * dimensional array. If the number of components of the pixel has to be dynamic then the ITK image has to be a
+   * VectorImage object that can store multi-component pixels in an arbitrarily high dimensional array.
+   *
+   * NRRD file can contain domain (space, time) axes and range (list, RGB, covariant-vector, etc.) axes in an arbitrary
+   * order. Therefore, axes of the NRRD file may have to be reordered to be able to read the NRRD file into an ITK
+   * image. In general, NRRD domain axes are mapped to ITK image axes then all range axes are used as additional image
+   * dimensions. The only exception is that it makes sense to pick one a range axis and use that as pixel component.
+   *
+   * In early ITK versions, any kind of NRRD range axis was used as pixel component by default (this is the
+   * UseAnyRangeAxisAsPixel strategy). However, this resulted in more complicated application code (developers need to
+   * instantiate an Image or VectorImage depending on what kind of axes are in the NRRD file) and potential unnecessary
+   * reordering of axes (e.g., sequence of 3D images was always converted into a single multi-component 3D image).
+   *
+   * Therefore additional strategies were introduced to allow more control over how NRRD axes are mapped to ITK image
+   * axes and pixel component. UseNonListRangeAxisAsPixel to ensure that the image can be always read into an Image
+   * object. UseScalarPixel to ensure that the image can be always read into an Image object with a scalar pixel type.
+   *
+   * \li UseAnyRangeAxisAsPixel Use the first a range axis (i.e., non-domain axis) of the NRRD file as pixel component.
+   *     This is the default (for backward compatibility).
+   *     If the NRRD file contains a non-list type (e.g., RGB, point, covariant-vector) range axis then the first one
+   * will be used as pixel component; the image will be read into an Image object. If the NRRD file only contains
+   * list-type range axis then it will be used as pixel component; the image will be read into a VectorImage.
+   * \li UseNonListRangeAxisAsPixel Use the first non-list range axis (i.e., non-domain axis) of the NRRD file as pixel
+   * component. The image is always read into an Image object, with potentially a non-scalar pixel type.
+   * \li UseScalarPixel Pixel type is scalar. All range (i.e., non-domain) axes of the NRRD file are added as additional
+   * domain axes. The image is always read into an Image object with a scalar pixel type.
+   *
+   * Examples:
+   *
+   *   NRRD file: list domain domain domain RGB-Color (or: domain domain domain list RGB-Color)
+   *   - UseAnyRangeAxisAsPixel:     Image<RGBPixel<double>,4> with axes: domain domain domain list
+   *   - UseNonListRangeAxisAsPixel: Image<RGBPixel<double>,4> with axes: domain domain domain list
+   *   - UseScalarPixel:             Image<double,5>           with axes: domain domain domain list RGB-Color
+   *
+   *   NRRD file: list domain domain domain (or: domain domain domain list)
+   *   - UseAnyRangeAxisAsPixel:     VectorImage<double,3> with axes: domain domain domain
+   *   - UseNonListRangeAxisAsPixel: Image<double,4>       with axes: domain domain domain list
+   *   - UseScalarPixel:             Image<double,4>       with axes: domain domain domain list
+   */
+  itkSetMacro(AxesReorder, AxesReorderEnum);
+  itkGetConstMacro(AxesReorder, AxesReorderEnum);
+
+  /** Use the first a range axis (i.e., non-domain axis) of the NRRD file as pixel component. */
+  void
+  SetAxesReorderToUseAnyRangeAxisAsPixel()
+  {
+    this->SetAxesReorder(AxesReorderEnum::UseAnyRangeAxisAsPixel);
+  }
+
+  /** UseNonListRangeAxisAsPixel Use the first non-list range axis (i.e., non-domain axis) of the NRRD file as pixel
+   * component. */
+  void
+  SetAxesReorderToUseNonListRangeAxisAsPixel()
+  {
+    this->SetAxesReorder(AxesReorderEnum::UseNonListRangeAxisAsPixel);
+  }
+
+  /** UseScalarPixel Pixel type is scalar. All range (i.e., non-domain) axes of the NRRD file are added as additional
+   * domain axes. */
+  void
+  SetAxesReorderToUseScalarPixel()
+  {
+    this->SetAxesReorder(AxesReorderEnum::UseScalarPixel);
+  }
+
 protected:
   NrrdImageIO();
   ~NrrdImageIO() override;
@@ -127,6 +217,8 @@ class ITKIONRRD_EXPORT NrrdImageIO : public ImageIOBase
   NrrdToITKComponentType(const int) const;
 
   const NrrdEncoding_t * m_NrrdCompressionEncoding{ nullptr };
+
+  AxesReorderEnum m_AxesReorder{ AxesReorderEnum::UseAnyRangeAxisAsPixel };
 };
 } // end namespace itk
 

Modules/IO/NRRD/itk-module.cmake

@@ -13,6 +13,8 @@ itk_module(
   PRIVATE_DEPENDS
     ITKNrrdIO
   TEST_DEPENDS
+    ITKImageCompose
+    ITKImageGrid
     ITKTestKernel
   FACTORY_NAMES
     ImageIO::Nrrd