[**:house: Home**](https://hackmd.io/s/rkkDP_l4M) | [:boy: **About**](https://hackmd.io/s/B149Z8v7b) | [**:microscope: Researches**](https://hackmd.io/s/rJPFNKlVz) | [**:rocket: Side projects**](https://hackmd.io/s/H1aS2qe4G) | [**:airplane: Life gallery**](https://hackmd.io/s/HJN4JslNM)

---

# [Algorithm] Preshower Detector Alignment

*<div style="text-align: center;" markdown="1">`minimum chi-squared` `gradient decent` `rotation matrix` `statistics`</div>*

## Introduction
[Preshower detector (ES)](http://cms.web.cern.ch/news/ecal-preshower) detector is built in the endcaps of the [CMS](http://cms.web.cern.ch/news/what-cms) detector between [tracker](http://cms.web.cern.ch/news/tracker-detector) (close to collision spot) and [Electromagnetic calorimeter (ECAL)](http://cms.web.cern.ch/news/electromagnetic-calorimeter) for increasing the precision of position reconstruction of photon and electron. There are four donut-shape planes placed in each endcap, i.e. each two planes are in each endcap side, and two of same side planes are called the *front* and *rear* plane.

[<div style="text-align: center;" markdown="1"><img src="http://hep1.phys.ntu.edu.tw/~alpha78718/cv/esidea.png" height="400"></div>](http://hep1.phys.ntu.edu.tw/~alpha78718/cv/esidea.png)

Each ES plane is made by several silicon modules with the width of 6.3 cm. The module is consist of 32 silicon strips. The silicon strips in the front planes are placed with vertical direction, while they are placed with horizontal direction in the rear planes. The direction gives the particular dimension sensitivity, the front plane is sensitive to $x$ direction, and rear plane is sensitive to $y$ direction. Basing on the $z$ coordinate of the CMS from $-$ to $+$, the four planes are denoted as $-\, Rear$, $-\,Front$, $+\,Front$ and $+\,Rear$. Thus, when the charge particle or photon passes through the silicon sensor, it remains the signal, called ***hit***, and then the position can be recognized.

Since the detector is placed manually and staying in the changeable large magnetic field space, ES planes' positions can be misaligned. This affect the position reconstruction of the particle in the LHC data. The ES alignment thus plays important role for the data quality. However, it is impossible to fix the ES planes' positions mechanically, we apply the position correction to the particle by transformation matrix. I maintain and update this framework since 2012 to 2017, faced the challenge of long-shutdown, Run II of CERN-LHC and the problem of magnetic field in the CMS.

## Techniques
The alignment algorithm for each plane are applied independently. The algorithm is the supervised learning which compares the ***measured hits*** of *ES* and ***predicted hits*** predicted from particle *tracks* in the *tracker*. As illustrated in the figure of [Introduction](#Introduction), after the collision, the particle tracks (dot lines) from the tracker (light yellow) are extrapolated to ES planes (light green) be the *predicted hits*. However, the tracker is expected to be well aligned, the *measured hits* on ES planes recored in the LHC data has to be corrected with respect to the *predicted hits* (target).

### 1. Definition
By using ***minimum $\chi^2$*** method, the position of *measured hits* are updated iteratively until the ES plane is aligned within the its position resolution. The $\chi^2$ is defined as

$$
\begin{split}
\chi^2&=\sum_{i=1}^{tracks}\epsilon_i^{T}V_i^{-1}\epsilon_i \\
&=\sum_{i=1}^{tracks}\begin{bmatrix} \epsilon_x\,\,\,\epsilon_y\end{bmatrix}_i
\begin{bmatrix} \sigma_x^2 & \sigma_{xy} \\ \sigma_{yx} & \sigma_y^2\end{bmatrix}^{-1}_i
\begin{bmatrix} \epsilon_x \\ \epsilon_y\end{bmatrix}_i\ ,
\end{split}
$$

where the ES plane is expected to only sensitive to 2-dimension, i.e. x and y axis; $\epsilon_x$ ($\epsilon_y$) is the position difference between the closest *measured hit* and *predicted hit* on x-axis (y-axis) in ES-plane coordinate, which is also called ***residual***. $V_i$ is the covariance matrix (error matrix) for tracks extrapolation and intrinsic resolution. If the measurement of x and y is uncorrelated, the errors are $\sigma_{xy}=\sigma_{yx}=0$ in anti-diagonal, the $\chi^2$ is simplified to the usual case as

$$
\chi^2=\sum_{i=1}^{tracks}\left(\frac{\epsilon_x^2}{\sigma_x^2}+\frac{\epsilon_y^2}{\sigma_y^2}\right)_i\ .
$$

The definition of residual is as

$$
\begin{split}
\epsilon &= v'_{predicted} - v'_{measured} \\
&= R_t(v_{predicted}-O_t)-R(v_{measured}-O)\ ,
\end{split}
$$

where $v'$ is the cartesian coordinate in the ES-plane coordinate, while $v$ is in the CMS coordinate. $R_t$ and $O_t$ is the updated $3\times3$ rotation matrix and ES-plane origin with the CMS coordinate in $t$ iteration, i.e. the $\epsilon$ and $\chi^2$ are changed iteratively; the original rotation matrix is identity, $R=I$, and $O$ is the default origin. The ES-plane coordinate is as the figure:

[<div style="text-align: center;" markdown="1"><img src="https://i.imgur.com/TmgYYff.png" height="150"></div>](https://i.imgur.com/TmgYYff.png)

which includes 3 rotations of the ***Euler angles***, $\alpha$, $\beta$ and $\gamma$, along the 3 transitions, $x$, $y$ and $z$, respectively. The rotation matrix is defined as

$$
R=\begin{bmatrix} \cos{\beta}\cos{\gamma}-\sin{\alpha}\sin{\beta}\sin{\gamma} & \cos{\beta}\sin{\gamma}+\sin{\alpha}\sin{\beta}\cos{\gamma} & -\cos{\alpha}\sin{\beta} \\
-\cos{\alpha}\sin{\gamma} & \cos{\alpha}\cos{\gamma} & \sin{\alpha} \\
\sin{\beta}\cos{\gamma}+\sin{\alpha}\cos{\beta}\sin{\gamma} & \sin{\beta}\sin{\gamma} - \sin{\alpha}\cos{\beta}\cos{\gamma} & \cos{\alpha}\cos{\beta}
\end{bmatrix}\ ,
$$

which rotates from $\gamma$, $\alpha$ to $\beta$.

### 2. Iteration and parameter update
According to the residual $\epsilon$ in above, the updating parameter are the rotation matrix, $R'(\alpha,\,\beta,\,\gamma)$, and the origin of ES-plane, $O'(x,\,y,\,z)$. The coordinate of ES plane can be presented with 6 dimensional vector as

$$
p = \begin{bmatrix} x\\ y \\ z \\ \alpha \\ \beta \\ \gamma \end{bmatrix}\ ,
$$

and the variance of all parameters is denoted as $\delta p=[\Delta x,\,\Delta y,\,\Delta z,\,\Delta\alpha,\,\Delta\beta,\,\Delta\gamma]^T$. The residual thus is as function of $p$, i.e. $\epsilon(p)$.

The parameter updating is using ***gradient descent method***, which is finding the best parameters to reach the minimum of $\chi^2$. The iteration for the $\chi^2$ is expected to climb down the hill until the valley. According to this concept, the step to reach the minimum can be obtained by the *chain rule of differential* as like:

$$
\begin{split}
\Delta x &= \frac{\Delta f(x)}{\Delta x}\frac{\Delta x}{\Delta f(x)}\Delta x \\
&=\dot{f}\times\ddot{f}^{-1}\\
&=gc^{-1}\ ,
\end{split}
$$

where $f(x)$ is a certain function having the minimum point as function of $x$; $g$ is the first-order-differential term $\dot{f}$ w.r.t $x$, called ***gradient (slope)*** ; and $c$ is the second-order-differential term $\ddot{f}$ w.r.t $x$, called ***second derivative (curvature)*** . Since the $f(x)$ has the minimum point, e.g. $f(x)=ax^2+bx+k$, the direction of updating $x$ (movement) is opposite from the sign of $g$, e.g. going down hill by increasing $x$, the $g$ is supposed to be negative. In the meanwhile, the $c$ term can ensure the $x$ is in the valley instead of the summit. As the two cases shown in following figures, the algorithm can ensure the $x$ update is going to reach the minimum.

[<div style="text-align: center;" markdown="1"><img src="https://i.imgur.com/cLEUV9A.png" width="500"></div>](https://i.imgur.com/cLEUV9A.png)

Thus, the $x$ updating in $t$ iteration becomes as

$$
x_t=x_{t-1}-\Delta x_t\ .
$$

For the alignment case, by replace the $f(x)$ to $\chi^2$ as function of $\epsilon_i(p)$

$$
\begin{split}
\delta p &= (\nabla_p^2\chi^2)^{-1}\,\nabla_p\chi^2 \\
&=\left[\sum_{i=1}^{tracks}J^T_iV_i^{-1}J_i\right]^{-1}\left[\sum_{i=1}^{tracks}J^T_iV_i^{-1}\epsilon_i\right] \\
&=M^{-1}P\ ,
\end{split}
$$

where $J_i$ is a Jacobian matrix of $\epsilon_i$ as
$$
J_i=\nabla_p\epsilon_i\ ,
$$

where it is $6\times2$ matrix, since $\nabla_p$ is a 6 dimensional operator as $[\frac{\partial}{\partial x},\,\frac{\partial}{\partial y},\,\frac{\partial}{\partial z},\,\frac{\partial}{\partial\alpha},\,\frac{\partial}{\partial\beta},\,\frac{\partial}{\partial\gamma} ]$. The $M$ is the curvature term as a $6\times6$ squared matrix, and $P$ is the gradient term as a $6\times1$ vector. As mentioned in above, since the direction of updating $p$ is opposite from $\delta p$, the updated coordinate of ES-plane with $t$ iterations is

$$
p=p_O-\delta p_1-\dots-\delta p_t\ .
$$

As the following figures from one of ES plane, the alignment algorithm changes each parameters until reaching the stopping point, i.e. the step of each dimension is zero.

<div style="text-align: center;" markdown="1"><img src="http://hep1.phys.ntu.edu.tw/~alpha78718/cv/esiters.png" height="400"></div>

### 3. Resolution
The alignment is achieved by the hit measurement of the ES plane, the performance is limited by the resolution of the silicon sensor. As mentioned in the [Introduction](#Introduction), according to the different structure of the silicon strips, front and rear plane have different sensitivities for the $x$ and $y$ axis, i.e. the resolutions of the measurement are different as well. Taking the example of the front plane, every group of the 32 silicon strips is well placed vertically in the 6.3 $\times$ 6.3 $\text{cm}^2$ module. The resolutions of $x$ and $y$ are calculated as

$$
\begin{split}
\sigma_x&=\frac{width}{\sqrt{12}}=\frac{6.3/32}{\sqrt{12}}&=0.055\,\text{cm}\ ,\\
\sigma_y&=\frac{length}{\sqrt{12}}=\frac{6.3}{\sqrt{12}}&=1.761\,\text{cm}\ ,
\end{split}
$$

where $\frac{1}{\sqrt{12}}$ is obtained by assuming the hit probability of silicon sensor is uniform distribution in every dimension, detail see [here](https://hackmd.io/s/HJ6D65HHG#1-Uniform-distribution). In contrast, the resolutions of the rear plane are $\sigma_x=1.761$ and $\sigma_y=0.055$, since it sensitive to $y$ direction. These uncertainties are included in the covariance matrix for the minimum $\chi^2$ method.

## Results
The validation is done by checking the residual distribution. Since the sensitivity difference in front and rear planes, the residual distributions are shown in the particular direction, i.e. the residual of $x$ direction for front plane, and the residual of $y$ direction for rear plane. AS shown in following figures, the residual distribution of the misaligned ES-plane is away from zero (first figure), while the mean of the aligned ES-plane (second figure) is zero within the resolution $\sigma_x=0.055\,\text{cm}$.

<div style="text-align: center;" markdown="1"><img src="http://hep1.phys.ntu.edu.tw/~alpha78718/cv/esbefor.png" height="250"><img src="http://hep1.phys.ntu.edu.tw/~alpha78718/cv/esafter.png" height="250"></div>

## References
- Github : <https://github.com/juifa-tsai/AlignmentTool>
- [**Sensor Alignment by Tracks**, *arXiv:physics/0306034*](https://arxiv.org/abs/physics/0306034v2)
- Conference note : <https://doi.org/10.1016/j.nima.2012.04.007> by Dr. K.Y. Kao
- Related Machine learning course : https://www.youtube.com/watch?v=yKKNr-QKz2Q

<br>

---
[:ghost: Github](https://github.com/juifa-tsai) | [:busts_in_silhouette: Linkedin ](https://www.linkedin.com/in/jui-fa-tsai-08ba0a93)