Merge pull request scilus#1057 from AlexVCaron/feat/precision

[ENH] Add precision to FRF scripts

scilpy/io/utils.py

@@ -22,6 +22,8 @@
 from scilpy.utils.spatial import RAS_AXES_NAMES
 
 
+FLOATING_POINTS_PRECISION = 12
+
 eddy_options = ["mb", "mb_offs", "slspec", "mporder", "s2v_lambda", "field",
                 "field_mat", "flm", "slm", "fwhm", "niter", "s2v_niter",
                 "cnr_maps", "residuals", "fep", "interp", "s2v_interp",
@@ -276,6 +278,14 @@ def add_skip_b0_check_arg(parser, will_overwrite_with_min,
         '--skip_b0_check', action='store_true', help=msg)
 
 
+def add_precision_arg(parser):
+    parser.add_argument('--precision', type=ranged_type(int, 1),
+                        default=FLOATING_POINTS_PRECISION,
+                        help='Precision for floating point values. Numbers '
+                             'are rounded up to \nthe number of decimals '
+                             'provided. [Default: %(default)s]')
+
+
 def add_verbose_arg(parser):
     parser.add_argument('-v', default="WARNING", const='INFO', nargs='?',
                         choices=['DEBUG', 'INFO', 'WARNING'], dest='verbose',

scripts/scil_frf_mean.py

@@ -18,6 +18,7 @@
 import numpy as np
 
 from scilpy.io.utils import (add_overwrite_arg,
+                             add_precision_arg,
                              assert_inputs_exist,
                              add_verbose_arg,
                              assert_outputs_exist)
@@ -33,6 +34,7 @@ def _build_arg_parser():
     p.add_argument('mean_frf', metavar='file',
                    help='Path of the output mean FRF file.')
 
+    add_precision_arg(p)
     add_verbose_arg(p)
     add_overwrite_arg(p)
 
@@ -66,7 +68,7 @@ def main():
 
     final_frf = np.mean(all_frfs, axis=0)
 
-    np.savetxt(args.mean_frf, final_frf)
+    np.savetxt(args.mean_frf, final_frf, fmt=f"%.{args.precision}f")
 
 
 if __name__ == "__main__":

scripts/scil_frf_memsmt.py

@@ -50,7 +50,8 @@
 from scilpy.image.utils import extract_affine
 from scilpy.io.btensor import generate_btensor_input
 from scilpy.io.image import get_data_as_mask
-from scilpy.io.utils import (add_overwrite_arg, add_verbose_arg,
+from scilpy.io.utils import (add_overwrite_arg, add_precision_arg,
+                             add_verbose_arg,
                              assert_inputs_exist, assert_outputs_exist,
                              assert_roi_radii_format, add_skip_b0_check_arg,
                              add_tolerance_arg,
@@ -169,6 +170,7 @@ def _build_arg_parser():
                    help='Path to the output CSF frf mask file, the voxels '
                         'used to compute the CSF frf.')
 
+    add_precision_arg(p)
     add_verbose_arg(p)
     add_overwrite_arg(p)
 
@@ -264,7 +266,7 @@ def main():
     frf_out = [args.out_wm_frf, args.out_gm_frf, args.out_csf_frf]
 
     for frf, response in zip(frf_out, responses):
-        np.savetxt(frf, response)
+        np.savetxt(frf, response, fmt=f"%.{args.precision}f")
 
 
 if __name__ == "__main__":

scripts/scil_frf_msmt.py

@@ -36,7 +36,8 @@
 from scilpy.dwi.utils import extract_dwi_shell
 from scilpy.gradients.bvec_bval_tools import check_b0_threshold
 from scilpy.io.image import get_data_as_mask
-from scilpy.io.utils import (add_overwrite_arg, add_skip_b0_check_arg,
+from scilpy.io.utils import (add_overwrite_arg, add_precision_arg,
+                             add_skip_b0_check_arg,
                              add_verbose_arg, assert_inputs_exist,
                              assert_outputs_exist, assert_roi_radii_format,
                              assert_headers_compatible)
@@ -136,6 +137,7 @@ def _build_arg_parser():
                    help='Path to the output CSF frf mask file, the voxels '
                         'used to compute the CSF frf.')
 
+    add_precision_arg(p)
     add_verbose_arg(p)
     add_overwrite_arg(p)
 
@@ -216,7 +218,7 @@ def main():
     frf_out = [args.out_wm_frf, args.out_gm_frf, args.out_csf_frf]
 
     for frf, response in zip(frf_out, responses):
-        np.savetxt(frf, response)
+        np.savetxt(frf, response, fmt=f"%.{args.precision}f")
 
 
 if __name__ == "__main__":

scripts/scil_frf_set_diffusivities.py

@@ -17,6 +17,7 @@
 import numpy as np
 
 from scilpy.io.utils import (add_overwrite_arg,
+                             add_precision_arg,
                              assert_inputs_exist,
                              add_verbose_arg,
                              assert_outputs_exist)
@@ -43,6 +44,7 @@ def _build_arg_parser():
                         'evaluated without the x 10**-4 factor. [%(default)s].'
                    )
 
+    add_precision_arg(p)
     add_verbose_arg(p)
     add_overwrite_arg(p)
 
@@ -59,7 +61,7 @@ def main():
 
     frf_file = np.loadtxt(args.frf_file)
     response = replace_frf(frf_file, args.new_frf, args.no_factor)
-    np.savetxt(args.output_frf_file, response)
+    np.savetxt(args.output_frf_file, response, fmt=f"%.{args.precision}f")
 
 
 if __name__ == "__main__":

scripts/scil_frf_ssst.py

@@ -20,6 +20,7 @@
 from scilpy.gradients.bvec_bval_tools import check_b0_threshold
 from scilpy.io.image import get_data_as_mask
 from scilpy.io.utils import (add_b0_thresh_arg, add_overwrite_arg,
+                             add_precision_arg,
                              add_skip_b0_check_arg, add_verbose_arg,
                              assert_inputs_exist, assert_outputs_exist,
                              assert_roi_radii_format,
@@ -80,6 +81,7 @@ def _build_arg_parser():
 
     add_b0_thresh_arg(p)
     add_skip_b0_check_arg(p, will_overwrite_with_min=True)
+    add_precision_arg(p)
     add_verbose_arg(p)
     add_overwrite_arg(p)
 
@@ -117,7 +119,7 @@ def main():
         min_fa_thresh=args.min_fa_thresh, min_nvox=args.min_nvox,
         roi_radii=roi_radii, roi_center=args.roi_center)
 
-    np.savetxt(args.frf_file, full_response)
+    np.savetxt(args.frf_file, full_response, fmt=f"%.{args.precision}f")
 
 
 if __name__ == "__main__":

scripts/tests/test_frf_mean.py

@@ -32,6 +32,23 @@ def test_execution_processing_msmt(script_runner, monkeypatch):
     assert ret.success
 
 
+def test_outputs_precision(script_runner, monkeypatch):
+    monkeypatch.chdir(os.path.expanduser(tmp_dir.name))
+    in_frf = os.path.join(SCILPY_HOME, 'commit_amico', 'wm_frf.txt')
+    ret = script_runner.run('scil_frf_mean.py', in_frf, in_frf, 'mfrfp.txt',
+                            '--precision', '4')
+    assert ret.success
+
+    expected = [
+        "0.0016 0.0004 0.0004 3076.7249",
+        "0.0012 0.0003 0.0003 3076.7249",
+        "0.0009 0.0003 0.0003 3076.7249"
+    ]
+    with open('mfrfp.txt', 'r') as result:
+            for i, line in enumerate(result.readlines()):
+                assert line.strip("\n") == expected[i]          
+
+
 def test_execution_processing_bad_input(script_runner, monkeypatch):
     monkeypatch.chdir(os.path.expanduser(tmp_dir.name))
     in_wm_frf = os.path.join(SCILPY_HOME, 'commit_amico', 'wm_frf.txt')

scripts/tests/test_frf_memsmt.py

@@ -130,6 +130,42 @@ def test_inputs_check(script_runner, monkeypatch):
     assert (not ret.success)
 
 
+def test_outputs_precision(script_runner, monkeypatch):
+    monkeypatch.chdir(os.path.expanduser(tmp_dir.name))
+    in_dwi_lin = os.path.join(SCILPY_HOME, 'btensor_testdata',
+                              'dwi_linear.nii.gz')
+    in_bval_lin = os.path.join(SCILPY_HOME, 'btensor_testdata',
+                               'linear.bvals')
+    in_bvec_lin = os.path.join(SCILPY_HOME, 'btensor_testdata',
+                               'linear.bvecs')
+    in_dwi_plan = os.path.join(SCILPY_HOME, 'btensor_testdata',
+                               'dwi_planar.nii.gz')
+    in_bval_plan = os.path.join(SCILPY_HOME, 'btensor_testdata',
+                                'planar.bvals')
+    in_bvec_plan = os.path.join(SCILPY_HOME, 'btensor_testdata',
+                                'planar.bvecs')
+    in_dwi_sph = os.path.join(SCILPY_HOME, 'btensor_testdata',
+                              'dwi_spherical.nii.gz')
+    in_bval_sph = os.path.join(SCILPY_HOME, 'btensor_testdata',
+                               'spherical.bvals')
+    in_bvec_sph = os.path.join(SCILPY_HOME, 'btensor_testdata',
+                               'spherical.bvecs')
+    ret = script_runner.run('scil_frf_memsmt.py', 'wm_frf.txt',
+                            'gm_frf.txt', 'csf_frf.txt', '--in_dwis',
+                            in_dwi_lin, in_dwi_plan, in_dwi_sph, '--in_bvals',
+                            in_bval_lin, in_bval_plan, in_bval_sph,
+                            '--in_bvecs', in_bvec_lin, in_bvec_plan,
+                            in_bvec_sph, '--in_bdeltas', '1', '-0.5', '0',
+                            '--min_nvox', '1', '--precision', '4', '-f')
+
+    assert ret.success
+
+    for frf_file in ['wm_frf.txt', 'gm_frf.txt', 'csf_frf.txt']:
+        with open(frf_file, "r") as f:
+            for item in f.readline().strip("\n").split(" "):
+                assert len(item.split(".")[1]) == 4
+
+
 def test_execution_processing(script_runner, monkeypatch):
     monkeypatch.chdir(os.path.expanduser(tmp_dir.name))
     in_dwi_lin = os.path.join(SCILPY_HOME, 'btensor_testdata',

scripts/tests/test_frf_msmt.py

@@ -75,6 +75,27 @@ def test_roi_radii_shape_parameter2(script_runner, monkeypatch):
     assert (not ret.success)
 
 
+def test_outputs_precision(script_runner, monkeypatch):
+    monkeypatch.chdir(os.path.expanduser(tmp_dir.name))
+    in_dwi = os.path.join(SCILPY_HOME, 'commit_amico',
+                          'dwi.nii.gz')
+    in_bval = os.path.join(SCILPY_HOME, 'commit_amico',
+                           'dwi.bval')
+    in_bvec = os.path.join(SCILPY_HOME, 'commit_amico',
+                           'dwi.bvec')
+    mask = os.path.join(SCILPY_HOME, 'commit_amico', 'mask.nii.gz')
+    ret = script_runner.run('scil_frf_msmt.py', in_dwi,
+                            in_bval, in_bvec, 'wm_frf.txt', 'gm_frf.txt',
+                            'csf_frf.txt', '--mask', mask, '--min_nvox', '20',
+                            '--precision', '4', '-f')
+    assert ret.success
+
+    for frf_file in ['wm_frf.txt', 'gm_frf.txt', 'csf_frf.txt']:
+        with open(frf_file, "r") as f:
+            for item in f.readline().strip("\n").split(" "):
+                assert len(item.split(".")[1]) == 4
+
+
 def test_execution_processing(script_runner, monkeypatch):
     monkeypatch.chdir(os.path.expanduser(tmp_dir.name))
     in_dwi = os.path.join(SCILPY_HOME, 'commit_amico',

scripts/tests/test_frf_set_diffusivities.py

@@ -3,6 +3,7 @@
 
 import os
 import tempfile
+import numpy as np
 
 from scilpy import SCILPY_HOME
 from scilpy.io.fetcher import fetch_data, get_testing_files_dict
@@ -34,6 +35,24 @@ def test_execution_processing_msmt(script_runner, monkeypatch):
     assert ret.success
 
 
+def test_outputs_precision(script_runner, monkeypatch):
+    monkeypatch.chdir(os.path.expanduser(tmp_dir.name))
+    in_frf = os.path.join(SCILPY_HOME, 'commit_amico', 'wm_frf.txt')
+    ret = script_runner.run('scil_frf_set_diffusivities.py', in_frf,
+                            '15,4,4,13,4,4,12,5,5', 'new_frf.txt',
+                            '--precision', '4', '-f')
+    assert ret.success
+
+    expected = [
+        "0.0015 0.0004 0.0004 3076.7249",
+        "0.0013 0.0004 0.0004 3076.7249",
+        "0.0012 0.0005 0.0005 3076.7249"
+    ]
+    with open('new_frf.txt', 'r') as result:
+        for i, line in enumerate(result.readlines()):
+            assert line.strip("\n") == expected[i] 
+
+
 def test_execution_processing__wrong_input(script_runner, monkeypatch):
     monkeypatch.chdir(os.path.expanduser(tmp_dir.name))
     in_frf = os.path.join(SCILPY_HOME, 'commit_amico', 'wm_frf.txt')

scripts/tests/test_frf_ssst.py

@@ -60,6 +60,18 @@ def test_roi_radii_shape_parameter(script_runner, monkeypatch):
     assert (not ret.success)
 
 
+def test_outputs_precision(script_runner, monkeypatch):
+    monkeypatch.chdir(os.path.expanduser(tmp_dir.name))
+    ret = script_runner.run('scil_frf_ssst.py', in_dwi,
+                            in_bval, in_bvec, 'frf.txt',
+                            '--precision', '4', '-f')
+    assert ret.success
+
+    with open("frf.txt", "r") as f:
+        for item in f.readline().strip("\n").split(" "):
+            assert len(item.split(".")[1]) == 4
+
+
 def test_execution_processing(script_runner, monkeypatch):
     monkeypatch.chdir(os.path.expanduser(tmp_dir.name))
     ret = script_runner.run('scil_frf_ssst.py', in_dwi,