diff --git a/docs/api-reference/index.md b/docs/api-reference/index.md index 3b023e75..9f3bd24b 100644 --- a/docs/api-reference/index.md +++ b/docs/api-reference/index.md @@ -46,7 +46,6 @@ conversions data - instrument_view load orso resolution diff --git a/docs/user-guide/data-reduction.md b/docs/user-guide/data-reduction.md new file mode 100644 index 00000000..c3e3694a --- /dev/null +++ b/docs/user-guide/data-reduction.md @@ -0,0 +1,169 @@ +# Focused reflectometry data reduction + +The goal of the reflectometry data reduction is to compute the sample reflectivity $R(Q)$ as a function of the momentum transfer $Q$. + +Based on [J. Stahn, A. Glavic, Focusing neutron reflectometry: Implementation and experience on the TOF-reflectometer Amor](#reference). + + +## Preliminaries + +The detector data consists of a list $EV$ of neutron detector events. +For each neutron event in the list we know its wavelength $\lambda$ and the pixel number $j$ of the detector pixel that it hit. +The detector pixel positions are known and so is the position and the orientation of the sample. +From this information we can compute the reflection angle $\theta$, and the momentum transfer $Q$ caused by the interaction with the sample. + +The purpose of this text is not to describe how the event coordinates $Q$ and $\theta$ are derived from the raw event data and the geometry information, so for now just take those relations for given. +For more details see the implementations for the respective instruments [Amor] and [Estia]. + +To simplify the description it is assumed that the sample- and reference measurements were made over the same length of time, and it is assumed the neutron intensity from the source did not vary between the two measurements. + + +## Model of event intensity in the detector + +The reflectivity of the sample is related to the intensity of neutron counts in the detector by the model +$$ +I_{\text{sam}}(\lambda, j) = F(\theta(\lambda, j, \mu_{\text{sam}}), w_{\text{sam}}) \cdot R(Q(\lambda, \theta(\lambda, j, \mu_{\text{sam}}))) \cdot I_{\text{ideal}}(\lambda, j) +$$ (model) +where $I_{\text{sam}}(\lambda, j)$ represents the expected number of neutrons detected in the $j$ pixel of the detector per unit of wavelength at the wavelength value $\lambda$. $I_{\text{ideal}}$ represents the expected number of neutrons detected if the sample was a perfect reflector and large enough so that the footprint of the focused beam on the sample was small compared to the sample. $F(\theta, w)$ is the fraction of the beam that hits the sample. It depends on the incidence angle $\theta$ and on the size of the sample represented by $w,$ and $\mu_{\text{sam}}$ is the sample rotation. + +The model does not hold for any $\lambda,j$. To make this explicit, let $M_{sam}$ represent the region of $\lambda,j$ where the model is expected to hold. + +The ideal intensity is estimated from a reference measurement on a neutron supermirror. +How it is computed will be described later, for now assume it is known. + +## Estimating $R(Q)$ +Move $F$ to the left-hand-side of equation {eq}`model` and integrate over all $\lambda, j\in M$ contributing to the $Q$-bin $[q_{i}, q_{i+1}]$ +$$ +\int_{M \cap Q(\lambda, \theta(\lambda, j, \mu_{\text{sam}})) \in [q_{i}, q_{i+1}]} \frac{I_{\text{sam}}(\lambda, j)}{F(\theta(\lambda, j, \mu_{\text{sam}}), w_{\text{sam}})} d\lambda \ dj = \\ +\int_{M \cap Q(\lambda, \theta(\lambda, j, \mu_{\text{sam}})) \in [q_{i}, q_{i+1}]} I_{\text{ideal}}(\lambda, j) R(Q(\lambda, \theta(\lambda, j, \mu_{\text{sam}}))) d\lambda \ dj. +$$ (reflectivity) +Notice that if the $Q$ binning is sufficiently fine then $R(Q)$ is approximately constant in the integration region. +Assuming the binning is fine enough $R(Q)$ can be moved outside the integral and isolated so that + +$$ + R(Q_{i+\frac{1}{2}}) \approx \frac{\int_{M \cap Q(\lambda, j, \mu_{\text{sam}}) \in [q_{i}, q_{i+1}]} \frac{I_{\text{sam}}(\lambda, j)}{F(\theta(\lambda, j, \mu_{\text{sam}}), w_{\text{sam}})} d\lambda \ dj }{\int_{M \cap Q(\lambda, \theta(\lambda, j, \mu_{\text{sam}})) \in [q_{i}, q_{i+1}]} I_{\text{ideal}}(\lambda, j) d\lambda \ dj} =: \frac{I_{measured}(Q_{i+\frac{1}{2}})}{I_{\text{ideal}}(Q_{i+\frac{1}{2}})} +$$ +for $Q_{i+\frac{1}{2}} \in [q_{i}, q_{i+1}]$. + +For the integral to make sense $M\sub M_{sam}$, but there might be other constraints limiting $M$ more, for now we leave it undefined. + + +## The reference intensity $I_{\text{ideal}}$ +$I_{\text{ideal}}$ is estimated from a reference measurement on a neutron supermirror with known reflectivity curve. +The reference measurement intensity $I_{\text{ref}}$ is modeled the same way the sample measurement was + +$$ +I_{\text{ref}}(\lambda, j) = F(\theta(\lambda, j, \mu_{\text{ref}}), w_{\text{ref}}) \cdot R_{\text{supermirror}}(Q(\lambda, \theta(\lambda, j, \mu_{\text{ref}}))) \cdot I_{\text{ideal}}(\lambda, j) +$$ +but in this case $R_{\text{supermirror}}(Q)$ is known. + +As before, the model does not hold for any $\lambda,j$. Let $M_{ref}$ represent the region of $\lambda,j$ where the model is expected to hold. $M_{ref}$ is typically not the same as $M_{sam}$. For example, if the supermirror reflectivity is not known for all $Q$ we are not going to have a model for the reference intensity at the $\lambda,j$ corresponding to those $Q$ and that is reflected in $M_{ref}$ but not in $M_{sam}$. + +Using the definition in {eq}`reflectivity` +$$ +I_{\text{ideal}}(Q_{i+\frac{1}{2}}) = \int_{M \cap Q(\lambda, \theta(\lambda, j, \mu_{\text{sam}})) \in [q_{i}, q_{i+1}]} \frac{I_{\text{ref}}(\lambda, j)}{F(\theta(\lambda, j, \mu_{\text{ref}}), w_{\text{ref}}) R_{\text{supermirror}}(Q(\lambda, \theta(\lambda, j, \mu_{\text{ref}})))} + d\lambda \ dj. +$$ +For this integral to make sense $M\sub M_{ref}$, so now we have an additional constraint on $M$ to keep in mind. + + +## Estimating intensities from detector counts + +The number of neutron counts in the detector is a Poisson process where the expected number of neutrons per pixel and unit of wavelength are the measurement intensities $I_{sam}$ and $I_{ref}$ defined above. +The expected number of counts can be estimated by the empirically observed count: +$$ +I_{measured}(Q_{i+\frac{1}{2}}) = \int_{M\cap Q(\lambda, \theta(\lambda, j, \mu_{\text{sam}})) \in [q_{i}, q_{i+1}]} \frac{I_{\text{sam}}(\lambda, j)}{F(\theta(\lambda, j, \mu_{\text{sam}}), w_{\text{sam}})} d\lambda \ dj = \\ + \approx +\sum_{\substack{k \in EV_{\text{sam}} \\ Q(\lambda_{k}, \theta(\lambda_{k}, j_{k}, \mu_{\text{sam}})) \in [q_{i}, q_{i+1}]\\ (\lambda_k, j_k)\in M}} \frac{1}{F(\theta(\lambda_{k}, j_{k}, \mu_{\text{sam}}), w_{\text{sam}})} +$$ +where $EV_{\text{sam}}$ refers to the event list from the sample experiment. + +The sum is compound Poisson distributed and for such random variables the variance is well approximated by the sum of squared summands +$$ +V\bigg[ \sum_{\substack{k \in EV_{\text{sam}} \\ Q(\lambda_{k}, \theta(\lambda_{k}, j_{k}, \mu_{\text{sam}})) \in [q_{i}, q_{i+1}]\\ (\lambda_k, j_k)\in M}} \frac{1}{F(\theta(\lambda_{k}, j_{k}, \mu_{\text{sam}}), w_{\text{sam}})} \bigg] \approx +\sum_{\substack{k \in EV_{\text{sam}} \\ Q(\lambda_{k}, \theta(\lambda_{k}, j_{k}, \mu_{\text{sam}})) \in [q_{i}, q_{i+1}]\\ (\lambda_k, j_k)\in M}} \frac{1}{F(\theta(\lambda_{k}, j_{k}, \mu_{\text{sam}}), w_{\text{sam}})^2}. +$$ + +The same estimates are used to approximate the ideal intensity: +$$ +I_{\text{ideal}}(Q_{i+\frac{1}{2}}) \approx \sum_{\substack{k \in EV_{\text{ref}} \\ Q(\lambda_{k}, \theta(\lambda_{k}, j_{k}, \mu_{\text{sam}})) \in [q_{i}, q_{i+1}]\\ (\lambda_k, j_k)\in M}} \frac{1}{F(\theta(\lambda_{k}, j_{k}, \mu_{\text{ref}}), w_{\text{ref}}) R_{\text{supermirror}}(Q(\lambda_{k}, \theta(\lambda_{k}, j_{k}, \mu_{\text{ref}})))} +$$ +The above two expressions together with {eq}`reflectivity` lets us express $R(Q_{i+{\frac{1}{2}}})$ entirely in terms of the neutron counts in the detector. + +## More efficient evaluation of the reference intensity +The above expression for the reference intensity is cumbersome to compute because it is a sum over the reference measurement event list, and the reference measurement is large compared to the sample measurement. + +Therefore we back up a bit. Consider the expression for the reference intensity, replacing the integrand with a generic $I(\lambda, j)$ it looks something like: + +$$ +I_{\text{ideal}}(Q_{i+\frac{1}{2}}) = \int_{Q(\lambda, \theta(\lambda, j, \mu_{\text{sam}})) \in [q_{i}, q_{i+1}]} I(\lambda, j) \ d\lambda \ dj. +$$ + +In the previous section we approximated the integral by summing over all events in the reference measurement. + +Alternatively, we could define a $\lambda$ grid with edges $\lambda_{k}$ for $k=1\ldots N$ and approximate the integration region as the union of a subset of the grid cells: + +$$ +\int_{Q(\lambda, \theta(\lambda, j, \mu_{\text{sam}})) \in [q_{i}, q_{i+1}]} I(\lambda, j) \ d\lambda \ dj +\approx \sum_{Q(\bar{\lambda}_{k+\frac{1}{2}},\ \theta(\bar{\lambda}_{k+\frac{1}{2}}, j, \mu_{\text{sam}})) \in [q_{i}, q_{i+1}]} I_{k+\frac{1}{2},j} +$$ +where + +$$ +I_{k+\frac{1}{2},j} = \int_{\lambda \in [\lambda_{k}, \lambda_{k+1}]} I(\lambda, j) \ d\lambda +$$ +and $\bar{\lambda}_{k+\frac{1}{2}} = (\lambda_{k} + \lambda_{k+1}) / 2$. + +Why would this be more efficient than the original approach? Note that $I_{k+\frac{1}{2}, j}$ does not depend on $\mu_{\text{sam}}$, and that it can be computed once and reused for all sample measurements. +This allows us to save computing time for each new sample measurement, as long as $|EV_{\text{ref}}| >> NM$ where $M$ is the number of detector pixels and $N$ is the size of the $\lambda$ grid. + +However, ideally computing the reference intensity should be quick compared to reducing the sample measurement. And since a reasonable value for $N$ is approximately $500$, and $M\approx 30000$, and a sample measurement is likely less than $10$ million events, the cost of computing the reference measurement is still considerable compared to reducing the sample measurement. + +Therefore there's one more approximation that is used to further reduce the cost of computing the reference intensity. +The description of the final approximation is instrument specific. +In the next section it is described specifically for the Amor instrument. + + +### Evaluating the reference intensity for the Amor instrument + +The Amor detector has three logical dimensions, `blade`, `wire` and `stripe`. It happens to be the case that $\theta(\lambda, j)$ is almost the same for all $j$ belonging to the same `stripe` of the detector. +We can express this as + +$$ +\theta(\lambda, j, \mu_{\text{sam}}) \approx \bar{\theta}(\lambda, \mathrm{bladewire}(j), \mu_{\text{sam}}) +$$ +where $\bar{\theta}$ is an approximation for $\theta$ that only depends on the blade and the wire of the pixel where the neutron was detected. +Then the above expression for the reference intensity can be rewritten as + +$$ +\int_{Q(\lambda, \bar{\theta}(\lambda, z, \mu_{\text{sam}})) \in [q_{i}, q_{i+1}]} \int_{\mathrm{bladewire}(j) = z} I(\lambda, j) \ dj \ d\lambda \ dz +\approx \sum_{Q(\bar{\lambda}_{k+\frac{1}{2}}, \bar{\theta}(\bar{\lambda}_{k+\frac{1}{2}}, z, \mu_{\text{sam}})) \in [q_{i}, q_{i+1}]} I_{k+\frac{1}{2},z} +$$ +where + +$$ +I_{k+\frac{1}{2},z} = \int_{\lambda \in [\lambda_{k}, \lambda_{k+1}]} \int_{\mathrm{bladewire}(j) = z} I(\lambda, j) \ dj \ d\lambda . +$$ +Like before, the benefit of doing this is that + +$$ + \int_{\lambda \in [\lambda_{k}, \lambda_{k+1}]} \int_{\mathrm{bladewire}(j) = z} I(\lambda, j) \ dj \ d\lambda +$$ +can be pre-computed because it doesn't depend on $\mu_{\text{sam}}$. +But unlike before $I_{k+\frac{1}{2},z}$ now has a much more manageable size, about 64x smaller than the first attempt. +This makes it comfortably smaller than the sample measurement. + + +## Reference + +``` +J. Stahn, A. Glavic, +Focusing neutron reflectometry: Implementation and experience on the TOF-reflectometer Amor, +Nuclear Instruments and Methods in Physics Research Section A: Accelerators, Spectrometers, Detectors and Associated Equipment, +Volume 821, +2016, +Pages 44-54, +ISSN 0168-9002, +https://doi.org/10.1016/j.nima.2016.03.007. +(https://www.sciencedirect.com/science/article/pii/S0168900216300250) +``` \ No newline at end of file diff --git a/docs/user-guide/index.md b/docs/user-guide/index.md index a86b0b21..ee652d08 100644 --- a/docs/user-guide/index.md +++ b/docs/user-guide/index.md @@ -5,5 +5,6 @@ maxdepth: 1 --- amor/index +data-reduction ```