API Reference
The following section outlines the API of shimming-toolbox.
Fieldmapping
- shimmingtoolbox.prepare_fieldmap.complex_difference(phase1, phase2)
Calculates the complex difference between 2 phase arrays (phase2 - phase1)
- Parameters:
phase1 (numpy.ndarray) -- Array containing phase data in radians
phase2 (numpy.ndarray) -- Array containing phase data in radians. Must be the same shape as phase1.
- Returns:
The difference in phase between each voxels of phase2 and phase1 (phase2 - phase1)
- Return type:
numpy.ndarray
- shimmingtoolbox.prepare_fieldmap.correct_2pi_offset(unwrapped, mag, mask, validity_threshold)
- Removes 2*pi offsets from unwrapped for a time series. If there is no offset, it returns the same array. The
'correct' offset is assumed to be at time 0.
- Parameters:
unwrapped (numpy.ndarray) -- 4d array of the spatially unwrapped phase
mag (numpy.ndarray) -- 4d array containing the magnitude values of the phase
mask (numpy.ndarray) -- 4d mask of the unwrapped phase array
validity_threshold (float) -- Threshold to create a mask on each timepoints and assume as reliable phase data
- Returns:
4d array of the unwrapped phase corrected if there were n*2*pi offsets between time points
- Return type:
numpy.ndarray
- shimmingtoolbox.prepare_fieldmap.prepare_fieldmap(list_nii_phase, echo_times, mag, unwrapper='prelude', mask=None, threshold=0.05, gaussian_filter=False, sigma=1, fname_save_mask=None)
Creates fieldmap (in Hz) from phase images. This function accommodates multiple echoes (2 or more) and phase difference. This function also accommodates 4D phase inputs, where the 4th dimension represents the time, in case multiple field maps are acquired across time for the purpose of real-time shimming experiments.
- Parameters:
list_nii_phase (list) -- List of nib.Nifti1Image phase values. The array can be [x, y], [x, y, z] or [x, y, z, t]. The values must range from [-pi to pi].
echo_times (list) -- List of echo times in seconds for each echo. The number of echotimes must match the number of echoes. It input is a phasediff (1 phase), input 2 echotimes.
unwrapper (str) -- Unwrapper to use for phase unwrapping. Supported: prelude.
mag (numpy.ndarray) -- Array containing magnitude data relevant for
phase
input. Shape must match phase[echo].mask (numpy.ndarray) -- Mask for masking output fieldmap. Must match shape of phase[echo].
threshold -- Threshold for masking if no mask is provided. Allowed range: [0, 1] where all scaled values lower than the threshold are set to 0.
gaussian_filter (bool) -- Option of using a Gaussian filter to smooth the fieldmaps (boolean)
sigma (float) -- Standard deviation of gaussian filter.
fname_save_mask (str) -- Filename of the mask calculated by the unwrapper
- Returns
numpy.ndarray: Unwrapped fieldmap in Hz.
Wrapper to different unwrapping algorithms.
- shimmingtoolbox.unwrap.unwrap_phase.unwrap_phase(nii_phase_wrapped, unwrapper='prelude', mag=None, mask=None, threshold=None, fname_save_mask=None)
Calls different unwrapping algorithms according to the specified unwrapper parameter. The function also allows to call the different unwrappers with more flexibility regarding input shape.
- Parameters:
nii_phase_wrapped (nib.Nifti1Image) -- 2D, 3D or 4D radian values [-pi to pi] to perform phase unwrapping. Supported shapes: [x, y], [x, y, z] or [x, y, z, t].
unwrapper (str, optional) -- Unwrapper algorithm name. Possible values:
prelude
.mag (numpy.ndarray) -- 2D, 3D or 4D magnitude data corresponding to phase data. Shape must be the same as
phase
.mask (numpy.ndarray) -- numpy array of booleans with shape of
phase
to mask during phase unwrapping.threshold (float) -- Prelude parameter, see prelude for more detail.
fname_save_mask (str) -- Filename of the mask calculated by the unwrapper
- Returns:
Unwrapped phase image.
- Return type:
numpy.ndarray
Wrapper to FSL Prelude (https://fsl.fmrib.ox.ac.uk/fsl/fslwiki/FUGUE/Guide#PRELUDE_.28phase_unwrapping.29)
- shimmingtoolbox.unwrap.prelude.prelude(nii_wrapped_phase, mag=None, mask=None, threshold=None, is_unwrapping_in_2d=False, fname_save_mask=None)
wrapper to FSL prelude
This function enables phase unwrapping by calling FSL prelude on the command line. A mask can be provided to mask the phase image provided. 2D unwrapping can be turned off. The output path can be specified. The temporary niis can optionally be saved.
- Parameters:
nii_wrapped_phase (nib.Nifti1Image) -- 2D or 3D radian numpy array to perform phase unwrapping. (2 pi interval)
mag (numpy.ndarray) -- 2D or 3D magnitude numpy array corresponding to the phase array
mask (numpy.ndarray, optional) -- numpy array of booleans with shape of complex_array to mask during phase unwrapping
threshold -- Threshold value for automatic mask generation (Use either mask or threshold, not both)
is_unwrapping_in_2d (bool, optional) -- prelude parameter to unwrap slice by slice
fname_save_mask (str) -- Filename of the mask calculated by the unwrapper
- Returns:
3D array with the shape of complex_array of the unwrapped phase output from prelude
- Return type:
numpy.ndarray
Masking
Image mask with shape API
- shimmingtoolbox.masking.shapes.shape_cube(data, len_dim1, len_dim2, len_dim3, center_dim1=None, center_dim2=None, center_dim3=None)
Creates a cube mask. Returns mask with the same shape as data.
- Parameters:
data (numpy.ndarray) -- Data to mask, must be 3 dimensional array.
len_dim1 (int) -- Length of the side of the square along first dimension (in pixels).
len_dim2 (int) -- Length of the side of the square along second dimension (in pixels).
len_dim3 (int) -- Length of the side of the square along third dimension (in pixels).
center_dim1 (int) -- Center of the square along first dimension (in pixels). If no center is provided, the middle is used.
center_dim2 (int) -- Center of the square along second dimension (in pixels). If no center is provided, the middle is used.
center_dim3 (int) -- Center of the square along third dimension (in pixels). If no center is provided, the middle is used.
- Returns:
Mask with booleans. True where the cube is located and False in the background.
- Return type:
numpy.ndarray
- shimmingtoolbox.masking.shapes.shape_sphere(data, radius, center_dim1=None, center_dim2=None, center_dim3=None)
Creates a spherical mask. Returns mask with the same shape as data.
- Parameters:
data (numpy.ndarray) -- Data to mask, must be 3 dimensional array.
radius (int) -- Radius of the sphere (in pixels).
center_dim1 (int) -- Center of the sphere along the first dimension (in pixels). If no center is provided, the middle is used.
center_dim2 (int) -- Center of the sphere along the second dimension (in pixels). If no center is provided, the middle is used.
center_dim3 (int) -- Center of the sphere along the third dimension (in pixels). If no center is provided, the middle is used.
- Returns:
Mask with booleans. True where the sphere is located and False in the background.
- Return type:
numpy.ndarray
- shimmingtoolbox.masking.shapes.shape_square(data, len_dim1, len_dim2, center_dim1=None, center_dim2=None)
Creates a square mask. Returns mask with the same shape as data.
- Parameters:
data (numpy.ndarray) -- Data to mask, must be 2 dimensional array.
len_dim1 (int) -- Length of the side of the square along first dimension (in pixels).
len_dim2 (int) -- Length of the side of the square along second dimension (in pixels).
center_dim1 (int) -- Center of the square along first dimension (in pixels). If no center is provided, the middle is used.
center_dim2 (int) -- Center of the square along second dimension (in pixels). If no center is provided, the middle is used.
- Returns:
Mask with booleans. True where the square is located and False in the background.
- Return type:
numpy.ndarray
- shimmingtoolbox.masking.shapes.shapes(data, shape, **kargs)
Wrapper to different shape masking functions.
- Parameters:
data (numpy.ndarray) -- Data to mask.
shape (str) -- Shape to mask, implemented shapes include: {'square', 'cube', 'sphere'}.
**kargs -- Refer to the specific function in this file for the specific arguments for each shape. See example section for more details.
- Returns:
Mask with booleans. True where the shape is located and False in the background.
- Return type:
numpy.ndarray
Examples
>>> dummy_data = np.ones([4,3]) >>> dummy_mask = shapes(dummy_data, 'square', center_dim1=1, center_dim2=1, len_dim1=1, len_dim2=3)
Image thresholding API
- shimmingtoolbox.masking.threshold.threshold(data, thr=30, scaled_thr=False)
Threshold an image
- Parameters:
data (threshold. For complex) -- Data to be masked
thr (float) -- Value to threshold the data: voxels will be set to zero if their value is equal or less than this
data --
values. (threshold is applied on the absolute) --
scaled_thr (bool) -- Specifies if the threshold is absolute or scaled [0, 1].
- Returns:
Boolean mask with same dimensions as data
- Return type:
numpy.ndarray
- shimmingtoolbox.masking.mask_utils.dilate_binary_mask(mask, shape='sphere', size=3)
Dilates a binary mask according to different shapes and kernel size
- Parameters:
mask (numpy.ndarray) -- 3d array containing the binary mask.
shape (str) -- 3d kernel to perform the dilation. Allowed shapes are: 'sphere', 'cross', 'line', 'cube', 'None'. 'line' uses 3 line kernels to extend in each directions by "(size - 1) / 2" only if that direction is smaller than (size - 1) / 2
size (int) -- Length of a side of the 3d kernel. Must be odd.
- Returns:
Dilated mask.
- Return type:
numpy.ndarray
Notes
Kernels for
- 'cross' size 3:
np.array([[[0, 0, 0], [0, 1, 0], [0, 0, 0]], [[0, 1, 0], [1, 1, 1], [0, 1, 0]], [[0, 0, 0], [0, 1, 0], [0, 0, 0]]])
- 'sphere' size 5:
np.array([[[0 0 0 0 0], [0 0 0 0 0], [0 0 1 0 0], [0 0 0 0 0], [0 0 0 0 0]], [[0 0 0 0 0], [0 0 1 0 0], [0 1 1 1 0], [0 0 1 0 0], [0 0 0 0 0]], [[0 0 1 0 0], [0 1 1 1 0], [1 1 1 1 1], [0 1 1 1 0], [0 0 1 0 0]], [[0 0 0 0 0], [0 0 1 0 0], [0 1 1 1 0], [0 0 1 0 0], [0 0 0 0 0]], [[0 0 0 0 0], [0 0 0 0 0], [0 0 1 0 0], [0 0 0 0 0], [0 0 0 0 0]]]
- 'cube' size 3:
np.array([[[1, 1, 1], [1, 1, 1], [1, 1, 1]], [[1, 1, 1], [1, 1, 1], [1, 1, 1]], [[1, 1, 1], [1, 1, 1], [1, 1, 1]]])
- 'line' size 3:
np.array([[[0, 0, 0], [0, 1, 0], [0, 0, 0]], [[0, 0, 0], [0, 1, 0], [0, 0, 0]], [[0, 0, 0], [0, 1, 0], [0, 0, 0]]])
- shimmingtoolbox.masking.mask_utils.resample_mask(nii_mask_from, nii_target, from_slices=None, dilation_kernel='None', dilation_size=3, path_output=None)
Select the appropriate slices from
nii_mask_from
usingfrom_slices
and resample ontonii_target
- Parameters:
nii_mask_from (nib.Nifti1Image) -- Mask to resample from. False or 0 signifies not included.
nii_target (nib.Nifti1Image) -- Target image to resample onto.
from_slices (tuple) -- Tuple containing the slices to select from nii_mask_from. None selects all the slices.
dilation_kernel (str) -- kernel used to dilate the mask. Allowed shapes are: 'sphere', 'cross', 'line' 'cube'. See
dilate_binary_mask()
for more details.dilation_size (int) -- Length of a side of the 3d kernel to dilate the mask. Must be odd. For example, a kernel of size 3 will dilate the mask by 1 pixel.
path_output (str) -- Path to output debug artefacts.
- Returns:
Mask resampled with nii_target.shape and nii_target.affine.
- Return type:
nib.Nifti1Image
Coils
- class shimmingtoolbox.coils.coil.Coil(profile, affine, constraints)
Coil profile object that stores coil profiles and there constraints
- dim
Dimension along specific axis. dim: 0,1,2 are spatial axes, while dim: 3 corresponds to the coil channel.
- Type:
Tuple[int]
- profile
(dim1, dim2, dim3, channels) 4d array of N 3d coil profiles
- Type:
np.ndarray
- affine
4x4 array containing the affine transformation associated with the NIfTI file of the coil profile. This transformation relates to the physical coordinates of the scanner (qform).
- Type:
np.ndarray
- required_constraints
List containing the required keys for
constraints
- Type:
list
- coef_sum_max
Contains the maximum value for the sum of the coefficients
- Type:
float
- coef_channel_minmax
List of
(min, max)
pairs for each coil channels. (None, None) is used to specify no bounds.- Type:
list
- name
Name of the coil.
- Type:
str
- __init__(profile, affine, constraints)
Initialize Coil
- Parameters:
profile (np.ndarray) -- Coil profile (dim1, dim2, dim3, channels) 4d array of N 3d coil profiles
affine (np.ndarray) -- 4x4 array containing the qform affine transformation for the coil profiles
constraints (dict) --
dict containing the constraints for the coil profiles. Required keys:
name (str): Name of the coil.
coef_sum_max (float): Contains the maximum value for the sum of the coefficients. None is used to specify no bounds
coef_channel_max (list): List of
(min, max)
pairs for each coil channels. (None, None) is used to specify no bounds.
Examples
# Example of constraints constraints = { 'name': "dummy coil", 'coef_sum_max': 10, # 8 channel coil 'coef_channel_minmax': [(-2, 2), (-2, 2), (-2, 2), (-2, 2), (-3, 3), (-3, 3), (-3, 3), (-3, 3)], }
- load_constraints(constraints)
Loads the constraints named in required_constraints as attribute to this class
- class shimmingtoolbox.coils.coil.ScannerCoil(coord_system, dim_volume, affine, constraints, order)
Coil class for scanner coils as they require extra arguments
- __init__(coord_system, dim_volume, affine, constraints, order)
Initialize Coil
- Parameters:
profile (np.ndarray) -- Coil profile (dim1, dim2, dim3, channels) 4d array of N 3d coil profiles
affine (np.ndarray) -- 4x4 array containing the qform affine transformation for the coil profiles
constraints (dict) --
dict containing the constraints for the coil profiles. Required keys:
name (str): Name of the coil.
coef_sum_max (float): Contains the maximum value for the sum of the coefficients. None is used to specify no bounds
coef_channel_max (list): List of
(min, max)
pairs for each coil channels. (None, None) is used to specify no bounds.
Examples
# Example of constraints constraints = { 'name': "dummy coil", 'coef_sum_max': 10, # 8 channel coil 'coef_channel_minmax': [(-2, 2), (-2, 2), (-2, 2), (-2, 2), (-3, 3), (-3, 3), (-3, 3), (-3, 3)], }
- shimmingtoolbox.coils.coil.convert_to_mp(shim_setting, manufacturers_model_name)
- Converts the ShimSettings tag from the json BIDS sidecar to the scanner units.
(i.e. For the Prisma fit DAC --> uT/m, uT/m^2 (1st order, 2nd order))
- Parameters:
shim_setting (list) -- List of coefficients. Found in the json BIDS sidecar under 'ShimSetting'.
manufacturers_model_name (str) -- Name of the model of the scanner. Found in the json BIDS sidecar under ManufacturersModelName'. Supported names: 'Prisma_fit'.
- Returns:
Coefficients with units converted.
- Return type:
list
- shimmingtoolbox.coils.spherical_harmonics.spherical_harmonics(orders, x, y, z)
Returns an array of spherical harmonic basis fields with the order/degree index along the 4th dimension.
- Parameters:
orders (numpy.ndarray) -- Degrees of the desired terms in the series expansion, specified as a vector of non-negative integers (
np.array(range(0, 3))
yields harmonics up to (n-1)-th order). Must be non negative.x (numpy.ndarray) -- 3-D arrays of grid coordinates
y (numpy.ndarray) -- 3-D arrays of grid coordinates (same shape as x)
z (numpy.ndarray) -- 3-D arrays of grid coordinates (same shape as x)
- Returns:
4d basis set of spherical harmonics with order/degree ordered along 4th dimension
- Return type:
numpy.ndarray
Examples
Initialize grid positions
>>> [x, y, z] = np.meshgrid(np.array(range(-10, 11)), np.array(range(-10, 11)), np.array(range(-10, 11)), indexing='ij')
0th-to-2nd order terms inclusive
>>> orders = np.array(range(0, 3)) >>> basis = spherical_harmonics(orders, x, y, z)
Notes
- basis[:, :, :,0] corresponds to the 0th-order constant term (globally=unity)
0: c
- basis[:, :, :, 1:3] to 1st-order linear terms
1: y
2: z
3: x
- basis[:, :, :, 4:8] to 2nd-order terms
4: xy
5: zy
6: z2
7: zx
8: x2y2
- Based on
spherical_harmonics.m by topfer@ualberta.ca
calc_spherical_harmonics_arb_points_cz.m by jaystock@nmr.mgh.harvard.edu
- shimmingtoolbox.coils.siemens_basis.siemens_basis(x, y, z, orders=(1, 2))
The function first wraps
shimmingtoolbox.coils.spherical_harmonics
to generate 1st and 2nd order spherical harmonicbasis
fields at the grid positions given by arraysX,Y,Z
. Following Siemens convention,basis
is then:Reordered along the 4th dimension as X, Y, Z, Z2, ZX, ZY, X2-Y2, XY
Rescaled to Hz/unit-shim, where "unit-shim" refers to the measure displayed in the Adjustments card of the Syngo console UI, namely:
1 micro-T/m for X,Y,Z gradients (= 0.042576 Hz/mm)
1 micro-T/m^2 for 2nd order terms (= 0.000042576 Hz/mm^2)
The returned
basis
is thereby in the form of ideal "shim reference maps", ready for optimization.- Parameters:
x (numpy.ndarray) -- 3-D arrays of grid coordinates, "Left->Right" grid coordinates in the patient coordinate system (i.e. NIfTI reference (RAS), units of mm)
y (numpy.ndarray) -- 3-D arrays of grid coordinates (same shape as x). "Posterior->Anterior" grid coordinates in the patient coordinate system (i.e. NIfTI reference (RAS), units of mm)
z (numpy.ndarray) -- 3-D arrays of grid coordinates (same shape as x). "Inferior->Superior" grid coordinates in the patient coordinate system (i.e. NIfTI reference, units of mm)
orders (tuple) -- Degrees of the desired terms in the series expansion, specified as a vector of non-negative integers (
(0:1:n)
yields harmonics up to n-th order, implemented 1st and 2nd order)
- Returns:
4-D array of spherical harmonic basis fields
- Return type:
numpy.ndarray
Notes
For now,
orders
is, in fact, as default [1:2]—which is suitable for the Prisma (presumably other Siemens systems as well) however, the 3rd-order shims of the Terra should ultimately be accommodated too. (Requires checking the Adjustments/Shim card to see what the corresponding terms and values actually are). So,basis
will return with 8 terms along the 4th dim if using the 1st and 2nd order.
- shimmingtoolbox.coils.coordinates.generate_meshgrid(dim, affine)
Generate meshgrid of size dim, with coordinate system defined by affine. :param dim: x, y and z dimensions. :type dim: tuple :param affine: 4x4 affine matrix :type affine: numpy.ndarray
- Returns:
List of numpy.ndarray containing meshgrid of coordinates
- Return type:
list
- shimmingtoolbox.coils.coordinates.get_main_orientation(cosines: list)
Returns the orientation of the slice axis by looking at the ImageOrientationPatientDICOM JSON tag
- Parameters:
cosines (list) -- list of 6 elements. The first 3 represent the x, y, z cosines of the first row. The last 3
x (represent the) --
y --
it (z cosines of the first column. This can be found in ImageOrientationPatientDICOM so) --
coordinates. (should be LPS) --
- Returns:
'SAG', 'COR' or 'TRA'
- Return type:
str
- shimmingtoolbox.coils.coordinates.phys_gradient(data, affine)
Calculate the gradient of
data
along physical coordinates defined byaffine
- Parameters:
data (numpy.ndarray) -- 3d array containing data to apply gradient
affine (numpy.ndarray) -- 4x4 array containing affine transformation
- Returns
numpy.ndarray: 3D matrix containing the gradient along the x direction in the physical coordinate system numpy.ndarray: 3D matrix containing the gradient along the y direction in the physical coordinate system numpy.ndarray: 3D matrix containing the gradient along the z direction in the physical coordinate system
- shimmingtoolbox.coils.coordinates.phys_to_vox_coefs(gx, gy, gz, affine)
Calculate the vector sum along the image coordinates defined by
affine
with coefficients in the patient coordinate system.- Parameters:
gx (numpy.ndarray) -- 3D matrix containing the coefs along the x direction in the patient coordinate system
gy (numpy.ndarray) -- 3D matrix containing the coefs along the y direction in the patient coordinate system
gz (numpy.ndarray) -- 3D matrix containing the coefs along the z direction in the patient coordinate system
affine (numpy.ndarray) -- 4x4 array containing affine transformation
- Returns:
3D matrix containing the coefs along the x direction in the image coordinate system numpy.ndarray: 3D matrix containing the coefs along the y direction in the image coordinate system numpy.ndarray: 3D matrix containing the coefs along the z direction in the image coordinate system
- Return type:
numpy.ndarray
- shimmingtoolbox.coils.coordinates.resample_from_to(nii_from_img, nii_to_vox_map, order=2, mode='nearest', cval=0.0, out_class=<class 'nibabel.nifti1.Nifti1Image'>)
Wrapper to nibabel's
resample_from_to
function. Resample image from_img to mapped voxel space to_vox_map. The wrapper adds support for 2D input data (adds a singleton) and for 4D time series. For more info, refer to nibabel.processing.resample_from_to.- Parameters:
nii_from_img (nibabel.Nifti1Image) -- Nibabel object with 2D, 3D or 4D array. The 4d case will be treated as a timeseries.
nii_to_vox_map (nibabel.Nifti1Image) -- Nibabel object with
order (int) -- Refer to nibabel.processing.resample_from_to
mode (str) -- Refer to nibabel.processing.resample_from_to
cval (scalar) -- Refer to nibabel.processing.resample_from_to
out_class -- Refer to nibabel.processing.resample_from_to
- Returns:
- Return a Nibabel object with the resampled data. The 4d case will have an extra dimension
for the different time points.
- Return type:
nibabel.Nifti1Image
- shimmingtoolbox.coils.biot_savart.biot_savart(centers, normals, radii, segment_numbers, fov_min, fov_max, fov_n)
Creates coil profiles for arbitrary loops, for use in multichannel shim examples that do not match spherical harmonics :param centers: List of 3D float center points for each loop in mm :type centers: list :param normals: List of 3D float normal vectors for each loop in mm :type normals: list :param radii: List of float radii for each loop in mm :type radii: list :param segment_numbers: List of integer number of segments for each loop approximation :type segment_numbers: list :param fov_min: Low 3D float corner of coil profile field of view (x, y, z) in mm :type fov_min: tuple :param fov_max: Inclusive high 3D float corner of coil profile field of view (x, y, z) in mm :type fov_max: tuple :param fov_n: Integer number of points for each dimension (x, y, z) in mm :type fov_n: tuple
- Returns:
(X, Y, Z, centers) coil profiles of magnetic field z-component in Hz/A -- (X, Y, Z, Channel)
- Return type:
numpy.ndarray
Shim
- shimmingtoolbox.shim.realtime_shim.realtime_shim(nii_fieldmap, nii_anat, pmu, json_fmap, nii_mask_anat_riro=None, nii_mask_anat_static=None, path_output=None)
This function will generate static and dynamic (due to respiration) Gz components based on a fieldmap time series and respiratory trace information obtained from Siemens bellows (PMUresp_signal.resp). An additional multi-gradient echo (MGRE) magnitude image is used to generate an ROI and resample the static and dynamic Gz component maps to match the MGRE image. Lastly the mean Gz values within the ROI are computed for each slice.
- Parameters:
nii_fieldmap (nibabel.Nifti1Image) -- Nibabel object containing fieldmap data in 4d where the 4th dimension is the timeseries. Fieldmap should be in Hz.
nii_anat (nibabel.Nifti1Image) -- Nibabel object containing a 3d image of the target data to shim.
pmu (PmuResp) -- PmuResp object containing the respiratory trace information.
json_fmap (dict) -- dict of the json sidecar corresponding to the fieldmap data (Used to find the acquisition timestamps).
nii_mask_anat_static (nibabel.Nifti1Image) -- Nibabel object containing the mask to specify the shimming region for the static component.
nii_mask_anat_riro (nibabel.Nifti1Image) -- Nibabel object containing the mask to specify the shimming region for the riro component.
path_output (str) -- Path to output figures and temporary variables. If none is provided, no debug output is provided.
- Returns:
tuple containing:
numpy.ndarray: 1D array of the x static_correction. The correction is in mT/m for each slice.
numpy.ndarray: 1D array of the y static_correction. The correction is in mT/m for each slice.
numpy.ndarray: 1D array of the z static_correction. The correction is in mT/m for each slice.
numpy.ndarray: 1D array of the x dynamic riro_correction. The correction is in (mT/m)*rms_pressure for each slice.
numpy.ndarray: 1D array of the y dynamic riro_correction. The correction is in (mT/m)*rms_pressure for each slice.
numpy.ndarray: 1D array of the z dynamic riro_correction. The correction is in (mT/m)*rms_pressure for each slice.
float: Average pressure of the pmu
float: RMS of the pmu pressure
- Return type:
(tuple)
- class shimmingtoolbox.shim.sequencer.RealTimeSequencer(nii_fieldmap, json_fmap, nii_anat, nii_static_mask, nii_riro_mask, slices, pmu: PmuResp, coils, method='least_squares', opt_criteria='mse', mask_dilation_kernel='sphere', mask_dilation_kernel_size=3, reg_factor=0, path_output=None)
Sequencer object that stores different nibabel object, and parameters. It's also doing real time optimization of the currents, and the evaluation of the shimming
- nii_fieldmap
Nibabel object containing fieldmap data in 4d where the 4th dimension is the timeseries. Also contains an affine transformation.
- Type:
nib.Nifti1Image
- json_fmap
Dict of the json sidecar corresponding to the fieldmap data (Used to find the acquisition timestamps).
- Type:
dict
- nii_anat
Nibabel object containing anatomical data in 3d.
- Type:
nib.Nifti1Image
- nii_static_mask
3D anat mask used for the optimizer to shim the region for the static component.
- Type:
nib.Nifti1Image
- nii_riro_mask
3D anat mask used for the optimizer to shim the region for the riro component.
- Type:
nib.Nifti1Image
- slices
1D array containing tuples of dim3 slices to shim according to the anat where the shape of anat: (dim1, dim2, dim3). Refer to
shimmingtoolbox.shim.sequencer.define_slices()
.- Type:
list
- coils
List of Coils containing the coil profiles. The coil profiles and the fieldmaps must have matching units (if fmap is in Hz, the coil profiles must be in hz/unit_shim). Refer to
shimmingtoolbox.coils.coil.Coil
. Make sure the extent of the coil profiles are larger than the extent of the fieldmap. This is especially true for dimensions with only 1 voxel(e.g. (50x50x1x10). Refer toshimmingtoolbox.shim.sequencer.extend_slice()
/shimmingtoolbox.shim.sequencer.update_affine_for_ap_slices()
- Type:
ListCoil
- method
Supported optimizer: 'least_squares', 'pseudo_inverse'. Note: refer to their specific implementation to know limits of the methods in:
shimmingtoolbox.optimizer
- Type:
str
- opt_criteria
Criteria for the optimizer 'least_squares'. Supported: 'mse': mean squared error, 'mae': mean absolute error, 'std': standard deviation.
- Type:
str
- reg_factor
Regularization factor for the current when optimizing. A higher coefficient will penalize higher current values while a lower factor will lower the effect of the regularization. A negative value will favour high currents (not preferred). Only relevant for 'least_squares' opt_method.
- Type:
float
- mask_dilation_kernel
Kernel used to dilate the mask. Allowed shapes are: 'sphere', 'cross', 'line' 'cube'. See
shimmingtoolbox.masking.mask_utils.dilate_binary_mask()
for more details.- Type:
str
- mask_dilation_kernel_size
Length of a side of the 3d kernel to dilate the mask. Must be odd. For example, a kernel of size 3 will dilate the mask by 1 pixel.
- Type:
int
- path_output
Path to the directory to output figures. Set logging level to debug to output debug artefacts.
- Type:
str
- optimizer
Object that contains everything needed for the optimization created from shimmingtoolbox.optimizer init method
- Type:
object
- bounds
List of the bounds for the currents for the real time optimization
- Type:
list
- acq_pressures
1D array that contains the acquisitions pressures
- Type:
np.ndarray
- eval(coef_static, coef_riro, mean_p, pressure_rms)
Evaluate the real time shimming by plotting and saving results
- Parameters:
coef_static (np.ndarray) -- coefficients got during the static optimization
coef_riro (np.ndarray) -- coefficients got during the real time optimization
mean_p (float) -- mean of the acquisitions pressures
pressure_rms (float) -- rms of the acquisitions pressures
- get_acq_pressures()
Get the acquisition pressures at the times when the fieldmap was acquired.
- Returns:
1D array that contains the acquisitions pressures
- Return type:
np.ndarray
- get_fieldmap(nii_fieldmap)
Get the fieldmap for the RealTimeSequencer class
- Parameters:
nii_fieldmap (nib.Nifti1Image) -- Nibabel object containing fieldmap data in 4d where the 4th dimension is the timeseries.
- Returns:
- Nibabel object containing fieldmap data in 4d where the 4th dimension
is the timeseries.
- Return type:
nib.Nifti1Image
- get_mask(nii_static_mask, nii_riro_mask)
Get both masks for the RealTimeSequencer Class
- Parameters:
nii_static_mask (nib.Nifti1Image) -- 3D anat mask used for the optimizer to shim the region for the static component.
nii_riro_mask (nib.Nifti1Image) -- 3D anat mask used for the optimizer to shim the region for the riro component.
- Returns:
- tuple containing:
nib.Nifti1Image: 3D anat mask used for the optimizer to shim the region for the static component.
nib.Nifti1Image: 3D anat mask used for the optimizer to shim the region for the riro component.
- Return type:
(tuple)
- get_real_time_parameters()
Get real time parameters used for the shimming
- Returns:
- tuple containing:
np.ndarray: 3D array containing the static data for the optimization
np.ndarray: 3D array containing the real time data for the optimization
float: Mean pressure of the respiratory trace.
- float: Root mean squared of the pressure trace. This is provided to compare results between scans,
multiply the riro coefficients by rms of the pressure to do so.
- Return type:
(tuple)
- opt_riro(i, optimizer, nii_mask_anat, slices_anat, dilation_kernel, dilation_size, path_output, shimwise_bounds)
Make the optmization of the currents for each slice
- Parameters:
i (int) -- Index of the shim group of the optimization
optimizer (Optimizer) -- Initialized Optimizer object
nii_mask_anat (nib.Nifti1Image) -- anat mask on which the optimization will be make
slices_anat -- 1D array containing tuples of dim3 slices to shim according to the anat, where the shape of anat is: (dim1, dim2, dim3). Refer to
shimmingtoolbox.shim.sequencer.define_slices()
.dilation_kernel (str) -- Kernel used to dilate the mask. Allowed shapes are: 'sphere', 'cross', 'line' 'cube'. See
shimmingtoolbox.masking.mask_utils.dilate_binary_mask()
for more details.dilation_size (int) -- Length of a side of the 3d kernel to dilate the mask. Must be odd. For example, a kernel of size 3 will dilate the mask by 1 pixel.
path_output (str) -- Path to the directory to output figures. Set logging level to debug to output debug artefacts.
shimwise_bounds (list) -- Bounds of the currents for the optimization
- Returns:
Riro coefficients of the coil profiles to shim (channels) [Hz/unit_pressure]
- Return type:
np.ndarray
- optimize_riro(nii_mask_anat)
- Parameters:
nii_mask_anat (nib.Nifti1Image) -- anat mask on which the optimization will be made
- Returns:
Riro coefficients of the coil profiles to shim (len(slices) x channels) [Hz/unit_pressure]
- plot_currents(static, riro=None)
Plot evolution of currents through shim groups
- Parameters:
static (np.ndarray) -- Array with the static currents
riro (np.ndarray) -- Array with the riro currents
- plot_pressure_points(acq_pressures, ylim)
Plot respiratory trace pressure points
- Parameters:
acq_pressures (np.ndarray) -- acquisitions pressures
ylim (float) -- Limit of the y-axis
- plot_shimmed_trace(unshimmed_trace, shim_trace_static, shim_trace_riro, shim_trace_static_riro)
Plot shimmed and unshimmed sum over the roi for each shim
- Parameters:
unshimmed_trace (np.ndarray) -- array with the trace of the nii_fieldmap data
shim_trace_static (np.ndarray) -- array with the trace of the nii_fieldmap data after the static shimming
shim_trace_riro (np.ndarray) -- array with the trace of the nii_fieldmap data after the riro shimming
shim_trace_static_riro (np.ndarray) -- array with the trace of the nii_fieldmap data after both shimming
- plot_static_riro(masked_unshimmed, masked_shim_static, masked_shim_static_riro, unshimmed, shimmed_static, shimmed_static_riro, i_t=0, i_slice=0, i_shim=0)
Plot Static and RIRO maps for a particular fieldmap slice, anat shim and timepoint
- Parameters:
masked_unshimmed (np.ndarray) -- Fieldmap masked before the shimming
masked_shim_static (np.ndarray) -- Fieldmap masked after static shimming
masked_shim_static_riro (np.ndarray) -- Fieldmap masked after the static and riro shimming
unshimmed (np.ndarray) -- Fieldmap not shimmed
shimmed_static (np.ndarray) -- Data of the nii_fieldmap after the static shimming
shimmed_static_riro (np.ndarray) -- Data of the nii_fieldmap after static and riro shimming
i_shim -- (int): index of the anat shim, where we want to plot the static and riro maps
i_slice -- (int): index of the slice, where we want to plot the static and riro maps
i_t -- (int): Index of the time, where we want to plot the static and riro maps
- print_rt_metrics(unshimmed, shimmed_static, shimmed_static_riro, shimmed_riro, mask)
Print to the console metrics about the realtime and static shim. These metrics isolate temporal and static components
Temporal: Compute the STD across time pixelwise, and then compute the mean across pixels. Static: Compute the MEAN across time pixelwise, and then compute the STD across pixels.
- Parameters:
unshimmed (np.ndarray) -- Fieldmap not shimmed
shimmed_static (np.ndarray) -- Data of the nii_fieldmap after the static shimming
shimmed_static_riro (np.ndarray) -- Data of the nii_fieldmap after static and riro shimming
shimmed_riro (np.ndarray) -- Data of the nii_fieldmap after the riro shimming
mask (np.ndarray) -- Mask where the shimming was done
- select_optimizer(unshimmed, affine, pmu: PmuResp | None = None)
Select and initialize the optimizer
- Parameters:
unshimmed (np.ndarray) -- 3D B0 map
affine (np.ndarray) -- 4x4 array containing the affine transformation for the unshimmed array
pmu (PmuResp) -- PmuResp object containing the respiratory trace information. Required for method 'least_squares_rt'.
- Returns:
Initialized Optimizer object
- Return type:
- shim()
Performs realtime shimming using one of the supported optimizers and an external respiratory trace.
- Returns:
- tuple containing:
np.ndarray: Static coefficients of the coil profiles to shim (len(slices) x channels) e.g. [Hz]
- np.ndarray: Riro coefficients of the coil profiles to shim (len(slices) x channels)
e.g. [Hz/unit_pressure]
float: Mean pressure of the respiratory trace.
- float: Root mean squared of the pressure.
This is provided to compare results between scans, multiply the riro coefficients by rms of the pressure to do so.
- Return type:
(tuple)
- class shimmingtoolbox.shim.sequencer.Sequencer(slices, mask_dilation_kernel, mask_dilation_kernel_size, reg_factor, path_output)
General class for the sequencer
- slices
1D array containing tuples of dim3 slices to shim according to the anat, where the shape of anat is: (dim1, dim2, dim3). Refer to
shimmingtoolbox.shim.sequencer.define_slices()
.- Type:
list
- mask_dilation_kernel
Kernel used to dilate the mask. Allowed shapes are: 'sphere', 'cross', 'line' 'cube'. See
shimmingtoolbox.masking.mask_utils.dilate_binary_mask()
for more details.- Type:
str
- mask_dilation_kernel_size
Length of a side of the 3d kernel to dilate the mask. Must be odd. For example, a kernel of size 3 will dilate the mask by 1 pixel.
- Type:
int
- reg_factor
Regularization factor for the current when optimizing. A higher coefficient will penalize higher current values while a lower factor will lower the effect of the regularization. A negative value will favour high currents (not preferred). Only relevant for 'least_squares' opt_method.
- Type:
float
- path_output
Path to the directory to output figures. Set logging level to debug to output debug
- Type:
str
- opt(i, optimizer, nii_mask_anat, slices_anat, dilation_kernel, dilation_size, path_output)
Make the optimization of the currents for one shim group.
- Parameters:
i (int) -- Index of the shim group for the optimization.
optimizer (Optimizer) -- Initialized Optimizer object.
nii_mask_anat (nib.Nifti1Image) -- 3D anat mask used for the optimizer to shim in the region of interest. (only consider voxels with non-zero values)
slices_anat (list) -- 1D list containing tuples of dim3 slices to shim according to the anat, where the shape of anat is: (dim1, dim2, dim3). Refer to
shimmingtoolbox.shim.sequencer.define_slices()
.dilation_kernel (str) -- Kernel used to dilate the mask. Allowed shapes are: 'sphere', 'cross', 'line', 'cube'. See
shimmingtoolbox.masking.mask_utils.dilate_binary_mask()
for more details.dilation_size (int) -- Length of a side of the 3d kernel to dilate the mask. Must be odd. For example, a kernel of size 3 will dilate the mask by 1 pixel.
path_output (str) -- Path to the directory to output figures. Set logging level to debug to output debug artefacts.
- Returns:
Coefficients of the coil profiles to shim (n_channels)
- Return type:
np.ndarray
- optimize(nii_mask_anat)
Optimization of the currents for each shim group. Wraps
shimmingtoolbox.shim.sequencer.Sequencer.opt()
.- Parameters:
nii_mask_anat (nib.Nifti1Image) -- 3D anat mask used for the optimizer to shim in the region of interest. (only consider voxels with non-zero values)
- Returns:
Coefficients of the coil profiles to shim (len(slices) x n_channels)
- Return type:
np.ndarray
- class shimmingtoolbox.shim.sequencer.ShimSequencer(nii_fieldmap, nii_anat, nii_mask_anat, slices, coils, method='least_squares', opt_criteria='mse', mask_dilation_kernel='sphere', mask_dilation_kernel_size=3, reg_factor=0, path_output=None)
ShimSequencer object to perform optimization of shim parameters for static and dynamic shimming. This object can also evaluate the shimming performance.
- nii_fieldmap
Nibabel object containing fieldmap data in 3d.
- Type:
nib.Nifti1Image
- nii_anat
Nibabel object containing anatomical data in 3d.
- Type:
nib.Nifti1Image
- nii_mask_anat
3D anat mask used for the optimizer to shim in the region of interest. (only consider voxels with non-zero values)
- Type:
nib.Nifti1Image
- coils
List of Coils containing the coil profiles. The coil profiles and the fieldmaps must have matching units (if fmap is in Hz, the coil profiles must be in hz/unit_shim). Refer to
shimmingtoolbox.coils.coil.Coil
. Make sure the extent of the coil profiles are larger than the extent of the fieldmap. This is especially true for dimensions with only 1 voxel(e.g. (50x50x1). Refer toshimmingtoolbox.shim.sequencer.extend_slice()
/shimmingtoolbox.shim.sequencer.update_affine_for_ap_slices()
- Type:
ListCoil
- method
Supported optimizer: 'least_squares', 'pseudo_inverse'. Note: refer to their specific implementation to know limits of the methods in:
shimmingtoolbox.optimizer
- Type:
str
- opt_criteria
Criteria for the optimizer 'least_squares'. Supported: 'mse': mean squared error, 'mae': mean absolute error, 'std': standard deviation.
- Type:
str
- nii_fieldmap_orig
Nibabel object containing the copy of the original fieldmap data
- Type:
nib.Nifti1Image
- fmap_is_extended
Tells whether the fieldmap has been extended by the object.
- Type:
bool
- calc_shimmed_anat_orient(coefs, list_shim_slice)
Calculate and save the shimmed anat orient
- Parameters:
coefs (np.ndarray) -- Coefficients of the coil profiles to shim (len(slices) x n_channels)
list_shim_slice (list) -- list of the index where there was a correction
- calc_shimmed_full_mask(unshimmed, correction, masks_fmap)
Calculate the shimmed full mask
- Parameters:
unshimmed (np.ndarray) -- Original fieldmap not shimmed
correction (np.ndarray) -- Corrections to apply to the fieldmap
masks_fmap (np.ndarray) -- Resampled mask on the fieldmap
- Returns:
- tuple containing:
np.ndarray: Masked shimmed fieldmap
np.ndarray: Binary mask in the fieldmap space
- Return type:
(tuple)
- display_shimmed_results(shimmed, unshimmed, masks_fmap, coef)
Print the efficiency of the corrections according to the opt_criteria
- Parameters:
shimmed (np.ndarray) -- Shimmed fieldmap
unshimmed (np.ndarray) -- Original fieldmap not shimmed
masks_fmap (np.ndarray) -- Resampled mask on the fieldmap
coef (np.ndarray) -- Coefficients of the coil profiles to shim (len(slices) x n_channels)
- eval(coef)
Calculate theoretical shimmed map and output figures.
- Args :
coef (np.ndarray): Coefficients of the coil profiles to shim (len(slices) x n_channels)
- evaluate_shimming(unshimmed, coef, merged_coils)
Evaluate the shimming and print the efficiency of the corrections.
- Parameters:
unshimmed (np.ndarray) -- Original fieldmap not shimmed
coef (np.ndarray) -- Coefficients of the coil profiles to shim (len(slices) x n_channels)
merged_coils (np.ndarray) -- Coils resampled on the original fieldmap
- Returns:
- tuple containing:
np.ndarray: Shimmed fieldmap
np.ndarray: Corrections to apply to the fieldmap
np.ndarray: Resampled mask on the fieldmap
list: List containing the indexes of the shimmed slices
- Return type:
(tuple)
- get_anat(nii_anat)
Get the target image and perform error checking.
- Parameters:
nii_anat (nib.Nifti1Image) -- Nibabel object containing anatomical data in 3d.
- Returns:
Nibabel object containing anatomical data in 3d.
- Return type:
nib.Nifti1Image
- get_fieldmap(nii_fieldmap)
Get the fieldmap and perform error checking.
- Parameters:
nii_fieldmap (nib.Nifti1Image) -- Nibabel object containing fieldmap data in 3d.
- Returns:
tuple containing:
nib.Nifti1Image: Nibabel object containing fieldmap data in 3d.
nib.Nifti1Image: Nibabel object containing the copy of the initial fieldmap data in 3d.
bool: Boolean indicating if the initial fieldmap has been changed.
- Return type:
(tuple)
- get_mask(nii_mask_anat)
Get the mask and perform error checking.
- Parameters:
nii_mask_anat (nib.Nifti1Image) -- 3D anat mask used for the optimizer to shim in the region of interest.(only consider voxels with non-zero values)
- Returns:
- 3D anat mask used for the optimizer to shim in the region of interest.
(Only consider voxels with non-zero values)
- Return type:
nib.Nifti1Image
- plot_currents(static)
Plot evolution of currents through shim groups
- Parameters:
static (np.ndarray) -- Array with the static coefficients
- plot_full_mask(unshimmed, shimmed_masked, mask)
Plot and save the static full mask
- Parameters:
unshimmed (np.ndarray) -- Original fieldmap not shimmed
shimmed_masked (np.ndarray) -- Masked shimmed fieldmap
mask (np.ndarray) -- Binary mask in the fieldmap space
- plot_partial_mask(unshimmed, shimmed, masks)
This figure shows a single fieldmap slice for all shim groups. The shimmed and unshimmed fieldmaps are in the background and the correction is overlaid in color.
- Parameters:
unshimmed (np.ndarray) -- Original fieldmap not shimmed
shimmed (np.ndarray) -- Shimmed fieldmap
masks (np.ndarray) -- Resampled mask on the fieldmap
- select_optimizer()
Select and initialize the optimizer
- Returns:
Initialized Optimizer object
- Return type:
- shim()
Performs shimming according to slices using one of the supported optimizers and coil profiles.
- Returns:
Coefficients of the coil profiles to shim (len(slices) x n_channels)
- Return type:
np.ndarray
- shimmingtoolbox.shim.sequencer.define_slices(n_slices: int, factor=1, method='sequential')
Define the slices to shim according to the output convention. (list of tuples)
- Parameters:
n_slices (int) -- Number of total slices.
factor (int) -- Number of slices per shim.
method (str) -- Defines how the slices should be sorted, supported methods include: 'interleaved', 'sequential', 'volume'. See Examples for more details.
- Returns:
1D list containing tuples of dim3 slices to shim. (dim1, dim2, dim3)
- Return type:
list
Examples
- ::
slices = define_slices(10, 2, 'interleaved') print(slices) # [(0, 5), (1, 6), (2, 7), (3, 8), (4, 9)] slices = define_slices(20, 5, 'sequential') print(slices) # [(0, 1, 2, 3, 4), (5, 6, 7, 8, 9), (10, 11, 12, 13, 14), (15, 16, 17, 18, 19)] slices = define_slices(20, method='volume') # 'volume' ignores the 'factor' option print(slices) # [(0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19)]
- shimmingtoolbox.shim.sequencer.extend_fmap_to_kernel_size(nii_fmap_orig, dilation_kernel_size, path_output=None)
Load the fmap and expand its dimensions to the kernel size
- Parameters:
nii_fmap_orig (nib.Nifti1Image) -- 3d (dim1, dim2, dim3) or 4d (dim1, dim2, dim3, t) nii to be extended
dilation_kernel_size -- Size of the kernel
path_output (str) -- Path to save the debug output
- Returns:
Nibabel object of the loaded and extended fieldmap
- Return type:
nib.Nifti1Image
- shimmingtoolbox.shim.sequencer.extend_slice(nii_array, n_slices=1, axis=2)
Adds n_slices on each side of the selected axis. It uses the nearest slice and copies it to fill the values. Updates the affine of the matrix to keep the input array in the same location.
- Parameters:
nii_array (nib.Nifti1Image) -- 3d or 4d array to extend the dimensions along an axis.
n_slices (int) -- Number of slices to add on each side of the selected axis.
axis (int) -- Axis along which to insert the slice(s), Allowed axis: 0, 1, 2.
- Returns:
Array extended with the appropriate affine to conserve where the original pixels were located.
- Return type:
nib.Nifti1Image
Examples
- ::
print(nii_array.get_fdata().shape) # (50, 50, 1, 10) nii_out = extend_slice(nii_array, n_slices=1, axis=2) print(nii_out.get_fdata().shape) # (50, 50, 3, 10)
- shimmingtoolbox.shim.sequencer.new_bounds_from_currents(currents, old_bounds)
Uses the currents to determine the appropriate bounds for the next optimization. It assumes that "old_coef + next_bound < old_bound".
- Parameters:
currents (np.ndarray) -- 2D array (n_shims x n_channels). Direct output from
_optimize()
.old_bounds (list) -- 1d list (n_channels) of tuples (min, max) containing the merged bounds of the previous optimization.
- Returns:
2d list (n_shim_groups x n_channels) of bounds (min, max) corresponding to each shim group and channel.
- Return type:
list
- shimmingtoolbox.shim.sequencer.parse_slices(fname_nifti)
Parse the BIDS sidecar associated with the input nifti file.
- Parameters:
fname_nifti (str) -- Full path to a NIfTI file
- Returns:
1D list containing tuples of dim3 slices to shim. (dim1, dim2, dim3)
- Return type:
list
- shimmingtoolbox.shim.sequencer.shim_max_intensity(nii_input, nii_mask=None)
- Find indexes of the 4th dimension of the input volume that has the highest signal intensity for each slice.
Based on: https://onlinelibrary.wiley.com/doi/10.1002/hbm.26018
- Parameters:
nii_input (nib.Nifti1Image) -- 4d volume where 4th dimension was acquired with different shim values
nii_mask (nib.Nifti1Image) -- Mask defining the spatial region to shim. If None: consider all voxels of nii_input.
- Returns:
1d array containing the index of the volume that maximizes signal intensity for each slice
- Return type:
np.ndarray
- shimmingtoolbox.shim.sequencer.update_affine_for_ap_slices(affine, n_slices=1, axis=2)
Updates the input affine to reflect an insertion of n_slices on each side of the selected axis
- Parameters:
affine (np.ndarray) -- 4x4 qform affine matrix representing the coordinates
n_slices (int) -- Number of pixels to add on each side of the selected axis
axis (int) -- Axis along which to insert the slice(s)
- Returns:
4x4 updated affine matrix
- Return type:
np.ndarray
- shimmingtoolbox.shim.b1shim.b1shim(b1, mask=None, algorithm=1, target=None, q_matrix=None, sar_factor=1.5)
Computes static optimized shim weights that minimize the B1+ field coefficient of variation over the masked region.
- Parameters:
b1 (numpy.ndarray) -- 4D array corresponding to the measured B1+ field. (x, y, n_slices, n_channels)
mask (numpy.ndarray) -- 3D array corresponding to the region where shimming will be performed. (x, y, n_slices)
algorithm (int) -- Number from 1 to 4 specifying which algorithm to use for B1+ optimization: 1 - Reduce the coefficient of variation of the B1+ field. Favors high B1+ efficiency. 2 - Magnitude least square (MLS) optimization targeting a specific B1+ value. Target value required. 3 - Maximizes the SAR efficiency (B1+/sqrt(SAR)). Q matrices required. 4 - Phase-only shimming.
target (float) -- Target B1+ value used by algorithm 2 in nT/V.
q_matrix (numpy.ndarray) -- Matrix used to constrain local SAR. If no matrix is provided, unconstrained optimization is performed, which might result in SAR excess at the scanner (n_channels, n_channels, n_vop).
sar_factor (float) -- Factor (=> 1) to which the maximum local SAR after shimming can exceed the phase-only shimming maximum local SAR. Values between 1 and 1.5 should work with Siemens scanners. High factors allow more shimming liberty but are more likely to result in SAR excess at the scanner.
- Returns:
Optimized and normalized 1D vector of complex shimming weights of length n_channels.
- Return type:
numpy.ndarray
- shimmingtoolbox.shim.b1shim.combine_maps(b1_maps, weights)
Combines the B1 field distribution of several channels into one map representing the total B1 field magnitude.
- Parameters:
b1_maps (numpy.ndarray) -- Complex B1 field for different channels (x, y, n_slices, n_channels).
weights (numpy.ndarray) -- 1D complex array of length n_channels.
- Returns:
B1 field distribution obtained when applying the provided shim weights.
- Return type:
numpy.ndarray
- shimmingtoolbox.shim.b1shim.complex_to_vector(weights)
Separates the real and imaginary components of a complex vector into a twice as long vector.
- Parameters:
weights (numpy.ndarray) -- 1D complex array of length n_channels.
- Returns:
1D array of length 2*n_channels. First/second half: real/imaginary.
- Return type:
numpy.ndarray
- shimmingtoolbox.shim.b1shim.load_siemens_vop(fname_sar_file)
Reads in a Matlab file in which the VOP matrices are stored and returns them as a numpy array.
- Parameters:
fname_sar_file -- Path to the 'SarDataUser.mat' file containing the scanner's VOPs. This file should be available at the scanner in 'C:/Medcom/MriProduct/PhysConfig'.
- Returns:
VOP matrices (n_coils, n_coils, n_VOPs)
- Return type:
numpy.ndarray
- shimmingtoolbox.shim.b1shim.max_sar(weights, q_matrix)
Returns the maximum local SAR corresponding to a set of shim weight and a set of Q matrices.
- Parameters:
weights (numpy.ndarray) -- 1D vector of complex shim weights. (length: n_channel)
q_matrix (numpy.ndarray) -- Q matrices used to compute the local energy deposition in the tissues.
(n_channels --
n_channels --
n_voxel) --
- Returns:
maximum local SAR.
- Return type:
float
- shimmingtoolbox.shim.b1shim.phase_only_shimming(b1_maps, init_phases=None)
Performs a phase-only RF-shimming to find a set of phases that homogenizes the B1+ field.
- Parameters:
b1_maps (numpy.ndarray) -- 4D array corresponding to the measured B1 field. (x, y, n_slices, n_channels)
init_phases (numpy.ndarray) -- 1D array of initial phase values used for optimization.
- Returns:
Optimized and normalized 1D vector of complex shimming weights of length n_channels.
- Return type:
numpy.ndarray
- shimmingtoolbox.shim.b1shim.vector_to_complex(weights)
Combines real and imaginary values contained in a vector into a half long complex vector.
- Parameters:
weights (numpy.ndarray) -- 1D array of length 2*n_channels. First/second half: real/imaginary.
- Returns:
1D complex array of length n_channels.
- Return type:
numpy.ndarray
Optimizer
- class shimmingtoolbox.optimizer.basic_optimizer.Optimizer(coils: List[Coil], unshimmed, affine)
Optimizer object that stores coil profiles and optimizes an unshimmed volume given a mask. Use optimize(args) to optimize a given mask. For basic optimizer, uses unbounded pseudo-inverse.
- coils
List of Coil objects containing the coil profiles and related constraints
- Type:
ListCoil
- unshimmed
3d array of unshimmed volume
- Type:
numpy.ndarray
- unshimmed_affine
4x4 array containing the qform affine transformation for the unshimmed array
- Type:
numpy.ndarray
- merged_coils
4d array containing all coil profiles resampled onto the target unshimmed array concatenated on the 4th dimension. See self.merge_coils() for more details
- Type:
numpy.ndarray
- merged_bounds
list of bounds corresponding to each merged coils: merged_bounds[3] is the (min, max) bound for merged_coils[..., 3]
- Type:
list
- __init__(coils: List[Coil], unshimmed, affine)
Initializes coils according to input list of Coil
- Parameters:
coils (ListCoil) -- List of Coil objects containing the coil profiles and related constraints
unshimmed (numpy.ndarray) -- 3d array of unshimmed volume
affine (numpy.ndarray) -- 4x4 array containing the affine transformation for the unshimmed array
- merge_bounds()
Merge the coil profile bounds into a single array.
- Returns:
list of bounds corresponding to each merged coils
- Return type:
list
- merge_coils(unshimmed, affine)
Uses the list of coil profiles to return a resampled concatenated list of coil profiles matching the unshimmed image. Bounds are also concatenated and returned.
- Parameters:
unshimmed (numpy.ndarray) -- 3d array of unshimmed volume
affine (numpy.ndarray) -- 4x4 array containing the affine transformation for the unshimmed array
- optimize(mask)
Optimize unshimmed volume by varying current to each channel
- Parameters:
mask (numpy.ndarray) -- 3d array of integers marking volume for optimization. Must be the same shape as unshimmed
- Returns:
- Coefficients corresponding to the coil profiles that minimize the objective function.
The shape of the array returned has shape corresponding to the total number of channels
- Return type:
numpy.ndarray
- set_merged_bounds(merged_bounds)
Changes the default bounds set in the coil profile
- Parameters:
merged_bounds -- Concatenated coil profile bounds
- set_unshimmed(unshimmed, affine)
Set the unshimmed array to a new array. Resamples coil profiles accordingly.
- Parameters:
unshimmed (numpy.ndarray) -- 3d array of unshimmed volume
affine -- (numpy.ndarray): 4x4 array containing the qform affine transformation for the unshimmed array
- class shimmingtoolbox.optimizer.lsq_optimizer.LsqOptimizer(coils: List[Coil], unshimmed, affine, opt_criteria='mse', reg_factor=0)
Bases:
Optimizer
Optimizer object that stores coil profiles and optimizes an unshimmed volume given a mask. Use optimize(args) to optimize a given mask. The algorithm uses a least squares solver to find the best shim. It supports bounds for each channel as well as a bound for the absolute sum of the channels.
- __init__(coils: List[Coil], unshimmed, affine, opt_criteria='mse', reg_factor=0)
Initializes coils according to input list of Coil
- Parameters:
coils (ListCoil) -- List of Coil objects containing the coil profiles and related constraints
unshimmed (numpy.ndarray) -- 3d array of unshimmed volume
affine (numpy.ndarray) -- 4x4 array containing the affine transformation for the unshimmed array
opt_criteria (str) -- Criteria for the optimizer 'least_squares'. Supported: 'mse': mean squared error, 'mae': mean absolute error, 'std': standard deviation.
reg_factor (float) -- Regularization factor for the current when optimizing. A higher coefficient will penalize higher current values while a lower factor will lower the effect of the regularization. A negative value will favour high currents (not preferred).
- get_initial_guess()
Calculates the initial guess according to the self.initial_guess_method
- Returns:
1d array (n_channels) containing the initial guess for the optimization
- Return type:
np.ndarray
- merge_bounds()
Merge the coil profile bounds into a single array.
- Returns:
list of bounds corresponding to each merged coils
- Return type:
list
- merge_coils(unshimmed, affine)
Uses the list of coil profiles to return a resampled concatenated list of coil profiles matching the unshimmed image. Bounds are also concatenated and returned.
- Parameters:
unshimmed (numpy.ndarray) -- 3d array of unshimmed volume
affine (numpy.ndarray) -- 4x4 array containing the affine transformation for the unshimmed array
- optimize(mask)
Optimize unshimmed volume by varying current to each channel
- Parameters:
mask (numpy.ndarray) -- 3D integer mask used for the optimizer (only consider voxels with non-zero values).
- Returns:
- Coefficients corresponding to the coil profiles that minimize the objective function.
The shape of the array returned has shape corresponding to the total number of channels
- Return type:
numpy.ndarray
- set_merged_bounds(merged_bounds)
Changes the default bounds set in the coil profile
- Parameters:
merged_bounds -- Concatenated coil profile bounds
- set_unshimmed(unshimmed, affine)
Set the unshimmed array to a new array. Resamples coil profiles accordingly.
- Parameters:
unshimmed (numpy.ndarray) -- 3d array of unshimmed volume
affine -- (numpy.ndarray): 4x4 array containing the qform affine transformation for the unshimmed array
- class shimmingtoolbox.optimizer.lsq_optimizer.PmuLsqOptimizer(coils, unshimmed, affine, opt_criteria, pmu: PmuResp, reg_factor=0)
Bases:
LsqOptimizer
Optimizer for the realtime component (riro) for this optimization: field(i_vox) = riro(i_vox) * (acq_pressures - mean_p) + static(i_vox) Unshimmed must be in units: [unit_shim/unit_pressure], ex: [Hz/unit_pressure]
This optimizer bounds the riro results to the coil bounds by taking the range of pressure that can be reached by the PMU.
- __init__(coils, unshimmed, affine, opt_criteria, pmu: PmuResp, reg_factor=0)
Initializes coils according to input list of Coil
- Parameters:
coils (ListCoil) -- List of Coil objects containing the coil profiles and related constraints
unshimmed (numpy.ndarray) -- 3d array of unshimmed volume
affine (numpy.ndarray) -- 4x4 array containing the affine transformation for the unshimmed array
opt_criteria (str) -- Criteria for the optimizer 'least_squares'. Supported: 'mse': mean squared error, 'mae': mean absolute error, 'std': standard deviation.
pmu (PmuResp) -- PmuResp object containing the respiratory trace information.
- get_initial_guess()
Calculates the initial guess according to the self.initial_guess_method
- Returns:
1d array (n_channels) containing the initial guess for the optimization
- Return type:
np.ndarray
- merge_bounds()
Merge the coil profile bounds into a single array.
- Returns:
list of bounds corresponding to each merged coils
- Return type:
list
- merge_coils(unshimmed, affine)
Uses the list of coil profiles to return a resampled concatenated list of coil profiles matching the unshimmed image. Bounds are also concatenated and returned.
- Parameters:
unshimmed (numpy.ndarray) -- 3d array of unshimmed volume
affine (numpy.ndarray) -- 4x4 array containing the affine transformation for the unshimmed array
- optimize(mask)
Optimize unshimmed volume by varying current to each channel
- Parameters:
mask (numpy.ndarray) -- 3D integer mask used for the optimizer (only consider voxels with non-zero values).
- Returns:
- Coefficients corresponding to the coil profiles that minimize the objective function.
The shape of the array returned has shape corresponding to the total number of channels
- Return type:
numpy.ndarray
- set_merged_bounds(merged_bounds)
Changes the default bounds set in the coil profile
- Parameters:
merged_bounds -- Concatenated coil profile bounds
- set_unshimmed(unshimmed, affine)
Set the unshimmed array to a new array. Resamples coil profiles accordingly.
- Parameters:
unshimmed (numpy.ndarray) -- 3d array of unshimmed volume
affine -- (numpy.ndarray): 4x4 array containing the qform affine transformation for the unshimmed array
Image manipulation
- shimmingtoolbox.image.concat_data(list_nii: List[Nifti1Image], axis=3, pixdim=None)
Concatenate data
- Parameters:
list_nii -- list of Nifti1Image
axis -- axis: 0, 1, 2, 3, 4.
pixdim -- pixel resolution to join to image header
- Returns:
concatenated image
- Return type:
ListNii
Numerical model
Create numerical model data for multi-echo B0 field mapping data
This module is for numerically simulating multi-echo B0 field mapping data. It considers features like: background B0 field, flip angle, echo time, and noise.
Typical usage example:
from shimmingtoolbox.simulate import *
b0_sim = NumericalModel(model="shepp-logan")
# Generate a background B0
b0_field = 13 # (Hz)
b0_sim.generate_deltaB0("linear", [0.0, b0_field])
# Simulate the signal data
FA = 15 # (degrees)
TE = [0.003, 0.015] # (seconds)
SNR = 50
b0_sim.simulate_measurement(FA, TE, SNR)
# Save simulation as NIfTI file (JSON sidecar also exported with parameters)
b0_sim.save('Phase', 'b0_mapping_data.nii', format='nifti')
- class shimmingtoolbox.simulate.numerical_model.NumericalModel(model=None, num_vox=128)
Multi-echo B0 field mapping data numerical simulator.
Simulate multi-echo B0 field mapping data in the presence of a B0 field. Can simulate data under ideal conditions or with noise. Export simulations in a NIfTI or
.mat
file formats.- gamma
Gyromagnetic ratio in rad * Hz / Tesla.
- Type:
float
- field_strength
Static field strength in Tesla.
- Type:
float
- handedness
Orientation of the cross-product for the Larmor equation. The value of this attribute is MRI vendor-dependent.
- measurement
Simulated measurement data array.
- proton_density
Default assumed brain proton density in %.
- T2_star
Default assumed brain T2* values in seconds at 3T.
- generate_deltaB0(field_type, params)
Generates a background B0 field.
Defines the starting volume. Sets the background B0 field to zeros.
- Parameters:
field_type -- Type of field to be generated. Available implementations are:
'linear'
.params -- List of parameters defining the field for the selected field type. If
field_type = 'linear'
, thenparams
are[m b]
where m (Hz/pixel) is the slope and b is the floor field (Hz).
- save(data_type, file_name, format=None)
Exports simulated data to a file with a JSON sidecar.
Resets the measurement class attribute to zero before simulating. Simulates the signal for each echo-time provided. If defined, adds noise to the complex simulated signal measurements using an SNR value.
- Parameters:
data_type -- Export data type. "Magnitude", "Phase", "Real", or "Imaginary".
file_name -- Filename of exported file, with or without file extension.
format -- File format for exported data. If no value given, will attempt to extract format from filename file extension, otherwise default to NIfTI.
- simulate_measurement(FA, TE, SNR=None)
Simulates a multi-echo measurement for field mapping
Resets the measurement class attribute to zero before simulating. Simulates the signal for each echo-time provided. If defined, adds noise to the complex simulated signal measurements using an SNR value.
- Parameters:
FA -- Flip angle in degrees.
TE -- Echo-times in seconds. Can be either a single value, list, or array.
SNR -- Signal-to-noise ratio used to define noise. If not set, no noise is added to the measurements.
Miscellaneous
- shimmingtoolbox.dicom_to_nifti.dicom_to_nifti(path_dicom, path_nifti, subject_id='sub-01', fname_config_dcm2bids='/home/docs/checkouts/readthedocs.org/user_builds/shimming-toolbox-py/checkouts/latest/config/dcm2bids.json', remove_tmp=False)
Converts dicom files into nifti files by calling dcm2bids
- Parameters:
path_dicom (str) -- Path to the input DICOM folder.
path_nifti (str) -- Path to the output NIfTI folder.
subject_id (str) -- Name of the imaged subject.
fname_config_dcm2bids (str) -- Path to the dcm2bids config JSON file.
remove_tmp (bool) -- If True, removes the tmp folder containing the NIfTI files created by dcm2niix.
- shimmingtoolbox.dicom_to_nifti.fix_tfl_b1(nii_b1, json_data)
Un-shuffles and rescales the magnitude and phase of complex B1+ maps acquired with Siemens' standard B1+ mapping sequence. Also computes a corrected affine matrix allowing the B1+ maps to be visualized in FSLeyes. :param nii_b1: Array of dimension (x, y, n_slices, 2*n_channels) as created by dcm2niix. :type nii_b1: numpy.ndarray :param json_data: Contains the different fields present in the json file corresponding to the nifti file. :type json_data: dict
- Returns:
NIfTI object containing the complex rescaled B1+ maps (x, y, n_slices, n_channels).
- Return type:
nib.Nifti1Image
- shimmingtoolbox.load_nifti.get_acquisition_times(nii_data, json_data)
Return the acquisition timestamps from a json sidecar. This assumes BIDS convention.
- Parameters:
nii_data (nibabel.Nifti1Image) -- Nibabel object containing the image timeseries.
json_data (dict) -- Json dict corresponding to a nifti sidecar.
- Returns:
Acquisition timestamps in ms.
- Return type:
numpy.ndarray
- shimmingtoolbox.load_nifti.load_nifti(path_data, modality='phase')
Load data from a directory containing NIFTI type file with nibabel. :param path_data: Path to the directory containing the file(s) to load :type path_data: str :param modality: Modality to read nifti (can be phase or magnitude) :type modality: str
- Returns:
List containing headers for every Nifti file dict: List containing all information in JSON format from every Nifti image numpy.ndarray: 5D array of all acquisition in time (x, y, z, echo, volume)
- Return type:
nibabel.Nifti1Image.Header
Note
If 'path' is a folder containing niftis, directly output niftis. It 'path' is a folder containing acquisitions, ask the user for which acquisition to use.
- shimmingtoolbox.load_nifti.read_nii(fname_nifti, auto_scale=True)
Reads a nifti file and returns the corresponding image and info. Also returns the associated json data. :param fname_nifti: direct path to the .nii or .nii.gz file that is going to be read :type fname_nifti: str :param auto_scale: Tells if scaling is done before return :type auto_scale:
bool
, optional- Returns:
Objet containing various data about the nifti file (returned by nibabel.load) json_data (dict): Contains the different fields present in the json file corresponding to the nifti file image (numpy.ndarray): For B0-maps, image contained in the nifti. Siemens phase images are rescaled between 0 and 2pi.
- Return type:
info (Nifti1Image)
- shimmingtoolbox.download.download_data(urls)
Download the binaries from a URL and return the destination filename Retry downloading if either server or connection errors occur on a SSL connection
- Parameters:
urls -- list of several urls (mirror servers) or single url (string)
- shimmingtoolbox.download.install_data(url, dest_folder, keep=False)
Download a data bundle from a URL and install in the destination folder.
- Parameters:
url -- URL or sequence thereof (if mirrors).
dest_folder -- destination directory for the data (to be created).
keep -- whether to keep existing data in the destination folder.
- Returns:
NoneType
Note
The function tries to be smart about the data contents.
Examples:
If the archive only contains a
README.md
, and the destination folder is${dst}
,${dst}/README.md
will be created. Note: an archive not containing a single folder is commonly known as a "tarbomb" because it puts files anywhere in the current working directory.If the archive contains a
${dir}/README.md
, and the destination folder is${dst}
,${dst}/README.md
will be created. Note: typically the package will be called${basename}-${revision}.zip
and contain a root folder named${basename}-${revision}/
under which all the other files will be located. The right thing to do in this case is to take the files from there and install them in${dst}
.Uses
download_data()
to retrieve the data.Uses
unzip()
to extract the bundle.
- shimmingtoolbox.download.unzip(compressed, dest_folder)
Extract compressed file to the
dest_folder
. Can handle.zip
,.tar.gz
. If none of this extension is found, simply copy the file indest_folder
.- Parameters:
compressed -- the compressed
.zip
or.tar.gz
filedest_folder -- the destination dir that expanded files are written to
- class shimmingtoolbox.pmu.PmuResp(fname_pmu)
PMU object containing the pressure values of a Siemens .resp file
- fname
Filename of the Siemens .resp file
- Type:
str
- data
Pressure values ranging from 0 to 4095
- Type:
numpy.ndarray
- start_time_mdh
Start time in milliseconds past midnight (mdh clock is expected to be the closest to the image header)
- Type:
int
- stop_time_mdh
Stop time in milliseconds past midnight (mdh clock is expected to be the closest to the image header)
- Type:
int
- start_time_mpcu
Start time in milliseconds past midnight
- Type:
int
- stop_time_mpcu
Stop time in milliseconds past midnight
- Type:
int
- interp_resp_trace(acquisition_times)
Interpolates
data
to the specifiedacquisition_times
- Parameters:
acquisition_times (numpy.ndarray) -- 1D array of the times in milliseconds past midnight of the desired times to interpolate the resp_trace. Times must be within
self.start_time_mdh
andself.stop_time_mdh
- Returns:
1D array with interpolated times
- Return type:
numpy.ndarray
- read_resp(fname_pmu)
Read a Siemens Physiological Log file. Returns a tuple with the logging data as numpy integer array and times in the form of milliseconds past midnight.
- Parameters:
fname_pmu -- Filename of the Siemens .resp file
- Returns:
A dict containing the
fname_pmu
infos. Contains the following keys:fname
data
start_time_mdh
stop_time_mdh
start_time_mpcu
stop_time_mpcu
- Return type:
dict
- shimmingtoolbox.utils.add_suffix(fname, suffix)
Add suffix between end of file name and extension.
- Parameters:
fname -- absolute or relative file name. Example:
t2.nii
suffix -- suffix. Example:
_mean
- Return:
file name string with suffix. Example:
t2_mean.nii
Examples:
add_suffix(t2.nii, _mean)
->t2_mean.nii
add_suffix(t2.nii.gz, a)
->t2a.nii.gz
- shimmingtoolbox.utils.create_fname_from_path(path, file_default)
Given a path, make sure it is not a directory, if it is add the default filename, if not, return the path
- Parameters:
path (str) -- filename or path to add the file_default to.
file_default (str) -- Name of the file + ext (example.nii.gz) to add to the path if the path is a directory.
- Returns:
Absolute path of a file
- Return type:
str
- shimmingtoolbox.utils.create_output_dir(path_output, is_file=False, output_folder_name='output')
Given a path, create the directory if it doesn't exist.
- Parameters:
path_output (str) -- Full path to either a folder or a file.
is_file (bool) -- True if the
path_output
is for a file, else False.output_folder_name (str) -- Name of sub-folder.
- shimmingtoolbox.utils.fill(data, invalid=None)
Replace the value of invalid 'data' cells (indicated by 'invalid') by the value of the nearest valid data cell
- Parameters:
data (numpy.ndarray)) -- array of any dimension
invalid (numpy.ndarray) -- a binary array of same shape as 'data'. True cells set where data value should be replaced. If None (default), use: invalid = np.isnan(data)
- Returns:
Return a filled array.
- Return type:
numpy.ndarray
- shimmingtoolbox.utils.iso_times_to_ms(iso_times)
Convert dicom acquisition times to ms
- Parameters:
iso_times (numpy.ndarray) -- 1D array of time strings from dicoms. Suported formats: "HHMMSS.mmmmmm" or "HH:MM:SS.mmmmmm"
- Returns:
1D array of times in milliseconds
- Return type:
numpy.ndarray
- shimmingtoolbox.utils.montage(X)
Concatenates images stored in a 3D array :param X: 3D array with the last dimension being the one in which the images are concatenated. :type X: numpy.ndarray
- Returns:
2D array of concatenated images.
- Return type:
numpy.ndarray
- shimmingtoolbox.utils.run_subprocess(cmd)
Wrapper for
subprocess.run()
.- Parameters:
cmd (list) -- list of arguments to be passed to the command line
- shimmingtoolbox.utils.save_nii_json(nii, json_data, fname_output)
Save the nii to a nifti file and dict to a json file.
- Parameters:
nii (nib.Nifti1Image) -- Nibabel object containing data save.
json_data (dict) -- Dictionary containing the json sidecar associated with the nibabel object.
fname_output (str) -- Output filename, supported types : '.nii', '.nii.gz'
- shimmingtoolbox.utils.set_all_loggers(verbose, list_exclude=('matplotlib',))
Set all loggers in the root manager to the verbosity level. Exclude any logger with the name in list_exclude
- Parameters:
verbose (str) -- Verbosity level: 'info', 'debug', 'warning', 'critical', 'error'
list_exclude -- List of string to exclude from logging
- shimmingtoolbox.utils.splitext(fname)
Split a fname (folder/file + ext) into a folder/file and extension.
Note: for .nii.gz the extension is understandably .nii.gz, not .gz (
os.path.splitext()
would want to do the latter, hence the special case).
- shimmingtoolbox.utils.st_progress_bar(*args, **kwargs)
Thin wrapper around tqdm.tqdm which checks SCT_PROGRESS_BAR muffling the progress bar if the user sets it to no, off, or false (case insensitive).
- shimmingtoolbox.utils.timeit(func)
Decorator to time a function. Decorate a function: @timeit on top of the function definition. The elapsed time will output in debug mode