5TTgen MSMT

docs/reference/commands/5ttgen.rst

@@ -15,7 +15,7 @@ Usage
 
     5ttgen algorithm [ options ] ...
 
--  *algorithm*: Select the algorithm to be used; additional details and options become available once an algorithm is nominated. Options are: freesurfer, fsl, gif, hsvs
+-  *algorithm*: Select the algorithm to be used; additional details and options become available once an algorithm is nominated. Options are: freesurfer, fsl, gif, hsvs, msmt
 
 Description
 -----------
@@ -476,3 +476,101 @@ See the Mozilla Public License v. 2.0 for more details.
 
 For more details, see http://www.mrtrix.org/.
 
+.. _5ttgen_msmt:
+
+5ttgen msmt
+===========
+
+Synopsis
+--------
+
+Generate a 5TT image from ODF images
+
+Usage
+-----
+
+::
+
+    5ttgen msmt odf_wm odf_gm odf_csf output [ options ]
+
+-  *odf_wm*: The input white matter ODF
+-  *odf_gm*: The input grey matter ODF
+-  *odf_csf*: The input cerebrospinal fluid ODF
+-  *output*: The output 5TT image
+
+Description
+-----------
+
+If the user does not manually provide a brain mask using the -mask option, the command will automatically determine a mask for the output 5TT image; this is necessary to conform to the expectations of the 5TT format. This mask will be computed based on where the sum of the ODFs exceeds 50% of what is expected from a voxel that contains a DWI signal that exactly matches one of the response functions.
+
+Options
+-------
+
+Options specific to the 'msmt' algorithm
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+- **-mask image** An input binary brain mask image
+
+Options common to all 5ttgen algorithms
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+- **-nocrop** Do NOT crop the resulting 5TT image to reduce its size (keep the same dimensions as the input image)
+
+- **-sgm_amyg_hipp** Represent the amygdalae and hippocampi as sub-cortical grey matter in the 5TT image
+
+Additional standard options for Python scripts
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+- **-nocleanup** do not delete intermediate files during script execution, and do not delete scratch directory at script completion.
+
+- **-scratch /path/to/scratch/** manually specify an existing directory in which to generate the scratch directory.
+
+- **-continue ScratchDir LastFile** continue the script from a previous execution; must provide the scratch directory path, and the name of the last successfully-generated file.
+
+Standard options
+^^^^^^^^^^^^^^^^
+
+- **-info** display information messages.
+
+- **-quiet** do not display information messages or progress status. Alternatively, this can be achieved by setting the MRTRIX_QUIET environment variable to a non-empty string.
+
+- **-debug** display debugging messages.
+
+- **-force** force overwrite of output files.
+
+- **-nthreads number** use this number of threads in multi-threaded applications (set to 0 to disable multi-threading).
+
+- **-config key value**  *(multiple uses permitted)* temporarily set the value of an MRtrix config file entry.
+
+- **-help** display this information page and exit.
+
+- **-version** display version information and exit.
+
+References
+^^^^^^^^^^
+
+* Smith, R. E.; Tournier, J.-D.; Calamante, F. & Connelly, A. Anatomically-constrained tractography: Improved diffusion MRI streamlines tractography through effective use of anatomical information. NeuroImage, 2012, 62, 1924-1938
+
+Tournier, J.-D.; Smith, R. E.; Raffelt, D.; Tabbara, R.; Dhollander, T.; Pietsch, M.; Christiaens, D.; Jeurissen, B.; Yeh, C.-H. & Connelly, A. MRtrix3: A fast, flexible and open software framework for medical image processing and visualisation. NeuroImage, 2019, 202, 116137
+
+--------------
+
+
+
+**Author:** Arkiev D'Souza (arkiev.dsouza@sydney.edu.au) & Robert E. Smith (robert.smith@florey.edu.au)
+
+**Copyright:** Copyright (c) 2008-2024 the MRtrix3 contributors.
+
+This Source Code Form is subject to the terms of the Mozilla Public
+License, v. 2.0. If a copy of the MPL was not distributed with this
+file, You can obtain one at http://mozilla.org/MPL/2.0/.
+
+Covered Software is provided under this License on an "as is"
+basis, without warranty of any kind, either expressed, implied, or
+statutory, including, without limitation, warranties that the
+Covered Software is free of defects, merchantable, fit for a
+particular purpose or non-infringing.
+See the Mozilla Public License v. 2.0 for more details.
+
+For more details, see http://www.mrtrix.org/.
+

python/mrtrix3/commands/5ttgen/__init__.py

@@ -14,4 +14,4 @@
 # For more details, see http://www.mrtrix.org/.
 
 # pylint: disable=unused-variable
-ALGORITHMS = ['freesurfer', 'fsl', 'gif', 'hsvs']
+ALGORITHMS = ['freesurfer', 'fsl', 'gif', 'hsvs', 'msmt']

python/mrtrix3/commands/5ttgen/msmt.py

@@ -0,0 +1,127 @@
+# Copyright (c) 2008-2024 the MRtrix3 contributors.
+#
+# This Source Code Form is subject to the terms of the Mozilla Public
+# License, v. 2.0. If a copy of the MPL was not distributed with this
+# file, You can obtain one at http://mozilla.org/MPL/2.0/.
+#
+# Covered Software is provided under this License on an "as is"
+# basis, without warranty of any kind, either expressed, implied, or
+# statutory, including, without limitation, warranties that the
+# Covered Software is free of defects, merchantable, fit for a
+# particular purpose or non-infringing.
+# See the Mozilla Public License v. 2.0 for more details.
+#
+# For more details, see http://www.mrtrix.org/.
+
+from mrtrix3 import app, image, run
+
+
+def usage(base_parser, subparsers):  # pylint: disable=unused-variable
+  parser = subparsers.add_parser('msmt', parents=[base_parser])
+  parser.set_author('Arkiev D\'Souza (arkiev.dsouza@sydney.edu.au) & Robert E. Smith (robert.smith@florey.edu.au)')
+  parser.set_synopsis('Generate a 5TT image from ODF images')
+  parser.add_description('If the user does not manually provide a brain mask using the -mask option,'
+                         ' the command will automatically determine a mask for the output 5TT image;'
+                         ' this is necessary to conform to the expectations of the 5TT format.'
+                         ' This mask will be computed based on where the sum of the ODFs exceeds 50%'
+                         ' of what is expected from a voxel that contains a DWI signal'
+                         ' that exactly matches one of the response functions.')
+
+  parser.add_argument('odf_wm', type=app.Parser.ImageIn(), help='The input white matter ODF')
+  parser.add_argument('odf_gm', type=app.Parser.ImageIn(), help='The input grey matter ODF')
+  parser.add_argument('odf_csf', type=app.Parser.ImageIn(), help='The input cerebrospinal fluid ODF')
+  parser.add_argument('output', type=app.Parser.ImageOut(), help='The output 5TT image')
+
+  options = parser.add_argument_group('Options specific to the \'msmt\' algorithm')
+  options.add_argument('-mask', type=app.Parser.ImageIn(), help='An input binary brain mask image')
+
+
+def execute(): # pylint: disable=unused-variable
+
+  target_voxel_grid = image.Header(app.ARGS.odf_wm)
+
+  class Tissue:
+    def __init__(self, inpath, name, expect_anisotropic):
+      self.name = name
+      header = image.Header(inpath)
+      do_regrid = not image.match(header, target_voxel_grid, up_to_dim=3)
+      if do_regrid:
+        app.warn(f'ODF image "{inpath}" not defined on same voxel grid as WM ODF image;'
+                 ' will have to resample')
+      else:
+        app.debug(f'{name} ODF image already on target voxel grid')
+      is_anisotropic = len(header.size()) > 3 and header.size()[3] > 1
+      if is_anisotropic != expect_anisotropic:
+        app.warn(f'Received {"anisotropic" if is_anisotropic else "isotropic"} ODF for {name}'
+                 f' but expected {"anisotropic" if expect_anisotropic else "isotropic"};'
+                 ' check order of input ODF images if this was not intentional')
+      self.lzeropath = f'{name}_lzero.mif'
+      command = ['mrgrid', inpath, 'regrid', '-template', app.ARGS.odf_wm, '-', '|'] \
+                if do_regrid                                                         \
+                else []
+      if is_anisotropic:
+        command.extend(['mrconvert', '-' if command else inpath, '-coord', '3', '0', '-axes', '0,1,2', '-', '|'])
+      command.extend(['mrcalc', '-' if command else inpath, '0.0', '-max', self.lzeropath])
+      run.command(command)
+      self.normpath = None
+
+    def normalise(self, volsum_image):
+      assert self.normpath is None
+      self.normpath = f'{self.name}_fraction.mif'
+      run.command(f'mrcalc {self.lzeropath} {volsum_image} -divide {self.normpath}')
+
+  # End definitiion of class Tissues
+
+  tissues = [Tissue(app.ARGS.odf_wm,  'WM',  True),
+             Tissue(app.ARGS.odf_gm,  'GM',  False),
+             Tissue(app.ARGS.odf_csf, 'CSF', False)]
+
+  # Compute volume sum in each voxel
+  volsum_image = 'totalvol.mif'
+  run.command(['mrmath', [item.lzeropath for item in tissues], 'sum', volsum_image])
+  # Normalise to tissue fractions
+  for tissue in tissues:
+    tissue.normalise(volsum_image)
+
+  # Create empty volume for SGM and pathology
+  empty_volume = 'empty_vol.mif'
+  run.command(f'mrcalc {volsum_image} inf -gt {empty_volume}')
+
+  # Concatenate volumes
+  result_unmasked = '5TT_unmasked.mif'
+  ####################################################################################################################
+  # 5TT volumes:          Cortical GM      Sub-cortical GM          WM                   CSF             Pathology   #
+  ####################################################################################################################
+  run.command(f'mrcat {tissues[1].normpath} {empty_volume} {tissues[0].normpath} {tissues[2].normpath} {empty_volume}'
+              f' {result_unmasked} -datatype float32')
+
+  # Tidy image by masking
+  result_masked = '5TT_masked.mif'
+  # If mask is provided, use it; if not, create one
+  if app.ARGS.mask:
+    # Case 1: Using provided brainmask
+    # Check if brainmask is on same voxel grid; if not, regrid
+    mask_header = image.Header(app.ARGS.mask)
+    if image.match(target_voxel_grid, mask_header, up_to_dim=3):
+      app.debug('Mask matches ODF voxel grid; no regridding of mask required')
+      run.command(f'mrcalc {result_unmasked} {app.ARGS.mask} -mult {result_masked}')
+    else:
+      app.warn('Mask has different voxel grid to WM ODF image; regridding')
+      run.command(f'mrgrid {app.ARGS.mask} -template {empty_volume} regrid - |'
+                  ' mrcalc - 0.5 -gt - |'
+                  f' mrcalc - {result_unmasked} -mult {result_masked}')
+  else:
+    # Case 2: Generate a mask that contains only those voxels where the sum of ODF l=0 terms exceeds 0.5/sqrt(4pi);
+    #   then select the largest component and fill any holes.
+    # 0.5/sqrt(4pi) = 0.1410473959
+    app.console('No 5TT mask provided; generating from input ODFs')
+    run.command(f'mrcalc {volsum_image} 0.1410473959 -gt - |'
+                ' maskfilter - clean - |'
+                ' maskfilter - fill - |'
+                f' mrcalc - {result_unmasked} -mult {result_masked}')
+
+  run.command(['mrconvert', result_masked, app.ARGS.output],
+    force=app.FORCE_OVERWRITE,
+    preserve_pipes=True)
+
+  return result_masked

testing/CMakeLists.txt

@@ -27,7 +27,7 @@ ExternalProject_Add(BinariesTestData
 ExternalProject_Add(ScriptsTestData
     PREFIX ${CMAKE_CURRENT_BINARY_DIR}/scripts_data
     GIT_REPOSITORY ${mrtrix_scripts_data_url}
-    GIT_TAG 76f47633cd0a37e901c42320f4540ecaffd51367
+    GIT_TAG ed1d3c630db9eafd6ce2808360e9e99529abfabc
     GIT_PROGRESS TRUE
     CONFIGURE_COMMAND ""
     BUILD_COMMAND ""
@@ -47,4 +47,4 @@ add_subdirectory(binaries)
 add_subdirectory(lib)
 add_subdirectory(scripts)
 add_subdirectory(tools)
-add_subdirectory(unit_tests)
+add_subdirectory(unit_tests)

testing/scripts/CMakeLists.txt

@@ -47,6 +47,8 @@ add_bash_script_test(5ttgen/hsvs_piping)
 add_bash_script_test(5ttgen/hsvs_template)
 add_bash_script_test(5ttgen/hsvs_whitespace)
 add_bash_script_test(5ttgen/hsvs_whitestem)
+add_bash_script_test(5ttgen/msmt_default)
+add_bash_script_test(5ttgen/msmt_masked)
 
 add_bash_script_test(dwi2mask/3dautomask_default)
 add_bash_script_test(dwi2mask/3dautomask_options)

testing/scripts/tests/5ttgen/msmt_default

@@ -0,0 +1,19 @@
+#!/bin/bash
+# Verify "5ttgen msmt" default operation
+# Output is compared to a prior output generateed by the command
+
+# First we need to get some multi-tissue ODFs
+dwi2response dhollander BIDS/sub-01/dwi/sub-01_dwi.nii.gz tmp_wm.txt tmp_gm.txt tmp_csf.txt \
+-fslgrad BIDS/sub-01/dwi/sub-01_dwi.bvec BIDS/sub-01/dwi/sub-01_dwi.bval \
+-mask BIDS/sub-01/dwi/sub-01_brainmask.nii.gz \
+-force
+dwi2fod msmt_csd BIDS/sub-01/dwi/sub-01_dwi.nii.gz \
+-fslgrad BIDS/sub-01/dwi/sub-01_dwi.bvec BIDS/sub-01/dwi/sub-01_dwi.bval \
+tmp_wm.txt tmp_wm.mif \
+tmp_gm.txt tmp_gm.mif \
+tmp_csf.txt tmp_csf.mif \
+-force
+
+5ttgen msmt tmp_wm.mif tmp_gm.mif tmp_csf.mif tmp.mif -force
+
+testing_diff_image tmp.mif 5ttgen/msmt/default.mif.gz -abs 1e-5

testing/scripts/tests/5ttgen/msmt_masked

@@ -0,0 +1,21 @@
+#!/bin/bash
+# Verify "5ttgen msmt" operation where a user-specified mask is provided
+# Output is compared to a prior output generateed by the command
+
+# First we need to get some multi-tissue ODFs
+dwi2response dhollander BIDS/sub-01/dwi/sub-01_dwi.nii.gz tmp_wm.txt tmp_gm.txt tmp_csf.txt \
+-fslgrad BIDS/sub-01/dwi/sub-01_dwi.bvec BIDS/sub-01/dwi/sub-01_dwi.bval \
+-mask BIDS/sub-01/dwi/sub-01_brainmask.nii.gz \
+-force
+dwi2fod msmt_csd BIDS/sub-01/dwi/sub-01_dwi.nii.gz \
+-fslgrad BIDS/sub-01/dwi/sub-01_dwi.bvec BIDS/sub-01/dwi/sub-01_dwi.bval \
+tmp_wm.txt tmp_wm.mif \
+tmp_gm.txt tmp_gm.mif \
+tmp_csf.txt tmp_csf.mif \
+-force
+
+5ttgen msmt tmp_wm.mif tmp_gm.mif tmp_csf.mif tmp.mif \
+-mask BIDS/sub-01/dwi/sub-01_brainmask.nii.gz \
+-force
+
+testing_diff_image tmp.mif 5ttgen/msmt/masked.mif.gz -abs 1e-5