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

Class to analyse event shapes. More...

#include <EventShapes2.h>

Inheritance diagram for Analysis2::EventShapes2:

List of all members.

Public Member Functions

virtual void analyze (const tPVector &particles)
 Analyze the given vector of particles.
Standard constructors and destructors.
 EventShapes2 ()
 The default constructor.
virtual ~EventShapes2 ()
 The destructor.
Functions used by the persistent I/O system.
void persistentOutput (PersistentOStream &os) const
 Function used to write out object persistently.
void persistentInput (PersistentIStream &is, int version)
 Function used to read in object persistently.

Static Public Member Functions

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

Protected Member Functions

void reset (const vector< Lorentz5Momentum > &part)
 Member to reset the particles to be considered.
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.
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 getX (const Lorentz5Momentum &p, const Energy &Ebeam)
 Single particle variables which do not depend on event shapes axes.
double getXi (const Lorentz5Momentum &p, const Energy &Ebeam)
 The scaled momentum $\xi=-\log\left( p/E_{\rm beam}\right)$.
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.

Private Member Functions

EventShapes2operator= (const EventShapes2 &)
 The assignment operator is private and must never be called.
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.
string _omthr
 Output option strings for $1-T$ distribution.
string _maj
 Output option strings for the major distribution.
string _min
 Output option strings for the minor distribution.
string _obl
 Output option strings for the oblateness distribution.
string _sph
 Output option strings for the sphericity distribution.
string _apl
 Output option strings for the aplanarity distribution.
string _pla
 Output option strings for the planarity distribution.
string _c
 Output option strings for the C distribution.
string _d
 Output option strings for the D distribution.
string _mhi
 Output option strings for the $M_{\rm high}$ distribution.
string _mlo
 Output option strings for the $M_{\rm low}$ distribution.
string _mdiff
 Output option strings for the $M_{\rm high}-M_{\rm low}$ distribution.
string _bmax
 Output option strings for the $B_{\rm max}$ distribution.
string _bmin
 Output option strings for the $B_{\rm min}$ distribution.
string _bsum
 Output option strings for the $B_{\rm max}+B_{\rm min}$ distribution.
string _bdiff
 Output option strings for the $B_{\rm max}-B_{\rm min}$ distribution.
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 ClassDescription
< EventShapes2
initEventShapes2
 The static object used to initialize the description of this class.

Standard Interfaced functions.

virtual void doinit () throw (InitException)
 Initialize this object after the setup phase before saving an EventGenerator to disk.
virtual void dofinish ()
 Finalize this object.
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.

Detailed Description

Class to analyse event shapes.

See also:
The interfaces defined for EventShapes2.

Definition at line 29 of file EventShapes2.h.


Member Function Documentation

virtual void Analysis2::EventShapes2::analyze ( const tPVector particles) [virtual]

Analyze the given vector of particles.

The default version calls analyze(tPPtr) for each of the particles.

Parameters:
particlesthe vector of pointers to particles to be analyzed

Reimplemented from Analysis2::Analysis2Base.

double Analysis2::EventShapes2::Bmax ( ) [inline, protected]

Jet broadening related event shapes.

The wide jet broadening

Methods that actually calculate the event shapes.

Calculate the hemisphere masses and jet broadenings

void Analysis2::EventShapes2::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 Analysis2::EventShapes2::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 Analysis2::EventShapes2::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

virtual IBPtr Analysis2::EventShapes2::clone ( ) const [inline, protected, virtual]

Make a simple clone of this object.

Returns:
a pointer to the new object.

Reimplemented from ThePEG::AnalysisHandler.

double Analysis2::EventShapes2::CParameter ( ) [inline, protected]

Linear momentum tensor related event shapes.

The C parameter

void Analysis2::EventShapes2::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).
virtual void Analysis2::EventShapes2::dofinish ( ) [protected, virtual]

Finalize this object.

Called in the run phase just after a run has ended. Used eg. to write out statistics.

Reimplemented from Analysis2::Analysis2Base.

virtual void Analysis2::EventShapes2::doinit ( ) throw (InitException) [protected, virtual]

Initialize this object after the setup phase before saving an EventGenerator to disk.

Exceptions:
InitExceptionif object could not be initialized properly.

Reimplemented from ThePEG::InterfacedBase.

vector<double> Analysis2::EventShapes2::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 Analysis2::EventShapes2::eigenvector ( const double  T[3][3],
const double &  lam 
) [private]

The eigenvector of.

Parameters:
Tto a given eigenvalue
lam
vector<Axis> Analysis2::EventShapes2::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 Analysis2::EventShapes2::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::AnalysisHandler.

double Analysis2::EventShapes2::getX ( const Lorentz5Momentum &  p,
const Energy &  Ebeam 
) [inline, protected]

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

Ratio of momentum to beam momentum

static void Analysis2::EventShapes2::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 Analysis2::Analysis2Base.

double Analysis2::EventShapes2::Mhigh2 ( ) [inline, protected]

Jet mass related event shapes.

The high hemishpere mass squared divided by the visible energy squared

EventShapes2& Analysis2::EventShapes2::operator= ( const EventShapes2 ) [private]

The assignment operator is private and must never be called.

In fact, it should not even be implemented.

Function used to read in object persistently.

Parameters:
isthe persistent input stream read from.
versionthe version number of the object when written.

Reimplemented from Analysis2::Analysis2Base.

Function used to write out object persistently.

Parameters:
osthe persistent output stream written to.

Reimplemented from Analysis2::Analysis2Base.

Energy Analysis2::EventShapes2::ptInT ( const Lorentz5Momentum &  p) [inline, protected]

Single particle variables related to one of the shape axis.

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

double Analysis2::EventShapes2::sphericity ( ) [inline, protected]

Quadratic momentum tensor related variables.

The sphericity

double Analysis2::EventShapes2::thrust ( ) [inline, protected]

Member functions to return thrust related shapes.

The thrust


Member Data Documentation

The jet broadenings.

The wide jet broadening

Definition at line 523 of file EventShapes2.h.

Hemisphere masses.

The high hemisphere mass

Definition at line 508 of file EventShapes2.h.

vector<double> Analysis2::EventShapes2::_thrust [private]

Values of axis related event shapes.

Values of thrust related variables

Definition at line 458 of file EventShapes2.h.

vector<Axis> Analysis2::EventShapes2::_thrustAxis [private]

Various event shape axes.

The thrust related axes

Definition at line 438 of file EventShapes2.h.

Whether or not certain event axes have been calculated.

Whether or not the thrust is calculated

Definition at line 478 of file EventShapes2.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 617 of file EventShapes2.h.


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