Fix linear operator in regularized iterative SENSE

[ckolbPTB]

For the regularized iterative SENSE reconstruction we allow the following to be minimised:
1/2*||Ex - y||2^2 + reg_weight/2 * ||B x - x_reg||_2^2
with B a LinearOperator.

If we calculate the derivative and rearrange we get
(E^H E + reg_weight* B^H B) x = E^H y + B^H x_reg
(c) Andreas Kofler

In our current implementation we were missing the B^H part.

[github-actions]

Coverage Report

File	Stmts	Miss	Cover	Missing
src/mrpro/algorithms/csm
inati.py	24	1	96%	44
walsh.py	16	1	94%	34
src/mrpro/algorithms/dcf
dcf_voronoi.py	53	4	92%	15, 48–49, 76
src/mrpro/algorithms/optimizers
adam.py	20	1	95%	69
pdhg.py	79	3	96%	177–178, 184
src/mrpro/algorithms/reconstruction
DirectReconstruction.py	28	16	43%	51–71, 85
IterativeSENSEReconstruction.py	13	1	92%	76
Reconstruction.py	50	22	56%	42, 54–56, 80–87, 104–113
RegularizedIterativeSENSEReconstruction.py	41	17	59%	96–100, 114–139
src/mrpro/data
AcqInfo.py	128	3	98%	26, 169, 207
CsmData.py	29	3	90%	15, 82–84
DcfData.py	45	8	82%	18, 66, 78–83
IData.py	67	9	87%	119, 125, 129, 159–167
IHeader.py	75	7	91%	75, 109, 127–131
KData.py	214	25	88%	111–112, 127, 134, 144, 152, 206–207, 245, 250–251, 270–281, 435, 437, 500, 515, 552, 583, 592
KHeader.py	153	17	89%	25, 119–123, 150, 199, 210, 217–218, 221, 228, 260–271
KNoise.py	31	15	52%	39–52, 56–61
KTrajectory.py	81	12	85%	108–113, 116–118, 203–207
MoveDataMixin.py	140	18	87%	15, 113, 129, 143–145, 207, 323–325, 338, 417, 437–438, 440, 455–456, 458
QData.py	39	7	82%	42, 65–73
Rotation.py	674	35	95%	100, 198, 335, 433, 477, 495, 581, 583, 592, 626, 628, 691, 768, 773, 776, 791, 808, 813, 889, 1077, 1082, 1085, 1109, 1113, 1240, 1242, 1250–1251, 1315, 1397, 1690, 1846, 1881, 1885, 1996
SpatialDimension.py	232	21	91%	34, 104, 141, 148, 154, 274–276, 289–291, 325, 343, 356, 369, 382, 395, 404–405, 420, 429
acq_filters.py	12	1	92%	47
src/mrpro/data/traj_calculators
KTrajectoryCalculator.py	25	2	92%	23, 45
KTrajectoryIsmrmrd.py	13	2	85%	41, 50
KTrajectoryPulseq.py	23	1	96%	55
src/mrpro/operators
CartesianSamplingOp.py	89	3	97%	118, 157, 280
ConstraintsOp.py	60	2	97%	46, 48
EndomorphOperator.py	65	2	97%	228, 234
FiniteDifferenceOp.py	27	2	93%	40, 105
FourierOp.py	158	3	98%	263, 381, 386
Functional.py	71	5	93%	20–22, 117, 119
GridSamplingOp.py	136	9	93%	72–73, 82–83, 90–91, 94, 96, 98
LinearOperator.py	168	9	95%	55, 91, 190, 220, 261, 270, 278, 295, 320
LinearOperatorMatrix.py	158	16	90%	82, 119, 152, 161, 166, 175–178, 191–194, 203, 215, 304, 331, 359
MultiIdentityOp.py	13	2	85%	43, 48
Operator.py	78	2	97%	25, 74
ProximableFunctionalSeparableSum.py	39	3	92%	50, 103, 110
SliceProjectionOp.py	173	8	95%	44, 61, 63, 69, 206, 227, 260, 300
WaveletOp.py	119	5	96%	151, 169, 204, 209, 232
ZeroPadOp.py	16	1	94%	30
src/mrpro/utils
filters.py	62	2	97%	44, 49
reshape.py	60	1	98%	191
slice_profiles.py	46	6	87%	20, 36, 113–116, 149
sliding_window.py	34	1	97%	34
split_idx.py	10	2	80%	43, 47
summarize_tensorvalues.py	11	9	18%	20–29
typing.py	18	11	39%	8–23
zero_pad_or_crop.py	31	6	81%	26, 30, 54, 57, 60, 63
TOTAL	4969	362	93%

Tests	Skipped	Failures	Errors	Time
2300	0 💤	0 ❌	0 🔥	1m 9s ⏱️

[github-actions]

📚 Documentation

📁 Download as zip
🔍 View online

[koflera]

a few other comments:

• in the title

  mrpro/src/mrpro/algorithms/reconstruction/RegularizedIterativeSENSEReconstruction.py

  Line 1 in 4aa5d59

  """Regularized Iterative SENSE Reconstruction by adjoint Fourier transform."""

  , the "adjoint Fourier Transform" seems unnecessary?
• in the docstring

  mrpro/src/mrpro/algorithms/reconstruction/RegularizedIterativeSENSEReconstruction.py

  Line 27 in 4aa5d59

  :math:`H x = b` with :math:`H = A^H W A + L B` and :math:`b = A^H W y + L x_0` where :math:`A`

  it should be :math:H x = b with :math:H = A^H W A + L B^H B and :math:b = A^H W y + L B^H x_0.
• I would also maybe change the capital L to l, to better distinguish scalars from linear operators. currently, the code would not make sense anyway for an entire regularization parameter map. Also, if with capital L it was implicitly intended to be able to have an entire reg parameter map, the docstring is wrong, cause it then would need to be :math:min_x \frac{1}{2}||W^\frac{1}{2} (Ax - y)||_2^2 + \frac{1}{2}||L(Bx - x_0)||_2^2, for which then the solution of the problem is however slightly different than what is currently implemented.
• maybe rename precondition_op to dcf_op in

  mrpro/src/mrpro/algorithms/reconstruction/RegularizedIterativeSENSEReconstruction.py

  Line 120 in 4aa5d59

  precondition_op = self.dcf.as_operator() if self.dcf is not None else IdentityOp()

  ?

[ckolbPTB]

I would also maybe change the capital L to l, to better distinguish scalars from linear operators

@koflera should the regularization strength then actually always be a float and not float | torch.Tensor?

[fzimmermann89]

please keep it as tensor or float.

[koflera]

I would also maybe change the capital L to l, to better distinguish scalars from linear operators

@koflera should the regularization strength then actually always be a float and not float | torch.Tensor?

I think it doesn't hurt to keep the option to have both. If the intention was to make it possible to have an entire Lambda-Map, I would say we can still update the code later on.

[koflera]

very last comment :)

[fzimmermann89]

Why not \lambda?

It is a latex equation.
The parameter is named differently anyways..?