From 976f13e2c71e79b1977eccffd21b7848b1cd2117 Mon Sep 17 00:00:00 2001 From: Alexandre D'Astous Date: Tue, 3 Nov 2020 00:22:07 -0500 Subject: [PATCH] Add resampling of fieldmap, riro and static to target acquisition --- setup.py | 1 + shimmingtoolbox/cli/realtime_zshim.py | 40 +++++++++++++++++++++++++-- test/cli/test_cli_realtime_zshim.py | 6 +++- 3 files changed, 44 insertions(+), 3 deletions(-) diff --git a/setup.py b/setup.py index 202419ba..aba63801 100644 --- a/setup.py +++ b/setup.py @@ -39,6 +39,7 @@ "pytest~=4.6.3", "pytest-cov~=2.5.1", "sklearn~=0.0", + "nilearn~=0.6.2" ], extras_require={ 'docs': ["sphinx>=1.6", "sphinx_rtd_theme>=0.2.4"], diff --git a/shimmingtoolbox/cli/realtime_zshim.py b/shimmingtoolbox/cli/realtime_zshim.py index 2b0e8032..4550a570 100644 --- a/shimmingtoolbox/cli/realtime_zshim.py +++ b/shimmingtoolbox/cli/realtime_zshim.py @@ -7,6 +7,7 @@ import nibabel as nib import json from sklearn.linear_model import LinearRegression +from nilearn.image import resample_img # TODO: remove matplotlib and dirtesting import from matplotlib.figure import Figure from shimmingtoolbox import __dir_testing__ @@ -33,11 +34,13 @@ help="3D nifti file with voxels between 0 and 1 used to weight the spatial region to shim.") @click.option('-resp', 'fname_resp', type=click.Path(), help="Siemens respiratory file containing pressure data.") +@click.option('-anat', 'fname_anat', type=click.Path(), + help="Filename of the anatomical image to apply the correction.") # TODO: Remove json file as input @click.option('-json', 'fname_json', type=click.Path(), help="Filename of json corresponding BIDS sidecar.") @click.option("-verbose", is_flag=True, help="Be more verbose.") -def realtime_zshim(fname_coil, fname_fmap, fname_mask, fname_resp, fname_json, verbose=True): +def realtime_zshim(fname_coil, fname_fmap, fname_mask, fname_resp, fname_json, fname_anat, verbose=True): """ Args: @@ -62,7 +65,7 @@ def realtime_zshim(fname_coil, fname_fmap, fname_mask, fname_resp, fname_json, v # TODO: Error handling might move to API if fieldmap.ndim != 4: - raise RuntimeError('fmap must be 4d (x, y, z, t)') + raise RuntimeError("fmap must be 4d (x, y, z, t)") nx, ny, nz, nt = fieldmap.shape # Load mask @@ -72,6 +75,12 @@ def realtime_zshim(fname_coil, fname_fmap, fname_mask, fname_resp, fname_json, v else: mask = np.ones_like(fieldmap) + # Load anat + nii_anat = nib.load(fname_anat) + anat = nii_anat.get_fdata() + if anat.ndim != 3: + raise RuntimeError("Anatomical image must be in 3d") + # Shim using sequencer and optimizer n_coils = coil.shape[-1] currents = np.zeros([n_coils, nt]) @@ -118,6 +127,20 @@ def realtime_zshim(fname_coil, fname_fmap, fname_mask, fname_resp, fname_json, v riro[i_x, i_y, i_z] = reg.coef_ static[i_x, i_y, i_z] = reg.intercept_ + # Resample masked_fieldmaps, riro and static to target anatomical image + nii_masked_fmap = nib.Nifti1Image(masked_fieldmaps, nii_fmap.affine) + nii_riro = nib.Nifti1Image(riro, nii_fmap.affine) + nii_static = nib.Nifti1Image(static, nii_fmap.affine) + + target_affine = nii_anat.affine[:-1, :-1] + nii_resampled_fmap = resample_img(nii_masked_fmap, target_affine=target_affine, interpolation='nearest') + nii_resampled_riro = resample_img(nii_riro, target_affine=target_affine, interpolation='nearest') + nii_resampled_static = resample_img(nii_static, target_affine=target_affine, interpolation='nearest') + + nib.save(nii_resampled_fmap, os.path.join(__dir_shimmingtoolbox__, 'resampled_fmap.nii.gz')) + nib.save(nii_resampled_riro, os.path.join(__dir_shimmingtoolbox__, 'resampled_riro.nii.gz')) + nib.save(nii_resampled_static, os.path.join(__dir_shimmingtoolbox__, 'resampled_static.nii.gz')) + # ================ PLOTS ================ # Calculate masked shim for spherical harmonics plot @@ -218,6 +241,19 @@ def realtime_zshim(fname_coil, fname_fmap, fname_mask, fname_resp, fname_json, v fname_figure = os.path.join(__dir_shimmingtoolbox__, 'realtime_zshim_pmu_vs_B0.png') fig.savefig(fname_figure) + # Show anatomical image + fig = Figure(figsize=(10, 10)) + ax = fig.add_subplot(2, 1, 1) + im = ax.imshow(anat[:-1, :-1, 10]) + fig.colorbar(im) + ax.set_title("Anatomical image [:-1, :-1, 10]") + ax = fig.add_subplot(2, 1, 2) + im = ax.imshow(nii_resampled_fmap.get_fdata()[0, :-1, :-1, 0]) + fig.colorbar(im) + ax.set_title("Resampled fieldmap [0, :-1, :-1, 0]") + fname_figure = os.path.join(__dir_shimmingtoolbox__, 'reatime_zshime_anat.png') + fig.savefig(fname_figure) + return fname_figure # Debug diff --git a/test/cli/test_cli_realtime_zshim.py b/test/cli/test_cli_realtime_zshim.py index 1fdfc531..c8f79b15 100644 --- a/test/cli/test_cli_realtime_zshim.py +++ b/test/cli/test_cli_realtime_zshim.py @@ -58,8 +58,12 @@ def test_cli_realtime_zshim(): fname_json = os.path.join(__dir_testing__, 'realtime_zshimming_data', 'nifti', 'sub-example', 'fmap', 'sub-example_magnitude1.json') + # Path for mag anat image + fname_anat = os.path.join(__dir_testing__, 'realtime_zshimming_data', 'nifti', 'sub-example', 'anat', + 'sub-example_unshimmed_e1.nii.gz') + result = runner.invoke(realtime_zshim, ['-fmap', fname_fieldmap, '-coil', fname_coil, '-mask', fname_mask, - '-resp', fname_resp, '-json', fname_json], + '-resp', fname_resp, '-json', fname_json, '-anat', fname_anat], catch_exceptions=False) assert result.exit_code == 0