6.7.6. EventShape

Introduction

This page will illustrate how to use the eventShape module, a tool that allows to calculate quantities sensitive to the global shape of the event. The goal of these variables is to deduce the original partonic state from the geometrical correlations among the final state particles. This topic is related to the continuum suppression, that is nothing but an application of the event shape concepts specific to the problem of separating \(e^+e^- \to B\bar{B}\) from \(e^+e^- \to q\bar{q}\) events.

Event shapes are discussed in several papers. Two comprehensive introductions can be found in:

The Pythia 6 manual

G. Fox, S. Wolfram, Event shapes in e+e- annihilation

For an introduction on both the theory of the event shapes variables and their use for continuum suppression and event characterization, you can also check the slides of the starter kit workshop .

What is the difference between event shape and continuum suppression?

This section is inspired by a B2question . The event shape variables are closely related to the continuum suppression one, but return completely different values and have also a different meaning. The Belle II continuum suppression (in short CS) uses the classic event shape variables, but calculates them separately on a B candidate and its rest of the event (the story is actually much longer: please see the Continuum suppression section). This makes the CS variables highly optimized for the continuum suppression when a B meson is reconstructed in an exclusive channel, but completely useless to any other application.

Unlike the CS variables, the event shape variables are calculated starting from lists of particles (usually a list of all the good tracks and a list of the good photons) intended to describe the whole event. While the CS variables are candidate-based (in a given event each B candidate will have different daughters and a different Rest Of Event, and therefore different CS variables), the event shape variables are event-based.

In summary, the event shape variables are suitable for:
  • studying the geometrical shape for the events

  • separating qq, BB, tau tau and other components without reconstruction any B or tau candidate

  • providing continuum suppression to all the non-B analyses

How to use eventShape

To access the event shape variables you have to append the EventShapeCalculator module to your path, which will calculate all the quantities you need and store them in the datastore making them available to the VariableManager. You can do this simply adding the modularAnalysis.buildEventShape() function to your analysis path:

modularAnalysis.buildEventShape(inputListNames=[], default_cleanup=True, allMoments=False, cleoCones=True, collisionAxis=True, foxWolfram=True, harmonicMoments=True, jets=True, sphericity=True, thrust=True, checkForDuplicates=False, path=None)[source]

Calculates the event-level shape quantities (thrust, sphericity, Fox-Wolfram moments…) using the particles in the lists provided by the user. If no particle list is provided, the function will internally create a list of good tracks and a list of good photons with (optionally) minimal quality cuts.

The results of the calculation are then stored into the EventShapeContainer dataobject, and are accessible using the variables of the EventShape group.

The user can switch the calculation of certain quantities on or off to save computing time. By default the calculation of the high-order moments (5-8) is turned off. Switching off an option will make the corresponding variables not available.

Warning

The user can provide as many particle lists as needed, using also combined particles, but the function will always assume that the lists are independent. If the lists provided by the user contain several times the same track (either with different mass hypothesis, or once as an independent particle and once as daughter of a combined particle) the results won’t be reliable. A basic check for duplicates is available setting the checkForDuplicate flags, but is usually quite time consuming.

Parameters
  • inputListNames – List of ParticleLists used to calculate the event shape variables. If the list is empty the default particleLists pi+:evtshape and gamma:evtshape are filled.

  • default_cleanup – If True, applies standard cuts on pt and cosTheta when defining the internal lists. This option is ignored if the particleLists are provided by the user.

  • path – Path to append the eventShape modules to.

  • thrust – Enables the calculation of thrust-related quantities (CLEO cones, Harmonic moments, jets).

  • collisionAxis – Enables the calculation of the quantities related to the collision axis .

  • foxWolfram – Enables the calculation of the Fox-Wolfram moments.

  • harmonicMoments – Enables the calculation of the Harmonic moments with respect to both the thrust axis and, if collisionAxis = True, the collision axis.

  • allMoments – If True, calculates also the FW and harmonic moments from order 5 to 8 instead of the low-order ones only.

  • cleoCones – Enables the calculation of the CLEO cones with respect to both the thrust axis and, if collisionAxis = True, the collision axis.

  • jets – Enables the calculation of the hemisphere momenta and masses. Requires thrust = True.

  • sphericity – Enables the calculation of the sphericity-related quantities.

  • checkForDuplicates – Perform a check for duplicate particles before adding them. This option is quite time consuming, instead of using it consider sanitizing the lists you are passing to the function.

This is a very minimal implementation of the event shape in a steering file:

import basf2 as b2
import modularAnalysis as ma

# create path
my_path = b2.create_path()

# add the event shape module with the default settings
ma.buildEventShape(path=my_path)

# Now the event shape variables are accessible
# to the VariableManager

For a complete steering file, see the example B2A704-EventShape.py .

Event shape variables

Once modularAnalysis.buildEventShape() has been added to the analysis path, the event shape variables will be accessible via the VariableManager as any other event-level variable. Here you can find some deeper theoretical explanation about their meaning, and how they are constructed.

First introduced by J.D. Bjorken and S.J. Brodsky in Phys. Rev. D1 (1970) 1416, the sphericity is a linear combination of the second and third eigenvalue of the sphericity tensor. This tensor \(S\) is constructed from the 3-vectors of all the particles (both charged and neutral) in the event:

\(S^{\alpha, \beta} = \frac{\sum_i p_i^\alpha p_i^\beta}{\sum_i p_i^2}\),

where \(i\) runs over all the reconstructed particles in the event and \(\alpha, \beta\) are the cartesian components of the momenta. \(S\) is not Lorentz-invariant, and it’s meaningful only if calculated in the center of mass frame. The eigenvectors of \(S\) define the ellipsoid that best matches with the particle distribution in the event. The eigenvalues \(\lambda_1, \lambda_2, \lambda_3\), ordered from the largest to the smallest, are proportional to the lengths of the axes of this ellipsoid. Since \(S\) is unitary, only two independent quantities can be constructed out of its eigenvalues. The VariableManager provides the sphericity (\(3(\lambda_2+\lambda_3)/2\)) and the aplanarity (\(3 lambda_3/2\)).

Thrust and thrust axis

The thrust axis is the axis that maximizes the sum the projections of all the particles’ momenta onto it. It is the best proxy for the main symmetry axis of the event and tends to be aligned to the first eigenvector of the sphericity matrix. The Thrust is the normalized sum the projections of all the particles’ momenta onto the thrust axis. It is close to 1 for a jet-like event, and close to 1/3 for a spherical one The calculation of the thrust axis allows the event to be divided into two hemispheres; total energy, momentum, and mass can be calculated for each of them individually. The VariableManager provides the thrust value (thrust), the three cartesian components of the thrust axis (thrustAxisX, thrustAxisY, thrustAxisZ), the polar direction of the thrust axis (thrustAxisCosTheta) and the metavariable useThrustFrame that allows to calculate any variable in a rotated reference frame where the z axis coincides with the thrust axis. Finally, you can access the invariant mass (forwardHemisphereMass), the total energy (forwardHemisphereEnergy), the total momentum (forwardHemisphereMomentum) and the components of the total momentum (forwardHemisphereX, forwardHemisphereY, forwardHemisphereZ) of the particles belonging to the forward emisphere. Similar variables are also defined for the backward hemisphere.

Fox-Wolfram moments

The Fox-Wolfram moments are among the most powerful quantities to define the shape of an event. They are based on the idea of expanding the geometrical distribution of the particles in the event in a suitable basis of special function. Neglecting the math, described in the references given above, the result is a rather simple definition:

\(H_l = \sum_{i,j}p_i p_j P_l(cos\alpha_{i,j})\),

where \(i,j\) runs over all the particles in the event, \(P_l\) is the Legendre polynomial of order \(l\), and \(\alpha_{i,j}\) is the angle between the momentum of the i-th and j-th particle. Instead of the raw moments, the event shape is always studied using the normalized moments \(R_i = H_i/H_0\), defined between 0 and 1. The VariableManager provides access to any Fox-Wolfram moment up to order 8 with the metavariables foxWolframH (raw moments) and foxWolframR (normalized moments). For convenience the normalized moments from 1 to 4 are also accessible using the variables foxWolframR1 to foxWolframR4.

Harmonic moments

The harmonic moments are the coefficients of the expansion of the event in spherical harmonics around an arbitrary axis. Unlike the Fox-Wolfram moments, they are not rotationally invariant. Currently the event shape tools allow you to calculate them respect to either the thrust axis, or the collision axis or both. They are defined as:

\(B_l = \sum_i\frac{p_i}{\sqrt{s}}P_l(cos\alpha_i)\),

where \(P_l\) is the Legendre polynomial of order \(l\), and \(\alpha_i\) is the angle between the momentum of the i-th particle and the chosen expansion axis. The VariableManager provides access to any harmonic moment up to order 8, calculated either using the collision axis or the thrust axis, via the metavariable harmonicMoment. For convenience the harmonic moments from 0 to 4 calculated using the thrust axis as polar axis are also accessible using the variables harmonicMomentThrust0 to harmonicMomentThrust4.

Cleo cones

The CLEO cones are a set of concentric cones with different opening angle having the thrust axis as symmetry axis. In our current implementation we use 9 cones, of opening angle between 10 and 90 degrees, 10 degrees apart one from the other. The measurement of the energy flux in each cone provides a set of variables that have been proven to be very powerful for continuum suppression. Jet like events will be more likely to have all the energy flow concentrated in the narrower cones, while spherical events will have a more evenly-distributed flow. The VariableManager provides access to any CLEO cone value up to order 8, calculated either using the collision axis or the thrust axis, via the metavariable cleoCone. For convenience the CLEO cones variables calculated using using the thrust axis are also accessible using the variables cleoConeThrust0 to cleoConeThrust8.

Complete list of the eventShape variables:

aplanarity

[Eventbased] Event aplanarity, defined as the 3/2 of the third sphericity eigenvalue.

backwardHemisphereEnergy

[Eventbased] Total energy the particles flying in the direction opposite to the thrust axis.

backwardHemisphereMass

[Eventbased] Invariant mass of the particles flying in the direction opposite to the thrust axis.

backwardHemisphereMomentum

[Eventbased] Total momentum the particles flying in the direction opposite to the thrust axis.

backwardHemisphereX

[Eventbased] X component of the total momentum of the particles flying in the direciton opposite to the thrust axis

backwardHemisphereY

[Eventbased] Y component of the total momentum of the particles flying in the direction opposite to the thrust axis

backwardHemisphereZ

[Eventbased] Z component of the total momentum of the particles flying in the direction opposite to the thrust axis

cleoCone(i, axisName)

[Eventbased] i-th order cleoCone, calculated respect to the axis axisName. The order can go up to 9th., the axisName can be either ‘thrust’ or ‘collision’

cleoConeThrust0

[Eventbased] 0th Cleo cone calculated respect to the thrust axis.

cleoConeThrust1

[Eventbased] 1st Cleo cone calculated respect to the thrust axis.

cleoConeThrust2

[Eventbased] 2nd Cleo cone calculated respect to the thrust axis.

cleoConeThrust3

[Eventbased] 3rd Cleo cone calculated respect to the thrust axis.

cleoConeThrust4

[Eventbased] 4th Cleo cone calculated respect to the thrust axis.

cleoConeThrust5

[Eventbased] 5th Cleo cone calculated respect to the thrust axis.

cleoConeThrust6

[Eventbased] 6th Cleo cone calculated respect to the thrust axis.

cleoConeThrust7

[Eventbased] 7th Cleo cone calculated respect to the thrust axis.

cleoConeThrust8

[Eventbased] 8th Cleo cone calculated respect to the thrust axis.

forwardHemisphereEnergy

[Eventbased] Total energy the particles flying in the same direction of the thrust axis.

forwardHemisphereMass

[Eventbased] Invariant mass of the particles flying in the same direction of the thrust axis.

forwardHemisphereMomentum

[Eventbased] Total momentum the particles flying in the same direction of the thrust axis.

forwardHemisphereX

[Eventbased] X component of the total momentum of the particles flying in the same direction of the thrust axis

forwardHemisphereY

[Eventbased] Y component of the total momentum of the particles flying in the same direction of the thrust axis

forwardHemisphereZ

[Eventbased] Z component of the total momentum of the particles flying in the same direction of the thrust axis

foxWolframH(i)

[Eventbased] i-th order Fox Wolfram moment. The order can go up to 8th.

foxWolframR(i)

[Eventbased] ratio of the i-th to the 0-th order Fox Wolfram moments. The order can go up to 8th.

foxWolframR1

[Eventbased] ratio of the 1-st to the 0-th order Fox Wolfram moments. This is just an alias of foxWolframR(1) defined for the user’s covenience.

foxWolframR2

[Eventbased] ratio of the 2-nd to the 0-th order Fox Wolfram moments. This is just an alias of foxWolframR(2) defined for the user’s covenience.

foxWolframR3

[Eventbased] ratio of the 3-rd to the 0-th order Fox Wolfram moments. This is just an alias of foxWolframR(3) defined for the user’s covenience.

foxWolframR4

[Eventbased] ratio of the 4-th to the 0-th order Fox Wolfram moments. This is just an alias of foxWolframR(4) defined for the user’s covenience.

harmonicMoment(i, axisName)

[Eventbased] i-th order harmonic moment, calculated respect to the axis axisName. The order can go up to 8th., the axisName can be either ‘thrust’ or ‘collision’

harmonicMomentThrust0

[Eventbased] Harmonic moment of the 0th order calculated respect to the thrust axis.

harmonicMomentThrust1

[Eventbased] Harmonic moment of the 1st order calculated respect to the thrust axis.

harmonicMomentThrust2

[Eventbased] Harmonic moment of the 2nd order calculated respect to the thrust axis.

harmonicMomentThrust3

[Eventbased] Harmonic moment of the 3rd order calculated respect to the thrust axis.

harmonicMomentThrust4

[Eventbased] Harmonic moment of the 4th order calculated respect to the thrust axis.

sphericity

[Eventbased] Event sphericity, defined as the linear combination of the sphericity eigenvlaues S = (3/2)(lambda2+lambda3)

thrust

[Eventbased] Event thrust.

thrustAxisCosTheta

[Eventbased] Cosine of the polar angle component of the thrust axis.

thrustAxisX

[Eventbased] X component of the thrust axis.

thrustAxisY

[Eventbased] Y component of the thrust axis.

thrustAxisZ

[Eventbased] Z component of the thrust axis.

useThrustFrame(variable)

Evaluates a variable value in the thrust reference frame.