PTB-MR / mrpro

MR image reconstruction and processing.
https://ptb-mr.github.io/mrpro/
Apache License 2.0
18 stars 2 forks source link

FourierOp Norm #542

Open fzimmermann89 opened 1 week ago

fzimmermann89 commented 1 week ago

At least two potential issues with Fourier Op:

Currently, the answer two both is "No" This adds the tests showcasing the potential problems

github-actions[bot] commented 1 week ago

Coverage

Coverage Report
FileStmtsMissCoverMissing
src/mrpro/algorithms/csm
   inati.py24196%44
   walsh.py16194%34
src/mrpro/algorithms/dcf
   dcf_voronoi.py53492%15, 48–49, 76
src/mrpro/algorithms/optimizers
   adam.py20195%69
src/mrpro/algorithms/reconstruction
   DirectReconstruction.py281643%51–71, 85
   IterativeSENSEReconstruction.py13192%76
   Reconstruction.py502256%42, 54–56, 80–87, 104–113
   RegularizedIterativeSENSEReconstruction.py411759%96–100, 114–139
src/mrpro/data
   AcqInfo.py128398%26, 169, 207
   CsmData.py29390%15, 82–84
   DcfData.py45882%18, 66, 78–83
   IData.py67987%119, 125, 129, 159–167
   IHeader.py75791%75, 109, 127–131
   KHeader.py1531789%25, 119–123, 150, 199, 210, 217–218, 221, 228, 260–271
   KNoise.py311552%39–52, 56–61
   KTrajectory.py69593%178–182
   MoveDataMixin.py1371887%15, 113, 129, 143–145, 207, 305–307, 320, 399, 419–420, 422, 437–438, 440
   QData.py39782%42, 65–73
   Rotation.py6743595%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.py2322291%33, 103, 128, 135, 141, 261–263, 276–278, 312, 330, 343, 356, 369, 382, 391–392, 407, 416, 420
   acq_filters.py12192%47
src/mrpro/data/_kdata
   KData.py1341887%109–110, 125, 132, 142, 150, 204–205, 243, 248–249, 268–279
   KDataRemoveOsMixin.py29293%44, 46
   KDataSelectMixin.py19289%48, 63
   KDataSplitMixin.py48394%53, 84, 93
src/mrpro/data/traj_calculators
   KTrajectoryCalculator.py25292%23, 45
   KTrajectoryIsmrmrd.py13285%41, 50
   KTrajectoryPulseq.py29197%54
src/mrpro/operators
   CartesianSamplingOp.py89298%118, 280
   ConstraintsOp.py60297%46, 48
   DensityCompensationOp.py14379%32–34
   EndomorphOperator.py65297%228, 234
   FiniteDifferenceOp.py27293%40, 105
   FourierOp.py160398%266, 384, 389
   Functional.py71593%20–22, 117, 119
   GridSamplingOp.py136993%72–73, 82–83, 90–91, 94, 96, 98
   LinearOperator.py1711293%55, 91, 190, 220, 261, 270, 278, 287, 295, 320, 418, 423
   LinearOperatorMatrix.py1581690%82, 119, 152, 161, 166, 175–178, 191–194, 203, 215, 304, 331, 359
   MultiIdentityOp.py13285%43, 48
   Operator.py78297%25, 74
   ProximableFunctionalSeparableSum.py39392%50, 103, 110
   SliceProjectionOp.py173895%44, 61, 63, 69, 206, 227, 260, 300
   WaveletOp.py120596%152, 170, 205, 210, 233
   ZeroPadOp.py16194%30
src/mrpro/utils
   filters.py62297%44, 49
   slice_profiles.py46687%20, 36, 113–116, 149
   sliding_window.py34197%34
   split_idx.py10280%43, 47
   summarize_tensorvalues.py11918%20–29
   typing.py181139%8–23
   zero_pad_or_crop.py31681%26, 30, 54, 57, 60, 63
TOTAL487735793% 

Tests Skipped Failures Errors Time
2293 0 :zzz: 18 :x: 0 :fire: 2m 46s :stopwatch:
github-actions[bot] commented 1 week ago

:books: Documentation

:file_folder: Download as zip :mag: View online

fzimmermann89 commented 1 week ago

Note to myself: wondering if this could also cause issues with self-supervised training.

Should check if F(subset of trajectory)(image) = F(image)[subset] -- which might also be a good test.

fzimmermann89 commented 1 week ago

Also, what was the outcome of #381 ?

fzimmermann89 commented 1 week ago

@ckolbPTB @schuenke

ckolbPTB commented 6 days ago

I thought to norm of FourierOp should always be 1?

I don't think this holds for arbitrary FourierOps. If a FourierOp only yields an undersampled k-space then the Norm will not be 1. Also if padding or cropping of image-space/k-space is required (mismatch of recon_matrix and encoding_matrix) the norm will not be 1.

Here I verified that FFT and NUFFT lead to the same result: https://github.com/PTB-MR/mrpro/blob/3c6d873cd72afe1feda0abba93d677036d783b9f/tests/operators/test_non_uniform_fast_fourier_op.py#L51

for recon_matrix == encoding_matrix and fully-sampled k-space.

The oversampling-parameter of the NUFFT I am not sure about. All the packages I know do not ensure that the norm is independent of this parameter but this does not mean it is the correct thing to do.

fzimmermann89 commented 6 days ago

It should also match for all of the common MR trajectories. Otherwise, a tiny trajectory change would result in a scaled image -- which would be quite surprising. Ill try if I can make these tests pass..

koflera commented 4 days ago

I thought to norm of FourierOp should always be 1?

I don't think this holds for arbitrary FourierOps. If a FourierOp only yields an undersampled k-space then the Norm will not be 1. Also if padding or cropping of image-space/k-space is required (mismatch of recon_matrix and encoding_matrix) the norm will not be 1.

Here I verified that FFT and NUFFT lead to the same result:

https://github.com/PTB-MR/mrpro/blob/3c6d873cd72afe1feda0abba93d677036d783b9f/tests/operators/test_non_uniform_fast_fourier_op.py#L51

for recon_matrix == encoding_matrix and fully-sampled k-space.

The oversampling-parameter of the NUFFT I am not sure about. All the packages I know do not ensure that the norm is independent of this parameter but this does not mean it is the correct thing to do.

For the Fast Fourier operator, the operator norm should always be be one (if norm="ortho" is used, which we always have, right?), regardless of full or undersampling. For the general FourierOp, the norm will indeed not always be one.

fzimmermann89 commented 4 days ago

My current opinion:

As we calculate the operator norm using the gram, it should be independent of the norm setting in the fft as long as it's the same in fft and ifft (we only use fft.gram for our operator norm)

Nufft and fft should match for on grid trajectories.

Nufft scaling should be independent of grid size.

An undersampled fft should always have norm <=1: The eigenvektor with energy only at sampled k space points has the highest ew. If and only if this is also ev to ew 1 of the fully sampled fft, the norm will be one. (Imdersamplimg only removes eigenvalues)

A zero padded fft should have norm 1 It's eogencektor will be the zero padded eogencektor of the normal fft.

I think depending on how to normalize, we can either make the zero padded(oversampled) or undersampled fft have norm 1.

What can't be norm 1 are transforms with row sum>1 - repeated samples in the Nufft. Here, dcf@fourier should have norm 1