Reflectometry rewrite #61

nvaytet · 2021-11-08T20:37:02Z

This is a major re-write of the reflectometry code. This aims to be a first step in the rewrite, where some functionality has been removed/disabled (see illumination correction fx) for simplicity, and will be added back in as part of follow-up PRs.

The idea is to write the reflectometry module in a more functional way, rather than an object-oriented approach.
This way, we are passing Scipp data arrays and variables around, instead of custom classes, which makes functions more re-useable/portable, as well as more in-line with the rest of the ess codebase.

In addition, the gravity correction calculation has been replaced by a custom conversion graph that is then supplied to transform_coords to perform the calculation of the two_theta scattering angle and the Q vector.

The illumination correction has been removed from the reduction workflow (for now), but not from the codebase.

Many tests have been deleted, trying to keep relevant tests for the functions (the transform_coords is being tested in scipp and scippneutron, so we don't need to test again all aspects of it here). Some of these tests will come back in future PRs, when we add the full functionality back in.

A basic reduction notebook, based on the ones from the ess-notebooks repo has been added, and this should in the future become the reference reduction workflow, once we add all the missing bits.

…f detector_id x nevets

…w_reflectometry_data

arm61 · 2021-11-16T12:46:28Z

src/ess/reflectometry/conversions.py

+    Generate a coordinate transformation graph for reflectometry.
+    """
+    graph = {**conversions.beamline(scatter=True), **conversions.elastic("tof")}
+    graph["two_theta"] = reflectometry_theta


Can we also produce something called 'theta' alongside this?

yes, it seems we can add something that is called theta which is just a renaming of two_theta by doing

graph["theta"] = "two_theta"

Are you sure we need/want a double meaning of theta? This seems like it will nearly with certainty lead to bugs.

I want a dimension called theta in my final dataset, which is equal to a half of two_theta

Then @nvaytet must have misunderstood you?

I think so. For now, this is fine, I will open an issue for this to be included.

I think I am very wary of using theta because I have things going wrong with that in the past. If think reflectometry may be an exception, since a mixup with \theta_{spherical} is not likely (unlike other techniques, which define z along incident_beam.

I have added suggestions based on the nomenclature used in the Amor paper (fig 5 https://doi.org/10.1016/j.nima.2016.03.007).

SimonHeybrock · 2021-12-09T08:31:55Z

src/ess/reflectometry/conversions.py

+    # Note: make use of in-place operations for better performance
+    drop *= grav * (m_n**2 / (2 * h**2))
+    drop += y
+    return sc.asin(sc.abs(drop) / L2) - sc.to_unit(sample_rotation, 'rad')


Suggested change

return sc.asin(sc.abs(drop) / L2) - sc.to_unit(sample_rotation, 'rad')

out = sc.abs(drop, out=drop)

out /= L2

out = sc.asin(out, out=out)

out -=sc.to_unit(sample_rotation, 'rad')

return out

does no extra allocation, vs. FOUR large allocations, I think?

SimonHeybrock · 2021-12-09T08:32:44Z

src/ess/reflectometry/conversions.py



 def reflectometry_q(wavelength: sc.Variable, theta: sc.Variable) -> sc.Variable:
    """
    Compute the Q vector from the theta angle computed as the difference
    between gamma and omega.
+    Note that this is different from the 'normal' Q defined in scippneutron where


Suggested change

Note that this is different from the 'normal' Q defined in scippneutron where

Note that this is identical to the 'normal' Q defined in scippneutron where

...except that the input is given as theta instead of two_theta?

SimonHeybrock · 2021-12-09T08:35:13Z

src/ess/v20/imaging/beamline.py

-            kind=ChopperKind.FRAME_OVERLAP)
-    }
+    beamline["chopper_wfm_1"] = sc.scalar(
+        sc.Dataset(


See above, a chopper "constructor" (even if it returns a dataset) would probably be useful?

…e dense or the event coordinates

SimonHeybrock · 2021-12-10T10:58:50Z

src/ess/choppers/make_chopper.py

+        if not np.all(np.diff(utils.cutout_angles_begin(chopper).values) > 0):
+            raise ValueError("Chopper begin cutout angles are not monotonic.")
+        if not np.all(np.diff(utils.cutout_angles_end(chopper).values) > 0):
+            raise ValueError("Chopper end cutout angles are not monotonic.")


Could this use https://scipp.github.io/generated/functions/scipp.issorted.html?

SimonHeybrock

Nice work!

SimonHeybrock · 2021-12-10T15:15:38Z

src/ess/choppers/make_chopper.py

@@ -52,9 +52,9 @@ def make_chopper(frequency: sc.Variable,
        widths = utils.cutout_angles_width(chopper)
        if (sc.min(widths) < sc.scalar(0.0, unit=widths.unit)).value:
            raise ValueError("Negative window width found in chopper cutout angles.")
-        if not np.all(np.diff(utils.cutout_angles_begin(chopper).values) > 0):
+        if not sc.issorted(utils.cutout_angles_begin(chopper), dim=widths.dim).value:


https://scipp.github.io/release/0.8.4/generated/functions/scipp.allsorted.html is your friend, forget that we have it.

ok, and just to make sure: do we also have a function that doesn't require the dim?

I don't think we do, we have always avoided that (even though there are a couple of exceptions).

nvaytet added 17 commits October 1, 2021 14:23

first draft and changing data loader structure

6dee37d

load amor data and convert to wavelength

24ea600

Merge branch 'main' into new_reflectometry_data

dcecd1a

Merge branch 'main' into new_reflectometry_data

6dbb87e

Merge branch 'main' into new_reflectometry_data

9de8eaa

manually broadcast positions to events to avoid automatic broadcast o…

0056e76

…f detector_id x nevets

use binning in wavelength to simplify theta calculation

35a1eb9

add qz coord

b9d0133

Merge branch 'new_reflectometry_data' of github.com:scipp/ess into ne…

20239e2

…w_reflectometry_data

start adding conversion graphs

87dc4e7

add skeleton notebook for amor reduction

ee80db7

mroe cleanup

6eec445

remove old files, we'll add functionality back in from the diff

7925e26

docstrings

8e03086

remove utils, add file managing with pooch and update ntoebook

9a635b0

update env file

efdc75c

Merge branch 'main' into new_reflectometry_data

d31e7de

nvaytet marked this pull request as draft November 8, 2021 20:39

nvaytet added 8 commits November 8, 2021 21:42

flake8

425dc4e

fix python tests

7f01825

fix unreleased erf function error

05f8359

add pythreejs and add normalization by supermirror in notebook

08f99fe

add graphviz

96edd88

workaround astype until 0.9 is released

881e215

scippneutron version 0.4.1

1c78219

remove TODO after 0.9 release

5fddb8e

nvaytet marked this pull request as ready for review November 15, 2021 19:35

nvaytet requested review from arm61 and SimonHeybrock November 15, 2021 19:35

arm61 reviewed Nov 16, 2021

View reviewed changes

SimonHeybrock reviewed Dec 9, 2021

View reviewed changes

nvaytet added 14 commits December 9, 2021 12:50

replace old Chopper constructor with a free function make_chopper

e1acf3c

use make_chopper in bemalines

d8c44a6

better use of in-place ops in coord transforms

5e4615e

update comment in docstring

d4ed7c2

add make_chopper docstring

4015110

use unit of L2

1f0ecc3

add more docstrings in chopper.utils

b6d2816

flake8

aaba230

update docs notebook

a3fb10c

update wfm notebook

e4da96f

try to cleanup api docs

e1df1ab

add toctrees in reference docs

6092471

branching optimization in tramsform coords when considering either th…

36c0e7b

…e dense or the event coordinates

update release notes

e7889a5

nvaytet mentioned this pull request Dec 9, 2021

Docs book theme #74

Merged

SimonHeybrock reviewed Dec 10, 2021

View reviewed changes

SimonHeybrock approved these changes Dec 10, 2021

View reviewed changes

nvaytet added 2 commits December 10, 2021 16:09

use issorted

4310120

flake8

99438ce

SimonHeybrock reviewed Dec 10, 2021

View reviewed changes

issorted().value -> allsorted()

785dbc6

nvaytet merged commit 0fa97b4 into main Dec 10, 2021

nvaytet deleted the new_reflectometry_data branch December 10, 2021 15:31

-    return sc.asin(sc.abs(drop) / L2) - sc.to_unit(sample_rotation, 'rad')
+    out = sc.abs(drop, out=drop)
+    out /= L2
+    out = sc.asin(out, out=out)
+    out -=sc.to_unit(sample_rotation, 'rad')
+    return out

	Note that this is different from the 'normal' Q defined in scippneutron where
	Note that this is identical to the 'normal' Q defined in scippneutron where

Reflectometry rewrite #61

Reflectometry rewrite #61

Uh oh!

Conversation

nvaytet commented Nov 8, 2021 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

arm61 Nov 16, 2021 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

SimonHeybrock left a comment

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Uh oh!

nvaytet commented Nov 8, 2021 •

edited

Loading

arm61 Nov 16, 2021 •

edited

Loading