herwig is hosted by Hepforge, IPPP Durham
Herwig++  2.7.0
Herwig::EventShapes Class Reference

The EventShapes class is designed so that certain event shapes, such as the thrust are only calculated once per event given the speed of the calculation. More...

#include <EventShapes.h>

Inheritance diagram for Herwig::EventShapes:

List of all members.

Public Member Functions

 EventShapes ()
 The default constructor.
void reset (const tPVector &part)
 Member to reset the particles to be considered.
void bookEEC (vector< double > &hi)
 Energy-energy correlation (EEC)
void normalizeEEC (vector< double > &hi, long evts)
 Before writing the histogram it has to be normalized according to the number of events.
double AEEC (vector< double > &hi, double &coschi)
 The asymmetry of EEC is calculated from a given $\cos\chi$ and EEC histogram, which is a vector<double> as described above.
double thrust ()
 Member functions to return thrust related shapes.
double thrustMajor ()
 The major.
double thrustMinor ()
 The minor.
double oblateness ()
 The oblateness.
Axis thrustAxis ()
 The thrust axis.
Axis majorAxis ()
 The major axis.
Axis minorAxis ()
 The minor axis.
double CParameter ()
 Linear momentum tensor related event shapes.
double DParameter ()
 The D parameter.
vector< double > linTenEigenValues ()
 The eigenvalues in descending order.
vector< Axis > linTenEigenVectors ()
 The eigenvectors in order of descending eigenvalue.
double sphericity ()
 Quadratic momentum tensor related variables.
double aplanarity ()
 The aplanarity.
double planarity ()
 The planarity.
Axis sphericityAxis ()
 The sphericity axis.
vector< double > sphericityEigenValues ()
 The sphericity eigenvalues.
vector< Axis > sphericityEigenVectors ()
 The sphericity eigenvectors.
double Mhigh2 ()
 Jet mass related event shapes.
double Mlow2 ()
 The low hemishpere mass squared divided by the visible energy squared.
double Mdiff2 ()
 The difference between the hemishpere masses squared divided by the visible energy squared.
double Bmax ()
 Jet broadening related event shapes.
double Bmin ()
 The narrow jet broadening.
double Bsum ()
 The sum of the jet broadenings.
double Bdiff ()
 The difference of the jet broadenings.
double getXi (const Lorentz5Momentum &p, const Energy &Ebeam)
 Single particle variables which do not depend on event shapes axes.
Energy getPt (const Lorentz5Momentum &p)
 Transverse momentum with respect to the beam.
double getRapidity (const Lorentz5Momentum &p)
 Rapidity with respect to the beam direction.
Energy ptInT (const Lorentz5Momentum &p)
 Single particle variables related to one of the shape axis.
Energy ptOutT (const Lorentz5Momentum &p)
 Transverse momentum with respect to the thrust axis out of the event plane.
double yT (const Lorentz5Momentum &p)
 Rapidity with respect to the thrust axis.
Energy ptInS (const Lorentz5Momentum &p)
 Transverse momentum with respect to the sphericity axis in the event plane.
Energy ptOutS (const Lorentz5Momentum &p)
 Transverse momentum with respect to the sphericity axis out of the event plane.
double yS (const Lorentz5Momentum &p)
 Rapidity with respect to the sphericity axis.

Static Public Member Functions

static void Init ()
 The standard Init function used to initialize the interfaces.

Protected Member Functions

Clone Methods.
virtual IBPtr clone () const
 Make a simple clone of this object.
virtual IBPtr fullclone () const
 Make a clone of this object, possibly modifying the cloned object to make it sane.

Private Member Functions

EventShapesoperator= (const EventShapes &)
 The assignment operator is private and must never be called.
void checkThrust ()
 Check whether the initialization of a certain class of event shapes has been calculated and if not do so.
void checkLinTen ()
 Check if the linear tensor related variables have been calculated and if not do so.
void checkSphericity ()
 Check if the quadratic tensor related variables have been calculated and if not do so.
void checkHemispheres ()
 Check if the hemisphere mass variables and jet broadenings have been calculated and if not do so.
void calcHemisphereMasses ()
 Methods that actually calculate the event shapes.
void calculateThrust ()
 Calculate the thrust and related axes.
void diagonalizeTensors (bool linear, bool cmboost)
 Diagonalize the tensors.
vector< double > eigenvalues (const double T[3][3])
 Quite general diagonalization of a symmetric Matrix T, given as an array of doubles.
Axis eigenvector (const double T[3][3], const double &lam)
 The eigenvector of.
vector< Axis > eigenvectors (const double T[3][3], const vector< double > &lam)
 The eigenvectors of.
void calcT (const vector< Momentum3 > &p, Energy2 &t, Axis &taxis)
 Member to calculate the thrust.
void calcM (const vector< Momentum3 > &p, Energy2 &m, Axis &maxis)
 Member to calculate the major.

Private Attributes

vector< Lorentz5Momentum > _pv
 Vector of particle momenta to be analysed.
bool _useCmBoost
 Whether ot not to boost to the CMS frame for the tensor diagonalizations.
vector< Axis > _thrustAxis
 Various event shape axes.
vector< Axis > _spherAxis
 The sphericity related axes.
vector< Axis > _linTenAxis
 The linearised tensor axes.
vector< double > _thrust
 Values of axis related event shapes.
vector< double > _spher
 Values of sphericity related variables.
vector< double > _linTen
 Values of linearized tensor related variables.
bool _thrustDone
 Whether or not certain event axes have been calculated.
bool _spherDone
 Whether or not the sphericity is calculated.
bool _linTenDone
 Whether or not the linearizes tensor is calculated.
bool _hemDone
 Whether or not the hemisphere masses have been calculated.
double _mPlus
 Hemisphere masses.
double _mMinus
 The low hemisphere mass.
double _bPlus
 The jet broadenings.
double _bMinus
 The narrow jet broadening.

Static Private Attributes

static NoPIOClassDescription
< EventShapes
initEventShapes
 The static object used to initialize the description of this class.

Detailed Description

The EventShapes class is designed so that certain event shapes, such as the thrust are only calculated once per event given the speed of the calculation.

See also:
The interfaces defined for EventShapes.

Definition at line 36 of file EventShapes.h.


Member Function Documentation

double Herwig::EventShapes::Bmax ( ) [inline]

Jet broadening related event shapes.

The wide jet broadening

Definition at line 259 of file EventShapes.h.

void Herwig::EventShapes::bookEEC ( vector< double > &  hi)

Energy-energy correlation (EEC)

Parameters:
hiis the histogram and has to be provided externally It is understood that the range of the histogam is -1 < cos(chi) < 1. hi.front() contains the bin [-1 < cos(chi) < -1+delta] and hi.back() the bin [1-delta < cos(chi) < 1]. delta = 2/hi.size(). We use classical indices to access the vector.

Methods that actually calculate the event shapes.

Calculate the hemisphere masses and jet broadenings

void Herwig::EventShapes::calcM ( const vector< Momentum3 > &  p,
Energy2 &  m,
Axis &  maxis 
) [private]

Member to calculate the major.

Parameters:
pThe three vectors
mThe major-squared (up to an Energy scale factor)
maxisThe major axis
void Herwig::EventShapes::calcT ( const vector< Momentum3 > &  p,
Energy2 &  t,
Axis &  taxis 
) [private]

Member to calculate the thrust.

Parameters:
pThe three vectors
tThe thrust-squared (up to an Energy scale factor)
taxisThe thrust axis
void Herwig::EventShapes::checkThrust ( ) [inline, private]

Check whether the initialization of a certain class of event shapes has been calculated and if not do so.

Check if thrust related variables have been calculated and if not do so

Definition at line 447 of file EventShapes.h.

virtual IBPtr Herwig::EventShapes::clone ( ) const [inline, protected, virtual]

Make a simple clone of this object.

Returns:
a pointer to the new object.

Implements ThePEG::InterfacedBase.

Definition at line 428 of file EventShapes.h.

double Herwig::EventShapes::CParameter ( ) [inline]

Linear momentum tensor related event shapes.

The C parameter

Definition at line 132 of file EventShapes.h.

void Herwig::EventShapes::diagonalizeTensors ( bool  linear,
bool  cmboost 
) [private]

Diagonalize the tensors.

Parameters:
linearswitch between diagonalization of linear/quadratic tensor.
cmboosttells whether to boost into cm frame of all momenta first, or not (default off, and no interface to this).
vector<double> Herwig::EventShapes::eigenvalues ( const double  T[3][3]) [private]

Quite general diagonalization of a symmetric Matrix T, given as an array of doubles.

The symmetry is not checked explicitly as this is clear in the context. It uses an explicit generic solution of the eigenvalue problem and no numerical approximation, based on Cardano's formula.

Parameters:
TMatrix to be diagonalised
Axis Herwig::EventShapes::eigenvector ( const double  T[3][3],
const double &  lam 
) [private]

The eigenvector of.

Parameters:
Tto a given eigenvalue
lam
vector<Axis> Herwig::EventShapes::eigenvectors ( const double  T[3][3],
const vector< double > &  lam 
) [private]

The eigenvectors of.

Parameters:
Tcorresponding to the eigenvectors
lam. The ordering of the vectors corresponds to the ordering of the eigenvalues.
virtual IBPtr Herwig::EventShapes::fullclone ( ) const [inline, protected, virtual]

Make a clone of this object, possibly modifying the cloned object to make it sane.

Returns:
a pointer to the new object.

Reimplemented from ThePEG::InterfacedBase.

Definition at line 433 of file EventShapes.h.

double Herwig::EventShapes::getXi ( const Lorentz5Momentum &  p,
const Energy &  Ebeam 
) [inline]

Single particle variables which do not depend on event shapes axes.

The scaled momentum $\xi=-\log\left( p/E_{\rm beam}\right)$.

Definition at line 298 of file EventShapes.h.

static void Herwig::EventShapes::Init ( ) [static]

The standard Init function used to initialize the interfaces.

Called exactly once for each class by the class description system before the main function starts or when this class is dynamically loaded.

Reimplemented from ThePEG::Interfaced.

double Herwig::EventShapes::Mhigh2 ( ) [inline]

Jet mass related event shapes.

The high hemishpere mass squared divided by the visible energy squared

Definition at line 227 of file EventShapes.h.

EventShapes& Herwig::EventShapes::operator= ( const EventShapes ) [private]

The assignment operator is private and must never be called.

In fact, it should not even be implemented.

Energy Herwig::EventShapes::ptInT ( const Lorentz5Momentum &  p) [inline]

Single particle variables related to one of the shape axis.

Transverse momentum with respect to the thrust axis in the event plane

Definition at line 326 of file EventShapes.h.

double Herwig::EventShapes::sphericity ( ) [inline]

Quadratic momentum tensor related variables.

The sphericity

Definition at line 172 of file EventShapes.h.

double Herwig::EventShapes::thrust ( ) [inline]

Member functions to return thrust related shapes.

The thrust

Definition at line 71 of file EventShapes.h.


Member Data Documentation

double Herwig::EventShapes::_bPlus [private]

The jet broadenings.

The wide jet broadening

Definition at line 661 of file EventShapes.h.

double Herwig::EventShapes::_mPlus [private]

Hemisphere masses.

The high hemisphere mass

Definition at line 646 of file EventShapes.h.

vector<double> Herwig::EventShapes::_thrust [private]

Values of axis related event shapes.

Values of thrust related variables

Definition at line 596 of file EventShapes.h.

vector<Axis> Herwig::EventShapes::_thrustAxis [private]

Various event shape axes.

The thrust related axes

Definition at line 576 of file EventShapes.h.

Whether or not certain event axes have been calculated.

Whether or not the thrust is calculated

Definition at line 616 of file EventShapes.h.

The static object used to initialize the description of this class.

Indicates that this is a concrete class with persistent data.

Definition at line 554 of file EventShapes.h.


The documentation for this class was generated from the following file: