Adding comments (Laurent)
[u/mrichter/AliRoot.git] / PWG2 / FLOW / AliSelectorFoF.cxx
CommitLineData
0ee287d5 1/* AliSelectorFoF.cxx, v1.1 01/02/2007 esimili Exp */
2/* derived from AliSelector.cxx,v 1.17 2006/08/31 jgrosseo Exp */
3
4// The class definition in esdV0.h has been generated automatically
5// by the ROOT utility TTree::MakeSelector(). This class is derived
6// from the ROOT class TSelector. For more information on the TSelector
7// framework see $ROOTSYS/README/README.SELECTOR or the ROOT User Manual.
8
9// The following methods are defined in this file:
10// Begin(): called everytime a loop on the tree starts,
11// a convenient place to create your histograms.
12// SlaveBegin(): called after Begin(), when on PROOF called only on the
13// slave servers.
14// Process(): called for each event, in this function you decide what
15// to read and fill your histograms.
16// SlaveTerminate: called at the end of the loop on the tree, when on PROOF
17// called only on the slave servers.
18// Terminate(): called at the end of the loop on the tree,
19// a convenient place to draw/fit your histograms.
20//
21// To use this file, try the following session on your Tree T:
22//
23// Root > T->Process("AliSelector.C")
24// Root > T->Process("AliSelector.C","some options")
25// Root > T->Process("AliSelector.C+")
26//
27
28#include "AliSelectorFoF.h"
29
30#include <TStyle.h>
31#include <TSystem.h>
32#include <TCanvas.h>
33#include <TRegexp.h>
34#include <TTime.h>
35#include <TFriendElement.h>
36#include <TTree.h>
37#include <TChain.h>
38#include <TFile.h>
39#include <TTimeStamp.h>
40#include <TMath.h>
41#include <TObject.h>
42#include <TObjArray.h>
43#include <TSelector.h>
44
92016a03 45#include "AliLog.h"
46#include "AliESD.h"
0ee287d5 47#include "AliESDtrack.h"
92016a03 48#include "AliESDv0.h"
0ee287d5 49
50#include <stdio.h>
51#include <stdlib.h>
52#include <iostream>
53using namespace std; //required for resolving the 'cout' symbol
54
55#include "AliFlowEvent.h"
56#include "AliFlowTrack.h"
57#include "AliFlowV0.h"
58#include "AliFlowConstants.h"
59#include "AliFlowSelection.h"
60#include "AliFlowMaker.h"
61#include "AliFlowAnalyser.h"
62#include "AliFlowWeighter.h"
63
64ClassImp(AliSelectorFoF)
65
66//-----------------------------------------------------------------------
67
68AliSelectorFoF::AliSelectorFoF() :
69 TSelector(),
70 fTree(0),
71 fESD(0),
72 fCountFiles(0),
73 fKineFile(0)
74{
75 //
76 // Constructor. Initialization of pointers
77 //
78
79 fFlowEventFileName = "fof_flowEvts.root" ;
80 fFlowWgtFileName = "fof_flowWgts.root" ;
81 fFlowAnalysisFileName = "fof_flowAnal.root" ;
82}
83
84//-----------------------------------------------------------------------
85
86AliSelectorFoF::~AliSelectorFoF()
87{
88 //
89 // Destructor
90 //
91
92 if (fTree) { fTree->ResetBranchAddresses() ; }
93
94 if (fESD)
95 {
96 delete fESD;
97 fESD = 0;
98 }
99}
100
101//-----------------------------------------------------------------------
102
103void AliSelectorFoF::CheckOptions()
104{
105 // checks the option string for the debug flag
106
107 AliLog::SetClassDebugLevel(ClassName(), AliLog::kInfo);
108
109 TString option = GetOption();
110
111 if (option.Contains("moredebug"))
112 {
113 printf("Enabling verbose debug mode for %s\n", ClassName());
114 AliLog::SetClassDebugLevel(ClassName(), AliLog::kDebug+1);
115 AliInfo(Form("Called with option %s.", option.Data()));
116 }
117 else if (option.Contains("debug"))
118 {
119 printf("Enabling debug mode for %s\n", ClassName());
120 AliLog::SetClassDebugLevel(ClassName(), AliLog::kDebug);
121 AliInfo(Form("Called with option %s.", option.Data()));
122 }
123}
124
125//-----------------------------------------------------------------------
126
127void AliSelectorFoF::Begin(TTree*)
128{
129 // The Begin() function is called at the start of the query.
130 // When running with PROOF Begin() is only called on the client.
131 // The tree argument is deprecated (on PROOF 0 is passed).
132
133 cout << " HERE I begin !!! " << endl ; cout << endl ;
134
135 CheckOptions();
136
137 AliDebug(AliLog::kDebug, "============BEGIN===========");
138
139 // flags
140 fDoNothing = kFALSE ; // kFALSE ;
92016a03 141 fOnFlyAnalysis = kTRUE ; // kTRUE ;
0ee287d5 142 fOnFlyWeight = kFALSE ; // kTRUE ;
143 fSaveFlowEvents= kTRUE ; // kTRUE ;
144
145// Maker part :
146 fFlowMaker = new AliFlowMaker() ;
147 // ESD Cuts
148 fFlowMaker->SetNHitsCut(1) ;
149 fFlowMaker->SetECut(0.01,100.) ;
150 fFlowMaker->PrintCutList() ;
151
152// Opens flow event (output) file
153 if(fSaveFlowEvents)
154 {
155 // flow events file (output)
156 fFlowfile = new TFile(fFlowEventFileName.Data(),"RECREATE") ;
157 fFlowfile->cd() ;
158 cout << " . Writing AliFlowEvents on : " << fFlowEventFileName.Data() << " . " << endl ;
159 }
160
161// Analysis part :
162 if(fOnFlyWeight || fOnFlyAnalysis)
163 {
164 cout << " . Here the flow selection ... " << endl ;
165 cout << endl ;
166
167 // AliFlowSelection...
168 fFlowSelect = new AliFlowSelection() ;
169 // Event Cuts
170 fFlowSelect->SetCentralityCut(-1) ;
171 fFlowSelect->SetRunIdCut(-1) ;
172 // R.P. calculation cuts
327288af 173 for(int j=0;j<AliFlowConstants::kHars;j++)
0ee287d5 174 {
da5aa0a0 175 fFlowSelect->SetEtaCut(0., 0.9, j, 1) ;
0ee287d5 176 fFlowSelect->SetPtCut(0.1, 10. , j, 1);
177 }
178 fFlowSelect->SetConstrainCut(kTRUE) ;
179 fFlowSelect->SetDcaGlobalCut(0.,0.1);
180 // Correlation analysis cuts (not all of them)
da5aa0a0 181 fFlowSelect->SetEtaPart(-0.9,0.9);
0ee287d5 182 fFlowSelect->SetPtPart(0.1,10.);
183 fFlowSelect->SetConstrainablePart(kTRUE);
184 fFlowSelect->SetDcaGlobalPart(0.,0.1);
185 // V0 analysis cuts (not all of them ... they are useless anyway)
186 fFlowSelect->SetV0Mass(0.4875,0.5078) ; // Mk0 = 0.49765
da5aa0a0 187 fFlowSelect->SetV0SideBands(0.1) ;
0ee287d5 188 fFlowSelect->SetV0Pt(0.1,10.) ;
189 fFlowSelect->SetV0Eta(-2.1,2.1) ;
190 // print list :
191 cout << " . Selection for R.P. calculation: " << endl ;
192 fFlowSelect->PrintSelectionList() ;
193 cout << " . Selection for correlation analysis: " << endl ;
194 fFlowSelect->PrintList() ;
195 cout << " . Selection for V0 analysis: " << endl ;
196 fFlowSelect->PrintV0List() ;
197
198// Weight part :
199 if(fOnFlyWeight)
200 {
201 cout << " . Here the flow weighter ... " << endl ;
202 cout << endl ;
203
204 // AliFlowWeighter...
205 fFlowWeighter = new AliFlowWeighter(fFlowSelect) ;
206
207 // Wgt file (create)
208 fFlowWgtFile = new TFile(fFlowWgtFileName.Data(),"RECREATE") ;
209 fFlowWeighter->SetWgtFile(fFlowWgtFile) ;
210 cout << " . Writing Wgt Histograms on : " << fFlowWeighter->GetWgtFileName() << " . " << endl ;
211 fFlowWeighter->Init() ;
212 }
213
214 if(fOnFlyAnalysis)
215 {
216 cout << " . Here the flow analysis ... " << endl ;
217 cout << endl ;
218
219 // AliFlowAnalyser...
220 fFlowAnal = new AliFlowAnalyser(fFlowSelect) ;
221
222 // Wgt file (read)
92016a03 223 TString wgtFileName = "fof_flowWgts.root" ;
224 TFile* wgtFile = new TFile(wgtFileName.Data(),"READ");
225 if(!wgtFile || wgtFile->IsZombie()) { cout << " . NO phi Weights . " << endl ;}
226 else { fFlowAnal->FillWgtArrays(wgtFile) ; cout << " . Weights from : " << fFlowAnal->GetWgtFileName() << " . " << endl ; }
0ee287d5 227
228 // analysis file (output)
229 fFlowAnal->SetHistFileName(fFlowAnalysisFileName.Data()) ;
230 cout << " . Writing Analysis Histograms on : " << fFlowAnal->GetHistFileName() << " . " << endl ;
231
232 // Analysis settings
4e566f2f 233 fFlowAnal->SetFlowForV0() ; // default kTRUE.
da5aa0a0 234 fFlowAnal->SetSub(1) ; //eta // ->SetChrSub() ; //charge // ->SetRndSub() ; //random (default)
92016a03 235 //fFlowAnal->SetV1Ep1Ep2() ; // default kFALSE.
236 //fFlowAnal->SetShuffle() ; // default kFALSE. shuffles track array
237 //fFlowAnal->SetRedoWgt(); // default kFALSE. recalculates phiWgt (even if phiWgt file is already there)
238 fFlowAnal->SetUsePhiWgt(kFALSE) ; // default kTRUE if phiWgt file is there, kFALSE if not (& phiWgt file is created)
0ee287d5 239 fFlowAnal->SetUseOnePhiWgt() ; // or // fFlowAnal->SetUseFirstLastPhiWgt() ; // uses 1 or 3 wgt histograms (default is 1)
92016a03 240 fFlowAnal->SetUseBayWgt(kFALSE) ; // default kFALSE. uses bayesian weights in P.id.
0ee287d5 241 //fFlowAnal->SetUsePtWgt(); // default kFALSE. uses pT as a weight for RP determination
242 //fFlowAnal->SetUseEtaWgt(); // default kFALSE. uses eta as a weight for RP determination
4e566f2f 243 //fFlowAnal->SetCustomRespFunc(); // default kFALSE. uses the combined response function from the ESD
0ee287d5 244 fFlowAnal->Init() ;
245 }
246 }
247
248 cout << " ... init is done " << endl ;
249}
250
251//-----------------------------------------------------------------------
252
253void AliSelectorFoF::SlaveBegin(TTree* tree)
254{
255 // The SlaveBegin() function is called after the Begin() function.
256 // When running with PROOF SlaveBegin() is called on each slave server.
257 // The tree argument is deprecated (on PROOF 0 is passed).
258
259 CheckOptions();
260
261 AliDebug(AliLog::kDebug, "=======SLAVEBEGIN========");
262 AliDebug(AliLog::kDebug, Form("Hostname: %s", gSystem->HostName()));
263 AliDebug(AliLog::kDebug, Form("Time: %s", gSystem->Now().AsString()));
264
265 if (tree != 0) { Init(tree) ; }
266}
267
268//-----------------------------------------------------------------------
269
270void AliSelectorFoF::Init(TTree *tree)
271{
272 // The Init() function is called when the selector needs to initialize
273 // a new tree or chain. Typically here the branch addresses of the tree
274 // will be set. It is normaly not necessary to make changes to the
275 // generated code, but the routine can be extended by the user if needed.
276 // Init() will be called many times when running with PROOF.
277
278 AliDebug(AliLog::kDebug, "=========Init==========");
279
280 fTree = tree;
281
282 if (fTree == 0)
283 {
284 AliDebug(AliLog::kError, "ERROR: tree argument is 0.");
285 return;
286 }
287
288 // Set branch address
289 fTree->SetBranchAddress("ESD", &fESD);
290 if (fESD != 0) { AliDebug(AliLog::kInfo, "INFO: Found ESD branch in chain.") ; }
291}
292
293//-----------------------------------------------------------------------
294
295Bool_t AliSelectorFoF::Notify()
296{
297 // The Notify() function is called when a new file is opened. This
298 // can be either for a new TTree in a TChain or when when a new TTree
299 // is started when using PROOF. Typically here the branch pointers
300 // will be retrieved. It is normaly not necessary to make changes
301 // to the generated code, but the routine can be extended by the
302 // user if needed.
303
304 AliDebug(AliLog::kDebug, "=========NOTIFY==========");
305 AliDebug(AliLog::kDebug, Form("Hostname: %s", gSystem->HostName()));
306 AliDebug(AliLog::kDebug, Form("Time: %s", TTimeStamp(time(0)).AsString()));
307
308 ++fCountFiles;
309 if (fTree)
310 {
311 TFile *f = fTree->GetCurrentFile();
312 AliDebug(AliLog::kInfo, Form("Processing %d. file %s", fCountFiles, f->GetName()));
313 }
314 else
315 {
316 AliDebug(AliLog::kError, "fTree not available");
317 }
318
319 DeleteKinematicsFile();
320
321 return kTRUE;
322}
323
324//-----------------------------------------------------------------------
325
326Bool_t AliSelectorFoF::Process(Long64_t entry)
327{
328 // The Process() function is called for each entry in the tree (or possibly
329 // keyed object in the case of PROOF) to be processed. The entry argument
330 // specifies which entry in the currently loaded tree is to be processed.
331 // It can be passed to either TTree::GetEntry() or TBranch::GetEntry()
332 // to read either all or the required parts of the data. When processing
333 // keyed objects with PROOF, the object is already loaded and is available
334 // via the fObject pointer.
335 //
336 // This function should contain the "body" of the analysis. It can contain
337 // simple or elaborate selection criteria, run algorithms on the data
338 // of the event and typically fill histograms.
339
340 // WARNING when a selector is used with a TChain, you must use
341 // the pointer to the current TTree to call GetEntry(entry).
342 // The entry is always the local entry number in the current tree.
343 // Assuming that fTree is the pointer to the TChain being processed,
344 // use fTree->GetTree()->GetEntry(entry).
345
346 AliDebug(AliLog::kDebug, Form("=========PROCESS========== Entry %lld", entry));
347
348 if(!fTree)
349 {
350 AliDebug(AliLog::kError, "ERROR: fTree is 0.") ;
351 return kFALSE ;
352 }
353
354 fEventNumber = entry ;
355 fTree->GetTree()->GetEntry(fEventNumber) ;
356
357 if(fESD) { AliDebug(AliLog::kDebug, Form("ESD: We have %d tracks.", fESD->GetNumberOfTracks())); }
358 cout << " event !!! " << entry << endl ;
359
360 fRunID = fESD->GetRunNumber() ;
92016a03 361 // fEventNumber = fESD->GetEventNumber() ;
362 fEventNumber = -1 ;
0ee287d5 363 fNumberOfTracks = fESD->GetNumberOfTracks() ;
364 fNumberOfV0s = fESD->GetNumberOfV0s() ;
365
366 cout << " *evt n. " << fEventNumber << " (run " << fRunID << ") " << endl ;
367 cout << " tracks: " << fNumberOfTracks << " , v0s " << fNumberOfV0s << endl ;
368
369 // Dummy Loop
370 if(fDoNothing) { cout << " ... doing nothing ... " << endl ; return kTRUE; }
371
372 // Instantiate a new AliFlowEvent
373 cout << " filling the flow event :| " << endl ;
374 fFlowEvent = fFlowMaker->FillFlowEvent(fESD) ;
92016a03 375 if(!fFlowEvent) { cout << "! something bad occurred !" << endl ; return kFALSE ; }
376 else { cout << "# event done :) " << entry << " # ok ! # " << endl ; }
0ee287d5 377
378 // Saves the AliFlowEvent
379 if(fSaveFlowEvents)
380 {
381 cout << " saving flow event :| " << endl ;
382 TString strID = "" ; strID += entry ;
383 fFlowfile->cd() ; fFlowEvent->Write(strID.Data()) ;
92016a03 384 cout << "# event saved :) " << strID.Data() << " # ok ! # " << endl ; cout << endl ;
0ee287d5 385 }
386
387 // Weights
388 if(fOnFlyWeight)
389 {
390 cout << " weightening :| " << endl ;
391 bool done = fFlowWeighter->WeightEvent(fFlowEvent) ;
92016a03 392 if(done) { cout << "# weights done :) " << entry << " # ok ! # " << endl ; }
0ee287d5 393 }
394
395 // On-Fly Analysis
396 if(fOnFlyAnalysis)
397 {
398 cout << " doing analysis :| " << endl ;
4e566f2f 399 bool done = fFlowAnal->Analyze(fFlowEvent) ;
92016a03 400 fFlowAnal->PrintEventQuantities() ;
401 if(done) { cout << "# analysis done :) " << entry << " # ok ! # " << endl ; }
0ee287d5 402 }
403
404 return kTRUE;
405}
406
407//-----------------------------------------------------------------------
408
409void AliSelectorFoF::SlaveTerminate()
410{
411 // The SlaveTerminate() function is called after all entries or objects
412 // have been processed. When running with PROOF SlaveTerminate() is called
413 // on each slave server.
414
415 AliDebug(AliLog::kDebug, "=======SLAVETERMINATE=======");
416
417 DeleteKinematicsFile();
418}
419
420//-----------------------------------------------------------------------
421
422void AliSelectorFoF::Terminate()
423{
424 // The Terminate() function is the last function to be called during
425 // a query. It always runs on the client, it can be used to present
426 // the results graphically or save the results to file.
427
428 AliDebug(AliLog::kDebug, "=========TERMINATE==========");
429
430 cout << " Finished ... " << endl ;
431
432 if(fDoNothing) { cout << " ... & nothing was done ! " << endl ; cout << endl ; return ; }
433
92016a03 434 cout << endl ;
0ee287d5 435 cout << " nTracks: " << fFlowMaker->GetNgoodTracks() << endl ;
436 cout << " nV0s: " << fFlowMaker->GetNgoodV0s() << endl ;
437 cout << " nTracks (|eta|<0.5): " << fFlowMaker->GetNgoodTracksEta() << endl ;
438 cout << " nTracks+: " << fFlowMaker->GetNposiTracks() << endl ;
439 cout << " nTracks-: " << fFlowMaker->GetNnegaTracks() << endl ;
440 cout << " nTracks unconstrained: " << fFlowMaker->GetNunconstrained() << endl ;
92016a03 441 cout << " Bayesian (e,mu,pi,k,p) : " ;
0ee287d5 442 for(int ii=0;ii<5;ii++) { cout << fFlowMaker->GetBayesianNorm(ii) << " " ; }
92016a03 443 cout << " . " << endl ; cout << endl ;
0ee287d5 444
445 if(fOnFlyWeight) { fFlowWeighter->Finish() ; delete fFlowWeighter ; }
446
447 if(fOnFlyAnalysis)
448 {
0ee287d5 449 fFlowAnal->Resolution() ; cout << endl;
0ee287d5 450 fFlowAnal->Finish() ; cout << endl;
92016a03 451 //cout << endl ; cout << "Particles normalized abundance:" << endl; // shows the bayesian vector
452 //fFlowAnal->PrintRunBayesian() ; cout << endl;
0ee287d5 453 delete fFlowAnal ;
454 }
455
456 if(fSaveFlowEvents) { fFlowfile->Close() ; cout << " file closed . " << endl ; }
457 delete fFlowMaker ;
458
459 delete fFlowSelect ;
460
461 cout << endl ;
462 return ;
463}
464
465//-----------------------------------------------------------------------
466
467TTree* AliSelectorFoF::GetKinematics()
468{
469 // Returns kinematics tree corresponding to current ESD active in fTree
470 // Loads the kinematics from the kinematics file, the file is identified by replacing "AliESDs" to
471 // "Kinematics" in the file path of the ESD file. This is a hack, to be changed!
472
473 if (!fKineFile)
474 {
475 if(!fTree->GetCurrentFile()) { return 0 ; }
476
477 TString fileName(fTree->GetCurrentFile()->GetName());
478 fileName.ReplaceAll("AliESDs", "Kinematics");
479
480 // temporary workaround for PROOF bug #18505
481 fileName.ReplaceAll("#Kinematics.root#Kinematics.root", "#Kinematics.root");
482
483 AliDebug(AliLog::kInfo, Form("Opening %s", fileName.Data()));
484
485 fKineFile = TFile::Open(fileName);
486 if(!fKineFile) { return 0 ; }
487 }
488
489 return dynamic_cast<TTree*> (fKineFile->Get(Form("Event%d/TreeK", fTree->GetTree()->GetReadEntry())));
490}
491
492//-----------------------------------------------------------------------
493
494void AliSelectorFoF::DeleteKinematicsFile()
495{
496 //
497 // Closes the kinematics file and deletes the pointer.
498 //
499
500 if (fKineFile)
501 {
502 fKineFile->Close();
503 delete fKineFile;
504 fKineFile = 0;
505 }
506}
507
508//-----------------------------------------------------------------------
509
510