# Real World Time Series Benchmark Datasets with Distribution Shifts:
Global Crude Oil Price and Volatility

Pranay Pasula111Accepted to IJCAI 2023, the 32nd International Joint
Conference on Artificial Intelligence (AI4TS).222Selected as a Contributed
Talk.333Awarded Best Paper Runner-Up / Honorable Mention. UC Berkeley EECS
<EMAIL_ADDRESS>
###### Abstract

The scarcity of task-labeled time-series benchmarks in the financial domain
hinders progress in continual learning. Addressing this deficit would foster
innovation in this area. Therefore, we present COB, Crude Oil Benchmark
datasets. COB includes 30 years of asset prices that exhibit significant
distribution shifts and optimally generates corresponding task (i.e., regime)
labels based on these distribution shifts for the three most important crude
oils in the world. Our contributions include creating real-world benchmark
datasets by transforming asset price data into volatility proxies, fitting
models using expectation-maximization (EM), generating contextual task labels
that align with real-world events, and providing these labels as well as the
general algorithm to the public. We show that the inclusion of these task
labels universally improves performance on four continual learning algorithms,
some state-of-the-art, over multiple forecasting horizons. We hope these
benchmarks accelerate research in handling distribution shifts in real-world
data, especially due to the global importance of the assets considered. We’ve
made the (1) raw price data, (2) task labels generated by our approach, (3)
and code for our algorithm available at oilpricebenchmarks.github.io.

##  1 Introduction

Benchmarks serve as testbeds for artificial intelligence (AI) algorithms and
have promoted the rapid acceleration of algorithm performance that has been
seen over the past ten years. Notable examples from the beginning of this ten
year period include ImageNet (Deng et al., 2009), COCO (Lin et al., 2014),
CIFAR-10 (Krizhevsky et al., 2009), and the Arcade Learning Environment
(Bellemare et al., 2013), the first of which is often credited with initiating
the ongoing deep learning revolution.

These benchmarks all contain images in the form of raw pixel data, which is an
important medium, but more recently, benchmarks have expanded to include
different modalities of data that are aimed at evaluating progress in a
variety of AI subfields and problem domains.

In this work, we focus on _temporal distribution shifts_ within the financial
domain and introduce time-series benchmark datasets of the three primary
global crude oil price markers, _West Texas Intermediate (WTI)_ , _Brent
Blend_ , and _Dubai Crude_.

For each dataset we provide over 30 years of data—daily spot prices for WTI
and Brent and monthly average prices for Dubai Crude—and see striking
distribution shifts throughout this time period for each asset.

Problem settings with time-series or otherwise sequential data often pose
serious issues to predictive AI models in several ways that the standard
supervised learning setting using (independent and identically distributed)
IID non-sequential data does not. Time-series often embodies obstacles, such
as

  * •

Correlatedness. While training predictive models, series data points that
occur near the same time tend to take similar values and have similar targets.
This biases the models to this local region of the input and target variable
spaces.

  * •

Non-stationarity. Detailed in Section 1.1.

  * •

Missing data. Subsets of data are often missing from series data. For example,
inadequate record keeping may be why the WTI and Brent datasets we introduce
are missing 3 and 2 percent of daily spot price data, respectively, but we
resample the datasets in a way that results in new series without missing
data.

  * •

Spurious data. Outliers, erroneous, or infeasible data skew model learning in
poor directions.

###  1.1 Non-Stationarity

Non-stationarity, or distribution shifts, threaten predictive models that have
been trained on data and then deployed to make predictions online, or in real-
time on different data. As the data distribution shifts further from the data
distribution the model was trained on, the worse the model performs.
Furthermore, this non-stationarity can manifest in fundamentally different
ways.

Let $x$ be the independent variables, $y$ be the target variables, and $P$ be
an appropriate probability measure, then we have the following types of non-
stationarity:

  * •

Covariate shift: A shift in the distribution of independent variables $P(X)$.

  * •

Prior probability shift: A shift in the distribution of target variables
$P(Y)$.

  * •

Concept drift: A shift in the distribution of the relationship between the
independent and target variables $P(Y\mid X)$.

###  1.2 Out-of-Distribution Detection

It’s crucial for online algorithms to quickly identify distribution shifts
because predictive models perform worse as the data they are predicting on
deviates further from the data that they were trained on. _But how does an
algorithm decide that a distribution shift has occurred?_ There’s a trade-off
between speed and confidence in deciding whether a recent sequence of data
includes a distribution shift.

###  1.3 Continual (Lifelong) Learning

A rapidly growing field of artificial intelligence that often handles data
online is _continual learning_ , which is sometimes used interchangeably with
_lifelong learning_ or _incremental learning_. Continual learning algorithms
aim to develop models that accrue useful knowledge over time to accomplish
some set of tasks, or equivalently, perform well on data distributions
corresponding to some set of underlying tasks while avoiding negative
interference or catastrophic forgetting.

We refer to the typical continual learning problem setting as one in which the
tasks are known, regime shifts happen instantaneously during a discrete
timestep, and both the regime-generated data and the regime labels are given
to the model at any time $t$.

###  1.4 Our Contribution

Our main contribution is the creation of three new real-world benchmark
datasets on prices of WTI, Brent, and Dubai crude oil, which are of critical
importance worldwide (Backus and Crucini, 2000; Kilian and Park, 2009).

We have

  1. 1.

Transformed the data into a proxy for volatility,

  2. 2.

Fitted the data into tasks using expectation-maximization (EM), a well-
established algorithm, to optimize for measures of information criteria,

  3. 3.

Provided contextual labels of major real-world events that align with price
and regime shifts,

  4. 4.

Demonstrated that our algorithmically task-labeled dataset generated tasks
that align with recessions and other spurious real-world events.

We hope that this work will accelerate the development of continual learning,
out-of-distribution detection, and other algorithms in handling issues posed
by real-world data that contains distribution shifts. Furthermore, we believe
that the societal importance of the assets represented by these datasets makes
these benchmarks especially appealing.

##  2 Datasets

A description of WTI, Brent, and Dubai Crude oils can be found in Section 4.

The raw data for WTI and Brent contains prices at a daily frequency while that
for Dubai Crude contains average prices at a monthly frequency. Therefore the
resampling process described in Section 3.3.1 applies only to WTI and Brent.

##  3 Dataset Creation Process

This section elucidates the steps that we have taken to construct the real-
world benchmark datasets that we are providing.

###  3.1 Converting Daily Spot Price to Volatility

#### WTI and Brent

We use the daily spot prices of WTI (U.S. EIA, 2022) and Brent (U.S. EIA,
2022) to derive weekly percent changes, a proxy of volatility, which is
important in evaluating _value at risk_ (VaR), a statistic that quantifies the
extent of potential financial losses within a position over a certain period
of time.

#### Dubai Crude

We use the monthly average prices of Dubai Crude (IMF, 2022) to derive monthly
percent changes, another proxy of volatility, but one that is attenuated by
the averaging over the monthly prices, as seen by comparing Figure 3 against
Figure 1 and Figure 2.

#### Heterogenous Dataset Difficulty

The difference in pricing frequency between the first two, WTI and Brent, and
the last, Dubai, allow us to create benchmark datasets with different
difficulties. Since Dubai Crude prices are monthly averages, much of the
pricing spread is smoothed out relative to if daily spot prices were used.
Unfortunately, the most reliable and granular Dubai Crude prices that had
history comparable in length to the WTI and Brent prices we use are the
monthly average prices we source from (IMF, 2022).

However, we see this as an opportunity to present benchmark datasets with
varying degrees of difficulty on various tasks, most obviously to make
predictions on, as target variables with lower spread are easier for
predictive models to handle.

###  3.2 Selecting Model Class Used to Create Tasks

#### Markov Regime Switching Models

Markov switching models are used on sequential data that is expected to
transition through a finite set of latent states, or tasks. These states can
and often differ significantly, as evidenced in Figure 1, Figure 2, and Figure
3.

(Kim et al., 1998) describes a Markov switching model that is particularly
well suited to handling the heteroskedasticity and mean reversion properties
of the raw data. The algorithm used in (Kim et al., 1998) has additional nice
properties, such as the use of expectation-maximization (EM), a form of
(local) maximum a posteriori (MAP) estimation.

Therefore, we use (Kim et al., 1998) to construct the benchmark datasets and
describe these steps in Section 3.3.

###  3.3 Creating Tasks Through Fitting model class

####  3.3.1 Appropriate Resampling

The raw WTI and Brent spot price data was captured daily, and the day-to-day
differences oscillated between positive and negative values so frequently that
the method we used to generate tasks (Kim et al., 1998) oscillated between
tasks so frequently that the results were unrealistic. Resampling the data by
every two days before using (Kim et al., 1998) didn’t resolve the issue, as
evidenced by Figure 5 in Appendix A.

We found that resampling weekly was an excellent frequency for (Kim et al.,
1998) to generate tasks that were reasonable and well optimized measures of
likelihood and information criteria.

####  3.3.2 Evaluating Stationarity

Since we’re putting forth the WTI, Brent, and Dubai crude oil datasets as
benchmark datasets, we want them to be amenable to several types of
evaluations. The most common evaluation on these series data are statistical
inference, specifically forecasting future asset prices. (Nagabandi et al.,
2018; He et al., 2019; Caccia et al., 2020; He and Sick, 2021; Pasula et al.,
2023) have shown that the provision of task, or context, information to
predictive models boosts their ability to accurately predict future values on
a variety of problem domains with challenging properties.

Therefore, it’s important to ensure that the time-series benchmark datasets
we’re providing are useful for evaluating tasks such as forecasting or
statistical inference. A critical step towards this is verifying that this
series data are not characterized by unit root processes, which can lead to
issues in statistical inference on this data (Phillips and Perron, 1988).

To evaluate whether unit roots are present in any of the benchmark datasets,
we perform an augmented Dickey-Fuller (ADF) test on each of the datasets. The
ADF tests on the datasets reject the null hypotheses that unit roots are
present in the data with all p-values $<$110-23$$.

####  3.3.3 Choosing the Number of Tasks

The problem of choosing the number of tasks $k$ is one of choosing the model
itself: once $k$ is specified, (Kim et al., 1998) uses EM to fit the data to a
partition of $k$ tasks. A rigorous way to guide model selection, and hence
$k$, is to use statistical measures of fit known as information criteria
(Stoica and Selen, 2004). The Akaike, Bayesian, and Hannan-Quinn information
criteria (AIC, BIC, HQIC) are three of the most commonly measures. Each
penalizes model deviation the same way but they vary the degree that they
penalize the number of model parameters and the loss of degrees of freedom.

We choose $k$ so that $\text{AIC}+\text{BIC}+\text{HQIC}$ is minimized under
the constraint that each task is present for at least one timestep for every
dataset. This results in $k=3$. We provide more information on our decision
process for this in Appendix B.

####  3.3.4 Smoothing probabilities

At any time $T$ with $t<T$, (Kim et al., 1998) uses an algorithm introduced in
(Kim, 1994) to smooth conditional probabilities $P(\mathcal{T}_{t}=i\mid
Y_{t})$ for every $t$ based on all observations available up through time $T$.

####  3.3.5 Assigning Tasks

For every timestep $t$, we assign task $\mathcal{T}_{t}$ to the argmax of the
smoothed probabilities at $t$ over all tasks $\\{1,2,\dotsc,k\\}$.

##  4 Datasets

###  4.1 West Texas Intermediate Crude Oil Daily Spot Price

A light (low density) and sweet (low sulfur content) crude oil that is ideal
for gasoline refining. This oil is extracted from land wells in the United
States and is transported to a location in the center of the US, Cushing,
Oklahoma. The cost required to transport WTI from its land-locked location to
overseas areas of the globe is a major drawback to this otherwise high-quality
crude oil and contributes to the positive Brent-WTI spread. WTI is the main
benchmark for oil that is used within the United States. We depict our
benchmark dataset creation algorithm for WTI in Figure 1.

Figure 1: Gray areas are periods that the NBER has classified as a recession.
(Top) WTI daily spot price. (Middle top) Excess returns, a proxy of volatility
and variance. (Middle bottom) Fitted model regime probabilities over time.
(Bottom) Task labels for each point in time, computed via argmax probabilites.

###  4.2 Brent Crude Oil Daily Spot Price

The most popular oil contracted, making up approximately two-thirds of all
global crude contracts, "Brent" is a blend of oil from four fields in the
North Sea and, like its WTI counterpart, is light and sweet. Because it can be
distilled into high-value products, such as gasoline and diesel fuel, and due
to it being sea-sourced and thus having relatively low transportation costs
globally, it has a far greater global reach and generally higher value than
WTI. We depict our benchmark dataset creation algorithm for Brent in Figure 2

Figure 2: (Top) Brent crude oil daily spot price per barrel. (Middle top)
Percent change over the weekly end price, a proxy of volatility and variance.
(Middle bottom) Fitted model regime probabilities over time. (Bottom) Regime
labels for each point in time, computed via argmax probabilites.

###  4.3 Dubai Crude Oil Monthly Average Price

Also known as Fateh, Dubai Crude is the least popular of the three primary
global crude oil price benchmarks, but it is used as such because it is one of
the few crude oils originating from the Persian Gulf that is available
immediately. Dubai Crude has medium density and higher sulfur content than WTI
and Brent but is priced comparably to both because of its central role in
pricing crude oil exports from the Persian Gulf to Asian-Pacific markets. We
depict our benchmark dataset creation algorithm for Dubai Crude in Figure 3

Figure 3: (Top) Dubai crude oil monthly average price per barrel. (Middle top)
Percent change over the _monthly average price_ , a proxy of volatility and
variance, but not as granular as the weekly end price that the WTI and Brent
datasets are using. (Middle bottom) Fitted model regime probabilities over
time. (Bottom) Regime labels for each point in time, computed via argmax
probabilites.

##  5 Experimental Design and Results

We use the regime, or task, labels derived from our algorithm as context
labels for four continual learning algorithms: MOLe Nagabandi et al. (2018),
MoB Pasula et al. (2023), MAML $k$-shot, and MAML-continuous. We evaluate on
these algorithms because (1) each uses meta-learned priors $\theta^{*}$ for
few-shot adaption to distribution shifts and (2) to cover a comprehensive set
of continual learning algorithm /emphclasses: a selector of expert models, a
weighted mixture of expert and constructive non-expert models, a few-shot
model that adapts from a history of $k$ data points, and a few-shot model that
adapts continously from each new data point, respectively.

The first two, MOLe Nagabandi et al. (2018) and MoB Pasula et al. (2023), are
modular algorithms, instantiating models from meta-learned priors when a new
task is detected. In comparison to MOLe, a state-of-the-art approach for its
time, MoB superceded MOLe by attaining both better overall performance and
requiring fewer additional models to do so on multiple diverse domains while
controlling for model architecture.

The latter two, MAML $k$-shot and MAML-continuous use a single model adapted
from a meta-learned prior. MAML $k$-shot adapts using the $k$ latest data
points, and MAML-continuous updates the model at each time step using the most
recent observation. We chose $k=13$ for WTI and Brent and $k=3$ for Dubai
Crude based on prior experience with similar time-series data and to align
with the notion of quarterly duration, which is of key importance in financial
domains.

We evaluate the usefulness of our contextual task labels by comparing the
percent change in forecasting accuracy over financially and seasonally
meaningful time-horizons with and without the inclusion of these task labels
as model inputs. Since the original Dubai Crude time-series were average
monthly prices, the horizons for this and for WTI and Brent differ. Results of
the improvements after including the task labels are shown in Figure 4.

By including our EM-based algorithm’s automatically generated contextual task
labels, we see universal improvement over each benchmark, time-horizon, and
model considered. Since the raw data for Dubai Crude have already been
averaged over month-long periods, we expect the performance changes on this
benchmark to be relatively small compared to those of WTI and Brent, which is
indeed what we see in Figure 4.

Figure 4: The Percent Improvement in MSE after including our algorithm’s
generated contextual task labels as model inputs for each continual learning
algorithm on each benchmark dataset over three financially and seasonally
meaningful Forecasting Horizons. Each bar plot shows the mean and standard
deviation "error" bars of 10 different seeds, so a total of 360 trials were
conducted.

##  6 Limitations

While our work on time-series benchmarks in the financial domain provides
valuable contributions, there are certain limitations to be acknowledged.
These limitations include:

_Domain Specificity_ : The benchmark datasets we introduce focus solely on the
financial domain, specifically crude oil prices. While this is an important
domain, the applicability of our benchmarks to other domains may be limited.

_Dataset Granularity_ : The datasets we provide have varying levels of
granularity, with daily spot prices for WTI and Brent and monthly average
prices for Dubai Crude. This difference in frequency may impact the
generalizability of our benchmarks across different time-series analysis
tasks.

_External Factors_ : Our benchmarks consider major real-world events and
contextual labels aligned with price and regime shifts. However, there may be
additional external factors and events that could impact the performance of AI
algorithms on the datasets, which are not explicitly captured in our work.

Addressing these limitations and expanding the scope of our benchmarks to
encompass diverse domains and data granularities would be valuable directions
for future research.

##  7 Discussion

In conclusion, we have introduced three novel time-series benchmark datasets
that encompass the global crude oil market. These datasets were meticulously
processed to fit into tasks using an expectation-maximization algorithm,
providing valuable insight into temporal distribution shifts and their
significant impact on predictive models. The benchmarks come equipped with
labels aligning with significant real-world events, creating a rich context
for further exploration and understanding.

Furthermore, we empirically show that including our algorithm’s automatically
generated task labels as model inputs universally improves performance over
all benchmarks, time-horizons, and models considered.

Our goal is to fuel advancements in areas such as continual learning and out-
of-distribution detection, while addressing challenges posed by sequential
data like non-stationarity and correlatedness. We believe that these
benchmarks, due to the societal importance of the assets they represent, can
make a substantial contribution to AI research, particularly in handling real-
world data that contains distribution shifts. As AI continues to permeate
every sector, the importance of such representative, real-world benchmarks
will continue to grow.

## Acknowledgments

We would like to extend our sincerest gratitude to the IJCAI 2023 AI4TS
Organizers and Program Committee. Their insightful reviews, constructive
feedback, and recognition of the importance of this work have been
instrumental in shaping the final manuscript. Their dedication to fostering
innovation and excellence in AI is truly commendable, and we’re honored that
they have awarded this work in recognition of promoting these ideals.

## Ethical Statement

There are no ethical issues.

## References

  * Backus and Crucini [2000] David K Backus and Mario J Crucini.  Oil prices and the terms of trade.  Journal of international Economics, 50(1):185–213, 2000.
  * Bellemare et al. [2013] Marc G Bellemare, Yavar Naddaf, Joel Veness, and Michael Bowling.  The arcade learning environment: An evaluation platform for general agents.  Journal of Artificial Intelligence Research, 47:253–279, 2013.
  * Caccia et al. [2020] Massimo Caccia, Pau Rodriguez, Oleksiy Ostapenko, Fabrice Normandin, Min Lin, Lucas Caccia, Issam Laradji, Irina Rish, Alexandre Lacoste, David Vazquez, et al.  Online fast adaptation and knowledge accumulation: a new approach to continual learning.  arXiv preprint arXiv:2003.05856, 2020.
  * Deng et al. [2009] Jia Deng, Wei Dong, Richard Socher, Li-Jia Li, Kai Li, and Li Fei-Fei.  Imagenet: A large-scale hierarchical image database.  In 2009 IEEE Conference on Computer Vision and Pattern Recognition, pages 248–255, 2009.
  * EIA [2022a] U.S. EIA.  Crude oil prices: Brent - europe [dcoilbrenteu], Sep 2022.  Source: U.S. Energy Information Administration. Retrieved from FRED, Federal Reserve Bank of St. Louis.
  * EIA [2022b] U.S. EIA.  Crude oil prices: West texas intermediate (wti) - cushing, oklahoma, Sep 2022.  Source: U.S. Energy Information Administration. Retrieved from FRED, Federal Reserve Bank of St. Louis.
  * He and Sick [2021] Yujiang He and Bernhard Sick.  Clear: An adaptive continual learning framework for regression tasks.  AI Perspectives, 3(1):1–16, 2021.
  * He et al. [2019] Xu He, Jakub Sygnowski, Alexandre Galashov, Andrei A Rusu, Yee Whye Teh, and Razvan Pascanu.  Task agnostic continual learning via meta learning.  arXiv preprint arXiv:1906.05201, 2019.
  * IMF [2022] IMF.  Global price of dubai crude [poildubusdm], Aug 2022.  Source: International Monetary Fund. Retrieved from FRED, Federal Reserve Bank of St. Louis.
  * Kilian and Park [2009] Lutz Kilian and Cheolbeom Park.  The impact of oil price shocks on the us stock market.  International Economic Review, 50(4):1267–1287, 2009.
  * Kim et al. [1998] Chang-Jin Kim, Charles R Nelson, and Richard Startz.  Testing for mean reversion in heteroskedastic data based on gibbs-sampling-augmented randomization.  Journal of Empirical finance, 5(2):131–154, 1998.
  * Kim [1994] Chang-Jin Kim.  Dynamic linear models with markov-switching.  Journal of econometrics, 60(1-2):1–22, 1994.
  * Krizhevsky et al. [2009] Alex Krizhevsky, Geoffrey Hinton, et al.  Learning multiple layers of features from tiny images.  2009\.
  * Lin et al. [2014] Tsung-Yi Lin, Michael Maire, Serge J. Belongie, Lubomir D. Bourdev, Ross B. Girshick, James Hays, Pietro Perona, Deva Ramanan, Piotr Dollár, and C. Lawrence Zitnick.  Microsoft COCO: common objects in context.  CoRR, abs/1405.0312, 2014.
  * Nagabandi et al. [2018] Anusha Nagabandi, Chelsea Finn, and Sergey Levine.  Deep online learning via meta-learning: Continual adaptation for model-based rl.  arXiv preprint arXiv:1812.07671, 2018.
  * Pasula et al. [2023] Pranay Pasula, Mengda Xu, and Sumitra Ganesh.  Mixture of basis for interpretable continual learning with distribution shifts.  https://openreview.net/pdf?id=5iibKv7Wk8W, 2023.
  * Phillips and Perron [1988] Peter CB Phillips and Pierre Perron.  Testing for a unit root in time series regression.  Biometrika, 75(2):335–346, 1988.
  * Stoica and Selen [2004] Petre Stoica and Yngve Selen.  Model-order selection: a review of information criterion rules.  IEEE Signal Processing Magazine, 21(4):36–47, 2004.

##  Appendix A Resampling Data

For daily frequency data, [Kim et al., 1998] returned unrealistic results, as
mentioned in Section 3.3.1. We resampled data by every two days, which
produced similarly unrealistic results. For WTI, the results for this
resampling frequency are shown in Figure 5. We found that resampling weekly
was an excellent frequency for [Kim et al., 1998] to generate tasks that were
reasonable and well optimized measures of likelihood and information criteria.

Figure 5: Unrealistic task generation on WTI daily spot price data with a
resampling frequency of every two days.

##  Appendix B Information Criteria

We compute the AIC, BIC, HQIC, and the sum of these of the fitted model for
$k=2,3,\dotsc,10$ for each oil and show the results in Table 1. Results for
$k>5$ are omitted because every model fit to this task count resulted in
greater sum of information criteria than the minimum (in bold), had at least
one task that didn’t have argmax probability for any timesteps, and didn’t
converge. We ran the EM model fitting algorithm on each $k$ 200 times, and if
no trial converged, then we list N under the row header Converge.

Oil | k Tasks | AIC | BIC | HQIC | Info Criteria Sum | Converge
---|---|---|---|---|---|---
WTI | 2 | 11025 | 11047 | 11033 | 33104 | Y
| 3 | 10963 | 11013 | 10981 | 32957 | Y
| 4 | 10956 | 11044 | 10989 | 32989 | N
| 5 | 10977 | 11115 | 11028 | 33119 | N
Brent | 2 | 10927 | 10949 | 10935 | 32812 | Y
| 3 | 10850 | 10899 | 10868 | 32617 | Y
| 4 | 10860 | 10949 | 10893 | 32702 | N
| 5 | 10882 | 11020 | 10933 | 32834 | N
Dubai | 2 | 2758 | 2774 | 2765 | 8297 | Y
| 3 | 2759 | 2795 | 2773 | 8328 | Y
| 4 | 2773 | 2837 | 2799 | 8409 | N
| 5 | 2790 | 2889 | 2830 | 8509 | N
Table 1: Information criteria values for each oil for various $k$

Even though Dubai Crude’s fitted model information criteria sum was minimized
for $k=2$, since WTI and Brent fitted models are both minimized for $k=3$ and
every task has argmax probability for at least one timestep for every oil for
$k=3$, we use $k=3$ for every oil.

Since Dubai Crude is sampled monthly, it’s not surprising that a model with
lesser $k$ than those for WTI and Brent, which are both resampled weekly,
results in lesser information criteria sum.


# Mapping the material distribution of a complex structure in an electron beam

Luise Poley  Department of Physics, Simon Fraser University, University Drive,
Burnaby, Canada  TRIUMF, Wesbrook Mall, Vancouver, Canada  Ulf Stolzenberg
II. Physikalisches Institut, Georg-August-Universität Göttingen, Friedrich-
Hund-Platz, Göttingen, Germany  Benjamin Schwenker  II. Physikalisches
Institut, Georg-August-Universität Göttingen, Friedrich-Hund-Platz, Göttingen,
Germany  Ariane Frey  II. Physikalisches Institut, Georg-August-Universität
Göttingen, Friedrich-Hund-Platz, Göttingen, Germany  Peter Göttlicher
Deutsches Elektronen-Synchrotron, Notkestraße, Hamburg, Germany  Carlos
Marinas  Instituto de Física Corpuscular, CSIC-Universidad de Valencia, c/
Catedrático José Beltrán, Paterna, Spain  Marcel Stanitzki  Deutsches
Elektronen-Synchrotron, Notkestraße, Hamburg, Germany  Bernd Stelzer
Department of Physics, Simon Fraser University, University Drive, Burnaby,
Canada  TRIUMF, Wesbrook Mall, Vancouver, Canada

###### Abstract

The simulation and analysis of High Energy Physics experiments require a
realistic simulation of the detector material and its distribution. The
challenge is to describe all active and passive parts of large scale detectors
like ATLAS in terms of their size, position and material composition. The
common method for estimating the radiation length by weighing individual
components, adding up their contributions and averaging the resulting material
distribution over extended structures provides a good general estimate, but
can deviate significantly from the material actually present.

A method has been developed to assess its material distribution with high
spatial resolution using the reconstructed scattering angles and hit positions
of high energy electron tracks traversing an object under investigation. The
study presented here shows measurements for an extended structure with a
highly inhomogeneous material distribution. The structure under investigation
is an End-of-Substructure-card prototype designed for the ATLAS Inner Tracker
strip tracker - a PCB populated with components of a large range of material
budgets and sizes.

The measurements presented here summarise requirements for data samples and
reconstructed electron tracks for reliable image reconstruction of large
scale, inhomogeneous samples, choices of pixel sizes compared to the size of
features under investigation as well as a bremsstrahlung correction for high
material densities and thicknesses.

###### Contents

  1. 1 Introduction
  2. 2 Multiple scattering models
  3. 3 Radiation length measurements
  4. 4 Measurements with test structures
  5. 5 Material distribution maps of extended structures
    1. 5.1 Structure under investigation
    2. 5.2 Measurement setup
  6. 6 Results
    1. 6.1 Comparison with known material
    2. 6.2 Pixel size
  7. 7 Impact on Tracking
  8. 8 Conclusion and Outlook

##  1 Introduction

Physics requirements on the vertex and momentum resolution of particles set
limits on the mass of vertex detectors for the next generation of collider
experiments [1]. Localised interactions of particles in the detector material
perturb the particle momentum, produce hadronic vertices with new secondary
particles or lead to photon conversions [2, 3]. To meet the conflicting goals
of minimising detector mass and also incorporating cooling, power distribution
and readout electronics, an integrated system design approach is needed [4].

Multiple scattering leads to random deflections of charged particles
traversing the detector. The magnitude of multiple scattering depends on the
local material budget measured in units of the material constant $X_{0}$. In
order to study the material composition of objects, a method has been
developed to estimate high-resolution 2D images of the present material
$X/X_{0}$ using high resolution tracking telescopes in test beam facilities
[5]. The amount of material can be extracted from the multiple scattering
deflections by using an appropriate theoretical model. The good spatial
resolution of the method opens a way to measure the overall scattering effect
in small areas. It can consequently be used to map fully passive materials
very precisely and has previously been used for small areas of low mass pixel
and strip detectors for the Belle II vertex detector [6].

For a large-scale application of the developed methods, a prototype End-of-
Substructure (EoS) card [7], an electronics board developed for the new Inner
Tracker [8] of the ATLAS detector [9], was chosen. The EoS card is a
comparably large structure (about 12$\times$5 cm) located within the ATLAS
tracking volume, which was designed to accommodate the combined ATLAS ITk
requirements of high reliability, radiation hardness and minimum material.
However, the functionalities of the EoS card (described in more detail in
section 5.1) require the use of several components containing high-density
materials. For a reliable description of the material distribution in the
detector simulation, a measurement of the material budget of the EoS card
(results shown in section 6) was therefore necessary.

##  2 Multiple scattering models

Multiple scattering is an unavoidable source of uncertainty in tracking and
vertex detector systems. The thickness $X$ of the material divided by a
material dependent quantity $X_{0}$, the radiation length constant, determines
the magnitude of scattering angles which can be expected. The radiation length
constant is approximately given by [10]

$\displaystyle X_{0}$
$\displaystyle=\frac{716.4\,\mathrm{g/}\mathrm{cm}^{2}\,A}{Z\left(Z+1\right)\,\ln{\left(\nicefrac{{287}}{{\sqrt{Z}}}\right)}}$
(1)

The $X$/$X_{0}$ measurements presented here are based on the reconstruction of
multiple scattering angles on a target material. Multiple scattering is the
net effect of a large number of individual Coulomb single scatterings on
nuclei of a traversed material and causes a direction change of the particle
trajectory which can be expressed by two projected angles $\vartheta_{u}$ and
$\vartheta_{v}$. There are several models, such as the Moliere [11] and
Highland [12] model, which can be used to describe the multiple scattering
angle distribution.

Moliere introduced a fully analytical model that describes the Gaussian core
as well as the non Gaussian tails of the multiple scattering distribution
[11]. The Highland model describes the Gaussian core of the multiple
scattering angle distribution. The standard deviation, called the Highland
width $\sigma_{\mathrm{HL}}$ is given by [13, 14]:

$\displaystyle\sigma_{\mathrm{HL}}$
$\displaystyle=\frac{0.0136\,\mathrm{GeV}}{\beta
p}\,\sqrt{\frac{X}{X_{0}}}\left(1+0.038\,\ln\left(\frac{X}{X_{0}}\right)\right)$
(2)

##  3 Radiation length measurements

The radiation length measurement approach presented here has already been used
to estimate material profiles of modules of the vertex detector of Belle II
[6]. More details on the approach can be found in [5]. The results presented
here are based on measurements conducted at the DESY II testbeam facility [15]
with a 2 GeV electron beam. The data analysis was conducted with the TestBeam
SoftWare (TBSW) framework [16]. An EUDET like reference telescope [17],
consisting of six Mimosa26 (M26) sensors, was used to measure the particle
trajectories. A schematic drawing of a generic telescope setup is depicted in
figure 1.

Figure 1: Schematic drawing of the telescope setup. The scattering target, the
EoS card, is placed between two telescope arms which consist of arrays of
three position sensitive M26 sensors each.

The basic idea of the radiation length estimation is to reconstruct multiple
scattering angles of charged particle transitions on a target plane centered
in a high-resolution telescope. Details about data analysis steps such as
clustering, tracking and alignment can be found in [18, 5]. The particle hits
on the sensors of the upstream telescope arm are used in a forward Kalman
filter to determine in-states $\eta_{\mathrm{in}}$ and their covariance matrix
$C_{\mathrm{in}}$ of the particle tracks on the target plane. The forward
filter is constructed such that the multiple scattering in air between
telescope planes is taken into account by adding virtual scattering planes
halfway between the telescope planes. The material of a telescope plane111The
overall thickness of each plane consists of 50 $\mu$m silicon and 25 $\mu$m
carbon is treated in an asymmetric way and contributes to the extrapolation to
the next downstream plane. In other words, the in-state at the scattering
target is possible without making an assumption on the material budget at that
plane. Hits in the downstream telescope arm are used in a time-reversed Kalman
filter pass starting from the last telescope plane. The unknown target
material does therefore not enter in the prediction of the track-states
$\eta_{\mathrm{in}}$ and $\eta_{\mathrm{out}}$ or their covariance matrices.

Without a magnetic field the tracks can be expressed by straight lines with
kinks on every material surface. The track states on the target plane with a
Cartesian coordinate system ($u$,$v$,$w$)222$u$ and $v$ are arranged in the
target plane, while $w$ is perpendicular to the surface can be therefore
expressed by two intersection coordinates $u$ and $v$, two slopes
$m_{u}=\nicefrac{{\mathrm{d}u}}{{\mathrm{d}w}}$ and
$m_{v}=\nicefrac{{\mathrm{d}v}}{{\mathrm{d}w}}$ and the particle charge $q$
divided by the particle momentum $p$, for example:

$\displaystyle\eta$
$\displaystyle=\left(\begin{array}[]{ccccc}m_{u}&m_{v}&u&v&q/p\end{array}\right)$
(4)

In- and out-states are matched according to a distance criterion. Calculating
the difference between the slope components of $\eta_{\mathrm{out}}$ and
$\eta_{\mathrm{in}}$ yields two projected scattering angles $\vartheta_{u}$
and $\vartheta_{v}$. The corresponding projected angle uncertainties can be
determined from error propagation when using the appropriate covariance matrix
entries. These uncertainties correspond to the telescope angle resolution
$\sigma_{\mathrm{reso}}$, which is an important parameter during the
calibration procedure. In order to perform spatially resolved measurements
each reconstructed angle must be associated with a intersection position. This
position is calculated as a weighted333Using the inverse intersection
variances as weights mean of the intersection entries in the track states.

The target plane is divided into small image pixels and the scattering angles
are filled into per pixel histograms using the estimated track intersection
coordinates. The reconstructed scattering angle distribution of the pixel
centered on position ($u$,$v$) can be described by a model function consisting
of convolution between a multiple scattering distribution
$f^{uv}_{\mathrm{msc}}$ and $f_{\mathrm{err}}$ for the angle reconstruction
error.

$\displaystyle
f^{uv}_{\mathrm{reco}}\left(\vartheta|\nicefrac{{X}}{{X_{0}}},p\right)=f^{uv}_{\mathrm{msc}}\left(\vartheta|\nicefrac{{X}}{{X_{0}}},p\right)\otimes
f_{\mathrm{err}}\left(\vartheta\right)\quad\mathrm{,}$ (5)

The scattering angle distribution $f^{uv}_{\mathrm{msc}}$ depends on the local
radiation length fraction $X$/$X_{0}$ of the traversed material and the beam
momentum. The density $f^{uv}_{\mathrm{msc}}$ can be expressed either by the
Highland model (eq. 2) or the Moliere model [11]. The angular resolution
function of the telescope $f_{\mathrm{err}}$ is approximated by a normal
distribution with zero mean and width $\sigma_{err}$. To first order, the
angle resolution is calculated from the uncertainties of the projected
scattering angles for in and out states. More precisely, it can be measured
from data during calibration as the width of the reconstructed angle
distribution when no material is inserted between the telescope arms. The
calibration factor $\lambda$ is the ratio between the calculated angle
resolution and the calibrated one.

In the Highland theory, we have an explicit formula for the fit model

$\displaystyle f^{uv}_{\mathrm{reco}}$
$\displaystyle=\sum_{i}w_{i,u,v}\mathcal{N}(\mu_{i,u,v},\sigma_{i,u,v})\quad\mathrm{,}$
(6)

where $\mathcal{N}(0,1)$ is a standard normal distribution. The sigma for
component $i$ in pixel $(u,v)$ is parametrised as

$\displaystyle\sigma_{i,u,v}=\sigma_{HL}(p(u,v),X/X_{0,i,u,v})^{2}+\lambda^{2}\sigma_{\mathrm{err}}^{2}\quad\mathrm{.}$
(7)

A parametrisation for the beam momentum $p$ is provided later in this section.
The weights $w_{i}$ can be interpreted as the fractional areas of a material
with thickness $X/X_{0,i}$ covering part of the pixel area. For homogeneous
pixels, a single component suffices to describe the central part of the
reconstructed angle distribution. For non homogeneous pixel, a mixture with
two components can improve the results. The total material length per pixel is
found as a weighed mean of the individual components.
As the radiation length estimation depends on accurate measurements of small
scattering angles 444$\mathcal{O}(100\,$µrad$)$ for thin materials such as 100
µm silicon, the position of the telescope planes and the intrinsic resolution
of the M26 sensors must be known precisely. To achieve this requirement a
telescope calibration step, which includes telescope alignment, measurements
of the cluster resolution and hot pixel masking must be conducted. Details on
the telescope calibration can be found in [18, 5]. Additionally, in order to
correct for small uncertainties in the telescope alignment, cluster
calibration and assumption of the nominal beam energy and telescope length a
calibration measurement with a well known material profile must be conducted.

In the frame of this calibration several measurement regions with a known
radiation length $X$/$X_{0}$ are defined. The reconstructed angle
distributions in these measurement regions can be described by an individual
model function in eq. 5. For each of these function the radiation length
parameter $X$/$X_{0}$ is fixed to the corresponding value. Calibration
parameters, such as $\lambda$, $\kappa$ and the exact linear dependence of the
beam energy $p(u,v)$ are shared by the model functions of all measurement
regions and can be determined by performing a simultaneous fit of the
reconstructed multiple scattering distributions.

The beam energy used in the fit model (see eq. 5) allows linear momentum
gradients and a scale correction $\kappa$ from the nominal momentum
$p_{\mathrm{nom}}$ set by the test beam facility.

$\displaystyle p\left(u,v\right)$
$\displaystyle=\kappa\left(p_{\mathrm{nom}}+\nabla_{u}\,p\cdot
u+\nabla_{v}\,p\cdot v\right)\quad\,\mathrm{.}$ (8)

In addition to the beam energy gradients, energy losses of electrons due to
bremsstrahlung in the scattering target decreases the momentum depending in
the local material distribution. The average energy loss can be modeled by the
Bethe Heitler model [19] and depends only on the local radiation length
$X$/$X_{0}$. Guided by experimental results, we used the following
parametrisation to correct for the effect of bremsstrahlung in our method.

$\displaystyle p_{\mathrm{wm}}^{\epsilon}(\nicefrac{{X}}{{X_{0}}},u,v)$
$\displaystyle=\left(1-\epsilon\right)\,p\left(u,v\right)+\epsilon\,p\left(u,v\right)\,\exp{\left(\nicefrac{{-X}}{{X_{0}}}(u,v)\right)}$
(9)

The beam momentum $p_{\mathrm{wm}}$ is a weighed mean of the momentum before
and after traversing the scattering target. The weight parameter $\epsilon$ is
a measure for the impact of the overall bremsstrahlung energy loss on the
multiple scattering distribution and has to be determined from measurements on
high $X$/$X_{0}$ materials. However, it is not a calibration parameter and the
task of determining $\epsilon$ is separate from the radiation length
calibration. Example measurements on aluminium and copper wedges with large
radiation lengths values will be presented in section 4.

For the data set presented here a calibration measurement was conducted with a
set of aluminium targets with a constant thickness. The seven aluminium layers
had thicknesses ranging from 0.5 to 6 mm of aluminium. Additionally,
scattering data without a scattering target were recorded. For each material
five million six hit tracks were reconstructed. The number of scattering
angles per measurement area was limited to 100 000 angles to ensure an
identical statistical weight of each angle distribution during the
simultaneous fit. The measurement with an air target is especially sensitive
on the calibration parameter $\lambda$ which is connected to the telescope
angle resolution. The calculated telescope angle resolution was found to be
$\sigma_{\mathrm{reso}}=\,$179 µrad. In order to determine the linear
dependence of the beam energy on the position on the target plane, measurement
areas should cover the whole range of the beam spot. Here measurement areas at
the four corners of the beam spot and a single measurement area near the
center of the beam spot are used. Figure 2 depicts one of the calibration
targets, aluminium with a thickness of 6 mm, with the corresponding
measurement areas marked by black and white rectangles.

Figure 2: Radiation length image of the 6 mm aluminium target. The measurement
areas providing scattering angle distributions for the simultaneous fit to
determine the calibration factors are shown as black rectangles. The
measurement areas are located in the central part of the image, where constant
radiation length values are observed.

The calibration parameters determined for $\epsilon$ values of 0.24, 0.30 and
0.36 are listed in table 1. According to these calibration parameters, the
calibrated telescope angle resolution $\lambda\,\sigma_{\mathrm{reso}}\approx
181\pm 1\,$µrad deviates from the calculated telescope angle resolution
$\sigma_{\mathrm{reso}}$ on the target plane by approximately 1$\,\%$.
$\kappa$ was found to be approximately 0.99, which indicates that there was a
small error in the assumed beam energy or telescope length. The beam energy
gradient in $v$ direction is compatible with
0$\,\nicefrac{{\mathrm{MeV}}}{{\mathrm{mm}}}$, while the gradient in $u$
direction is $\nabla_{u}\,p=-10\pm 1\nicefrac{{\mathrm{MeV}}}{{\mathrm{mm}}}$.

$\epsilon$ | $\lambda$ | $\kappa$ |  $\Delta_{u}p$ [MeV/mm] |  $\Delta_{v}p$ [MeV/mm]
---|---|---|---|---
0.24 |  1.011$\pm$0.005 |  0.988$\pm$0.003 |  -9$\pm$2 |  -1$\pm$3
0.30 |  1.011$\pm$0.005 |  0.987$\pm$0.003 |  -10$\pm$2 |  -1$\pm$3
0.36 |  1.012$\pm$0.005 |  0.985$\pm$0.003 |  -10$\pm$2 |  -1$\pm$3
Table 1: Calibration results for $\epsilon$ values of 0.24, 0.30 and 0.36.
$\lambda$ is the calibration factor of the telescope resolution. $\kappa$ is
the calibration factor of the scattering angle width and $\Delta_{u}p$ and
$\Delta_{v}p$ correspond to beam energy gradients in $u$ and $v$ direction.
The given errors correspond to the statistical uncertainties of the fit.

Beam energy gradients in this direction have been observed frequently in our
previous measurements [5]. The non-zero energy gradient in $u$ direction is
caused by the beam energy selection process at DESY, which includes the
horizontal spread of the multi-energetic particle beam via dipole magnet[15].
Two of the 40 scattering angle distributions used during the calibration
process together with their fit functions are depicted in figure 3.

Figure 3: Global fit to obtain calibration constants $\kappa$, $\lambda$,
$\Delta_{u}$ and $\Delta_{v}$. The fit minimises the joint $\chi^{2}$ for
fitting the scattering angles in 30 measurement areas on aluminium plates with
known thickness. The thickness is fixed and only the calibration constants are
floated. The data and fit is is shown for two example areas and show good
agreement.

The presented $X$/$X_{0}$ calibration procedure contains a validation step.
Using the determined calibration factors ($\lambda$, $\kappa$ and the beam
energy gradients) the angle distributions, which were used during the
calibration, are fitted to estimate the radiation length or analogously the
aluminium thickness. The determined values should match the thickness
assumptions used during the $X$/$X_{0}$ calibration. These estimated
thicknesses together with the mechanically measured aluminium thickness values
are depicted in figure 4. The results of the estimation show a good agreement
with the mechanically measured values. As can be seen, the selected small
differences in the $\epsilon$ values only have a minor effect on the thickness
measurements. The deviations between measurements with different $\epsilon$
values are smaller than the statistical uncertainties. Only for a $\epsilon$
value of zero a small deviation at the largest aluminium thickness is
observed. This is not surprising as the largest radiation length value
employed during the calibration corresponds to 6.7 % and the bremsstrahlung
correction only has a large impact for $X$/$X_{0}$ values beyond 10 %.

Figure 4: Comparison between mechanically measured thickness values of the
calibration targets with constant thicknesses (ranging from 0.5 to 6 mm of
aluminium) and the estimated thickness after the calibration for four
different $\epsilon$ values. As can be seen the estimated thicknesses are
compatible with the mechanically measured thicknesses within the limits of the
statistical uncertainty for nearly all thicknesses. Only for a thickness of 6
mm the $\epsilon=0$ estimated values show small deviations from the mechanical
measurement.

##  4 Measurements with test structures

Test structures of aluminium and copper with a wedge like profile provide a
test bench for the presented imaging algorithm. The samples are homogeneous
and the thickness profile is accurately known. The analysis focuses on
extracting quantitative radiation length profiles from scattering angles along
the wedge and experimentally finding good values for the weight parameter
$\epsilon$.The measurements where conducted in 2017 and 2018 at the DESY test
beam facility.

A schematic drawing and a photograph of the wedge are depicted in figure 5.
The wedge ranges from 6 mm at its thinnest point up to 30 mm of aluminium,
which corresponds to a radiation length $X$/$X_{0}$ of 33.7 % 555According to
[13] the radiation length constant of aluminium is $X_{0}$ = 88.97 mm. The
wedge is a rather large structure with a total length of 78 mm. The slope
angle of the wedge was determined to be $24.296^{\circ}$ with an optical
coordinate measurement machine (optical CMM). Due to the limited size of the
M26 sensitive area of approximately 10$\,\times\,$20 $\mathrm{mm}^{2}$ several
measurements at different positions along the wedge had to be conducted. The
distance between two measurement positions was 7 mm and a total of eleven
wedge positions were measured. Before the measurements calibration steps were
conducted as described in the last section.

Figure 5: Schematic drawing and photograph of the measured aluminium wedge.
The wedge has a continuous thickness profile ranging from approximately 5.4 mm
to 30 mm of aluminium which corresponds to a $X$/$X_{0}$ value of roughly 34
%. Images from [5]

The $X$/$X_{0}$ image of the wedge is depicted in figure 6. The image pixel
inside the black rectangle are used to determine the radiation length along
the slope of the wedge and compare them to values determined from the
precision measurement of the slope angle. The material profiles for six
different $\epsilon$ values are depicted in figure 7. The alternating grey and
white areas indicate the eleven different measurement regions along the wedge.
As can be seen the $X$/$X_{0}$ estimation with epsilon values between 0.26 and
0.34 match the expected profile along the slope of the wedge quite well. The
estimated radiation length values without bremsstrahlung correction, i.e. an
epsilon of zero, show large deviations from the mechanically measured profile.

Figure 6: Radiation length $X$/$X_{0}$ image of the aluminium wedge. Several
features of the wedge and the image pixels used to generate the profile in
figure 7 are indicated. Images from [5] Figure 7: Material profile of the
aluminium wedge. The different measurement positions along the wedge are
indicated by the alternating grey and white areas. The estimated $X$/$X_{0}$
values with epsilon 0.30 and 0.32 match the mechanically measured profile
(indicated by the black dashed line) the best. The estimated $X$/$X_{0}$
values without the bremsstrahlung correction ($\epsilon=0$) show large
deviations from the mechanical measurement, especially for large $X$/$X_{0}$
values.

Figure 8 displays the difference
$\Delta\nicefrac{{X}}{{X_{0}}}=\left(\nicefrac{{X}}{{X_{0}}}\right)^{\mathrm{meas}}-\left(\nicefrac{{X}}{{X_{0}}}\right)^{\mathrm{mech}}$
between the estimated radiation length values
$\left(\nicefrac{{X}}{{X_{0}}}\right)^{\mathrm{meas}}$ and the expected values
based on mechanical measurements
$\left(\nicefrac{{X}}{{X_{0}}}\right)^{\mathrm{mech}}$ as a function of the
corresponding aluminium thickness of the wedge for five different $\epsilon$
values. As can be seen an $\epsilon$ value between 0.3 and 0.32 shows the best
agreement with the expected values.

Figure 8: Relative difference between the estimated radiation length values at
specific points along the aluminium wedge and the corresponding values of the
mechanical measurement. Estimated values with a weight parameter $\epsilon$ =
0.30 and 0.32 match the mechanically measured values the best. Accordingly,
the optimal weight parameter lies between 0.30 and 0.32.

In order to confirm this weight factor and test its validity, a complementary
measurement with a copper wedge was conducted [5]. The measurements confirmed
our previous findings that for materials with large radiation length in the
order of 10% $X$/$X_{0}$, energy losses due to bremsstrahlung have to be
included in radiation length estimation to ensure unbiased $X$/$X_{0}$
measurements.

##  5 Material distribution maps of extended structures

###  5.1 Structure under investigation

The End-of-Substructure (EoS) card is a component designed for the Upgrade of
the ATLAS Inner Detector: the ITk strip tracker [8]. The ITk strip tracker
consists of four barrel layers in the central region and six end-cap discs in
the forward region, each populated with silicon sensor modules mounted on both
sides. Each barrel layer and end-cap disc will comprise identical
substructures: staves in the central region, holding 14 modules on either
side, and petals in the forward region, holding nine modules on either side.
EoS cards are mounted on both sides of each substructure to act as the
interface between modules and off-detector electronics. For the foreseen
detector geometry, particle tracks with a pseudo-rapidity $\eta$ between about
1.2 and 2.0, tracks are likely to traverse at least one EoS card.

The ITk is designed to have a material budget of less than 2 % $X/X_{0}$ [8]
along a particle track. However, locally, a higher material budget can not be
avoided due to the need for certain components, such as copper coils and heat
sinks on an EoS card. While locally high material budgets can not always be
avoided, it is important to model their material budget correctly in tracking
and tracking simulations.

A prototype EoS card was therefore chosen as a test object for these material
budget measurement studies: as an extended structure (about 12$\times$5 cm,
see figure 9) with a large range of components it provides an opportunity to
study material budget imaging for a range of structures and feature sizes. In
addition, the determination of its overall material budget provides a valuable
input for ITk tracker simulations.

Figure 9: Photograph of the EoS card under investigation. In order to protect
wire bonds on the chips mounted on the EoS, sensitive areas are protected by a
cover (shown in red) on the circuit board.

The EoS card studied here is a populated printed circuit board (PCB). It is
designed to distribute the control signals and power to the front-end modules
and collects, serialises and transmits the up-stream data to an off-detector
data acquisition system. For the prototype shown here, the active high data
rate part is performed within an ASIC designed for the LHC experiments, the
GBTX [20]. Another ASIC, the GBT-SCA [21], collects slowly changing
environmental information. The communication to the off-detector electronics
is performed with a bi-directional 5 Gbit/s optical link. For the prototype,
the electrical-optical conversion is done using commercial SFP$+$ plugins,
which are planned to be replaced by low weight radiation hard optical modules
for the final design. In addition to the active components, passive C-L-C pi-
filters are used for different functional blocks.

The PCB itself is a multi-layer stack of interleaving glass-epoxy (FR4) and 10
copper foils (see figure 10). It has a total thickness of 1.2 mm ($\pm$10 %).
The central layer is a 50 $\upmu$m polyimide layer. Most of the homogeneous
material is made up by six copper layers - nearly full planes of individual
thicknesses from 9 $\upmu$m to 35 $\upmu$m, adding up to a combined thickness
of 120 $\upmu$m.

In the transverse structure, the PCB layers are interrupted by vias, vertical
connections between layers. Most of them are designed as vias connecting
through the whole PCB (through vias). They have a central hole with a radius
of 0.15 mm, visible as low material areas of circular shape. Corresponding
drills have a depth of 1.2 mm and their cylindrical surfaces are plated with
$\approx 20$ $\upmu$m of copper, visible as areas with more material. An
additional contribution to areas with visibly more material comes from copper
rings with a radius of 0.3 mm surrounding the drill holes in each of the ten
layers (polyimide and copper), which add up to a total thickness of 0.2 mm of
copper.

A different set of holes in the PCB corresponds to vias from the outer layer
to layer 3, with a depth of only 0.23 mm and the diameters of both copper
rings and drill holes reduced by 1 mm each.

Figure 10: Layer stack of the EoS card and the structure of a through going
interconnect (via).

Known features were used to compare the estimated material with the known
material budget in the corresponding area (see section 6.1).

###  5.2 Measurement setup

Measurements were conducted at Beamline 21 at the DESY II testbeam facility
using electrons using the DATURA beam telescope [22]. While other types of
particle beams, e.g. muon beams, would be better suited for multiple
scattering measurements due to lower energy losses and no bremsstrahlung,
electron beams provide a good compromise between availability with high beam
currents and traceable range of scattering angles.

Several data runs were performed for the mapping and data analysis of the
structure under investigation. It is important to note that the measurement
setup was unchanged between different measurements: a constant beam energy of
2 GeV was used in combination of a telescope sensor threshold setting of 5
(corresponding to a threshold 5 $\sigma$ above the sensor noise). For all
measurements related to the EoS card, the telescope planes were positioned as
listed in table 2.

Plane | Position $[$mm$]$
---|---
0 | 0
1 | 147.5
2 | 303.5
3 | 523.5
4 | 673.5
5 | 828.5
Table 2: Positions of telescope planes (along the beam axis, with the first
traversed plane 0 positioned at $z=0$) for all material distribution
measurements. The structure under investigation was placed between planes 2
and 3.

Data runs with the following targets were performed (in accordance with the
measurement routine explained in section 3):

  * •

no target (”air run”)
data taking without any target for telescope calibration

  * •

calibration runs
data taking with aluminium plates of known, precise thicknesses to calibrate
calculated material using known samples, about 1 million telescope events per
calibration target. The obtained reconstruction constants are summarised in
table 3

  * •

material distribution mapping
Mapping the full EoS card by systematically moving the structures, mounted on
precision translation stages, with respect to the beam. 5-6 million telescope
events were collected per position.

  * •

high statistics runs
for two structures of particular interest (readout chips with a heat sink),
large data samples were collected (about 50 million telescope events for a
single stage position) in order to produce high-resolution maps

Constant | Value
---|---
$\lambda$ | 1.044
$\kappa$ | 0.95
$\Delta_{u}p$ | 7.3 MeV/mm
$\Delta_{v}p$ | 0.8 MeV/mm
Table 3: Calibration constants calculated based on calibration measurements as
defined in section 3. All measured constants are in good agreement with the
values used to determine a reasonable bremsstrahlung correction.

The size of an individual material distribution map is determined by the size
of the beam telescope’s MIMOSA-26 sensors [23] ($1\times
2\,\mathrm{\text{cm}^{2}}$ [17]), since only tracks passing through the active
sensor area contribute to the reconstructed image (see figure 11). In order to
use the full telescope sensor area, a beam collimator with the corresponding
size ($1\times 2\,\mathrm{\text{cm}^{2}}$) was chosen for these measurements.

Figure 11: Number of reconstructed tracks per pixel in a material distribution
map from 5 million events: 400 tracks per pixel were required for a reliable
fit. The white frame shows the area assumed to be usable per stage position as
a basis of estimate of the stepping size.

In order to reconstruct the material budget using the fit described in section
3, a minimum number of 400 reconstructed tracks per pixel was required in
order to consider the resulting result reliable. As a result, pixels at the
edges of the imaging window are only reconstructed for large data samples,
where sufficient tracks are found even in the edge regions of the image. The
stepping size between subsequent stage positions was therefore chosen to
minimise overlap between adjacent images (and therefore minimise the number of
stage positions where data had to be collected) while providing hermetic
coverage between adjacent images. Stage positions were therefore changed by 19
mm in x and 9 mm in y (see figure 11), with 5-6 million telescope events per
stage position being collected. For each beam position, data was collected for
45-55 minutes (for 4-5 million telescope events), at a trigger rate of about
1.5 kHz. Mapping the full structure (excluding detail images mapped with
larger data samples) required a total of 39 hours.

This size was chosen as a compromise to maximise the reconstructed area per
position (see figure 12(a)) and precision of the reconstructed value (see
figure 12(b)), while minimising the beam time required per position.

(a) Illustration of the impact of the number of reconstructed events on the
usable image area: for a minimum track requirement of 400 tracks, a minimum of
4 million events was found to result in a reasonably large reconstructed image
area

(b) Correlation of the reconstructed material budget for an image
reconstructed using 4 million and 66 million events. While low material
budgets showed a good correlation for both data sets (correlation factors r =
0.93 for values $<4.5\,\mathrm{\%\,X/X_{0}}$ and r = 0.89 for values
$>4.5\,\mathrm{\%\,X/X_{0}}$), large discrepancies were found for large
material budgets, since small variations in the histogram can have large
impacts on the fitted distribution width.

Figure 12: Impact of the size of the data set used for reconstruction on the
reconstructed image area and pixel value.

Images were reconstructed requiring a minimum of 400 tracks per pixel, which
was found to provide reliable results while leading to a sufficiently large
reconstructed image area per position (see figure 13.

Figure 13: Image area that can be reconstructed from a 4 million data sample
using different required minimum track cuts.

##  6 Results

Due to the size of the EoS card under investigation with respect to the
support structure used to mount it in the electron beam, the EoS was mounted
on the support structure in two different positions and mapped in two series
of movements (see figures 14(a) and 14(b)).

(a) End of Structure card (top half) mapped in electron beam

(b) End of Structure card (bottom half) mapped in electron beam

Figure 14: Both sections of an EoS card mapped in electron beam in two
sessions

Both sections of the EoS card were afterwards joined using an algorithm that
varied the planar rotation angle and displacement to find a maximum overlap
between the common EoS regions on both maps (see figure 15).

Figure 15: Full map of $X/X_{0}$ for an EoS card, calculated bremsstrahlung
correction factor of $\epsilon=0.3$, computed from two separate measurements.
Combining both individual measurements into one overall EoS map showed a good
agreement between both separate measurements.

An improvement of the imaging algorithm was the implementation of a
bremsstrahlung correction described above. By accounting for energy losses
caused by bremsstrahlung, the resulting overestimation of the reconstructed
material budget was compensated.

The obtained map was used for a comparison with the common approach of
averaging all material over the full component area: calculating an average
for the full EoS card results in 3.7 % $X/X_{0}$, assuming that no extra
volumes are added for any of the components with more material, e.g. the
copper cooling plate (see figure 16).

Figure 16: Distribution of reconstructed material per $100\times
100\,\mathrm{\upmu\text{(}m)^{2}}$ pixel for the full EoS card area (with and
without correction). While the average amount of material is small, individual
pixels show up to more than 20 times the estimated average (for values with
bremsstrahlung correction).

This average can be used to calculate the discrepancy between average and
actual amount of material per bin (see figure 17).

Figure 17: Map showing the local difference between an overall average and the
actual material

The map shows several areas where material is underestimated, in addition to
several components with significantly more material than estimated from the
average.

###  6.1 Comparison with known material

In order to judge the quality of the reconstruction method, known features of
the EoS card were used for comparison:

  * •

two pin headers with a diameter of 0.64$\times$0.64 mm2, made out of phosphor
bronze, were mounted to extend over the edge of the EoS card (see figure
18(a)) and could be used for a direct material comparison.

(a) Region of EoS card with two pins extending over the edge of the EoS (edge
of the PCB marked in white)

(b) Average material per vertical scan line per pin, calculated with
bremsstrahlung correction

Figure 18: Analysis of reconstructed material for two pin headers extending
over the edge of an EoS card, i.e. with no background material except for air.
The pin headers were treated as homogeneous (except for some solder residue on
the right pin), which agreed well with the estimated values.

While the exact composition of the bronze pins is unknown, it can be estimated
as a standard bronze alloy consisting of 60 % copper
($X_{0}=14.36\,\mathrm{\text{mm}})$ and 40 % tin
($X_{0}=12.06\,\mathrm{\text{mm}})$, leading to a combined radiation length of
13.44 mm. With a thickness of 0.64 mm, each pin can be expected to correspond
to a material contribution of 4.8 % $X/X_{0}$, which agrees well with the
estimated values within uncertainties. Additionally, both pins are consistent
with each other (except for an area with solder residue on one pin) and no
gradient is observed over the length of the pin, which indicates a correctly
determined beam energy calibration factor.

  * •

Two wires with a cross-section of 1 mm2 can be seen leading away from the EoS
card (see figure 19).

Figure 19: Mapped image area showing two wires with insulation leading away
from the EoS card.

Due to the round cross-section of the cable, a fit function was used to
calculate the radiation length of both the metal core of the cable and the
insulation material (see figure 20).

Figure 20: Reconstructed material along one line across an individual cable
and applied fit function: the low material regions correspond to image areas
with insulation only, the centre shows the metal core of the cable.

The fit results were used to calculate the radiation length X0 of both the
metal and the insulation comprising both cables (see figures 21(a) and 21(b)).
The determined radiation lengths show good agreement along each cable and
between the two cables, which were of the same type. While the exact material
contributions of both cables are unknown, the determined values are in good
agreement with a metal on the inside and a plastic on the outside of the
cable.

(a) Radiation length of metal cable core, determined from fit, for left and
right cable (determined with applied bremsstrahlung correction)

(b) Radiation length of cable insulation material, determined from fit, for
left and right cable (determined with applied bremsstrahlung correction)

Figure 21: Radiation lengths of cable materials determined per image line. The
results are consistent within uncertainties within uncertainties for both
cables along the cable length.

  * •

drill holes with vias can be seen in several areas of the EoS card: they
traverse all layers of the EoS, with a 20 $\upmu$m copper layer on the inside
of the drill hole for electrical contact (see figure 10). For a comparison
with the expected material, a drill hole was chosen which was not covered by
the protective layer on the EoS card (see figure 9). A picture of the hole
used for analysis is shown in figure 22.

Figure 22: Detail area of an EoS imaging scan showing a drill hole with a vias
in an unobstructed area of the EoS card. Dashed lines indicate positions used
to analyse the material distribution of the drill hole (see figure 23).

As shown in figure 10, an EoS FR4 stackup without any copper amounts to
approximately 0.75 % $X/X_{0}$, which is expected for the material surrounding
the drill hole. Since there is no material within the hole, the reconstructed
material is expected to be zero (within uncertainties) in that area. Figure 23
shows the obtained material distribution.

Figure 23: Material profile along horizontal and vertical profile lines
positioned around a drill hole with vias.

It was found that the calculated values agree well with the expectations: in
the absence of copper traces in the PCB, the reconstructed material reaches a
value of 0.75 % $X/X_{0}$, as well as reconstructed material consistent with
zero (within uncertainties) at the centre of the drill hole. The ring
surrounding the drill hole was found to be significantly lower than expected
from the copper wall inside the vias, which corresponds to 1.2 mm of copper
(8.4 % $X/X_{0}$). It has a thickness of about 20 $\upmu$m or less than half a
pixel length, which leads to a combined value for PCB and vias areas.

Figure 24: Scattering angle distribution from a pixel at the edge of a vias
showing the contributions from areas with little material and significantly
more scattering material. While the standard approach of a single Gaussian fit
(dashed red line) does provide a reasonable approximation of the narrow peak
(corresponding to air), it underestimates the impact of the wider peak. In a
second approach, the two peaks were first fitted separately (black line
fitting the narrow peak (air), cyan line corresponding to the wide peak
(copper-plated wall of the vias)) and the obtained fit parameters used as
input parameters to a double Gaussian fit (violet line), which shows a better
approximation of the scattering angle distribution.

In figure 24, this effect is demonstrated: a single Gaussian fit results in an
overall material estimate of 0.5 % $X/X_{0}$. Fitting both peaks separately
yields two distinct Gaussians: a narrow peak (0.01 % $X/X_{0}$, consistent
with air) and a wide peak (5.7 % $X/X_{0}$, which is lower than expected for
1.2 mm of copper, but can be expected due to the thin wall thickness and large
scattering angles, which could scatter the particles out of the copper wall).
Weighing them based on the areas under the respective curves leads to a
combined value of 2.5 % $X/X_{0}$, which is distinctly larger than the
estimate based on a single Gaussian fit and, judging from the distributions, a
more reasonable approximation.

Studies of dedicated areas on the EoS card found a good agreement between
reconstructed and known material, where available.

###  6.2 Pixel size

The resolution that can be achieved from electrons is determined by several
factors:

  * •

the intrinsic resolution of the beam telescope [17]

  * •

offsets in the reconstructed tracks which can occur for thick layers of
materials with high density

  * •

the minimum pixel size that can be used for the image reconstruction process.

The image pixel size is usually determined by the size of the data sample:
since a minimum number of tracks is required to perform a reliable fit, the
overall image needs to be divided into pixels containing the required minimum
of tracks. Scattering angles from different materials within the same pixel
lead to a smeared Gaussian distribution. Applying a single Gaussian fit
function to the resulting angle distribution can naively be assumed to result
in a reconstructed scattering angle width corresponding to an average of the
different material contributions.

This assumption was tested in these measurements by collecting a large data
set for an area under investigation that was known to have an uneven material
distribution: a readout chip with a heat sink. A data set of 56 million events
was obtained to provide a sufficient number of tracks per pixel to ensure a
reliable fit and the same data sample was reconstructed using pixel sizes of
50, 100, 200 and 400 $\upmu$m (see figure 25)

Figure 25: Image reconstructed for a readout chip with a cooling plate for
different pixel sizes. A correlation plot for the same area reconstructed
using different pixel sizes shows that, for this area, a larger pixel size can
cause a systematic shift towards lower values of material.

Comparing the reconstructed images showed that up to a pixel size of 200
$\upmu$m, the maximum and minimum reconstructed material budget values were
comparable for all images. For a pixel size of 400 $\upmu$m the reconstructed
image was found to lead to an overall reduction of the average material
budget. Additionally, a pixel size of 400 $\upmu$m was found to result in
pixels where the obtained track angle distribution could not be described by
the applied fit function anymore and the fit failed (see figure 26).

Figure 26: Example of a scattering angle distribution that can not be
described by a single Gaussian fit: due to an inhomogeneous material
distribution within a pixel, the scattering angle distribution consists of an
overlay of Gaussian distributions and can not be approximated by a Gaussian
fit function. Using a combination of two independent Gaussian distributions
(cyan line) improves the description of the observed scattering angle
distribution.

For a quantifiable comparison of the impact of larger pixel sizes, images with
retroactively increased pixel sizes were constructed for comparison with their
reconstructed equivalent: matrices of smaller pixels were averages into larger
pixels and the resulting image was compared with an image reconstructed with
that same pixel size (see figure 27).

Figure 27: Deviation maps of images constructed by averaging smaller pixels
with the corresponding image reconstructed with a larger pixel size

The resulting maps showed that the pixel size used for reconstruction has
little impact on the reconstructed material budget for homogeneous material,
but does affect the reconstructed value in very inhomogeneous regions. For the
region under investigation, the size of individual features was about 200
$\upmu$m, and the obtained maps showed only localised discrepancies for pixel
sizes up to that size, with only small and localised discrepancies between
maps reconstructed using pixel sizes of 50 $\upmu$m and 100 $\upmu$m. However,
larger pixel sizes show systematically lower material budgets than smaller
pixel sizes, which affect the overall reconstructed average and lead to a
systematically lower reconstructed value.

Studies with random number generated track angles show that in cases of pixels
with different fractions of material within the pixel, fitting the resulting
superposition of Gaussian distribution leads to a calculated distribution
width that is dominated by the biggest fraction of material contained in the
pixel rather than the expected combination of both values.

Since the EoS region shown here has a mostly low material budget with pillars
of more material occupying a smaller fraction of space, pixels can be assumed
to frequently have smaller fractions of large material areas and therefore to
disproportionately reflect the lower material budget areas, which leads to a
general underestimation. Studying this area led to the conclusions that in
general, pixel sizes for image reconstruction should be chosen to be smaller
than features on the imaged object in order to avoid the surrounding area to
dominate the obtained material budget.

Additionally, work was started on a method to account for uneven material
distributions within one pixel by fitting track angle distributions with a
mixture of two Gaussian densities and using their widths to calculate a
weighted mean (see figure 26) rather than a single Gaussian fit (as defined in
the algorithm used here, see section 3), which is planned to be included in
future versions of the reconstruction algorithm.

##  7 Impact on Tracking

After reconstructing the material distribution for the EoS card under
investigation, the obtained map was used to study the impact of the material
distribution on the tracking performance. For this study, a subset of 1
million events in one area of the EoS card with a range of material
thicknesses was used (image area with pin headers, see figure 18(a)). Tracks
were reconstructed for three scenarios:

  * •

air, where no material was estimated for the EoS card layer

  * •

constant, where a constant, flat material density of 3.7 % was assumed for the
whole EoS area, corresponding to the calculated average over the full EoS area

  * •

map, where the estimated material distribution was used for the reconstruction

For each scenario, a straight line fit of tracks was applied, where both
multiple scattering and energy loss were taken into account. The straight line
fit was implemented using a Kalman filter starting from a hit on the first
telescope plane. In the absence of a magnetic field, the momentum affects the
track fit only via the multiple scattering variance from the Highland model.
The initial momentum at the first plane was set to the beam energy and
subsequently decreased by the mean energy loss caused by bremsstrahlung
depending on the traversed material.

A first comparison of the material description was conducted by extrapolating
tracks found in the first three telescope planes onto the first telescope
sensor behind the EoS card. The closest hit within a window of $\pm 0.8$ mm
from the extrapolated track intersection was matched. The map of track
intersections in figure 28 shows a smooth beam profile. This map is
independent of the applied material description and therefore identical for
all investigated cases.

Figure 28: Number of track intersections per pixel reconstructed based on a
triplet reconstructed in the first three telescope planes and matched to a hit
in the first telescope plane behind the EoS card.

The standard deviation of position residuals
$x_{\mathrm{fit}}-x_{\mathrm{hit}}$ has a multiple scattering component that
scales with the square root of the traversed material at the EoS card. The non
homogeneous material distribution at the EoS card is indeed visible when
mapping the standard deviation of residuals, see figure 29.

Figure 29: Map of the RMS of residuals per pixel for tracks reconstructed
based on a triplet in the first three telescope planes and matched to a hit in
the first telescope plane behind the EoS card.

The improved material description by a detailed $X/X_{0}$ map leads to an
improved estimation of the Kalman filter for the variance of extrapolated
track intersections $x_{\mathrm{fit}}$. Figure 30 shows the distribution of
standardised residuals (pulls) for all three cases (material descriptions as
air, constant and map).

Figure 30: Comparison of pull values for hits in the EoS plane using three
different material models for the EoS card. The width of the distributions
illustrates that the air model underestimates the scattering material, while
the constant distribution overestimates it in areas. The map model provides
the most realistic estimate.

A Gaussian fit to the standardised residual histogram shows that the map
description leads to the best fit to a normal distribution and the fitted
standard deviation is closest to one, i.e. provides the best estimate for the
uncertainty of extrapolated hit positions. In all three cases non Gaussian
tails are observed in the standardised pulls. These tails are already present
in the residual histogram (see figure 31) and likely originate from non
Gaussian multiple scattering or wrongly matched hits during track finding.

Figure 31: Residual distribution for all hits in a track sample reconstructed
from the same data set using different material descriptions. As expected, the
residuals are not affected by the material description of the EoS card.

Only tracks with hits on all six sensors and a total track $\chi^{2}$-value
below 40 where accepted. In addition, new hits were added to the track if the
$\chi^{2}$-value of the predicted position residual was $<10$. It should be
noted that the material description at the central EoS card determined if a
hit from a downstream sensor was added to the track. As a consequence, the
tracking was repeated for all three material descriptions, which led to three
different track samples:

  * •

air: 2.4 million tracks

  * •

constant: 3.5 million tracks

  * •

map: 3.4 million tracks

A material description as air was found to result in many missed tracks,
because large scatterings at the EoS card resulted in high $\chi^{2}$ values
due to the underestimated material at the EoS card. The $\chi^{2}$ probability
distribution of the track fit is shown in figure 32.

Figure 32: Frequency distributions of the $\chi^{2}$-probabilities for tracks
reconstructed from the same data set using different material descriptions.
Due to the different numbers of reconstructed tracks, the distributions were
normalised for comparability. As expected, the air sample shows a large
fraction of tracks with an underestimated scattering angle, while the flat
estimate of the constant distribution, based on the full EoS card with several
areas with high material density, overestimates the material for a large
fraction of tracks in this subset of data. The Kolmogorov-Smirnov test was
applied to quantify the difference between a map description with and without
inclusion of a bremsstrahlung correction (see table 4).

For a material description as a map, the distribution is closest to a uniform
distribution in $[0,1]$ (see table 4). While the map description was found to
be an improvement, multiple scattering and energy loss by bremsstrahlung are
non Gaussian physical processes and provide only an approximation to actual
scattering processes, resulting in deviations from a uniform distribution.

Model | P-value
---|---
Air | 0.06
Constant | 0.38
Map, without bremsstrahlung | 0.60
Map, with bremsstrahlung | 0.68
Table 4: P-values of the Kolmogorov-Smirnov test comparing the observed
distribution of the $\chi^{2}$-probabilities of the track fit with the ideal
flat distribution. As expected, the P-value improves with improved material
description.

This study confirmed that implementing the reconstructed material distribution
for tracking reduces over- or underestimating track uncertainties. It also
helps to avoid the loss of tracks due to an underestimate of scattering
material. The study confirmed that the explicit description of non homogeneous
materials provides an advantage over using averaged material.

##  8 Conclusion and Outlook

The studies presented here were conducted to test the radiation length imaging
algorithm, developed for electron beams, for structures of different
complexities. Test structures with a range of material thicknesses were used
to expand the existing algorithm for measurements of objects with thicknesses
of up to 40 % $X/X_{0}$. The measurements on wedges motivated the
implementation of an empirical correction for the energy loss of electrons by
bremsstrahlung when passing the device under study.

The correction was tested further in measurements of an EoS card for the ATLAS
ITk as a structure with highly inhomogeneous material distribution. The
measurements found a good agreement between the present material, where known,
and the reconstructed values.

The measurements investigated the impact of different parameters on the
reconstructed material distribution and found reliable working points for
parameters such as stepping size between stage positions for adjacent images,
minimum number of required tracks per pixel and minimum number of events to be
used for analysis per image. The measurements confirmed that the
reconstruction of scattering angles for an object consisting of different
layers, e.g. the wires surrounded by isolation material and the PCB with
embedded copper layers, results in the expected overall material budget. In
contrast, areas with highly inhomogeneous material distributions within one
reconstructed image pixel showed a systematic shift in the reconstructed
material. A single Gaussian distribution was found insufficient to describe
more complex material distributions within one imaged pixel, which can be
avoided by choosing a pixel size below the feature size of the imaged object.

Future modifications of the analysis code used here are therefore planned to
incorporate a double Gaussian fit which does take mixed contributions within
the same pixel into account: instead of implementing a single fit to
approximate mixed contributions, using a fit function that contains Gaussian
contributions with different distributions widths has been found to produce
more reliable results for pixels with mixed contributions.

## Acknowledgements

The work presented here was partially funded by the Canada Foundation for
Innovation (CFI) and the Natural Sciences and Engineering Research Council
(NSERC) of Canada as well as the Alexander von Humboldt Foundation. We
acknowledge the financial support by the Federal Ministry of Education and
Research of Germany.

The measurements leading to these results have been performed at the Test Beam
Facility at DESY Hamburg (Germany), a member of the Helmholtz Association
(HGF). The authors would like to thank the DESY testbeam coordinators and beam
telescope support, especially Dr. Jan Dreyling-Eschweiler.

## References

  * [1] William E. Cooper.  Low-Mass Materials and Vertex Detector Systems.  PoS, Vertex2013:036, 2013.
  * [2] ATLAS Collaboration.  Study of the material of the ATLAS Inner Detector for Run 2 of the LHC.  JINST, 12(12):P12009, 2017.
  * [3] A.M. Sirunyan et. al.  Precision measurement of the structure of the CMS inner tracking system using nuclear interactions.  Journal of Instrumentation, 13(10):P10034–P10034, oct 2018.
  * [4] M. Demarteau et. al.  Planning the Future of U.S. Particle Physics (Snowmass 2013): Chapter 8: Instrumentation Frontier, 2014.
  * [5] U. Stolzenberg.  Radiation length measurements with high-resolution telescopes.  PhD thesis, Universität Göttingen, 2019.  II.Physik-UniGö-Diss-2019/06.
  * [6] U. Stolzenberg, A. Frey, B. Schwenker, P. Wieduwilt, C. Marinas, and F. Lütticke.  Radiation length imaging with high-resolution telescopes.  Nuclear Instruments and Methods in Physics Research Section A: Accelerators, Spectrometers, Detectors and Associated Equipment, 845:173 – 176, 2017.  Proceedings of the Vienna Conference on Instrumentation 2016.
  * [7] C. Wanotayaroj, H. Ceslik, H. Colbow, S. D. Cornell, P. Goettlicher, I. M. Gregor, A. Melnik, M. Stanitzki, and J. Wolff.  The End-of-Substructure (EoS) Card for the Strip Tracker Upgrade of the ATLAS experiment.  In 2018 IEEE Nuclear Science Symposium and Medical Imaging Conference Proceedings (NSS/MIC), pages 1–3, 2018.
  * [8] ATLAS Collaboration.  Technical Design Report for the ATLAS Inner Tracker Strip Detector.  Technical Report CERN-LHCC-2017-005. ATLAS-TDR-025, CERN, Geneva, Apr 2017\.
  * [9] ATLAS Collaboration.  The ATLAS Experiment at the CERN Large Hadron Collider.  Journal of Instrumentation, 3(08):S08003, 2008.
  * [10] C. Caso et al.  Review of Particle Physics.  Eur.Phys.J., C3, 1998.
  * [11] G. Molière.  Theorie der Streuung schneller geladener Teilchen II. Mehrfach- und Vielfachstreuung.  Zeitschrift Naturforschung Teil A, 3:78–97, 1947.
  * [12] V. Highland.  Some practical remarks on multiple scattering.  Nucl. Instr. and Meth., 129(2):497 – 499, 1975.
  * [13] M. Tanabashi, K. Hagiwara, and Hikasa.  Review of Particle Physics.  Phys.Rev.D, 98, 2018.
  * [14] G. R. Lynch and O. I. Dahl.  Approximations to multiple Coulomb scattering.  Nucl. Instr. Meth., B58:6 – 10, 1991.
  * [15] R. Diener, J. Dreyling-Eschweiler, H. Ehrlichmann, I.M. Gregor, U. Kötz, U. Krämer, N. Meyners, N. Potylitsina-Kube, A. Schütz, P. Schütze, and M. Stanitzki.  The DESY II test beam facility.  Nuclear Instruments and Methods in Physics Research Section A: Accelerators, Spectrometers, Detectors and Associated Equipment, 922:265 – 286, 2019.
  * [16] Testbeam Software Developers, 2016.  Test Beam Software Framework (TBSW), https://bitbucket.org/testbeam/tbsw.
  * [17] I. Rubinskiy.  An EUDET/AIDA Pixel Beam Telescope for Detector Development.  volume 37, pages 923 – 931, 2012.  Proceedings of the 2nd International Conference on Technology and Instrumentation in Particle Physics (TIPP 2011) .
  * [18] B. Schwenker.  Development and validation of a model for the response of the Belle II vertex detector.  PhD thesis, Universität Göttingen, 2014.  II.Physik-UniGö-Diss-2014/05, CERN-THESIS-2014-165.
  * [19] Heitler W. Bethe H. and Dirac Paul Adrien Maurice.  On the stopping of fast particles and on the creation of positive electrons.  Proc. R. Soc. Lond. A, 146:83–112, 1934.
  * [20] P Leitao, S Feger, et al.  Test bench development for the radiation hard GBTX ASIC.  Journal of Instrumentation, 10(01):C01038–C01038, Jan 2015\.
  * [21] A. Caratelli and S. Bonacini and K. Kloukinas and A. Marchioro and P. Moreira and R. De Oliveira and C. Paillard.  The GBT-SCA, a radiation tolerant ASIC for detector control and monitoring applications in HEP experiments.  Journal of Instrumentation, 10(03), mar 2015.
  * [22] Jansen, H. and Spannagel, S. and Behr, J. and others.  Performance of the EUDET-type beam telescopes.  EPJ Techniques and Instrumentation, 3(7), 2016.
  * [23] J. Baudot et al.  First test results of MIMOSA-26, a fast CMOS sensor with integrated zero suppression and digitized output.  In 2009 IEEE Nuclear Science Symposium Conference Record (NSS/MIC), pages 1169–1173, 10 2009.