Skip to content

Commit

Permalink
Merge branch 'cd_pipeline' of https://github.com/PTB-MR/mrpro into cd…
Browse files Browse the repository at this point in the history
…_pipeline
  • Loading branch information
schuenke committed Oct 11, 2024
2 parents 156e7de + d51e2aa commit e7d0b45
Show file tree
Hide file tree
Showing 33 changed files with 75 additions and 60 deletions.
17 changes: 13 additions & 4 deletions .pre-commit-config.yaml
Original file line number Diff line number Diff line change
Expand Up @@ -14,27 +14,27 @@ repos:
- id: mixed-line-ending

- repo: https://github.com/astral-sh/ruff-pre-commit
rev: v0.4.7
rev: v0.6.9
hooks:
- id: ruff # linter
args: [--fix]
- id: ruff-format # formatter

- repo: https://github.com/crate-ci/typos
rev: v1.21.0
rev: v1.25.0
hooks:
- id: typos

- repo: https://github.com/pre-commit/mirrors-mypy
rev: v1.10.0
rev: v1.11.2
hooks:
- id: mypy
pass_filenames: false
always_run: true
args: [src, tests, examples]
additional_dependencies:
- numpy
- torch>=2.3.0
- torch>=2.4.0
- types-requests
- einops
- pydicom
Expand All @@ -43,3 +43,12 @@ repos:
- xsdata
- "--index-url=https://download.pytorch.org/whl/cpu"
- "--extra-index-url=https://pypi.python.org/simple"
ci:
autofix_commit_msg: |
[pre-commit] auto fixes from pre-commit hooks
autofix_prs: false
autoupdate_branch: ''
autoupdate_commit_msg: '[pre-commit] pre-commit autoupdate'
autoupdate_schedule: monthly
skip: [mypy]
submodules: false
12 changes: 6 additions & 6 deletions README.md
Original file line number Diff line number Diff line change
Expand Up @@ -19,18 +19,18 @@ MR image reconstruction and processing package specifically developed for PyTorc

- **ISMRMRD support** MRpro supports [ismrmrd-format](https://ismrmrd.readthedocs.io/en/latest/) for MR raw data.
- **PyTorch** All data containers utilize PyTorch tensors to ensure easy integration in PyTorch-based network schemes.
- **Cartesian and non-Cartesian trajectories** MRpro can reconstruct data obtained with Cartesian and non-Cartesian (e.g. radial, spiral...) sapling schemes. MRpro automatically detects if FFT or nuFFT is required to reconstruction the k-space data.
- **Cartesian and non-Cartesian trajectories** MRpro can reconstruct data obtained with Cartesian and non-Cartesian (e.g. radial, spiral...) sapling schemes. MRpro automatically detects if FFT or nuFFT is required to reconstruct the k-space data.
- **Pulseq support** If the data acquisition was carried out using a [pulseq-based](http://pulseq.github.io/) sequence, the seq-file can be provided to MRpro and the used trajectory is automatically calculated.
- **Signal models** A range of different MR signal models are implemented (e.g. T1 recovery, WASABI).
- **Regularised image reconstruction** Regularised image reconstruction algorithms including Wavelet-based compressed sensing or total variation regularised image reconstruction are available.
- **Regularized image reconstruction** Regularized image reconstruction algorithms including Wavelet-based compressed sensing or total variation regularized image reconstruction are available.

## Examples

In the following we show some code snippets to highlight the use of MRpro. Each code snippet only shows the main steps. A complete working notebook can be found in the provided link.
In the following, we show some code snippets to highlight the use of MRpro. Each code snippet only shows the main steps. A complete working notebook can be found in the provided link.

### Simple reconstruction

Read in the data and trajectoy and reconstruct an image by applying a density compensation function and then the adjoint of the Fourier operator and the adjoint of the coil sensitivity operator.
Read the data and trajectory and reconstruct an image by applying a density compensation function and then the adjoint of the Fourier operator and the adjoint of the coil sensitivity operator.

```python
# Read the trajectory from the ISMRMRD file
Expand All @@ -46,7 +46,7 @@ Full example: <https://github.com/PTB-MR/mrpro/blob/main/examples/direct_reconst

### Estimate quantitative parameters

Quantitative parameter maps can be obtained by creating a functional to be minimised and calling a non-linear solver such as ADAM.
Quantitative parameter maps can be obtained by creating a functional to be minimized and calling a non-linear solver such as ADAM.

```python
# Define signal model
Expand Down Expand Up @@ -87,4 +87,4 @@ We recommend to use [Microsoft Visual Studio Code](https://code.visualstudio.com

### Style

Please have a look at our [contributor guide](https://ptb-mr.github.io/mrpro/contributor_guide.html) for more information on the structure of the repository, naming conventions and other useful information.
Please look at our [contributor guide](https://ptb-mr.github.io/mrpro/contributor_guide.html) for more information on the repository structure, naming conventions, and other useful information.
2 changes: 1 addition & 1 deletion examples/iterative_sense_reconstruction.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -59,7 +59,7 @@
"\n",
"$ y = Ax + n $\n",
"\n",
"where $n$ describes complex Gaussian noise. The image $x$ can be obtained by minimising the functionl $F$\n",
"where $n$ describes complex Gaussian noise. The image $x$ can be obtained by minimizing the functional $F$\n",
"\n",
"$ F(x) = ||W^{\\frac{1}{2}}(Ax - y)||_2^2 $\n",
"\n",
Expand Down
2 changes: 1 addition & 1 deletion examples/iterative_sense_reconstruction.py
Original file line number Diff line number Diff line change
Expand Up @@ -26,7 +26,7 @@
#
# $ y = Ax + n $
#
# where $n$ describes complex Gaussian noise. The image $x$ can be obtained by minimising the functionl $F$
# where $n$ describes complex Gaussian noise. The image $x$ can be obtained by minimizing the functional $F$
#
# $ F(x) = ||W^{\frac{1}{2}}(Ax - y)||_2^2 $
#
Expand Down
2 changes: 1 addition & 1 deletion src/mrpro/algorithms/prewhiten_kspace.py
Original file line number Diff line number Diff line change
Expand Up @@ -21,7 +21,7 @@ def prewhiten_kspace(kdata: KData, knoise: KNoise, scale_factor: float | torch.T
More information can be found in [ISMa]_ [HAN2014]_ [ROE1990]_.
If the the data has more samples in the 'other'-dimensions (batch/slice/...),
If the data has more samples in the 'other'-dimensions (batch/slice/...),
the noise covariance matrix is calculated jointly over all samples.
Thus, if the noise is not stationary, the noise covariance matrix is not accurate.
In this case, the function should be called for each sample separately.
Expand Down
2 changes: 1 addition & 1 deletion src/mrpro/data/Rotation.py
Original file line number Diff line number Diff line change
Expand Up @@ -413,7 +413,7 @@ def from_euler(cls, seq: str, angles: torch.Tensor | NestedSequence[float] | flo
chosen to be the basis vectors.
The three rotations can either be in a global frame of reference
(extrinsic) or in a body centred frame of reference (intrinsic), which
(extrinsic) or in a body centered frame of reference (intrinsic), which
is attached to, and moves with, the object under rotation [EULa]_.
Parameters
Expand Down
1 change: 1 addition & 0 deletions src/mrpro/data/_kdata/KData.py
Original file line number Diff line number Diff line change
Expand Up @@ -12,6 +12,7 @@
import numpy as np
import torch
from einops import rearrange

from mrpro.data._kdata.KDataRearrangeMixin import KDataRearrangeMixin
from mrpro.data._kdata.KDataRemoveOsMixin import KDataRemoveOsMixin
from mrpro.data._kdata.KDataSelectMixin import KDataSelectMixin
Expand Down
1 change: 1 addition & 0 deletions src/mrpro/data/_kdata/KDataProtocol.py
Original file line number Diff line number Diff line change
Expand Up @@ -3,6 +3,7 @@
from typing import Literal, Protocol, Self

import torch

from mrpro.data.KHeader import KHeader
from mrpro.data.KTrajectory import KTrajectory

Expand Down
1 change: 1 addition & 0 deletions src/mrpro/data/_kdata/KDataRearrangeMixin.py
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,7 @@
from typing import Self

from einops import rearrange

from mrpro.data._kdata.KDataProtocol import _KDataProtocol
from mrpro.data.AcqInfo import AcqInfo
from mrpro.utils import modify_acq_info
Expand Down
1 change: 1 addition & 0 deletions src/mrpro/data/_kdata/KDataRemoveOsMixin.py
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,7 @@
from typing import Self

import torch

from mrpro.data._kdata.KDataProtocol import _KDataProtocol
from mrpro.data.KTrajectory import KTrajectory

Expand Down
1 change: 1 addition & 0 deletions src/mrpro/data/_kdata/KDataSelectMixin.py
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,7 @@
from typing import Literal, Self

import torch

from mrpro.data._kdata.KDataProtocol import _KDataProtocol
from mrpro.utils import modify_acq_info

Expand Down
1 change: 1 addition & 0 deletions src/mrpro/data/_kdata/KDataSplitMixin.py
Original file line number Diff line number Diff line change
Expand Up @@ -5,6 +5,7 @@

import torch
from einops import rearrange, repeat

from mrpro.data._kdata.KDataProtocol import _KDataProtocol
from mrpro.data.EncodingLimits import Limits
from mrpro.utils import modify_acq_info
Expand Down
2 changes: 1 addition & 1 deletion src/mrpro/operators/CartesianSamplingOp.py
Original file line number Diff line number Diff line change
Expand Up @@ -12,7 +12,7 @@
class CartesianSamplingOp(LinearOperator):
"""Cartesian Sampling Operator.
Selects points on a Cartisian grid based on the the k-space trajectory.
Selects points on a Cartesian grid based on the k-space trajectory.
Thus, the adjoint sorts the data into regular Cartesian sampled grid based on the k-space
trajectory. Non-acquired points are zero-filled.
"""
Expand Down
2 changes: 1 addition & 1 deletion src/mrpro/operators/EinsumOp.py
Original file line number Diff line number Diff line change
Expand Up @@ -37,7 +37,7 @@ class EinsumOp(LinearOperator):
of :math:`N` vectors :math:`x1, x2, ..., xN`. Then, the operation defined by
:math:`A @ x: = diag(A1, A2,..., AN) * [x1, x2, ..., xN]^T`
can be implemented by the einsum rule ``"... i j, ... j -> ... i"``.
This is the default behaviour of the operator.
This is the default behavior of the operator.
"""

def __init__(self, matrix: torch.Tensor, einsum_rule: str = '... i j, ... j -> ... i') -> None:
Expand Down
10 changes: 5 additions & 5 deletions src/mrpro/operators/LinearOperator.py
Original file line number Diff line number Diff line change
Expand Up @@ -55,7 +55,7 @@ def operator_norm(
dim
the dimensions of the tensors on which the operator operates.
For example, for a matrix-vector multiplication example, a batched matrix tensor with shape (4,30,80,160),
input tensors of shape (4,30,160) to be multiplied, and dim = None, it is understood that the the
input tensors of shape (4,30,160) to be multiplied, and dim = None, it is understood that the
matrix representation of the operator corresponds to a block diagonal operator (with 4*30 matrices)
and thus the algorithm returns a tensor of shape (1,1,1) containing one single value.
In contrast, if for example, dim=(-1,), the algorithm computes a batched operator
Expand Down Expand Up @@ -132,7 +132,7 @@ def operator_norm(

return op_norm

@overload # type: ignore[override]
@overload
def __matmul__(self, other: LinearOperator) -> LinearOperator: ...

@overload
Expand Down Expand Up @@ -179,7 +179,7 @@ def __mul__(self, other: torch.Tensor) -> LinearOperator:
"""
return LinearOperatorElementwiseProductLeft(self, other)

def __rmul__(self, other: torch.Tensor) -> LinearOperator: # type: ignore[misc]
def __rmul__(self, other: torch.Tensor) -> LinearOperator:
"""Operator elementwise right multiplication with tensor.
Returns lambda x: other*self(x)
Expand Down Expand Up @@ -213,7 +213,7 @@ class LinearOperatorElementwiseProductRight(
):
"""Operator elementwise right multiplication with a tensor.
Peforms Tensor*LinearOperator(x)
Performs Tensor*LinearOperator(x)
"""

def adjoint(self, x: torch.Tensor) -> tuple[torch.Tensor,]:
Expand All @@ -226,7 +226,7 @@ class LinearOperatorElementwiseProductLeft(
):
"""Operator elementwise left multiplication with a tensor.
Peforms LinearOperator(Tensor*x)
Performs LinearOperator(Tensor*x)
"""

def adjoint(self, x: torch.Tensor) -> tuple[torch.Tensor,]:
Expand Down
4 changes: 2 additions & 2 deletions src/mrpro/operators/Operator.py
Original file line number Diff line number Diff line change
Expand Up @@ -45,7 +45,7 @@ def __mul__(self, other: torch.Tensor) -> Operator[*Tin, Tout]:
"""
return OperatorElementwiseProductLeft(self, other)

def __rmul__(self, other: torch.Tensor) -> Operator[*Tin, Tout]: # type: ignore[misc]
def __rmul__(self, other: torch.Tensor) -> Operator[*Tin, Tout]:
"""Operator multiplication with tensor.
Returns lambda x: other*self(x)
Expand Down Expand Up @@ -94,7 +94,7 @@ def forward(self, *args: *Tin) -> Tout:
class OperatorElementwiseProductRight(Operator[*Tin, Tout]):
"""Operator elementwise right multiplication with a tensor.
Peforms Tensor*Operator(x)
Performs Tensor*Operator(x)
"""

def __init__(self, operator: Operator[*Tin, Tout], tensor: torch.Tensor):
Expand Down
8 changes: 4 additions & 4 deletions src/mrpro/operators/SliceProjectionOp.py
Original file line number Diff line number Diff line change
Expand Up @@ -396,12 +396,12 @@ def projection_matrix(
# by adding all possible combinations of 0 and 1 to the point and flooring
offsets = torch.tensor(list(itertools.product([0, 1], repeat=3)))
# all points that influence a pixel
# x,y,8-neighbours,(2*w+1)-raylength,3-dimensions input_shape.xinput_shape.yinput_shape.z)
# x,y,8-neighbors,(2*w+1)-raylength,3-dimensions input_shape.xinput_shape.yinput_shape.z)
points_influencing_pixel = (
einops.rearrange(pixel_rotated_y_x_zyx, ' y x zyxdim -> y x 1 1 zyxdim')
+ einops.rearrange(ray, ' ray zyxdim -> 1 1 1 ray zyxdim')
+ einops.rearrange(offsets, ' neighbours zyxdim -> 1 1 neighbours 1 zyxdim')
).floor() # y x neighbours ray zyx
+ einops.rearrange(offsets, ' neighbors zyxdim -> 1 1 neighbors 1 zyxdim')
).floor() # y x neighbors ray zyx
# directional distance in source volume coordinate system
distance = pixel_rotated_y_x_zyx[:, :, None, None, :] - points_influencing_pixel
# Inverse rotation projects this back to the original coordinate system, i.e
Expand All @@ -415,7 +415,7 @@ def projection_matrix(

source = einops.rearrange(
points_influencing_pixel,
'y x neighbours raylength zyxdim -> (y x) (neighbours raylength) zyxdim',
'y x neighbors raylength zyxdim -> (y x) (neighbors raylength) zyxdim',
).int()

# mask of only potential source points inside the source volume
Expand Down
2 changes: 1 addition & 1 deletion src/mrpro/operators/WaveletOp.py
Original file line number Diff line number Diff line change
Expand Up @@ -52,7 +52,7 @@ def __init__(
For complex images the wavelet coefficients are calculated for real and imaginary part separately.
For a 2D image, the coefficients are labelled [aa, (ad_n, da_n, dd_n), ..., (ad_1, da_1, dd_1)] where a refers
For a 2D image, the coefficients are labeled [aa, (ad_n, da_n, dd_n), ..., (ad_1, da_1, dd_1)] where a refers
to the approximation coefficients and d to the detail coefficients. The index indicates the level.
Parameters
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -9,7 +9,7 @@ class TransientSteadyStateWithPreparation(SignalModel[torch.Tensor, torch.Tensor
r"""Signal model for transient steady state.
This signal model describes the behavior of the longitudinal magnetization during continuous acquisition after
a preparation pulse. The effect of the preparation pulse is modelled by a scaling factor applied to the
a preparation pulse. The effect of the preparation pulse is modeled by a scaling factor applied to the
equilibrium magnetization. A delay after the preparation pulse can be defined. During this time T1 relaxation to M0
occurs. Data acquisition starts after this delay. Perfect spoiling is assumed and hence T2 effects are not
considered in the model. In addition this model assumes :math:`TR << T1` and :math:`TR << T1^*`
Expand All @@ -18,20 +18,20 @@ class TransientSteadyStateWithPreparation(SignalModel[torch.Tensor, torch.Tensor
Let's assume we want to describe a continuous acquisition after an inversion pulse, then we have three parts:
[Part A: 180° inversion pulse][Part B: spoiler gradient][Part C: Continuous data acquisition]
- Part A: The 180° pulse leads to an inversion of the equilibrium magnetisation (:math:`M_0`) to :math:`-M_0`.
- Part A: The 180° pulse leads to an inversion of the equilibrium magnetization (:math:`M_0`) to :math:`-M_0`.
This can be described by setting the scaling factor ``m0_scaling_preparation`` to -1.
- Part B: Commonly after an inversion pulse a strong spoiler gradient is played out to compensate for non-perfect
inversion. During this time the magnetisation :math:`M_z(t)` follows the signal model:
inversion. During this time the magnetization :math:`M_z(t)` follows the signal model:
:math:`M_z(t) = M_0 + (s * M_0 - M_0)e^{(-t / T1)}` where :math:`s` is ``m0_scaling_preparation``.
- Part C: After the spoiler gradient the data acquisition starts and the magnetisation :math:`M_z(t)` can be
- Part C: After the spoiler gradient the data acquisition starts and the magnetization :math:`M_z(t)` can be
described by the signal model: :math:`M_z(t) = M_0^* + (M_{init} - M_0^*)e^{(-t / T1^*)}` where the initial
magnetisation is :math:`M_{init} = M_0 + (s*M_0 - M_0)e^{(-\Delta t / T1)}` where :math:`s` is
magnetization is :math:`M_{init} = M_0 + (s*M_0 - M_0)e^{(-\Delta t / T1)}` where :math:`s` is
``m0_scaling_preparation`` and :math:`\Delta t` is ``delay_after_preparation``. The effective longitudinal
relaxation time is :math:`T1^* = 1/(1/T1 - ln(cos(\alpha)/TR)`
where :math:`TR` is ``repetition_time`` and :math:`\alpha` is ``flip_angle``.
The steady-state magnetisation is :math:`M_0^* = M_0 T1^* / T1`.
The steady-state magnetization is :math:`M_0^* = M_0 T1^* / T1`.
References
----------
Expand Down
2 changes: 1 addition & 1 deletion tests/algorithms/test_prewhiten_kspace.py
Original file line number Diff line number Diff line change
Expand Up @@ -44,7 +44,7 @@ def test_prewhiten_kspace_cpu(random_kheader):
_test_prewhiten_kspace(random_kheader, 'cpu')


@pytest.mark.cuda()
@pytest.mark.cuda
def test_prewhiten_kspace_cuda(random_kheader):
"""Prewhitening of k-space data on GPU."""
_test_prewhiten_kspace(random_kheader, 'cuda')
4 changes: 2 additions & 2 deletions tests/conftest.py
Original file line number Diff line number Diff line change
Expand Up @@ -158,7 +158,7 @@ def random_full_ismrmrd_header(request) -> xsd.ismrmrdschema.ismrmrdHeader:
)


@pytest.fixture()
@pytest.fixture
def random_ismrmrd_file(random_acquisition, random_noise_acquisition, full_header):
with tempfile.NamedTemporaryFile(suffix='.h5') as file:
dataset = ismrmrd.Dataset(file.name)
Expand All @@ -184,7 +184,7 @@ def random_kheader(request, random_full_ismrmrd_header, random_acq_info):
return kheader


@pytest.fixture()
@pytest.fixture
def random_acq_info(random_acquisition):
"""Random (not necessarily valid) AcqInfo."""
acq_info = AcqInfo.from_ismrmrd_acquisitions([random_acquisition])
Expand Down
2 changes: 1 addition & 1 deletion tests/data/test_csm_data.py
Original file line number Diff line number Diff line change
Expand Up @@ -33,7 +33,7 @@ def test_CsmData_smoothing_width(csm_method, ellipse_phantom, random_kheader):
assert torch.equal(csm_using_spatial_dimension.data, csm_using_int.data)


@pytest.mark.cuda()
@pytest.mark.cuda
@pytest.mark.parametrize('csm_method', [CsmData.from_idata_walsh, CsmData.from_idata_inati])
def test_CsmData_cuda(csm_method, ellipse_phantom, random_kheader):
"""CsmData obtained on GPU in CUDA memory."""
Expand Down
2 changes: 1 addition & 1 deletion tests/data/test_dcf_data.py
Original file line number Diff line number Diff line change
Expand Up @@ -75,7 +75,7 @@ def test_dcf_spiral_traj_voronoi_singlespiral():
torch.testing.assert_close(dcf_three_broadcast.data[..., 1, :-ignore_last], dcf_single.data[..., 0, :-ignore_last])


@pytest.mark.cuda()
@pytest.mark.cuda
@pytest.mark.parametrize(('n_kr', 'n_ka', 'n_k0'), [(10, 6, 20)])
def test_dcf_rpe_traj_voronoi_cuda(n_kr, n_ka, n_k0):
"""Voronoi-based dcf calculation for RPE trajectory in CUDA memory."""
Expand Down
4 changes: 2 additions & 2 deletions tests/data/test_idata.py
Original file line number Diff line number Diff line change
Expand Up @@ -76,15 +76,15 @@ def test_IData_to_complex128(random_kheader, random_test_data):
assert idata_complex128.data.dtype == torch.complex128


@pytest.mark.cuda()
@pytest.mark.cuda
def test_IData_cuda(random_kheader, random_test_data):
"""Move IData object to CUDA memory."""
idata = IData.from_tensor_and_kheader(data=random_test_data, kheader=random_kheader)
idata_cuda = idata.cuda()
assert idata_cuda.data.is_cuda


@pytest.mark.cuda()
@pytest.mark.cuda
def test_IData_cpu(random_kheader, random_test_data):
"""Move IData object to CUDA memory and back to CPU memory."""
idata = IData.from_tensor_and_kheader(data=random_test_data, kheader=random_kheader)
Expand Down
Loading

0 comments on commit e7d0b45

Please sign in to comment.