KaliVeda
Toolkit for HIC analysis
KVReconDataAnalyser.cpp
1 //Created by KVClassFactory on Thu Jul 19 2018
2 //Author: John Frankland
3 
4 #include "KVBase.h"
5 #include "KVReconDataAnalyser.h"
6 #include "KVDBRun.h"
7 #include "KVMultiDetArray.h"
8 #include "KVDataAnalysisTask.h"
9 #include "KVDataSet.h"
10 #include "TChain.h"
11 #include "TObjString.h"
12 #include "TChain.h"
13 #include "KVAvailableRunsFile.h"
14 
15 #include <KVClassFactory.h>
16 #include <KVTriggerConditions.h>
17 
18 using namespace std;
19 
21 
22 
23 
26  : fSelector(nullptr), theChain(nullptr)
27 {
28 }
29 
30 
31 
34 
36 {
37  //Reset task variables
39  theChain = nullptr;
40  fSelector = nullptr;
41  TotalEntriesToRead = 0;
42 }
43 
44 
45 
46 
47 
48 
51 
53 {
54  // Checks the task variables
55 
57 
58  cout << "============> Analysis summary <=============" << endl;
59  cout << "Analysis of runs " << GetRunList().GetList() << " with the class ";
60  cout << "\"" << GetUserClass().name << "\"." << endl;
61  if (GetNbEventToRead()) {
62  cout << GetNbEventToRead() << " events will be processed." << endl;
63  }
64  else {
65  cout << "All events will be processed." << endl;
66  }
67  cout << "=============================================" << endl;
68 
69  return kTRUE;
70 }
71 
72 
73 
74 
79 
81 {
82  //Run the interactive analysis
83 
84  //make the chosen dataset the active dataset ( = gDataSet; note this also opens database
85  //and positions gDataBase & gExpDB).
86  GetDataSet()->cd();
87  fSelector = nullptr;
88 
89  theChain = new TChain("ReconEvents");
90  theChain->SetDirectory(0); // we handle delete
91 
92  // open and add to TChain all required files
93  // we force the opening of the files to avoid problems with xrootd which sometimes
94  // seems to have a little difficulty
95  for (auto& run : GetRunList()) {
96  TString fullPathToRunfile = gDataSet->GetFullPathToRunfile(GetDataType(), run);
97  cout << "Opening file " << fullPathToRunfile << endl;
98  TFile* f = gDataSet->OpenRunfile<TFile>(GetDataType(), run);
99  cout << "Adding file " << fullPathToRunfile;
100  cout << " to the TChain." << endl;
101  dynamic_cast<TChain*>(theChain)->Add(fullPathToRunfile);
102  if (f && !f->IsZombie()) {
103  // update run infos in available runs file if necessary
105  if (ARF->InfosNeedUpdate(run, gSystem->BaseName(fullPathToRunfile))) {
106  if (!((TTree*)f->Get("ReconEvents"))) {
107  Error("SubmitTask", "No tree named ReconEvents is present in the current file");
108  delete theChain;
109  return;
110  }
111  TEnv* treeInfos = (TEnv*)((TTree*)f->Get("ReconEvents"))->GetUserInfo()->FindObject("TEnv");
112  if (treeInfos) {
113  TString kvversion = treeInfos->GetValue("KVBase::GetKVVersion()", "");
114  TString username = treeInfos->GetValue("gSystem->GetUserInfo()->fUser", "");
115  if (kvversion != "") ARF->UpdateInfos(run, gSystem->BaseName(fullPathToRunfile), kvversion, username);
116  }
117  else {
118  Info("SubmitTask", "No TEnv object associated to the tree");
119  }
120  }
121  }
122  }
123  TotalEntriesToRead = theChain->GetEntries();
124  TString option = Form("EventsReadInterval=%lld,", GetAnalysisTask()->GetStatusUpdateInterval());
125  option += Form("FullRunList=%s", GetFullRunList().GetList().Data());
126 
127  // Add any user-defined options
128  if (GetUserClassOptions() != "") {
129  option += ",";
130  option += GetUserClassOptions();
131  }
132 
133  TObject* new_selector = GetInstanceOfUserClass();
134 
135  if (!new_selector || !new_selector->InheritsFrom("TSelector")) {
136  cout << "The selector \"" << GetUserClassName() << "\" is not valid." << endl;
137  cout << "Process aborted." << endl;
138  SafeDelete(new_selector);
139  }
140  else {
141  SafeDelete(new_selector);
142  Info("SubmitTask", "Beginning TChain::Process...");
143 #ifdef WITH_CPP11
144  if (GetProofMode() != KVDataAnalyser::EProofMode::None) dynamic_cast<TChain*>(theChain)->SetProof(kTRUE);
145 #else
146  if (GetProofMode() != KVDataAnalyser::None) dynamic_cast<TChain*>(theChain)->SetProof(kTRUE);
147 #endif
148  TString analysis_class;
149  if (GetAnalysisTask()->WithUserClass())
150  analysis_class.Form("%s%s", GetUserClass().full_path_imp().Data(), GetACliCMode());
151  else analysis_class = GetUserClassName();
152 
153  if (GetNbEventToRead()) {
154  theChain->Process(analysis_class, option.Data(), GetNbEventToRead());
155  }
156  else {
157  theChain->Process(analysis_class, option.Data());
158  }
159  }
160  delete theChain;
161  fSelector = nullptr; //deleted by TChain/TTreePlayer
162 }
163 
164 
165 
181 
183 {
184  // Called by currently-processed TSelector when a new file in the TChain is opened.
185  //
186  // We call gMultiDetArray->SetParameters for the current run.
187  // Whether or not only physics parameters are set, or the full set of calibrations and identifications
188  // for each detector/identification telescope is determined by the environment variable
189  //~~~
190  // [dataset].ReconAnalysis.WithCalibInfos: [yes/no]
191  //~~~
192  // This can be overridden for any individual analysis by setting the analysis class option
193  //~~~~
194  // WithCalibInfos=yes
195  //~~~~
196  //
197  // Infos on currently read file/tree are printed.
198 
199  Bool_t physics_parameters_only = !gDataSet->GetDataSetEnv("ReconAnalysis.WithCalibInfos", kTRUE);
200  if (fSelector->IsOptGiven("WithCalibInfos"))
201  physics_parameters_only = (fSelector->GetOpt("WithCalibInfos") != "yes");
202 
203  Int_t run = GetRunNumberFromFileName(theChain->GetCurrentFile()->GetName());
205  KVMultiDetArray::MakeMultiDetector(GetDataSet()->GetName(), run);
206 
207  KVDBRun* CurrentRun = gExpDB->GetDBRun(run);
208  SetCurrentRun(CurrentRun);
209  fSelector->SetCurrentRun(CurrentRun);
210 
211  cout << endl << " =================== New Run =================== " <<
212  endl << endl;
213 
214  CurrentRun->Print();
215  if (CurrentRun->GetSystem()) {
216  if (CurrentRun->GetSystem()->GetKinematics())
217  CurrentRun->GetSystem()->GetKinematics()->Print();
218  }
219 
220  cout << endl << " ================================================= " <<
221  endl << endl;
222 
223  PrintTreeInfos();
224  Info("preInitRun", "Data written with series %s, release %d", GetDataSeries().Data(),
225  GetDataReleaseNumber());
226  fRustines.InitializePatchList(GetDataSet()->GetName(), GetDataType(), run, GetDataSeries(),
227  GetDataReleaseNumber(), theChain->GetCurrentFile()->GetStreamerInfoCache());
228  fRustines.Print();
229 }
230 
231 
232 
235 
237 {
238  // apply any required patches to data
239  if (fRustines.HasActivePatches()) fRustines.Apply(fSelector->GetEvent());
240 }
241 
242 
243 
245 
247 {
248  return (TEnv*)theChain->GetTree()->GetUserInfo()->FindObject("TEnv");
249 }
250 
251 
252 
268 
270 {
271  // When called from the InitRun() method of a user's analysis class, this method will ensure that only data
272  // compatible with the experimental trigger will be provided for analysis in the user's Analysis() method.
273  //
274  // This will be done by searching for a KVTriggerConditions plugin class defined for the currently-analysed
275  // dataset, defined like so:
276  //
277  //~~~~
278  //+Plugin.KVTriggerConditions: [dataset] [classname] [libname] "[default constructor]()"
279  //~~~~
280  //
281  // An object of the plugin class will be instantiated, and then its overridden
282  // KVTriggerConditions::SetTriggerConditionsForRun() method will be called with 2 arguments:
283  // - a pointer to the user's analysis class derived from KVReconEventSelector;
284  // - the number of the run currently being analysed
285 
286  TPluginHandler* ph = KVBase::LoadPlugin("KVTriggerConditions", GetDataSet()->GetName());
287  if (!ph) {
288  Info("SetTriggerConditionsForRun",
289  "No definition of trigger conditions available for dataset %s",
290  GetDataSet()->GetName());
291  return;
292  }
293  std::unique_ptr<KVTriggerConditions> trig((KVTriggerConditions*) ph->ExecPlugin(0));
294  trig->SetTriggerConditionsForRun(fSelector, run);
295 }
296 
297 
298 
301 
303 {
304  // Print informations on currently analysed TTree
305  TEnv* treeInfos = GetReconDataTreeInfos();
306  if (!treeInfos) return;
307  cout << endl << "----------------------------------------------------------------------------------------------------" << endl;
308  cout << "INFORMATIONS ON VERSION OF KALIVEDA USED TO GENERATE FILE:" << endl << endl;
309  fDataVersion = treeInfos->GetValue("KVBase::GetKVVersion()", "(unknown)");
310  cout << "version = " << fDataVersion << endl ;
311  cout << "build date = " << treeInfos->GetValue("KVBase::GetKVBuildDate()", "(unknown)") << endl ;
312  cout << "source directory = " << treeInfos->GetValue("KVBase::GetKVSourceDir()", "(unknown)") << endl ;
313  cout << "KVROOT = " << treeInfos->GetValue("KVBase::GetKVRoot()", "(unknown)") << endl ;
314  if (strcmp(treeInfos->GetValue("KVBase::bzrBranchNick()", "(unknown)"), "(unknown)")) {
315  cout << "BZR branch name = " << treeInfos->GetValue("KVBase::bzrBranchNick()", "(unknown)") << endl ;
316  cout << "BZR revision #" << treeInfos->GetValue("KVBase::bzrRevisionNumber()", "(unknown)") << endl ;
317  cout << "BZR revision ID = " << treeInfos->GetValue("KVBase::bzrRevisionId()", "(unknown)") << endl ;
318  cout << "BZR revision date = " << treeInfos->GetValue("KVBase::bzrRevisionDate()", "(unknown)") << endl ;
319  }
320  else if (strcmp(treeInfos->GetValue("KVBase::gitBranch()", "(unknown)"), "(unknown)")) {
321  cout << "git branch = " << treeInfos->GetValue("KVBase::gitBranch()", "(unknown)") << endl;
322  cout << "git commit = " << treeInfos->GetValue("KVBase::gitCommit()", "(unknown)") << endl;
323  }
324  cout << endl << "INFORMATIONS ON GENERATION OF FILE:" << endl << endl;
325  cout << "Generated by : " << treeInfos->GetValue("gSystem->GetUserInfo()->fUser", "(unknown)") << endl ;
326  cout << "Analysis task : " << treeInfos->GetValue("AnalysisTask", "(unknown)") << endl ;
327  cout << "Job name : " << treeInfos->GetValue("BatchSystem.JobName", "(unknown)") << endl ;
328  cout << "Job submitted from : " << treeInfos->GetValue("LaunchDirectory", "(unknown)") << endl ;
329  cout << "Runs : " << treeInfos->GetValue("Runs", "(unknown)") << endl ;
330  cout << "Number of events requested : " << treeInfos->GetValue("NbToRead", "(unknown)") << endl ;
331  cout << endl << "----------------------------------------------------------------------------------------------------" << endl;
332 
333  // if possible, parse fDataVersion into series and release number
334  // e.g. if fDataVersion = "1.8.10":
335  // => fDataSeries = "1.8" fDataReleaseNum = 10
336  // e.g. if fDataVersion = "1.11/00":
337  // => fDataSeries = "1.11" fDataReleaseNum = 0
338  Int_t a, b, c;
339  if (fDataVersion != "(unknown)") {
340  if (fDataVersion.Contains("/")) {
341  if (sscanf(fDataVersion.Data(), "%d.%d/%d", &a, &b, &c) == 3) {
342  fDataSeries.Form("%d.%d", a, b);
343  fDataReleaseNum = c;
344  }
345  }
346  else {
347  if (sscanf(fDataVersion.Data(), "%d.%d.%d", &a, &b, &c) == 3) {
348  fDataSeries.Form("%d.%d", a, b);
349  fDataReleaseNum = c;
350  }
351  }
352  }
353  else {
354  fDataSeries = "";
355  fDataReleaseNum = -1;
356  }
357 }
358 
359 
int Int_t
#define SafeDelete(p)
#define f(i)
#define c(i)
bool Bool_t
constexpr Bool_t kFALSE
constexpr Bool_t kTRUE
Option_t Option_t option
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t b
char * Form(const char *fmt,...)
R__EXTERN TSystem * gSystem
void Print(Option_t *opt="") const override
Definition: KV2Body.cpp:813
Handles lists of available runs for different datasets and types of data.
Bool_t InfosNeedUpdate(const run_index_t &run, const KVString &filename)
void UpdateInfos(const run_index_t &run, const KVString &filename, const KVString &kvversion, const KVString &username)
static TPluginHandler * LoadPlugin(const Char_t *base, const Char_t *uri="0")
Definition: KVBase.cpp:772
Description of an experimental run in database ,,.
Definition: KVDBRun.h:40
void Print(Option_t *option="") const override
Definition: KVDBRun.cpp:61
KVDBSystem * GetSystem() const
Definition: KVDBRun.cpp:240
KV2Body * GetKinematics()
Definition: KVDBSystem.cpp:80
virtual Bool_t CheckTaskVariables()
void Reset() override
TString GetFullPathToRunfile(const KVString &type, const run_index_t &run) const
Definition: KVDataSet.cpp:929
KVAvailableRunsFile * GetAvailableRunsFile(const Char_t *type) const
Definition: KVDataSet.cpp:50
KVString GetDataSetEnv(const Char_t *type, const Char_t *defval="") const
Definition: KVDataSet.cpp:782
FileType * OpenRunfile(const KVString &type, const run_index_t &run)
Definition: KVDataSet.h:164
KVDBRun * GetDBRun(Int_t number) const
Definition: KVExpDB.h:89
static KVMultiDetArray * MakeMultiDetector(const Char_t *dataset_name, Int_t run=-1, TString classname="KVMultiDetArray", KVExpDB *db=nullptr)
static void SetMakeMultiDetectorPhysicsParametersOnly(Bool_t on=kTRUE)
Manages user analysis of reconstructed experimental data.
void preAnalysis() override
apply any required patches to data
TEnv * GetReconDataTreeInfos() const
void PrintTreeInfos()
Print informations on currently analysed TTree.
Bool_t CheckTaskVariables(void) override
Checks the task variables.
virtual void SetTriggerConditionsForRun(int) override
void SubmitTask() override
void Reset() override
Reset task variables.
Set trigger conditions for analysis of reconstructed data.
virtual void SetProof(Bool_t on=kTRUE, Bool_t refresh=kFALSE, Bool_t gettreeheader=kFALSE)
virtual const char * GetValue(const char *name, const char *dflt) const
virtual TObject * FindObject(const char *name) const
virtual Bool_t InheritsFrom(const char *classname) const
Longptr_t ExecPlugin(int nargs)
void Form(const char *fmt,...)
virtual const char * BaseName(const char *pathname)
BinData::ErrorType GetDataType(const TGraph *gr, DataOptions &fitOpt)
void Error(const char *location, const char *fmt,...)
void Info(const char *location, const char *fmt,...)
Add
TArc a
ClassImp(TPyArg)