KaliVeda
1.13/01
Heavy-Ion Analysis Toolkit
|
Base class for relativistic kinematics of massive particles.
Implements all kinematical manipulations necessary for obtaining angles, kinetic energies, transverse energies etc. etc., defining and changing relativistic reference frames.
Unless otherwise stated,
This class derives from TLorentzVector - a particle is basically a Lorentz 4-vector with added attributes - and therefore any methods which are not documented here will be found in that class
This also means that all particles and kinematics are relativistic in KaliVeda.
Methods defined/redefined in this class:
Particle properties can be defined either using one of the constructors :
or with the usual 'Set' methods:
or, for the following method:
Set momentum 3-vector, leaving rest mass unchanged, with :
opt = "cart"
or "cartesian"
: using cartesian components (px,py,pz)
opt = "spher"
or "spherical"
: set momentum 3-vector using spherical coordinatespx
= magnitude, py
= theta [degrees], pz
= phi [degrees].Other methods:
GetEtran() and GetTransverseEnergy() return the non-relativistic transverse energy defined as the kinetic energy multiplied by the squared sinus of the polar angle: i.e. \(E_{tran} = E\sin^2(\theta)\) where \(\theta\) is the polar angle of the particle measured with respect to the beam (z-) axis.
GetREtran() and GetRTransverseEnergy() return the relativistic transverse energy which is frame invariant.
Particle kinematics in different frames can be accessed with method InFrame(). This method does not modify the particle's kinematics or default reference frame (see 2 below). If you want to access lots of information from this frame, it is probably more efficient to define and store it in the particle's list of reference frames (see 3 below).
Example: inspect kinematics of particle (accessed through pointer KVParticle* p
)
In this case we need to pass two arguments to the KVFrameTransform constructor: before C++11, this required an explicit call to the constructor (top line); since C++11, uniform initialization allows to just give the expected arguments, using {}
to enclose the arguments for any aggregate (class) types.
Particle kinematics can be modified using method ChangeFrame():
Example: change kinematics of particle (accessed through pointer KVParticle* p
)
Rather than changing the reference frame of the particle, you can define and use several different reference frames while keeping the original kinematics as the default. Each frame can be used independently, and new frames can be defined based on any of the existing frames:
Example: (for a particle accessed through pointer KVParticle* p
{.cpp}):
Note that the same frame can be defined directly from the original particle by using a combined boost-then-rotation transform:
Note that the frame "rotated_moving_frame"
is directly accessible even if it is defined in two steps as a rotation of the "moving_frame"
.
The total number of defined frames (including the default frame) is given by GetNumberOfDefinedFrames(). Calling Print() will show all reference frames defined for the particle:
Indentation indicates the relationships between frames: "rotated_moving_frame"
is a child frame of "moving_frame"
. The first line is the default kinematics, which by default has no name (a name can be set using SetFrameName()).
The current default kinematics can always be accessed using GetCurrentDefaultKinematics(): this method returns the default kinematical frame for the particle when accessed from any kinematical frame:
Let us consider a particle for which the different reference frames in the previous paragraph have been defined. For an example, imagine that the default kinematics are that of particle with rest mass 939 MeV moving with a momentum of 250 MeV/c at an angle of \(45^o\) to the +ve z-direction:
Now if we want to change the default kinematical frame for this particle:
The relationships between frames are preserved, i.e. if we present the frames as graphs:
with "lab" as default frame:
with "rotated_moving_frame" as default frame:
If you call SetFrame() several times with the same frame name [note that frame names are case insensitive], the existing reference frame will be updated to use the new transformation, which will be applied to the kinematics of the particle in the 'parent' frame used to define the frame. Any frames which were defined based on the frame will be updated too.
Example: for the previous particle & frame definitions, after resetting the default kinematics:
Note that in this case, the "rotated_moving_frame" is updated automatically to take account of the new velocity of "moving_frame".
However, if you change the kinematics of the particle in its original (default) frame, you have to update the kinematics in all defined frames by hand, it does not occur automatically:
AddGroup() and BelongsToGroup() methods allow to sort particles in an event into subsets based on various criteria such as "emitted by QP", "backwards emitted", or "include in calorimetry". The group attribution can be done 'by hand' or using a KVParticleCondition object to define a selection. The number of such groups for a given particle is unlimited.
Group names are case insensitive. Example:
Groups are valid for all kinematical frames, and can be defined using a pointer to the nucleus in any of its defined kinematical frames:
Definition at line 398 of file KVParticle.h.
Classes | |
class | FrameList |
Public Types | |
enum | { kIsOK = BIT(14) , kIsOKSet = BIT(15) , kIsDetected = BIT(16) } |
Public Types inherited from TLorentzVector | |
typedef Double_t | Scalar |
Public Types inherited from TObject | |
enum | EDeprecatedStatusBits |
enum | EStatusBits |
Static Public Member Functions | |
static Double_t | C () |
Static Public Member Functions inherited from TObject | |
static Longptr_t | GetDtorOnly () |
static Bool_t | GetObjectStat () |
static void | SetDtorOnly (void *obj) |
static void | SetObjectStat (Bool_t stat) |
Protected Member Functions | |
void | AddGroups (KVUniqueNameList *un) |
list of groups added to the current one More... | |
void | SetGroups (KVUniqueNameList *un) |
Define for the particle a new list of groups. More... | |
Protected Member Functions inherited from TObject | |
virtual void | DoError (int level, const char *location, const char *fmt, va_list va) const |
void | MakeZombie () |
Protected Attributes | |
TVector3 * | fE0 |
the momentum of the particle before it is slowed/stopped by an absorber More... | |
KVNameValueList | fParameters |
a general-purpose list of parameters associated with this particle More... | |
Protected Attributes inherited from TObject | |
kOnlyPrepStep | |
Private Attributes | |
FrameList | fBoosted {this} |
Bool_t | fFrameCopyOnly {kFALSE} |
if != nullptr, this object is a representation of the particle in a kinematical frame More... | |
TString | fFrameName |
non-persistent frame name field, sets when calling SetFrame method More... | |
KVUniqueNameList | fGroups |
list of momenta of the particle in different Lorentz-boosted frames More... | |
TString | fName |
non-persistent name field - Is useful More... | |
KVParticle * | fOriginal {nullptr} |
parent kinematical frame More... | |
KVParticle * | fParentFrame {nullptr} |
Static Private Attributes | |
static Double_t | kSpeedOfLight = TMath::C() * 1.e-07 |
speed of light in cm/ns More... | |
Friends | |
class | KVKinematicalFrame |
Additional Inherited Members | |
Public Attributes inherited from TLorentzVector | |
kNUM_COORDINATES | |
kSIZE | |
kT | |
kX | |
kY | |
kZ | |
Public Attributes inherited from TObject | |
kBitMask | |
kCanDelete | |
kCannotPick | |
kHasUUID | |
kInconsistent | |
kInvalidObject | |
kIsOnHeap | |
kIsReferenced | |
kMustCleanup | |
kNoContextMenu | |
kNotDeleted | |
kObjInCanvas | |
kOverwrite | |
kSingleKey | |
kWriteDelete | |
kZombie | |
#include <KVParticle.h>
anonymous enum |
Enumerator | |
---|---|
kIsOK | |
kIsOKSet | |
kIsDetected |
Definition at line 551 of file KVParticle.h.
KVParticle::KVParticle | ( | ) |
Definition at line 40 of file KVParticle.cpp.
create particle with given mass and momentum vector
Definition at line 79 of file KVParticle.cpp.
create particle with given mass and momentum vector
Definition at line 92 of file KVParticle.cpp.
KVParticle::KVParticle | ( | const KVParticle & | obj | ) |
copy ctor
Definition at line 63 of file KVParticle.cpp.
|
virtual |
|
private |
used to inhibit full copying of particles in different kinematical frames
Definition at line 487 of file KVParticle.cpp.
Associate this particle with the given named group. Optional argument "from" allows to put a condition on the already stored group list, is set to "" by default
Definition at line 436 of file KVParticle.cpp.
|
protected |
list of groups added to the current one
Definition at line 469 of file KVParticle.cpp.
Check if particle belong to a given group return kTRUE if groupname="". return kFALSE if no group has be defined
Definition at line 509 of file KVParticle.cpp.
|
static |
Static function. Returns speed of light in cm/ns units.
Definition at line 117 of file KVParticle.cpp.
Make existing reference frame 'newdef' the new default frame for particle kinematics. The current default frame will then be accessible from the list of frames using its name (if set with SetFrameName).
You can change/set the name of the previous default frame with 'defname'. If no name was set and none given, it will be renamed "default" by default.
Definition at line 622 of file KVParticle.cpp.
void KVParticle::ChangeFrame | ( | const KVFrameTransform & | t, |
const KVString & | name = "" |
||
) |
Permanently modify kinematics of particle according to the given transformation. You can optionally set the name of this new default kinematics. NB the current kinematics will be lost. If you want to keep it after changing the default kinematics, define the new frame with SetFrame and then use ChangeDefaultFrame.
Definition at line 601 of file KVParticle.cpp.
Reset particle properties i.e. before creating/reading a new event.
Reimplemented from TObject.
Reimplemented in KVINDRAReconNuc, KVFAZIAReconNuc, KVReconstructedNucleus, and KVNucleus.
Definition at line 295 of file KVParticle.cpp.
Copy this to obj Particle kinematics are copied using operator=(const KVParticle&) List of particle's groups is copied The particle's name is copied The list of parameters associated with the particle is copied
Reimplemented from TObject.
Reimplemented in KVSimNucleus, KVNucleus, KVReconstructedNucleus, KVINDRAReconNuc, and KVFAZIAReconNuc.
Definition at line 269 of file KVParticle.cpp.
|
private |
[in] | frame | name of a kinematical frame |
Definition at line 937 of file KVParticle.cpp.
|
private |
PRIVATE method for internal use only Returns pointer to parent frame of 'f'
Definition at line 965 of file KVParticle.cpp.
|
inline |
Definition at line 686 of file KVParticle.h.
|
inline |
Definition at line 772 of file KVParticle.h.
|
inline |
Definition at line 654 of file KVParticle.h.
|
inline |
Definition at line 640 of file KVParticle.h.
|
inline |
Definition at line 623 of file KVParticle.h.
|
inline |
Definition at line 631 of file KVParticle.h.
KVParticle const * KVParticle::GetFrame | ( | const Char_t * | frame, |
Bool_t | warn_and_return_null_if_unknown = kTRUE |
||
) | const |
Return the momentum of the particle in the Lorentz-boosted frame corresponding to the name "frame" given as argument (see SetFrame() for definition of different frames). If the default frame name has been set (see KVEvent::SetFrameName) and 'frame' is the name of this default frame (KVParticle::fFrameName), we return the address of the particle itself.
This frame may have been defined by a direct transformation of the original kinematics of the particle (using SetFrame(newframe,...)) or by a transformation of the kinematics in another user-defined frame (using SetFrame(newframe,oldframe,...)).
Note that frames are not "dynamic": if any changes are made to the original particle's kinematics after definition of a frame, if you want these changes to affect also the other frames you need to update them by hand by calling KVParticle::UpdateAllFrames().
Frame names are case insensitive: "CM" or "cm" or "Cm" are all good...
By default, if no frame with the given name is found, we return nullptr and print a warning. If warn_and_return_null_if_unknown=kFALSE
, we return the address of the particle itself, i.e. the original/default kinematics. [Note that this is an inversion of the previous default behaviour]
Note that the properties of the particle returned by this method can not be modified: this is deliberate, as any modifications e.g. to kinematics will have no effect in any other frames.
The returned pointer corresponds to a "pseudoparticle" in the desired frame, therefore you can use any KVParticle method in order to access the kinematics of the particle in the boosted frame, e.g.
Definition at line 855 of file KVParticle.cpp.
Definition at line 801 of file KVParticle.h.
|
inline |
Definition at line 758 of file KVParticle.h.
|
inline |
Definition at line 616 of file KVParticle.h.
|
inline |
Definition at line 711 of file KVParticle.h.
|
inline |
Definition at line 636 of file KVParticle.h.
|
inline |
Definition at line 572 of file KVParticle.h.
|
inline |
Definition at line 606 of file KVParticle.h.
|
virtual |
|
inline |
Definition at line 794 of file KVParticle.h.
|
inline |
Definition at line 753 of file KVParticle.h.
|
inlineprivate |
Definition at line 536 of file KVParticle.h.
|
inlineprivate |
Definition at line 532 of file KVParticle.h.
|
inline |
Definition at line 816 of file KVParticle.h.
|
inlineprivate |
Definition at line 513 of file KVParticle.h.
|
inline |
Definition at line 690 of file KVParticle.h.
|
inline |
Definition at line 728 of file KVParticle.h.
|
inline |
Definition at line 650 of file KVParticle.h.
|
inline |
Definition at line 645 of file KVParticle.h.
should be in fm
Definition at line 666 of file KVParticle.h.
|
inline |
Definition at line 682 of file KVParticle.h.
|
inlineprivate |
Definition at line 521 of file KVParticle.h.
|
inline |
Definition at line 627 of file KVParticle.h.
|
inline |
Definition at line 610 of file KVParticle.h.
|
inline |
Definition at line 673 of file KVParticle.h.
TVector3 KVParticle::GetVelocity | ( | ) | const |
returns velocity vector in cm/ns units
Definition at line 1013 of file KVParticle.cpp.
|
inline |
Definition at line 677 of file KVParticle.h.
Double_t KVParticle::GetVperp | ( | ) | const |
returns transverse velocity in cm/ns units sign is +ve if py>=0, -ve if py<0
Definition at line 1033 of file KVParticle.cpp.
|
inline |
should be in fm
Definition at line 658 of file KVParticle.h.
[in] | frame | name of a previously defined kinematical frame |
Definition at line 783 of file KVParticle.h.
KVParticle KVParticle::InFrame | ( | const KVFrameTransform & | t | ) |
Use this method to obtain 'on-the-fly' some information on particle kinematics in a different reference frame. The default kinematics of the particle are unchanged.
Definition at line 582 of file KVParticle.cpp.
default initialisation
Definition at line 50 of file KVParticle.cpp.
|
inline |
Definition at line 777 of file KVParticle.h.
|
inline |
Definition at line 736 of file KVParticle.h.
Bool_t KVParticle::IsOK | ( | ) | const |
Determine whether this particle is considered "good" or not for analysis, depending on a previous call to SetIsOK(Bool_t flag). If SetIsOK has not been called, IsOK returns kTRUE by default.
Definition at line 318 of file KVParticle.cpp.
List all stored groups.
Definition at line 558 of file KVParticle.cpp.
Reimplemented from TObject.
Definition at line 353 of file KVParticle.cpp.
KVParticle & KVParticle::operator= | ( | const KVParticle & | rhs | ) |
KVParticle assignment operator.
Definition at line 371 of file KVParticle.cpp.
print out characteristics of particle
Reimplemented from TLorentzVector.
Reimplemented in KVSimNucleus, KVNucleus, KVReconstructedNucleus, KVINDRAReconNuc, and KVFAZIAReconNuc.
Definition at line 212 of file KVParticle.cpp.
recursive print out of all defined kinematical frames
Definition at line 189 of file KVParticle.cpp.
void KVParticle::RemoveAllGroups | ( | ) |
Remove all groups.
Definition at line 546 of file KVParticle.cpp.
Remove group from list of groups.
Definition at line 530 of file KVParticle.cpp.
void KVParticle::ResetEnergy | ( | ) |
Used for simulated particles after passage through some absorber/detector. The passage of a particle through the different absorbers modifies its kinetic energies, indeed the particle may be stopped in the detector. Calling this method will reset the particle's momentum to its initial value i.e. before it entered the first absorber. Particles which have not encountered any absorbers/detectors are left as they are.
Definition at line 393 of file KVParticle.cpp.
|
inlineprivate |
Definition at line 509 of file KVParticle.h.
|
inline |
Definition at line 706 of file KVParticle.h.
|
inline |
Definition at line 591 of file KVParticle.h.
Definition at line 596 of file KVParticle.h.
Definition at line 717 of file KVParticle.h.
Definition at line 601 of file KVParticle.h.
void KVParticle::SetFrame | ( | const Char_t * | frame, |
const KVFrameTransform & | ft | ||
) |
Define a Lorentz-boosted and/or rotated frame in which to calculate this particle's momentum and energy.
The new frame will have the name given in the string "frame", which can then be used to access the kinematics of the particle in different frames using GetFrame() (frame names are case-insensitive). Calling this method with the name of an existing frame will update the kinematics of the particle in that frame using the given transform (note that kinematics in all frames defined as 'subframes' of this frame will also be updated).
USING BOOSTS The boost velocity vector is that of the boosted frame with respect to the original frame of the particles in the event. The velocity vector can be given either in cm/ns units (default) or in units of 'c' (beta=kTRUE).
E.g. to define a frame moving at 0.1c in the +ve z-direction with respect to the original event frame:
(...supposing a valid pointer KVParticle* my_part...) TVector3 vframe(0,0,0.1); my_part->SetFrame("my_frame", KVFrameTransform(vframe, kTRUE));
or with velocity in cm/ns units (default): TVector3 vframe(0,0,3); my_part->SetFrame("my_frame", KVFrameTransform(vframe)); OR my_part->SetFrame("my_frame", vframe);
USING ROTATIONS According to the conventions adopted for the TRotation and TLorentzRotation classes, we actually use the inverse of the TLorentzRotation to make the transformation, to get the coordinates corresponding to a rotated coordinate system, not the coordinates of a rotated vector in the same coordinate system => you do not need to invert the transformation matrix
E.g. if you want to define a new frame whose coordinate axes are rotated with respect to the original axes, you can set up a TRotation like so:
TRotation rot; TVector3 newX, newY, newZ; // the new coordinate axes rot.RotateAxes(newX, newY, newZ);
If you are using one of the two global variables which calculate the event tensor (KVTensP and KVTensPCM) you can obtain the transformation to the tensor frame using:
TRotation rot; KVTensP* tens_gv;// pointer to tensor global variable tens_gv->GetTensor()->GetRotation(rot);// see KVTenseur3::GetRotation
Then the new frame can be defined by my_part->SetFrame("my_frame", KVFrameTransform(rot)); OR my_part->SetFrame("my_frame", rot);
USING COMBINED BOOST AND ROTATION You can define a frame using both a boost and a rotation like so: my_part->SetFrame("my_frame", KVFrameTransform(vframe,rot,kTRUE));
ACCESSING KINEMATICS IN NEW FRAMES In order to access the kinematics of the particle in the new frame:
my_part->GetFrame("my_frame")->GetTransverseEnergy();// transverse energy in "my_frame"
Definition at line 741 of file KVParticle.cpp.
void KVParticle::SetFrame | ( | const Char_t * | newframe, |
const Char_t * | oldframe, | ||
const KVFrameTransform & | ft | ||
) |
Define new kinematical frame by transformation from existing frame See SetFrame(const Char_t*,const KVFrameTransform&) for details on defining kinematically-transformed frames.
Definition at line 993 of file KVParticle.cpp.
|
inlineprivate |
Definition at line 505 of file KVParticle.h.
Set the (non-persistent) name of the reference frame for this particle's kinematics. Also sets a (persistent) parameter "frameName"
Definition at line 806 of file KVParticle.h.
|
protected |
Define for the particle a new list of groups.
Definition at line 457 of file KVParticle.cpp.
|
inline |
Definition at line 732 of file KVParticle.h.
Set acceptation/rejection status for the particle.
In order to 'forget' this status (accept all particles) use ResetIsOK()
Definition at line 339 of file KVParticle.cpp.
Change particle KE, keeping momentum direction constant If momentum is zero (i.e. no direction defined) the particle will be given a velocity in the positive z-direction
Definition at line 230 of file KVParticle.cpp.
Definition at line 568 of file KVParticle.h.
Definition at line 576 of file KVParticle.h.
Definition at line 580 of file KVParticle.h.
Set Momentum components (in MeV/c) if option is "cart" or "cartesian" we give cartesian components (x,y,z) if option is "spher" or "spherical" we give components (rho,theta,phi) in spherical system with theta, phi in DEGREES
Definition at line 1049 of file KVParticle.cpp.
set momentum with kinetic energy t and unit direction vector d (d is normalised first in case it is not a unit vector)
Definition at line 169 of file KVParticle.cpp.
Set Name of the particle.
Definition at line 412 of file KVParticle.cpp.
|
inlineprivate |
Definition at line 540 of file KVParticle.h.
|
inline |
Definition at line 820 of file KVParticle.h.
|
inlineprivate |
Definition at line 517 of file KVParticle.h.
Definition at line 437 of file KVParticle.h.
Definition at line 699 of file KVParticle.h.
Definition at line 441 of file KVParticle.h.
Definition at line 445 of file KVParticle.h.
Definition at line 449 of file KVParticle.h.
Definition at line 461 of file KVParticle.h.
Definition at line 453 of file KVParticle.h.
Definition at line 457 of file KVParticle.h.
void KVParticle::SetRandomMomentum | ( | Double_t | T, |
Double_t | thmin, | ||
Double_t | thmax, | ||
Double_t | phmin, | ||
Double_t | phmax, | ||
Option_t * | opt = "isotropic" |
||
) |
Give randomly directed momentum to particle with kinetic energy T Direction will be between (thmin,thmax) [degrees] limits in polar angle, and (phmin,phmax) [degrees] limits in azimuthal angle.
If opt = "" or "isotropic" (default) : direction is isotropically distributed over the solid angle If opt = "random" : direction is randomly distributed over solid angle
Based on KVPosition::GetRandomDirection().
Definition at line 136 of file KVParticle.cpp.
Definition at line 465 of file KVParticle.h.
Definition at line 469 of file KVParticle.h.
Definition at line 695 of file KVParticle.h.
TLorentzVector setters should not be used.
Definition at line 425 of file KVParticle.h.
Definition at line 429 of file KVParticle.h.
Definition at line 433 of file KVParticle.h.
Set velocity of particle (in cm/ns units)
Definition at line 1080 of file KVParticle.cpp.
Definition at line 473 of file KVParticle.h.
Definition at line 485 of file KVParticle.h.
Definition at line 489 of file KVParticle.h.
Definition at line 477 of file KVParticle.h.
Definition at line 481 of file KVParticle.h.
void KVParticle::UpdateAllFrames | ( | ) |
Call this method to update particle kinematics in all defined frames if you change the kinematics of the particle in its original/default frame.
Definition at line 913 of file KVParticle.cpp.
|
friend |
Definition at line 400 of file KVParticle.h.
|
mutableprivate |
Definition at line 420 of file KVParticle.h.
|
protected |
the momentum of the particle before it is slowed/stopped by an absorber
Definition at line 496 of file KVParticle.h.
if != nullptr, this object is a representation of the particle in a kinematical frame
Definition at line 502 of file KVParticle.h.
|
private |
non-persistent frame name field, sets when calling SetFrame method
Definition at line 408 of file KVParticle.h.
|
mutableprivate |
list of momenta of the particle in different Lorentz-boosted frames
list of TObjString for manage different group name
Definition at line 421 of file KVParticle.h.
|
private |
non-persistent name field - Is useful
Definition at line 407 of file KVParticle.h.
|
private |
parent kinematical frame
Definition at line 501 of file KVParticle.h.
|
protected |
a general-purpose list of parameters associated with this particle
Definition at line 497 of file KVParticle.h.
|
private |
Definition at line 500 of file KVParticle.h.
speed of light in cm/ns
Definition at line 422 of file KVParticle.h.