Thrust

The quantity thrust $T$ is defined by [Bra64]

$$T = \max_{\vert\mathbf{n}\vert = 1} \, \frac{\displaystyle \... ...bf{p}_i\vert} {\displaystyle \sum_i \vert\mathbf{p}_i\vert} ~,$$ (273)

and the thrust axis $\mathbf{v}_1$ is given by the $\mathbf{n}$ vector for which maximum is attained. The allowed range is $1/2 \leq T \leq 1$, with a 2-jet event corresponding to $T \approx 1$ and an isotropic event to $T \approx 1/2$.

In passing, we note that this is not the only definition found in the literature. The definitions agree for events studied in the c.m. frame and where all particles are detected. However, a definition like

$$T = 2 \, \max_{\vert\mathbf{n}\vert = 1} \, \frac{\displayst... ..._i \right\vert } {\displaystyle \sum_i \vert\mathbf{p}_i\vert}$$ (274)

(where $\theta(x)$ is the step function, $\theta(x) = 1$ if $x > 0$, else $\theta(x) = 0$) gives different results than the one above if e.g. only charged particles are detected. It would even be possible to have $T > 1$; to avoid such problems, often an extra fictitious particle is introduced to balance the total momentum [Bra79].

Eq. () may be rewritten as

$$T = \max_{\epsilon_i = \pm 1} \, \frac{\displaystyle \left\v... ...\right\vert } {\displaystyle \sum_i \vert\mathbf{p}_i\vert} ~.$$ (275)

(This may also be viewed as applying eq. () to an event with $2 n$ particles, $n$ carrying the momenta $\mathbf{p}_i$ and $n$ the momenta $- \mathbf{p}_i$, thus automatically balancing the momentum.) To find the thrust value and axis this way, $2^{n-1}$ different possibilities would have to be tested. The reduction by a factor of 2 comes from $T$ being unchanged when all $\epsilon_i \to - \epsilon_i$. Therefore this approach rapidly becomes prohibitive. Other exact methods exist, which `only' require about $4n^2$ combinations to be tried.

In the implementation in PYTHIA, a faster alternative method is used, in which the thrust axis is iterated from a starting direction $\mathbf{n}^{(0)}$ according to

$$\mathbf{n}^{(j+1)} = \frac{\displaystyle \sum_i \epsilon( \m... ...hbf{n}^{(j)} \cdot \mathbf{p}_i) \, \mathbf{p}_i \right\vert }$$ (276)

(where $\epsilon(x) = 1$ for $x > 0$ and $\epsilon(x) = -1$ for $x < 0$). It is easy to show that the related thrust value will never decrease, $T^{(j+1)} \geq T^{(j)}$. In fact, the method normally converges in 2-4 iterations. Unfortunately, this convergence need not be towards the correct thrust axis but is occasionally only towards a local maximum of the thrust function [Bra79]. We know of no foolproof way around this complication, but the danger of an error may be lowered if several different starting axes $\mathbf{n}^{(0)}$ are tried and found to agree. These $\mathbf{n}^{(0)}$ are suitably constructed from the $n'$ (by default 4) particles with the largest momenta in the event, and the $2^{n' -1}$ starting directions $\sum_i \epsilon_i \, \mathbf{p}_i$ constructed from these are tried in falling order of the corresponding absolute momentum values. When a predetermined number of the starting axes have given convergence towards the same (best) thrust axis this one is accepted.

In the plane perpendicular to the thrust axis, a major [MAR79] axis and value may be defined in just the same fashion as thrust, i.e.

$$M_a = \max_{\vert\mathbf{n}\vert = 1, \, \mathbf{n} \cdot \m... ...f{p}_i\vert} {\displaystyle \sum_i \vert\mathbf{p}_i\vert } ~.$$ (277)

In a plane more efficient methods can be used to find an axis than in three dimensions [Wu79], but for simplicity we use the same method as above. Finally, a third axis, the minor axis, is defined perpendicular to the thrust and major ones, and a minor value $M_i$ is calculated just as thrust and major. The difference between major and minor is called oblateness, $O = M_a -M_i$. The upper limit on oblateness depends on the thrust value in a not so simple way. In general $O \approx 0$ corresponds to an event symmetrical around the thrust axis and high $O$ to a planar event.

As in the case of sphericity, a generalization to arbitrary momentum dependence may easily be obtained, here by replacing the $\mathbf{p}_i$ in the formulae above by $\vert\mathbf{p}_i\vert^{r-1} \, \mathbf{p}_i$. This possibility is included, although so far it has not found any experimental use.