A global variable is an analysis tool for condensing the information in a multibody event into one or a few characteristic values. A simple example is the event multiplicity (the number of particles in each event), which can be used to characterize heavy-ion collision events in terms of violence or centrality.
In KaliVeda, the base class for a multibody event is KVEvent, which is basically a collection of nuclei (base class KVNucleus) (see Multi-particle Events chapter). Therefore the global variable classes below can be used with any event described by a class derived from KVEvent, containing particles described by a class which inherits from KVNucleus.
There are several base classes used for implementing global variables:
Many specific implementations are provided - see the list in the Global Variables module. Each global variable class defines how to calculate a given variable, but not for which particles or in which kinematical reference frame. This is done by calling the SetSelection()
and SetFrame()
methods of KVVarGlob.
Typical usage:
"var1"); // daughter class implementing 1-body global variable
SomeVarGlob VG(
VG.SetSelection( [particle selection criteria] );for kinematics] );
VG.SetFrame( [reference frame
VG.Init();
while( [loop over events] )
{while( [loop over particles in event] )
{// calculate contribution of particle to variable
VG.Fill( [particle] );
}
// perform any necessary calculations
VG.Calculate();
// accessing the result
auto vg = VG.GetValue(); // retrieve unique or principal value of variable
}
Note that although the Fill()
method is called for all particles, only those which satsify the conditions given to SetSelection()
will be used to calculate the variable.
All user analysis classes for analysing multibody event data inherit from KVEventSelector, which manages a list of global variables (a KVGVList). Users can add global variables to this list in their InitAnalysis()
or InitRun()
method (see Data Analysis) by calling
const Char_t* class_name, const Char_t* name) KVVarGlob* KVEventSelector::AddGV(
which returns a pointer to the new global variable, allowing to call directly its SetSelection()
or SetFrame()
method, as necessary.
For every event which is analysed, the list of global variables defined by the user is calculated automatically by KVEventSelector just before calling the user’s Analysis()
method, using all particles which satisfy the user’s particle selection criteria for the analysis, plus any additional selections applied to each global variable through its SetSelection()
method.
Therefore, in Analysis()
, one may retrieve the value(s) of each global variable for the event by using method
const Char_t* name) KVVarGlob* KVEventSelector::GetGV(
and the appropriate methods of KVVarGlob.
The typical usage case shown above, in the context of a KVEventSelector analysis class, would look something like this:
void SomeAnalysis::InitAnalysis()
{auto vg = AddGV("SomeVarGlob", "var1");
vg->SetSelection( [particle selection criteria] );for kinematics] );
vg->SetFrame( [reference frame
}
Bool_t SomeAnalysis::Analysis()
{// accessing the result
auto vg_val = GetGV("var1")->GetValue();
}
TTree
branches with global variable valuesThe KVGVList class provides methods to automatically create a branch for the value of each global variable in the list. In your analysis class, you can get a pointer to the internal list of global variables using method
KVGVList* KVEventSelector::GetGVList()
Then in your InitAnalysis()
or InitRun()
method, after creating a TTree
and filling the list of global variables required for your analysis with AddGV()
(see above), you can create the necessary branches in your TTree
like so:
new TTree(...);
TTree* myTree =
AddGV(...);
AddGV(...); GetGVList()->MakeBranches(myTree);
In your Analysis()
method in order to fill the branches with the values of all the variables for the current event you just need to call
GetGVList()->FillBranches(); myTree->Fill();
For multi-valued global variables (see for example KVFlowTensor), the default behaviour is to create a branch for each of the values calculated by the class. Each branch will have a name of the form toto.VARn
where toto
is the name of the global variable (i.e. the name given to AddGV
) and VARn
is the name of the n-th value calculated.
If you don’t want all values to be stored in your tree, you can choose to keep only the first N values like so:
"KVFlowTensor", "tensor")->SetMaxNumBranches(3); AddGV(
The resulting tree will contain the following branches:
tensor.FlowAngle
tensor.Sphericity
tensor.Coplanarity
corresponding to the first three values of this global variable.
Calling KVVarGlob::SetMaxNumBranches
with argument 0 (zero) will effectively prevent the creation of any branches for a given global variable. This can be used both for multi- and single-valued global variables.
Starting from version 1.12 of KaliVeda, global variables used in data analysis can be used to define new kinematical frames which can in turn be used by any variables which occur after it in the list.
In order to do so, call method KVVarGlob::SetNewFrameDefinition() with a lambda function having the following signature:
const KVVarGlob* vg){ e->SetFrame("_frame_name_", ...); } [](KVEvent* e,
When called (e.g. by KVGVList), the KVVarGlob pointer gives access to the global variable.
As an example of use, imagine that KVZmax is used to find the heaviest (largest Z) fragment in the forward CM hemisphere, then the velocity of this fragment is used to define a “QP_FRAME” in order to calculate the KVFlowTensor in this frame:
auto vg = AddGV("KVZmax", "zmax");
"CM");
vg->SetFrame("V>0", [](const KVNucleus* n){ return n->GetVpar()>0; }} );
vg->SetSelection( {
vg->SetNewFrameDefinition(const KVVarGlob* v){
[](KVEvent* e, "QP_FRAME", static_cast<const KVZmax*>(v)->GetZmax(0)->GetVelocity());
e->SetFrame(
});"KVFlowTensor", "qp_tensor");
vg = AddGV("QP_FRAME"); // frame will have been defined before tensor is filled vg->SetFrame(
Starting from version 1.12 of KaliVeda, conditions (‘cuts’) can be set on each variable which decide whether or not to retain an event for analysis. If any variable in the list fails the test, processing of the list is abandoned.
Selection criteria are set using lambda expressions. In this example, the variable “mtot” must have a value of at least 4 for the event to be retained:
auto mtot = AddGV("KVMult","mtot");
const KVVarGlob* v){ return v->GetValue()>=4; }); mtot->SetEventSelection([](
Any event selection criterion is tested as soon as each variable has been calculated. If the test fails, no further variables in the list will be calculated, and the current event will not be processed further: the user’s Analysis() method will NOT be called.
If you want to implement a new global variable, there are two possibilities. First, check that what you want to do is not possible with the KVVGSum class, or by changing the selection criteria and/or kinematical reference frame of an existing global variable: there is no need to write a new class if you want to calculate the total charge of fragments with 5 ≤ Z ≤ 15 emitted in polar angle range 10o ≤ θ ≤ 35o in the centre-of-mass frame: just use a KVZtot
auto vg = AddGV("KVZtot","zsum_z5to15_t10to35");
"CM");
vg->SetFrame("5<=Z<=15", [](const KVNucleus* n){ return (n->GetZ()>4 && n->GetZ()<16); });
vg->SetSelection("10<=theta<=35", [](const KVNucleus* n){ return (n->GetTheta()>=10 && n->GetTheta()<=35); }); vg->AddSelection(
Note that in the second selection, no reference is made to the "CM"
reference frame: the n
pointer in the lambda function will receive a pointer to the nucleus in the "CM"
frame, as defined by the call to SetFrame()
.
If not, you can write your own global variable class. You need to decide first of all which base class to use (looking at the inheritance of existing global variables may help you make your choice). If your global variable will calculate the mean (and possibly also the variance) of some quantity, the base class to use is KVVarGlobMean
. If your global variable will calculate a single value, use KVVarGlob1
. If neither of the previous two cases apply, use the most general base class, KVVarGlob
.
Next, decide if your global variable is of type 1-body (calculation based on properties of individual particles of the event, independently of the others), 2-body (calculation based on correlations between pairs of particles), or N-body (multibody correlations) (looking at existing classes may help to decide).
Then, use one of the following methods on the command line in order to generate skeleton ‘.h’ and ‘.cpp’ files for your class:
0] KVVarGlob::MakeClass("MyClassName", "The description of my class", type)
kaliveda[0] KVVarGlob1::MakeClass("MyClassName", "The description of my class", type)
kaliveda[0] KVVarGlobMean::MakeClass("MyClassName", "The description of my class", type) kaliveda[
with ‘type’ equal to one of KVVarGlob::kOneBody
, KVVarGlob::kTwoBody
, KVVarGlob::kNBody
.
When you have modified the skeleton according to your needs (see comments in generated code for help), you can test the compilation of your class in an interactive session:
0] .L MyClassName.cpp+ kaliveda[
You need to set up your environment correctly in order to tell ROOT where to find the sources for your global variables (this is essential for analysis tasks to run in batch mode), and in order to add your variables to the list available for KaliVeda analysis.
First modify (or create if it doesn’t exist) your .rootrc
file, adding/modifying the following line:
+Unix.*.Root.MacroPath: /full/path/to/directory/with/source/files +ACLiC.IncludePaths: -I/full/path/to/directory/with/source/files
You can use environment variables in this definition, as long as you enclose them in '$()'
, e.g. '$(HOME)'
. If you have several source directories, you can put them together:
+Unix.*.Root.MacroPath: /first/path:/second/path:/another/path +ACLiC.IncludePaths: -I/first/path -I/second/path -I/another/path
After relaunching ROOT/KaliVeda, you will now be able to compile your class(es) even if you launch in a different directory.
Next, modify (or create if it doesn’t exist) your .kvrootrc
file, adding for each of your global variables a line like this:
"MyClassName()" +Plugin.KVVarGlob: MyClassName MyClassName MyClassName.cpp+
Then you can use your global variable in an analysis in exactly the same way as other variables:
[in my_analyis.cpp]
auto myvg = AddGV("MyClassName", "myVarGlob");
You should try to avoid adding any new methods to your class (i.e. not defined by the interface of KVVarGlob). If you need additional parameters (named double
values) and/or options (named TString
values) for your global variables, use the methods provided by KVVarGlob (see “Options, parameters, reference frames, particle selection, etc.” in the class documentation).
If there is no alternative but to add a new method, you can get away with it:
[in my_analyis.cpp]#include "MyClassName.h"
...auto myvg = dynamical_cast<MyClassName*>(AddGV("MyClassName", "myVarGlob"));