Add anisotropic regularization for FILD measurements showcaseAdd notebook with external file references and LICENSE

[bosschmidt]

Describe your contribution

CIL user showcase 13

Checklist when you are ready to request a review

• I have performed a self-review of my code
• I have created a new folder, containing my contributions which are in the form of jupyter notebooks, with any necessary supporting python files, and a LICENSE file.
• I have added a description of my contribution(s) to the top of my file(s)
• If publicly available, I have added a link to the dataset used near the top of my file(s)
• I have added the CIL version I ran with near the top of my file(s)
• The content of this Pull Request (the Contribution) is intentionally submitted for inclusion in CIL-User-Showcase.
• I confirm that the contribution does not violate any intellectual property rights of third parties.
• I confirm that I have added license headers to all of the files I am contributing (with a license of my choice)
• Change pull request label to 'Waiting for review'

Note: for an example of a contribution, where a license header, description, data link and CIL version has been added, please
see: example_contribution

[MargaretDuff]

Thanks @bosschmidt - looking forward to seeing you next week! We will take a look at this, probably after the user meeting

[MargaretDuff]

Hi @bosschmidt! Thanks again for your submission and good to see you at the user meeting last week.

I have just taken a quick look. The notebook lists a set of custom files needed to run the notebook:

• my_custom_operator.py: Custom FILD operator implementation
• manolo_colormap.py: Custom colormap for visualization
• WaveletOperator.py: CIL's wavelet transform implementation
• b_2D.npy: Experimental FILD data
• 34614_Bo_v2023_4_weight_function.nc: Simulated system matrix

It is nice if these notebooks can be run by the user but understand if data cannot be shared. It especially would be good to share the "my_custom_operator.py" file, so people can see the custom CIL operator.

[MargaretDuff]

Hi @bosschmidt! Thanks again for your submission and good to see you at the user meeting last week.

I have just taken a quick look. The notebook lists a set of custom files needed to run the notebook:

• my_custom_operator.py: Custom FILD operator implementation
• manolo_colormap.py: Custom colormap for visualization
• WaveletOperator.py: CIL's wavelet transform implementation
• b_2D.npy: Experimental FILD data
• 34614_Bo_v2023_4_weight_function.nc: Simulated system matrix

It is nice if these notebooks can be run by the user but understand if data cannot be shared. It especially would be good to share the "my_custom_operator.py" file, so people can see the custom CIL operator.

Discussed this with Bo offline and there is now a link in the notebook to a github where these are stored.

[jakobsj]

Hi @bosschmidt,
Thank you so much for this excellent and comprehensive notebook. It provides a super nice example with detailed guidance on setting up a custom non-CT inverse problem, solving it a number of combinations of regularizers (thereby illustrating the optimization prototyping capabilities of CIL) and providing instructive commentary and conclusions. Also thank you @MargaretDuff for nicely handling and editing.

I had a few mostly minor I think editorial comments:

1. "Physical meaning: Each row represents the detector response to a specific velocity-space point" maybe I misunderstood but shouldn't it say each column instead of row being a detector response (since column j of the matrix is what you obtain if you multiply the matrix onto the unit vector with a 1 at index j and all other elements 0).
2. "Interpolate weight functions to higher different dimension on domain geometry to generate the forward operator:" Could you elaborate a bit on the motivation for this interpolation (inverse crime related?) and the sizes of grid before and after interpolation and which grid is used for what (simulation and reconstruction respectively)?
3. "First-Order Tikhonov Regularization" Since a nonnegativity constraint is also enforced, and Tikhonov is the name of the unconstrained formulation (not just two-norm regularization), it would be more correct to refer to this problem as nonnegativity constrained Tikhonov.
4. "Includes 10% measurement noise and 1% background noise" Could you elaborate what each of measurement noise and background noise are, perhaps include a mathematical formula (in addition to the code) for how the simulated observed data is computed?
5. In the simulated data besov section "interpolator = RegularGridInterpolator((y, x), result_array_2D)" could you elaborate why an interpolation is used on the reconstruced image?
6. Then in the real data results, no RegularGridInterpolater call, but instead imshow with bicubic interpolation - why the difference?

[bosschmidt]

Thanks for the detailed review and thoughtful comments @jakobsj ! Let me address each point:

1. The matrix-vector product interpretation is correct, but we describe the model in terms of weight functions (rows) as this better represents the physical measurement process: each row shows how any velocity-space distribution contributes (or "is weighted") to the signal at a specific detector pixel--Ref. W W Heidbrink et al 2007 Plasma Phys. Control. Fusion 49 1457 explains this concept in detail.
2. This is indeed related to avoiding inverse crime. We use different grid dimensions for simulation (92x52 = 4784 points) and reconstruction (91x51 = 4641 points).
3. You're right. The terminology should be "non-negativity constrained first-order Tikhonov regularization" throughout the notebook.
4. The noise model follows equation (4) in Schmidt et al 2024 Nucl. Fusion 64 076009:

• Measurement noise (10%) represents detector-specific noise in the FILD measurements.
• Background noise (1%) accounts for other factors like neutron radiation and gamma rays.

5. You're right - the interpolation was only used for comparison with ground truth and isn't essential. These three lines can be removed.
6. The difference in interpolation was due to different needs: the real data section doesn't require comparison with ground truth, so no interpolation to a common grid is needed. The imshow with bicubic interpolation is purely for visualization purposes.

Would you prefer I implement these clarifications directly in the notebook?

[jakobsj]

Hi @bosschmidt, yes please, that would be great! For bullet 1., I think what you say makes sense. My concern was/is about the phrase "detector response" which in my understanding (but I may be wrong or different usage in different areas) is the resulting measurement vector if applying the forward a canonical basis vector. Could a solution here be to avoid the phrase "detector response" and explain it in words and with a reference to the source you cited?

[bosschmidt]

@jakobsj Yes, sounds good. I will attend to this by the middle of next week.

[MargaretDuff]

Hi @jakobsj - ready for one (last) review