Skip to content
Snippets Groups Projects
Commit 65249c41 authored by Julius Steiglechner's avatar Julius Steiglechner
Browse files

Decide whether to use absolute or complex values.

parent 7f3f61d2
No related branches found
No related tags found
No related merge requests found
"""This module provides methods to denoise mp2rage MRI.
"""
Provide methods to denoise mp2rage MRI.
adapted from RobustCombination function from Jose Marques,
adapted from RobustCombination function from Jose Marques,
https://github.com/JosePMarques/MP2RAGE-related-scripts
Applications:
Applications:
https://journals.plos.org/plosone/article?id=10.1371/journal.pone.0099676
https://journals.plos.org/plosone/article?id=10.1371/journal.pone.0210803
"""
......@@ -18,45 +19,74 @@ from .unite_mp2rage import unite_mp2rage_robust
def quadratic_formula_positive(a, b, c):
"""Give positive solution of quadratic formula."""
return (-b + np.sqrt(b**2 - 4 * a * c)) / (2 * a)
def quadratic_formula_negative(a, b, c):
"""Give negative solution of quadratic formula."""
return (-b - np.sqrt(b**2 - 4 * a * c)) / (2 * a)
def combine_mp2rage_robust(uniform_filepath: Union[str, Path],
inv1_fielpath: Union[str, Path],
inv2_filepath: Union[str, Path],
denoised_output_filepath: Union[str, Path],
regularization_factor: float = 6) -> None:
"""Suppress background noise inside and outside the skull.
Args:
uniform_filepath (Union[str, Path]): [description]
inv1_fielpath (Union[str, Path]): [description]
inv2_filepath (Union[str, Path]): [description]
denoised_output_filepath (Union[str, Path]): [description]
regularization_factor (float, optional): [description]. Defaults to 6.
Note:
usually the regularization factor shouldn't be greater then 10, but that
is not the ase when the image is bias field corrected, in which case the
noise estimated at the edge of the image might not be such a good measure
def combine_mp2rage_robust(
uniform_filepath: Union[str, Path],
inv1_fielpath: Union[str, Path],
inv2_filepath: Union[str, Path],
denoised_output_filepath: Union[str, Path],
regularization_factor: float = 6,
value_mode: str = "absolute",
) -> None:
"""
Suppress background noise inside and outside the skull.
Note
----
usually the regularization factor shouldn't be greater then 10, but
that is not the case when the image is bias field corrected, in which
case the noise estimated at the edge of the image might not be such a
good measure.
Parameters
----------
uniform_filepath : Union[str, Path]
DESCRIPTION.
inv1_fielpath : Union[str, Path]
DESCRIPTION.
inv2_filepath : Union[str, Path]
DESCRIPTION.
denoised_output_filepath : Union[str, Path]
DESCRIPTION.
regularization_factor : float, optional
DESCRIPTION. The default is 6.
value_mode : str, optional
DESCRIPTION. The default is 'absolute'.
Returns
-------
None
DESCRIPTION.
"""
if value_mode not in ["absolute", "complex"]:
raise ValueError("Value mode has to be absolute or complex.")
mp2rage_nii = nib.load(uniform_filepath)
inv1_nii = nib.load(inv1_fielpath)
inv2_nii = nib.load(inv2_filepath)
mp2rage_img = mp2rage_nii.get_fdata()
inv1_img = inv1_nii.get_fdata()
inv2_img = inv2_nii.get_fdata()
mp2rage_img = mp2rage_nii.dataobj[:]
if value_mode == "complex":
inv1_img = inv1_nii.dataobj[:]
inv2_img = inv2_nii.dataobj[:]
elif value_mode == "absolute":
inv1_img = np.abs(inv1_nii.dataobj[:])
inv2_img = np.abs(inv2_nii.dataobj[:])
# rescale the range of features to scale the range[-.5, .5]
if mp2rage_img.min() >= 0 and mp2rage_img.max() >= 0.51:
# assumes that it is getting only positive values
mp2rage_img = (mp2rage_img - mp2rage_img.max() / 2) / mp2rage_img.max()
integerformat = 1
integerformat = False
elif mp2rage_img.min() < 0:
# rescale min max
min_mp2rage = np.min(mp2rage_img)
......@@ -65,9 +95,9 @@ def combine_mp2rage_robust(uniform_filepath: Union[str, Path],
min_mp2rage) - .5
inv1_img = inv1_img / (max_mp2rage - min_mp2rage) * 2
inv2_img = inv2_img / (max_mp2rage - min_mp2rage) * 2
integerformat = 1
integerformat = False
else:
integerformat = 0
integerformat = True
# computes correct polarity to inv1 image
inv1_img = np.sign(mp2rage_img) * inv1_img
......@@ -101,7 +131,7 @@ def combine_mp2rage_robust(uniform_filepath: Union[str, Path],
mp2rage_robust = unite_mp2rage_robust(inv1_estimated, inv2_img,
noiselevel**2)
if integerformat == 0:
if integerformat:
denoised_mp2rage_img = mp2rage_robust
else:
denoised_mp2rage_img = np.round(4095 * (mp2rage_robust + 0.5))
......
......@@ -3,11 +3,13 @@
import numpy as np
def unite_mp2rage_robust(TI1_img: np.ndarray,
TI2_img: np.ndarray,
beta: float = 1e-5) -> np.ndarray:
def unite_mp2rage_robust(
TI1_img: np.ndarray,
TI2_img: np.ndarray,
beta: float = 1e-5,
) -> np.ndarray:
"""Combine the two images in a robust way.
See also:
https://archive.ismrm.org/2013/0269.html
https://journals.plos.org/plosone/article?id=10.1371/journal.pone.0099676#pone.0099676.e004
......@@ -25,4 +27,4 @@ def unite_mp2rage(TI1_img: np.ndarray,
mp2rage = np.real(np.conj(TI1_img) * TI2_img) / np.clip(
sum_imgs, a_min=beta, a_max=np.max(sum_imgs))
return mp2rage
\ No newline at end of file
return mp2rage
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment