updating physics selection for heavy-ion data
[u/mrichter/AliRoot.git] / ANALYSIS / AliPhysicsSelection.cxx
CommitLineData
61899827 1/* $Id: AliPhysicsSelection.cxx 35782 2009-10-22 11:54:31Z jgrosseo $ */
2
3/**************************************************************************
4 * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
5 * *
6 * Author: The ALICE Off-line Project. *
7 * Contributors are mentioned in the code where appropriate. *
8 * *
9 * Permission to use, copy, modify and distribute this software and its *
10 * documentation strictly for non-commercial purposes is hereby granted *
11 * without fee, provided that the above copyright notice appears in all *
12 * copies and that both the copyright notice and this permission notice *
13 * appear in the supporting documentation. The authors make no claims *
14 * about the suitability of this software for any purpose. It is *
15 * provided "as is" without express or implied warranty. *
16 **************************************************************************/
17
18//-------------------------------------------------------------------------
19// Implementation of Class AliPhysicsSelection
528640ed 20// This class selects collision candidates from data runs, applying selection cuts on triggers
21// and background rejection based on the content of the ESD
22//
23// Usage:
24//
a2ce3799 25// Create the object:
528640ed 26// fPhysicsSelection = new AliPhysicsSelection;
a2ce3799 27//
28// For MC data, call
29// fPhysicsSelection->SetAnalyzeMC()
528640ed 30//
31// To check if an event is a collision candidate, use:
32// fPhysicsSelection->IsCollisionCandidate(fESD)
33//
34// After processing save the resulting histograms to a file with (a folder physics_selection
35// will be created that contains the histograms):
36// fPhysicsSelection->SaveHistograms("physics_selection")
37//
38// To print statistics after processing use:
39// fPhysicsSelection->Print();
40//
73cc8654 41// The BX ids corresponding to real bunches crossings p2 are
42// automatically selected. You cannot process runs with different
43// filling schemes if this option is set. If you want to disable this,
44// use:
45// fPhysicsSelection->SetUseBXNumbers(0);
46//
47//
48// If you are analizing muons and you want to keep the muon triggers
49// besides the CINT1B you can set:
50// fPhysicsSelection->SetUseMuonTriggers();
51//
85c71ba7 52//
53// To compute the Background automatically using the control triggers
54// use:
55// fPhysicsSelection->SetComputeBG();
56// this will show the value of the Beam Gas, accidentals and good
57// events as additional rows in the statistic tables, but it will NOT
58// subtract the background automatically.
59// This option will only work for runs taken with the CINT1
60// suite. This options enables automatically also the usage of BX
61// numbers. You can only process one run at a time if you require this
62// options, because it uses the bunch intensity estimated run by run.
63//
64// The BG will usually be more important in the so-called "bin 0": the
65// class can also compute the statistics table for events in this
66// bin. Since the definition of bin 0 may in general change from
67// analysis to analysis, the user needs to provide a callback
68// implementing the definition of bin zero. The callback should be
69// implemented as a method in the analysis task and should override
70// the IsEventInBinZero method of AliAnalysisTaskSE, and should thus
71// have the the following prototype:
72// Bool_t IsEventInBinZero();
73// It should return true if the event is in the bin 0 and it is set by
74// passing to the physics selection the NAME of the task where the
75// callback is implemented:
76// fPhysicsSelection->SetBin0Callback("MyTask").
77//
78//
decf6fd4 79// Usually the class selects the trigger scheme by itself depending on the run number.
80// Nevertheless, you can do that manually by calling AddCollisionTriggerClass() and AddBGTriggerClass()
81// Example:
82// To define the class CINT1B-ABCE-NOPF-ALL as collision trigger (those will be accepted as
83// collision candidates when they pass the selection):
84// AddCollisionTriggerClass("+CINT1B-ABCE-NOPF-ALL #769 #3119");
85// To select on bunch crossing IDs in addition, use:
86// AddCollisionTriggerClass("+CINT1B-ABCE-NOPF-ALL #769 #3119");
87// To define the class CINT1A-ABCE-NOPF-ALL as a background trigger (those will only be counted
88// for the control histograms):
89// AddBGTriggerClass("+CINT1A-ABCE-NOPF-ALL");
90// You can also specify more than one trigger class in a string or you can require that some are *not*
91// present. The following line would require CSMBA-ABCE-NOPF-ALL, but CSMBB-ABCE-NOPF-ALL is not allowed
92// to be present:
93// AddBGTriggerClass("+CSMBA-ABCE-NOPF-ALL -CSMBB-ABCE-NOPF-ALL");
94//
85c71ba7 95// Origin: Jan Fiete Grosse-Oetringhaus, CERN
96// Michele Floris, CERN
61899827 97//-------------------------------------------------------------------------
98
99#include <Riostream.h>
100#include <TH1F.h>
101#include <TH2F.h>
102#include <TList.h>
103#include <TIterator.h>
104#include <TDirectory.h>
296dd262 105#include <TObjArray.h>
61899827 106
107#include <AliPhysicsSelection.h>
108
109#include <AliTriggerAnalysis.h>
110#include <AliLog.h>
111
112#include <AliESDEvent.h>
85c71ba7 113#include <AliAnalysisTaskSE.h>
114#include "AliAnalysisManager.h"
141265a2 115#include "TPRegexp.h"
61899827 116
117ClassImp(AliPhysicsSelection)
118
119AliPhysicsSelection::AliPhysicsSelection() :
296dd262 120 AliAnalysisCuts("AliPhysicsSelection", "AliPhysicsSelection"),
121 fCurrentRun(-1),
a2ce3799 122 fMC(kFALSE),
296dd262 123 fCollTrigClasses(),
124 fBGTrigClasses(),
125 fTriggerAnalysis(),
126 fBackgroundIdentification(0),
91bea6e7 127 fHistBunchCrossing(0),
73cc8654 128 fHistTriggerPattern(0),
91bea6e7 129 fSkipTriggerClassSelection(0),
85c71ba7 130 fUsingCustomClasses(0),
131 fSkipV0(0),
132 fBIFactorA(1),
133 fBIFactorC(1),
141265a2 134 fBIFactorAC(1),
85c71ba7 135 fComputeBG(0),
141265a2 136 fBGStatOffset(0),
73cc8654 137 fUseBXNumbers(1),
138 fUseMuonTriggers(0),
85c71ba7 139 fFillingScheme(""),
ff097e3f 140 fBin0CallBack(""),
141 fBin0CallBackPointer(0)
61899827 142{
143 // constructor
144
296dd262 145 fCollTrigClasses.SetOwner(1);
146 fBGTrigClasses.SetOwner(1);
147 fTriggerAnalysis.SetOwner(1);
85c71ba7 148 fHistStatistics[0] = 0;
149 fHistStatistics[1] = 0;
296dd262 150
61899827 151 AliLog::SetClassDebugLevel("AliPhysicsSelection", AliLog::kWarning);
152}
153
154AliPhysicsSelection::~AliPhysicsSelection()
155{
156 // destructor
85c71ba7 157
296dd262 158 fCollTrigClasses.Delete();
159 fBGTrigClasses.Delete();
160 fTriggerAnalysis.Delete();
61899827 161
85c71ba7 162 if (fHistStatistics[0])
163 {
164 delete fHistStatistics[0];
165 fHistStatistics[0] = 0;
166 }
167 if (fHistStatistics[1])
61899827 168 {
85c71ba7 169 delete fHistStatistics[1];
170 fHistStatistics[1] = 0;
61899827 171 }
172
173 if (fHistBunchCrossing)
174 {
175 delete fHistBunchCrossing;
176 fHistBunchCrossing = 0;
177 }
73cc8654 178 if (fHistTriggerPattern)
179 {
180 delete fHistTriggerPattern;
181 fHistTriggerPattern = 0;
182 }
716f045f 183 if (fBackgroundIdentification)
184 {
185 delete fBackgroundIdentification;
186 fBackgroundIdentification = 0;
187 }
61899827 188}
296dd262 189
c2ba5a61 190UInt_t AliPhysicsSelection::CheckTriggerClass(const AliESDEvent* aEsd, const char* trigger, Int_t& triggerLogic) const
296dd262 191{
192 // checks if the given trigger class(es) are found for the current event
c2ba5a61 193 // format of trigger: +TRIGGER1 -TRIGGER2 [#XXX] [&YY] [*ZZ]
296dd262 194 // requires TRIGGER1 and rejects TRIGGER2
0c6c629b 195 // in bunch crossing XXX
1209509c 196 // if successful, YY is returned (for association between entry in fCollTrigClasses and AliVEvent::EOfflineTriggerTypes)
c2ba5a61 197 // triggerLogic is filled with ZZ, defaults to kCINT1
296dd262 198
decf6fd4 199 Bool_t foundBCRequirement = kFALSE;
200 Bool_t foundCorrectBC = kFALSE;
201
0c6c629b 202 UInt_t returnCode = AliVEvent::kUserDefined;
c2ba5a61 203 triggerLogic = kCINT1;
0c6c629b 204
296dd262 205 TString str(trigger);
206 TObjArray* tokens = str.Tokenize(" ");
207
208 for (Int_t i=0; i < tokens->GetEntries(); i++)
209 {
210 TString str2(((TObjString*) tokens->At(i))->String());
211
decf6fd4 212 if (str2[0] == '+' || str2[0] == '-')
296dd262 213 {
decf6fd4 214 Bool_t flag = (str2[0] == '+');
215
216 str2.Remove(0, 1);
217
218 if (flag && !aEsd->IsTriggerClassFired(str2))
219 {
220 AliDebug(AliLog::kDebug, Form("Rejecting event because trigger class %s is not present", str2.Data()));
221 delete tokens;
222 return kFALSE;
223 }
224 if (!flag && aEsd->IsTriggerClassFired(str2))
225 {
226 AliDebug(AliLog::kDebug, Form("Rejecting event because trigger class %s is present", str2.Data()));
227 delete tokens;
228 return kFALSE;
229 }
296dd262 230 }
decf6fd4 231 else if (str2[0] == '#')
296dd262 232 {
decf6fd4 233 foundBCRequirement = kTRUE;
234
235 str2.Remove(0, 1);
236
237 Int_t bcNumber = str2.Atoi();
238 AliDebug(AliLog::kDebug, Form("Checking for bunch crossing number %d", bcNumber));
239
240 if (aEsd->GetBunchCrossNumber() == bcNumber)
241 {
242 foundCorrectBC = kTRUE;
243 AliDebug(AliLog::kDebug, Form("Found correct bunch crossing %d", bcNumber));
244 }
296dd262 245 }
0c6c629b 246 else if (str2[0] == '&' && !fUsingCustomClasses)
247 {
248 str2.Remove(0, 1);
249
1209509c 250 returnCode = str2.Atoll();
0c6c629b 251 }
c2ba5a61 252 else if (str2[0] == '*')
253 {
254 str2.Remove(0, 1);
255
256 triggerLogic = str2.Atoi();
257 }
decf6fd4 258 else
259 AliFatal(Form("Invalid trigger syntax: %s", trigger));
296dd262 260 }
261
262 delete tokens;
decf6fd4 263
264 if (foundBCRequirement && !foundCorrectBC)
265 return kFALSE;
266
0c6c629b 267 return returnCode;
296dd262 268}
2d45a1a7 269
270//______________________________________________________________________________
271TObject *AliPhysicsSelection::GetStatistics(Option_t *option) const
272{
273// Get the statistics histograms ("ALL" and "BIN0")
274 TString opt(option);
275 opt.ToUpper();
276 Int_t ihist = 0;
277 if (opt == "ALL") ihist = kStatIdxAll;
278 if (opt == "BIN0") ihist = kStatIdxBin0;
279 return fHistStatistics[ihist];
280}
281
282//______________________________________________________________________________
0c6c629b 283UInt_t AliPhysicsSelection::IsCollisionCandidate(const AliESDEvent* aEsd)
61899827 284{
285 // checks if the given event is a collision candidate
0c6c629b 286 //
287 // returns a bit word describing the fired offline triggers (see AliVEvent::EOfflineTriggerTypes)
61899827 288
141265a2 289 if (fCurrentRun != aEsd->GetRunNumber()) {
95e6f4ec 290 if (!Initialize(aEsd->GetRunNumber()))
291 AliFatal(Form("Could not initialize for run %d", aEsd->GetRunNumber()));
141265a2 292 if(fComputeBG) SetBIFactors(aEsd); // is this safe here?
293 }
61899827 294 const AliESDHeader* esdHeader = aEsd->GetHeader();
295 if (!esdHeader)
296 {
297 AliError("ESD Header could not be retrieved");
298 return kFALSE;
299 }
300
a2ce3799 301 // check event type; should be PHYSICS = 7 for data and 0 for MC
302 if (!fMC)
303 {
304 if (esdHeader->GetEventType() != 7)
305 return kFALSE;
306 }
307 else
308 {
309 if (esdHeader->GetEventType() != 0)
310 AliFatal(Form("Invalid event type for MC: %d", esdHeader->GetEventType()));
311 }
758941d4 312
0c6c629b 313 UInt_t accept = 0;
296dd262 314
315 Int_t count = fCollTrigClasses.GetEntries() + fBGTrigClasses.GetEntries();
316 for (Int_t i=0; i < count; i++)
61899827 317 {
296dd262 318 const char* triggerClass = 0;
319 if (i < fCollTrigClasses.GetEntries())
320 triggerClass = ((TObjString*) fCollTrigClasses.At(i))->String();
321 else
322 triggerClass = ((TObjString*) fBGTrigClasses.At(i - fCollTrigClasses.GetEntries()))->String();
61899827 323
296dd262 324 AliDebug(AliLog::kDebug, Form("Processing trigger class %s", triggerClass));
61899827 325
296dd262 326 AliTriggerAnalysis* triggerAnalysis = static_cast<AliTriggerAnalysis*> (fTriggerAnalysis.At(i));
61899827 327
296dd262 328 triggerAnalysis->FillTriggerClasses(aEsd);
61899827 329
c2ba5a61 330 Int_t triggerLogic = 0;
331 UInt_t singleTriggerResult = CheckTriggerClass(aEsd, triggerClass, triggerLogic);
332
0c6c629b 333 if (singleTriggerResult)
296dd262 334 {
335 triggerAnalysis->FillHistograms(aEsd);
336
85c71ba7 337 Bool_t isBin0 = kFALSE;
338 if (fBin0CallBack != "") {
339 AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
340 if (!mgr) {
341 AliError("Cannot get the analysis manager");
342 }
343 else {
344 isBin0 = ((AliAnalysisTaskSE*)mgr->GetTask(fBin0CallBack.Data()))->IsEventInBinZero();
345 }
ff097e3f 346 } else if (fBin0CallBackPointer) {
347 isBin0 = (*fBin0CallBackPointer)(aEsd);
85c71ba7 348 }
ff097e3f 349
c2ba5a61 350 // hardware trigger
cc9d9320 351 Int_t fastORHW = triggerAnalysis->SPDFiredChips(aEsd, 1); // SPD number of chips from trigger bits (!)
c2ba5a61 352 Int_t fastORHWL1 = triggerAnalysis->SPDFiredChips(aEsd, 1, kFALSE, 2); // SPD number of chips from trigger bits in second layer (!)
85c71ba7 353 Bool_t v0AHW = fSkipV0 ? 0 :(triggerAnalysis->V0Trigger(aEsd, AliTriggerAnalysis::kASide, kTRUE) == AliTriggerAnalysis::kV0BB);// should replay hw trigger
354 Bool_t v0CHW = fSkipV0 ? 0 :(triggerAnalysis->V0Trigger(aEsd, AliTriggerAnalysis::kCSide, kTRUE) == AliTriggerAnalysis::kV0BB);// should replay hw trigger
c2ba5a61 355
cc9d9320 356 // offline trigger
357 Int_t fastOROffline = triggerAnalysis->SPDFiredChips(aEsd, 0); // SPD number of chips from clusters (!)
c2ba5a61 358 Int_t fastOROfflineL1 = triggerAnalysis->SPDFiredChips(aEsd, 0, kFALSE, 2); // SPD number of chips from clusters in second layer (!)
359 Bool_t v0A = fSkipV0 ? 0 :triggerAnalysis->IsOfflineTriggerFired(aEsd, AliTriggerAnalysis::kV0A);
360 Bool_t v0C = fSkipV0 ? 0 :triggerAnalysis->IsOfflineTriggerFired(aEsd, AliTriggerAnalysis::kV0C);
85c71ba7 361 Bool_t v0ABG = fSkipV0 ? 0 :triggerAnalysis->IsOfflineTriggerFired(aEsd, AliTriggerAnalysis::kV0ABG);
362 Bool_t v0CBG = fSkipV0 ? 0 :triggerAnalysis->IsOfflineTriggerFired(aEsd, AliTriggerAnalysis::kV0CBG);
f4ca8f20 363 Bool_t v0BG = v0ABG || v0CBG;
85c71ba7 364
365 // fmd
366 Bool_t fmdA = triggerAnalysis->IsOfflineTriggerFired(aEsd, AliTriggerAnalysis::kFMDA);
367 Bool_t fmdC = triggerAnalysis->IsOfflineTriggerFired(aEsd, AliTriggerAnalysis::kFMDC);
368 Bool_t fmd = fmdA || fmdC;
61899827 369
85c71ba7 370 // SSD
371 Int_t ssdClusters = triggerAnalysis->SSDClusters(aEsd);
c2ba5a61 372
373 // ZDC
374 Bool_t zdcA = triggerAnalysis->IsOfflineTriggerFired(aEsd, AliTriggerAnalysis::kZDCA);
375 Bool_t zdcC = triggerAnalysis->IsOfflineTriggerFired(aEsd, AliTriggerAnalysis::kZDCC);
376
85c71ba7 377
378 // Some "macros"
379 Bool_t mb1 = (fastOROffline > 0 || v0A || v0C) && (!v0BG);
380 Bool_t mb1prime = (fastOROffline > 1 || (fastOROffline > 0 && (v0A || v0C)) || (v0A && v0C) ) && (!v0BG);
381
382 // Background rejection
8663e267 383 Bool_t bgID = kFALSE;
384 if (fBackgroundIdentification)
385 bgID = ! fBackgroundIdentification->IsSelected(const_cast<AliESDEvent*> (aEsd));
386
c2ba5a61 387 /*Int_t ntrig = fastOROffline; // any 2 hits
85c71ba7 388 if(v0A) ntrig += 1;
389 if(v0C) ntrig += 1; //v0C alone is enough
390 if(fmd) ntrig += 1;
c2ba5a61 391 if(ssdClusters>1) ntrig += 1;*/
85c71ba7 392
c2ba5a61 393 // replay hardware trigger (should only remove events for MC)
394
395 Bool_t hwTrig = kFALSE;
396 switch (triggerLogic)
397 {
398 case kCINT1: hwTrig = fastORHW > 0 || v0AHW || v0CHW; break;
399 case kCMBS2A: hwTrig = fastORHWL1 > 1 && v0AHW; break;
400 case kCMBS2C: hwTrig = fastORHWL1 > 1 && v0CHW; break;
401 case kCMBAC: hwTrig = v0AHW && v0CHW; break;
402 default: AliFatal(Form("Undefined trigger logic %d", triggerLogic)); break;
403 }
85c71ba7 404
73cc8654 405 // Fill trigger pattern histo
406 Int_t tpatt = 0;
407 if (fastORHW>0) tpatt+=1;
408 if (v0AHW) tpatt+=2;
409 if (v0CHW) tpatt+=4;
410 fHistTriggerPattern->Fill( tpatt );
411
85c71ba7 412 // fill statistics and return decision
413 const Int_t nHistStat = 2;
414 for(Int_t iHistStat = 0; iHistStat < nHistStat; iHistStat++){
415 if (iHistStat == kStatIdxBin0 && !isBin0) continue; // skip the filling of bin0 stats if the event is not in the bin0
416
417 fHistStatistics[iHistStat]->Fill(kStatTriggerClass, i);
418
85c71ba7 419 // We fill the rest only if hw trigger is ok
420 if (!hwTrig)
421 {
422 AliDebug(AliLog::kDebug, "Rejecting event because hardware trigger is not fired");
423 continue;
424 } else {
425 fHistStatistics[iHistStat]->Fill(kStatHWTrig, i);
426 }
733f0542 427
85c71ba7 428
429 // v0 BG stats
430 if (v0ABG)
431 fHistStatistics[iHistStat]->Fill(kStatV0ABG, i);
432 if (v0CBG)
433 fHistStatistics[iHistStat]->Fill(kStatV0CBG, i);
434
435 // We fill the rest only if mb1 && ! v0BG
436 if (mb1)
437 fHistStatistics[iHistStat]->Fill(kStatMB1, i);
438 else continue;
439
440 if (mb1prime)
441 fHistStatistics[iHistStat]->Fill(kStatMB1Prime, i);
442
443 if (fmd)
444 fHistStatistics[iHistStat]->Fill(kStatFMD, i);
445
c2ba5a61 446 //if(ntrig >= 2 && !v0BG)
447 // fHistStatistics[iHistStat]->Fill(kStatAny2Hits, i);
85c71ba7 448
449 if (fastOROffline > 0)
450 fHistStatistics[iHistStat]->Fill(kStatFO1, i);
451 if (fastOROffline > 1)
452 fHistStatistics[iHistStat]->Fill(kStatFO2, i);
c2ba5a61 453 if (fastOROfflineL1 > 1)
454 fHistStatistics[iHistStat]->Fill(kStatFO2L1, i);
cc9d9320 455
85c71ba7 456 if (v0A)
457 fHistStatistics[iHistStat]->Fill(kStatV0A, i);
458 if (v0C)
459 fHistStatistics[iHistStat]->Fill(kStatV0C, i);
c2ba5a61 460
461 if (zdcA)
462 fHistStatistics[iHistStat]->Fill(kStatZDCA, i);
463 if (zdcC)
464 fHistStatistics[iHistStat]->Fill(kStatZDCC, i);
465 if (zdcA && zdcC)
466 fHistStatistics[iHistStat]->Fill(kStatZDCAC, i);
85c71ba7 467
468 // if (fastOROffline > 1 && !v0BG)
469 // fHistStatistics[iHistStat]->Fill(kStatFO2NoBG, i);
470
c2ba5a61 471 //if (fastOROffline > 0 && (v0A || v0C) && !v0BG)
472 // fHistStatistics[iHistStat]->Fill(kStatFO1AndV0, i);
296dd262 473
85c71ba7 474 if (v0A && v0C && !v0BG && !bgID)
475 fHistStatistics[iHistStat]->Fill(kStatV0, i);
476
c2ba5a61 477 Bool_t offlineAccepted = kFALSE;
478
479 switch (triggerLogic)
480 {
481 case kCINT1: offlineAccepted = mb1; break;
482 case kCMBS2A: offlineAccepted = fastOROfflineL1 > 1 && v0A; break;
483 case kCMBS2C: offlineAccepted = fastOROfflineL1 > 1 && v0C; break;
484 case kCMBAC: offlineAccepted = v0A && v0C; break;
485 default: AliFatal(Form("Undefined trigger logic %d", triggerLogic)); break;
486 }
85c71ba7 487
c2ba5a61 488 if ( offlineAccepted )
85c71ba7 489 {
490 if (!v0BG || fSkipV0)
491 {
492 if (!v0BG) fHistStatistics[iHistStat]->Fill(kStatOffline, i);
61899827 493
85c71ba7 494 if (fBackgroundIdentification && bgID)
495 {
496 AliDebug(AliLog::kDebug, "Rejecting event because of background identification");
497 fHistStatistics[iHistStat]->Fill(kStatBG, i);
498 }
499 else
500 {
c2ba5a61 501 AliDebug(AliLog::kDebug, Form("Accepted event for histograms with trigger logic %d", triggerLogic));
296dd262 502
85c71ba7 503 fHistStatistics[iHistStat]->Fill(kStatAccepted, i);
504 if(iHistStat == kStatIdxAll) fHistBunchCrossing->Fill(aEsd->GetBunchCrossNumber(), i); // Fill only for all (avoid double counting)
505 if((i < fCollTrigClasses.GetEntries() || fSkipTriggerClassSelection) && (iHistStat==kStatIdxAll))
0c6c629b 506 accept |= singleTriggerResult; // only set for "all" (should not really matter)
85c71ba7 507 }
508 }
509 else
510 AliDebug(AliLog::kDebug, "Rejecting event because of V0 BG flag");
511 }
512 else
c2ba5a61 513 AliDebug(AliLog::kDebug, Form("Rejecting event because trigger logic %d is not fulfilled", triggerLogic));
296dd262 514 }
296dd262 515 }
516 }
517
518 if (accept)
0c6c629b 519 AliDebug(AliLog::kDebug, Form("Accepted event as collision candidate with bit mask %d", accept));
61899827 520
296dd262 521 return accept;
61899827 522}
a2ce3799 523
85c71ba7 524Int_t AliPhysicsSelection::GetTriggerScheme(UInt_t runNumber) const
a2ce3799 525{
526 // returns the current trigger scheme (classes that are accepted/rejected)
527
c2ba5a61 528 Bool_t pp = kTRUE;
529 if (runNumber >= 136849) // last pp run 2010
530 pp = kFALSE;
531
a2ce3799 532 if (fMC)
c2ba5a61 533 {
534 if (pp)
535 return 0; // pp
536 else
537 return -1; // HI
538 }
a2ce3799 539
540 // TODO dependent on run number
4b7e8f3b 541
542 switch (runNumber)
543 {
544 // CSMBB triggers
545 case 104044:
546 case 105054:
547 case 105057:
548 return 2;
549 }
550
c2ba5a61 551 // defaults
552 if (pp)
553 return 1;
554 else
555 return 3;
a2ce3799 556}
85c71ba7 557
558const char * AliPhysicsSelection::GetFillingScheme(UInt_t runNumber) {
559
560 if(fMC) return "MC";
561
562 if (runNumber >= 104065 && runNumber <= 104160) {
563 return "4x4a";
564 }
565 else if (runNumber >= 104315 && runNumber <= 104321) {
566 return "4x4a*";
567 }
568 else if (runNumber >= 104792 && runNumber <= 104803) {
569 return "4x4b";
570 }
571 else if (runNumber >= 104824 && runNumber <= 104892) {
572 return "4x4c";
573 }
574 else if (runNumber == 105143 || runNumber == 105160) {
575 return "16x16a";
576 }
577 else if (runNumber >= 105256 && runNumber <= 105268) {
578 return "4x4c";
5f337f3d 579 }
580 else if (runNumber >= 114786 && runNumber <= 116684) {
ca234d7a 581 return "Single_2b_1_1_1";
85c71ba7 582 }
5f337f3d 583 else if (runNumber >= 117048 && runNumber <= 117120) {
584 return "Single_3b_2_2_2";
585 }
586 else if (runNumber >= 117220 && runNumber <= 119163) {
587 return "Single_2b_1_1_1";
588 }
589 else if (runNumber >= 119837 && runNumber <= 119862) {
590 return "Single_4b_2_2_2";
591 }
592 else if (runNumber >= 119902 && runNumber <= 120691) {
593 return "Single_6b_3_3_3";
594 }
595 else if (runNumber >= 120741 && runNumber <= 122375) {
596 return "Single_13b_8_8_8";
597 }
141265a2 598 else if (runNumber >= 130148 && runNumber <= 130375) {
599 return "125n_48b_36_16_36";
600 }
8fc17e87 601 else if (runNumber >= 130601 && runNumber <= 130640) {
602 return "1000ns_50b_35_14_35";
603 }
85c71ba7 604 else {
605 AliError(Form("Unknown filling scheme (run %d)", runNumber));
606 }
607
608 return "Unknown";
609}
610
85c71ba7 611const char * AliPhysicsSelection::GetBXIDs(UInt_t runNumber, const char * trigger) {
612
613 if (!fUseBXNumbers || fMC) return "";
614
615 if (runNumber >= 104065 && runNumber <= 104160) {
616 if (!strcmp("CINT1B-ABCE-NOPF-ALL",trigger)) return " #2128 #3019";
617 else if(!strcmp("CINT1A-ABCE-NOPF-ALL",trigger)) return " #346 #3465";
618 else if(!strcmp("CINT1C-ABCE-NOPF-ALL",trigger)) return " #1234 #1680";
619 else if(!strcmp("CINT1-E-NOPF-ALL",trigger)) return " #790";
e0bf2892 620 // else AliError(Form("Unknown trigger: %s", trigger));
85c71ba7 621 }
622 else if (runNumber >= 104315 && runNumber <= 104321) {
623 if (!strcmp("CINT1B-ABCE-NOPF-ALL",trigger)) return " #2000 #2891";
624 else if(!strcmp("CINT1A-ABCE-NOPF-ALL",trigger)) return " #218 #3337";
625 else if(!strcmp("CINT1C-ABCE-NOPF-ALL",trigger)) return " #1106 #1552";
626 else if(!strcmp("CINT1-E-NOPF-ALL",trigger)) return " #790";
e0bf2892 627 // else AliError(Form("Unknown trigger: %s", trigger));
85c71ba7 628 }
629 else if (runNumber >= 104792 && runNumber <= 104803) {
630 if (!strcmp("CINT1B-ABCE-NOPF-ALL",trigger)) return " #2228 #3119";
631 else if(!strcmp("CINT1A-ABCE-NOPF-ALL",trigger)) return " #2554 #446";
632 else if(!strcmp("CINT1C-ABCE-NOPF-ALL",trigger)) return " #1334 #769";
633 else if(!strcmp("CINT1-E-NOPF-ALL",trigger)) return " #790";
e0bf2892 634 // else AliError(Form("Unknown trigger: %s", trigger));
85c71ba7 635 }
636 else if (runNumber >= 104824 && runNumber <= 104892) {
637 if (!strcmp("CINT1B-ABCE-NOPF-ALL",trigger)) return " #3119 #769";
638 else if(!strcmp("CINT1A-ABCE-NOPF-ALL",trigger)) return " #2554 #446";
639 else if(!strcmp("CINT1C-ABCE-NOPF-ALL",trigger)) return " #1334 #2228";
640 else if(!strcmp("CINT1-E-NOPF-ALL",trigger)) return " #790";
e0bf2892 641 // else AliError(Form("Unknown trigger: %s", trigger));
85c71ba7 642 }
643 else if (runNumber == 105143 || runNumber == 105160) {
85c71ba7 644 if (!strcmp("CINT1B-ABCE-NOPF-ALL",trigger)) return " #1337 #1418 #2228 #2309 #3119 #3200 #446 #527";
645 else if(!strcmp("CINT1A-ABCE-NOPF-ALL",trigger)) return " #1580 #1742 #1904 #2066 #2630 #2792 #2954 #3362";
646 else if(!strcmp("CINT1C-ABCE-NOPF-ALL",trigger)) return " #845 #1007 #1169 #1577 #3359 #3521 #119 #281 ";
647 else if(!strcmp("CINT1-E-NOPF-ALL",trigger)) return " #790";
e0bf2892 648 // else AliError(Form("Unknown trigger: %s", trigger));
85c71ba7 649 }
650 else if (runNumber >= 105256 && runNumber <= 105268) {
651 if (!strcmp("CINT1B-ABCE-NOPF-ALL",trigger)) return " #3019 #669";
652 else if(!strcmp("CINT1A-ABCE-NOPF-ALL",trigger)) return " #2454 #346";
653 else if(!strcmp("CINT1C-ABCE-NOPF-ALL",trigger)) return " #1234 #2128";
654 else if(!strcmp("CINT1-E-NOPF-ALL",trigger)) return " #790";
e0bf2892 655 // else AliError(Form("Unknown trigger: %s", trigger));
73cc8654 656 } else if (runNumber >= 114786 && runNumber <= 116684) { // 7 TeV 2010, assume always the same filling scheme
ca234d7a 657 if (!strcmp("CINT1B-ABCE-NOPF-ALL",trigger)) return " #346";
658 else if(!strcmp("CINT1A-ABCE-NOPF-ALL",trigger)) return " #2131";
659 else if(!strcmp("CINT1C-ABCE-NOPF-ALL",trigger)) return " #3019";
660 else if(!strcmp("CINT1-E-NOPF-ALL",trigger)) return " #1238";
e0bf2892 661 // else AliError(Form("Unknown trigger: %s", trigger));
85c71ba7 662 }
5f337f3d 663 else if (runNumber >= 117048 && runNumber <= 117120) {
664 // return "Single_3b_2_2_2";
665 if (!strcmp("CINT1B-ABCE-NOPF-ALL",trigger)) return " #346 #1240 ";
666 else if (!strcmp("CINT1A-ABCE-NOPF-ALL",trigger)) return " #2131 ";
667 else if (!strcmp("CINT1C-ABCE-NOPF-ALL",trigger)) return " #3019 ";
668 else if (!strcmp("CINT1-E-NOPF-ALL",trigger)) return " #1238";
e0bf2892 669 // else AliError(Form("Unknown trigger: %s", trigger));
5f337f3d 670
671 }
c2ba5a61 672 else if ((runNumber >= 117220 && runNumber <= 118555) || (runNumber >= 118784 && runNumber <= 119163))
673 {
5f337f3d 674 // return "Single_2b_1_1_1";
675 if (!strcmp("CINT1B-ABCE-NOPF-ALL",trigger)) return " #346 ";
676 else if (!strcmp("CINT1A-ABCE-NOPF-ALL",trigger)) return " #2131 ";
677 else if (!strcmp("CINT1C-ABCE-NOPF-ALL",trigger)) return " #3019 ";
678 else if (!strcmp("CINT1-E-NOPF-ALL",trigger)) return " #1238 ";
e0bf2892 679 // else AliError(Form("Unknown trigger: %s", trigger));
5f337f3d 680 }
c2ba5a61 681 else if (runNumber >= 118556 && runNumber <= 118783) {
682 // return "Single_2b_1_1_1";
683 // same as previous but was misaligned by 1 BX in fill 1069
684 if (!strcmp("CINT1B-ABCE-NOPF-ALL",trigger)) return " #345 ";
685 else if (!strcmp("CINT1A-ABCE-NOPF-ALL",trigger)) return " #2130 ";
686 else if (!strcmp("CINT1C-ABCE-NOPF-ALL",trigger)) return " #3018 ";
687 else if (!strcmp("CINT1-E-NOPF-ALL",trigger)) return " #1238 ";
688 // else AliError(Form("Unknown trigger: %s", trigger));
689 }
5f337f3d 690 else if (runNumber >= 119837 && runNumber <= 119862) {
691 // return "Single_4b_2_2_2";
692 if (!strcmp("CINT1B-ABCE-NOPF-ALL",trigger)) return " #669 #3019 ";
693 else if (!strcmp("CINT1A-ABCE-NOPF-ALL",trigger)) return " #346 #2454 ";
694 else if (!strcmp("CINT1C-ABCE-NOPF-ALL",trigger)) return " #1234 #2128 ";
695 else if (!strcmp("CINT1-E-NOPF-ALL",trigger)) return " #1681 #3463";
e0bf2892 696 // else AliError(Form("Unknown trigger: %s", trigger));
5f337f3d 697
698 }
699 else if (runNumber >= 119902 && runNumber <= 120691) {
700 // return "Single_6b_3_3_3";
701 if (!strcmp("CINT1B-ABCE-NOPF-ALL",trigger)) return " #346 #546 #746 ";
702 else if (!strcmp("CINT1A-ABCE-NOPF-ALL",trigger)) return " #2131 #2331 #2531 ";
703 else if (!strcmp("CINT1C-ABCE-NOPF-ALL",trigger)) return " #3019 #3219 #3419 ";
704 else if (!strcmp("CINT1-E-NOPF-ALL",trigger)) return " #1296 #1670";
e0bf2892 705 // else AliError(Form("Unknown trigger: %s", trigger));
5f337f3d 706 }
707 else if (runNumber >= 120741 && runNumber <= 122375) {
708 // return "Single_13b_8_8_8";
709 if (!strcmp("CINT1B-ABCE-NOPF-ALL",trigger)) return " #346 #446 #546 #646 #1240 #1340 #1440 #1540 ";
710 else if (!strcmp("CINT1A-ABCE-NOPF-ALL",trigger)) return " #946 #2131 #2231 #2331 #2431 ";
711 else if (!strcmp("CINT1C-ABCE-NOPF-ALL",trigger)) return " #3019 #3119 #3219 #3319 #3519 ";
712 else if (!strcmp("CINT1-E-NOPF-ALL",trigger)) return " #1835 #2726";
e0bf2892 713 // else AliError(Form("Unknown trigger: %s", trigger));
5f337f3d 714
bc18a83e 715 }
716 else if (runNumber >= 130148 && runNumber <= 130375) {
717 TString triggerString = trigger;
718 static TString returnString = " ";
719 returnString = "";
720 if (triggerString.Contains("B")) returnString += " #346 #396 #446 #496 #546 #596 #646 #696 #1240 #1290 #1340 #1390 #1440 #1490 #1540 #1590 ";
721 if (triggerString.Contains("A")) returnString += " #755 #805 #855 #905 #955 #1005 #1799 #1849 #1899 #2131 #2181 #2231 #2281 #2331 #2381 #2431 #2481 #2531 #2581 #2631 #2846 #3016 #3066 #3116 #3166 #3216 #3266 #3316 #3366 #3425 #3475 #3525 ";
722 if (triggerString.Contains("C")) returnString += " #3019 #3069 #3119 #3169 #3219 #3269 #3319 #3369 #14 #64 #114 #746 #796 #846 #908 #958 #1008 #1640 #1690 #1740 #2055 #2125 #2175 #2225 #2275 #2325 #2375 #2425 #2475 #2534 #2584 #2634 ";
723 // Printf("0x%x",returnString.Data());
724 // Printf("%s",returnString.Data());
725 return returnString.Data();
8fc17e87 726 }
727 else if (runNumber >= 130601 && runNumber <= 130640) {
728 TString triggerString = trigger;
729 static TString returnString = " ";
730 returnString = "";
731 if (triggerString.Contains("B")) returnString += " #346 #386 #426 #466 #506 #546 #586 #1240 #1280 #1320 #1360 #1400 #1440 #1480 ";
732 if (triggerString.Contains("A")) returnString += " #626 #666 #706 #746 #786 #826 #866 #1520 #1560 #1600 #1640 #1680 #1720 #1760 #2076 #2131 #2171 #2211 #2251 #2291 #2331 #2371 #2414 #2454 #2494 #2534 #2574 #2614 #2654 #2694 #2734 #2774 #2814 "; //#2854 #2894 #2934 not present in this run
733 if (triggerString.Contains("C")) returnString += " #3019 #3059 #3099 #3139 #3179 #3219 #3259 #3299 #3339 #3379 #3419 #3459 #3499 #3539 #115 #629 #669 #709 #749 #789 #829 #869 #909 #949 #989 #1029 #1069 #1109 #1149 #1523 #1563 #1603 #1643 "; //#1683 #1723 #1763 not present in this run
734 return returnString.Data();
5f337f3d 735 }
736
85c71ba7 737 else {
8fc17e87 738 AliWarning(Form("Unknown run %d, using all BXs!",runNumber));
85c71ba7 739 }
740
741 return "";
742}
61899827 743
85c71ba7 744Bool_t AliPhysicsSelection::Initialize(Int_t runNumber)
61899827 745{
746 // initializes the object for the given run
747 // TODO having the run number here and parameters hardcoded is clearly temporary, a way needs to be found to have a CDB-like configuration also for analysis
748
29e8486e 749 Bool_t oldStatus = TH1::AddDirectoryStatus();
750 TH1::AddDirectory(kFALSE);
751
85c71ba7 752 if(!fBin0CallBack)
753 AliError("Bin0 Callback not set: will not fill the statistics for the bin 0");
754
755 if (fMC) {
0c6c629b 756 // override BX and bg options in case of MC
85c71ba7 757 fComputeBG = kFALSE;
758 fUseBXNumbers = kFALSE;
759 }
760
a2ce3799 761 Int_t triggerScheme = GetTriggerScheme(runNumber);
91bea6e7 762 if (!fUsingCustomClasses && fCurrentRun != -1 && triggerScheme != GetTriggerScheme(fCurrentRun))
a2ce3799 763 AliFatal("Processing several runs with different trigger schemes is not supported");
296dd262 764
85c71ba7 765 if(fComputeBG && fCurrentRun != -1 && fCurrentRun != runNumber)
766 AliFatal("Cannot process several runs because BG computation is requested");
767
768 if(fComputeBG && !fUseBXNumbers)
e0bf2892 769 AliFatal("Cannot compute BG if BX numbers are not used");
85c71ba7 770
771 if(fUseBXNumbers && fFillingScheme != "" && fFillingScheme != GetFillingScheme(runNumber))
772 AliFatal("Cannot process runs with different filling scheme if usage of BX numbers is requested");
773
85c71ba7 774 fFillingScheme = GetFillingScheme(runNumber);
775
61899827 776 AliInfo(Form("Initializing for run %d", runNumber));
777
758941d4 778 // initialize first time?
a2ce3799 779 if (fCurrentRun == -1)
a2ce3799 780 {
91bea6e7 781 if (fUsingCustomClasses) {
782 AliInfo("Using user-provided trigger classes");
783 } else {
784 switch (triggerScheme)
785 {
c2ba5a61 786 case -1:
787 // MC Heavy-Ion
788
789 fCollTrigClasses.Add(new TObjString(Form("&%u *%d", (UInt_t) AliVEvent::kMB, (Int_t) kCMBS2A)));
790 fCollTrigClasses.Add(new TObjString(Form("&%u *%d", (UInt_t) AliVEvent::kMB, (Int_t) kCMBS2C)));
791 fCollTrigClasses.Add(new TObjString(Form("&%u *%d", (UInt_t) AliVEvent::kMB, (Int_t) kCMBAC)));
792 break;
793
a2ce3799 794 case 0:
c2ba5a61 795 // MC Proton-Proton
796
797 fCollTrigClasses.Add(new TObjString(Form("&%u *%d", (UInt_t) AliVEvent::kMB, (Int_t) kCINT1)));
a2ce3799 798 break;
91bea6e7 799
a2ce3799 800 case 1:
c2ba5a61 801 // Proton-Proton
802
4a6e052f 803 // trigger classes used before August 2010
1209509c 804 fCollTrigClasses.Add(new TObjString(Form("%s%s &%u","+CINT1B-ABCE-NOPF-ALL", GetBXIDs(runNumber,"CINT1B-ABCE-NOPF-ALL"), (UInt_t) AliVEvent::kMB)));
805 fBGTrigClasses.Add (new TObjString(Form("%s%s &%u","+CINT1A-ABCE-NOPF-ALL", GetBXIDs(runNumber,"CINT1A-ABCE-NOPF-ALL"), (UInt_t) AliVEvent::kMB)));
806 fBGTrigClasses.Add (new TObjString(Form("%s%s &%u","+CINT1C-ABCE-NOPF-ALL", GetBXIDs(runNumber,"CINT1C-ABCE-NOPF-ALL"), (UInt_t) AliVEvent::kMB)));
807 fBGTrigClasses.Add (new TObjString(Form("%s%s &%u","+CINT1-E-NOPF-ALL", GetBXIDs(runNumber,"CINT1-E-NOPF-ALL"), (UInt_t) AliVEvent::kMB)));
0c6c629b 808
809 // Muon trigger have the same BXIDs of the corresponding CINT triggers
1209509c 810 fCollTrigClasses.Add(new TObjString(Form("%s%s &%u","+CMUS1B-ABCE-NOPF-MUON", GetBXIDs(runNumber,"CINT1B-ABCE-NOPF-ALL"), (UInt_t) AliVEvent::kMUON)));
811 fBGTrigClasses.Add (new TObjString(Form("%s%s &%u","+CMUS1A-ABCE-NOPF-MUON", GetBXIDs(runNumber,"CINT1A-ABCE-NOPF-ALL"), (UInt_t) AliVEvent::kMUON)));
812 fBGTrigClasses.Add (new TObjString(Form("%s%s &%u","+CMUS1C-ABCE-NOPF-MUON", GetBXIDs(runNumber,"CINT1C-ABCE-NOPF-ALL"), (UInt_t) AliVEvent::kMUON)));
813 fBGTrigClasses.Add (new TObjString(Form("%s%s &%u","+CMUS1-E-NOPF-MUON" , GetBXIDs(runNumber,"CINT1-E-NOPF-ALL"), (UInt_t) AliVEvent::kMUON)));
4a6e052f 814
815 // triggers classes used from August 2010
816 // MB
bc18a83e 817 fCollTrigClasses.Add(new TObjString(Form("%s%s &%u", "+CINT1-B-NOPF-ALLNOTRD" , GetBXIDs(runNumber, "B"), (UInt_t) AliVEvent::kMB)));
818 fBGTrigClasses.Add (new TObjString(Form("%s%s &%u", "+CINT1-AC-NOPF-ALLNOTRD", GetBXIDs(runNumber, "AC"), (UInt_t) AliVEvent::kMB)));
819 fBGTrigClasses.Add (new TObjString(Form("%s%s &%u", "+CINT1-E-NOPF-ALLNOTRD" , GetBXIDs(runNumber, "E"), (UInt_t) AliVEvent::kMB)));
820
821 // MUON
822 fCollTrigClasses.Add(new TObjString(Form("%s%s &%u", "+CMUS1-B-NOPF-ALLNOTRD" , GetBXIDs(runNumber, "B"), (UInt_t) AliVEvent::kMUON)));
823 fBGTrigClasses.Add (new TObjString(Form("%s%s &%u", "+CMUS1-AC-NOPF-ALLNOTRD", GetBXIDs(runNumber, "AC"), (UInt_t) AliVEvent::kMUON)));
824 fBGTrigClasses.Add (new TObjString(Form("%s%s &%u", "+CMUS1-E-NOPF-ALLNOTRD" , GetBXIDs(runNumber, "E"), (UInt_t) AliVEvent::kMUON)));
825
826 // High Multiplicity
827 fCollTrigClasses.Add(new TObjString(Form("%s%s &%u", "+CSH1-B-NOPF-ALLNOTRD" , GetBXIDs(runNumber, "B"), (UInt_t) AliVEvent::kHighMult)));
828 fBGTrigClasses.Add (new TObjString(Form("%s%s &%u", "+CSH1-AC-NOPF-ALLNOTRD" , GetBXIDs(runNumber, "AC"), (UInt_t) AliVEvent::kHighMult)));
829 fBGTrigClasses.Add (new TObjString(Form("%s%s &%u", "+CSH1-E-NOPF-ALLNOTRD" , GetBXIDs(runNumber, "E"), (UInt_t) AliVEvent::kHighMult)));
4a6e052f 830
141265a2 831 // WARNING: IF YOU ADD MORE TRIGGER CLASSES, PLEASE CHECK THAT THE REGULAR EXPRESSION IN GetStatRow IS STILL VALID
832
a2ce3799 833 break;
834
4b7e8f3b 835 case 2:
c2ba5a61 836 // Proton-Proton
837
1209509c 838 fCollTrigClasses.Add(new TObjString(Form("+CSMBB-ABCE-NOPF-ALL &%u", (UInt_t) AliVEvent::kMB)));
839 fBGTrigClasses.Add(new TObjString(Form("+CSMBA-ABCE-NOPF-ALL -CSMBB-ABCE-NOPF-ALL &%u", (UInt_t) AliVEvent::kMB)));
840 fBGTrigClasses.Add(new TObjString(Form("+CSMBC-ABCE-NOPF-ALL -CSMBB-ABCE-NOPF-ALL &%u", (UInt_t) AliVEvent::kMB)));
0c6c629b 841 break;
842
843 case 3:
c2ba5a61 844 // Heavy-ion
845
846 fCollTrigClasses.Add(new TObjString(Form("+CMBS2A-B-NOPF-ALL &%u *%d", (UInt_t) AliVEvent::kMB, (Int_t) kCMBS2A)));
847 fBGTrigClasses.Add (new TObjString(Form("+CMBS2A-A-NOPF-ALL &%u *%d", (UInt_t) AliVEvent::kMB, (Int_t) kCMBS2A)));
848 fBGTrigClasses.Add (new TObjString(Form("+CMBS2A-C-NOPF-ALL &%u *%d", (UInt_t) AliVEvent::kMB, (Int_t) kCMBS2A)));
849 fBGTrigClasses.Add (new TObjString(Form("+CMBS2A-E-NOPF-ALL &%u *%d", (UInt_t) AliVEvent::kMB, (Int_t) kCMBS2A)));
850
851 fCollTrigClasses.Add(new TObjString(Form("+CMBS2C-B-NOPF-ALL &%u *%d", (UInt_t) AliVEvent::kMB, (Int_t) kCMBS2C)));
852 fBGTrigClasses.Add (new TObjString(Form("+CMBS2C-A-NOPF-ALL &%u *%d", (UInt_t) AliVEvent::kMB, (Int_t) kCMBS2C)));
853 fBGTrigClasses.Add (new TObjString(Form("+CMBS2C-C-NOPF-ALL &%u *%d", (UInt_t) AliVEvent::kMB, (Int_t) kCMBS2C)));
854 fBGTrigClasses.Add (new TObjString(Form("+CMBS2C-E-NOPF-ALL &%u *%d", (UInt_t) AliVEvent::kMB, (Int_t) kCMBS2C)));
855
856 fCollTrigClasses.Add(new TObjString(Form("+CMBAC-B-NOPF-ALL &%u *%d", (UInt_t) AliVEvent::kMB, (Int_t) kCMBAC)));
857 fBGTrigClasses.Add (new TObjString(Form("+CMBAC-A-NOPF-ALL &%u *%d", (UInt_t) AliVEvent::kMB, (Int_t) kCMBAC)));
858 fBGTrigClasses.Add (new TObjString(Form("+CMBAC-C-NOPF-ALL &%u *%d", (UInt_t) AliVEvent::kMB, (Int_t) kCMBAC)));
859 fBGTrigClasses.Add (new TObjString(Form("+CMBAC-E-NOPF-ALL &%u *%d", (UInt_t) AliVEvent::kMB, (Int_t) kCMBAC)));
860
4b7e8f3b 861 break;
91bea6e7 862
a2ce3799 863 default:
864 AliFatal(Form("Unsupported trigger scheme %d", triggerScheme));
91bea6e7 865 }
a2ce3799 866 }
91bea6e7 867
a2ce3799 868 Int_t count = fCollTrigClasses.GetEntries() + fBGTrigClasses.GetEntries();
869
870 for (Int_t i=0; i<count; i++)
871 {
872 AliTriggerAnalysis* triggerAnalysis = new AliTriggerAnalysis;
873 triggerAnalysis->SetAnalyzeMC(fMC);
874 triggerAnalysis->EnableHistograms();
875 triggerAnalysis->SetSPDGFOThreshhold(1);
4011b280 876 triggerAnalysis->SetDoFMD(kFALSE);
a2ce3799 877 fTriggerAnalysis.Add(triggerAnalysis);
878 }
879
85c71ba7 880 // TODO: shall I really delete this?
881 if (fHistStatistics[0])
882 delete fHistStatistics[0];
883 if (fHistStatistics[1])
884 delete fHistStatistics[1];
296dd262 885
85c71ba7 886 fHistStatistics[kStatIdxBin0] = BookHistStatistics("_Bin0");
887 fHistStatistics[kStatIdxAll] = BookHistStatistics("");
a2ce3799 888
889 if (fHistBunchCrossing)
890 delete fHistBunchCrossing;
891
892 fHistBunchCrossing = new TH2F("fHistBunchCrossing", "fHistBunchCrossing;bunch crossing number;", 4000, -0.5, 3999.5, count, -0.5, -0.5 + count);
73cc8654 893
894 if (fHistTriggerPattern)
895 delete fHistTriggerPattern;
896
897 const int ntrig=3;
85c71ba7 898 Int_t n = 1;
73cc8654 899 const Int_t nbinTrig = TMath::Nint(TMath::Power(2,ntrig));
900
901 fHistTriggerPattern = new TH1F("fHistTriggerPattern", "Trigger pattern: FO + 2*v0A + 4*v0C",
902 nbinTrig, -0.5, nbinTrig-0.5);
903 fHistTriggerPattern->GetXaxis()->SetBinLabel(1,"NO TRIG");
904 fHistTriggerPattern->GetXaxis()->SetBinLabel(2,"FO");
905 fHistTriggerPattern->GetXaxis()->SetBinLabel(3,"v0A");
906 fHistTriggerPattern->GetXaxis()->SetBinLabel(4,"FO & v0A");
907 fHistTriggerPattern->GetXaxis()->SetBinLabel(5,"v0C");
908 fHistTriggerPattern->GetXaxis()->SetBinLabel(6,"FO & v0C");
909 fHistTriggerPattern->GetXaxis()->SetBinLabel(7,"v0A & v0C");
910 fHistTriggerPattern->GetXaxis()->SetBinLabel(8,"FO & v0A & v0C");
911
912
913 n = 1;
a2ce3799 914 for (Int_t i=0; i < fCollTrigClasses.GetEntries(); i++)
915 {
a2ce3799 916 fHistBunchCrossing->GetYaxis()->SetBinLabel(n, ((TObjString*) fCollTrigClasses.At(i))->String());
917 n++;
918 }
919 for (Int_t i=0; i < fBGTrigClasses.GetEntries(); i++)
920 {
a2ce3799 921 fHistBunchCrossing->GetYaxis()->SetBinLabel(n, ((TObjString*) fBGTrigClasses.At(i))->String());
922 n++;
923 }
85c71ba7 924
925
926
a2ce3799 927 }
928
929 Int_t count = fCollTrigClasses.GetEntries() + fBGTrigClasses.GetEntries();
296dd262 930 for (Int_t i=0; i<count; i++)
61899827 931 {
a2ce3799 932 AliTriggerAnalysis* triggerAnalysis = static_cast<AliTriggerAnalysis*> (fTriggerAnalysis.At(i));
933
296dd262 934 switch (runNumber)
935 {
a2ce3799 936 case 104315:
296dd262 937 case 104316:
938 case 104320:
a2ce3799 939 case 104321:
296dd262 940 case 104439:
941 triggerAnalysis->SetV0TimeOffset(7.5);
942 break;
a2ce3799 943 default:
944 triggerAnalysis->SetV0TimeOffset(0);
296dd262 945 }
f4ca8f20 946 }
947
758941d4 948 fCurrentRun = runNumber;
949
29e8486e 950 TH1::AddDirectory(oldStatus);
951
61899827 952 return kTRUE;
953}
954
85c71ba7 955TH2F * AliPhysicsSelection::BookHistStatistics(const char * tag) {
bc18a83e 956 // add 6 rows to count for the estimate of good, accidentals and
957 // BG and the ratio of BG and accidentals to total +ratio goot to
958 // first col + 2 for error on good.
959 // TODO: Remember the the indexes of rows for the BG selection. Add new member fBGRows[] and use kStat as indexes
85c71ba7 960
961 Int_t count = fCollTrigClasses.GetEntries() + fBGTrigClasses.GetEntries();
962#ifdef VERBOSE_STAT
141265a2 963 Int_t extrarows = fComputeBG != 0 ? 8 : 0;
85c71ba7 964#else
141265a2 965 Int_t extrarows = fComputeBG != 0 ? 3 : 0;
85c71ba7 966#endif
967 TH2F * h = new TH2F(Form("fHistStatistics%s",tag), Form("fHistStatistics - %s ;;",tag), kStatAccepted, 0.5, kStatAccepted+0.5, count+extrarows, -0.5, -0.5 + count+extrarows);
968
969 h->GetXaxis()->SetBinLabel(kStatTriggerClass, "Trigger class");
970 h->GetXaxis()->SetBinLabel(kStatHWTrig, "Hardware trigger");
971 h->GetXaxis()->SetBinLabel(kStatFO1, "FO >= 1");
972 h->GetXaxis()->SetBinLabel(kStatFO2, "FO >= 2");
c2ba5a61 973 h->GetXaxis()->SetBinLabel(kStatFO2L1, "FO (L1) >= 2");
85c71ba7 974 h->GetXaxis()->SetBinLabel(kStatV0A, "V0A");
975 h->GetXaxis()->SetBinLabel(kStatV0C, "V0C");
976 h->GetXaxis()->SetBinLabel(kStatFMD, "FMD");
85c71ba7 977 h->GetXaxis()->SetBinLabel(kStatV0ABG, "V0A BG");
978 h->GetXaxis()->SetBinLabel(kStatV0CBG, "V0C BG");
c2ba5a61 979 h->GetXaxis()->SetBinLabel(kStatZDCA, "ZDCA");
980 h->GetXaxis()->SetBinLabel(kStatZDCC, "ZDCC");
981 h->GetXaxis()->SetBinLabel(kStatZDCAC, "ZDCA & ZDCC");
85c71ba7 982 h->GetXaxis()->SetBinLabel(kStatMB1, "(FO >= 1 | V0A | V0C) & !V0 BG");
983 h->GetXaxis()->SetBinLabel(kStatMB1Prime, "(FO >= 2 | (FO >= 1 & (V0A | V0C)) | (V0A &v0C) ) & !V0 BG");
c2ba5a61 984 //h->GetXaxis()->SetBinLabel(kStatFO1AndV0, "FO >= 1 & (V0A | V0C) & !V0 BG");
85c71ba7 985 h->GetXaxis()->SetBinLabel(kStatV0, "V0A & V0C & !V0 BG & !BG ID");
986 h->GetXaxis()->SetBinLabel(kStatOffline, "Offline Trigger");
c2ba5a61 987 //h->GetXaxis()->SetBinLabel(kStatAny2Hits, "2 Hits & !V0 BG");
85c71ba7 988 h->GetXaxis()->SetBinLabel(kStatBG, "Background identification");
989 h->GetXaxis()->SetBinLabel(kStatAccepted, "Accepted");
990
991 Int_t n = 1;
992 for (Int_t i=0; i < fCollTrigClasses.GetEntries(); i++)
993 {
994 h->GetYaxis()->SetBinLabel(n, ((TObjString*) fCollTrigClasses.At(i))->String());
995 n++;
996 }
997 for (Int_t i=0; i < fBGTrigClasses.GetEntries(); i++)
998 {
999 h->GetYaxis()->SetBinLabel(n, ((TObjString*) fBGTrigClasses.At(i))->String());
1000 n++;
1001 }
1002
1003 if(fComputeBG) {
141265a2 1004 fBGStatOffset = n;
1005 h->GetYaxis()->SetBinLabel(n++, Form("BG (A+C) (Mask [0x%x])", fComputeBG));
85c71ba7 1006 h->GetYaxis()->SetBinLabel(n++, "ACC");
1007#ifdef VERBOSE_STAT
1008 h->GetYaxis()->SetBinLabel(n++, "BG (A+C) % (rel. to CINT1B)");
1009 h->GetYaxis()->SetBinLabel(n++, "ACC % (rel. to CINT1B)");
1010 h->GetYaxis()->SetBinLabel(n++, "ERR GOOD %");
1011 h->GetYaxis()->SetBinLabel(n++, "GOOD % (rel. to 1st col)");
1012 h->GetYaxis()->SetBinLabel(n++, "ERR GOOD");
1013#endif
1014 h->GetYaxis()->SetBinLabel(n++, "GOOD");
1015 }
1016
1017 return h;
1018}
1019
e9247450 1020void AliPhysicsSelection::Print(Option_t *option) const
61899827 1021{
1022 // print the configuration
e9247450 1023 TString msg;
a2ce3799 1024 Printf("Configuration initialized for run %d (MC: %d):", fCurrentRun, fMC);
e9247450 1025 msg += Form("Configuration initialized for run %d (MC: %d):\n", fCurrentRun, fMC);
61899827 1026
296dd262 1027 Printf("Collision trigger classes:");
1028 for (Int_t i=0; i < fCollTrigClasses.GetEntries(); i++)
1029 Printf("%s", ((TObjString*) fCollTrigClasses.At(i))->String().Data());
61899827 1030
296dd262 1031 Printf("Background trigger classes:");
1032 for (Int_t i=0; i < fBGTrigClasses.GetEntries(); i++)
1033 Printf("%s", ((TObjString*) fBGTrigClasses.At(i))->String().Data());
1034
1035 AliTriggerAnalysis* triggerAnalysis = dynamic_cast<AliTriggerAnalysis*> (fTriggerAnalysis.At(0));
61899827 1036
296dd262 1037 if (triggerAnalysis)
1038 {
1039 if (triggerAnalysis->GetV0TimeOffset() > 0)
1040 Printf("V0 time offset active: %.2f ns", triggerAnalysis->GetV0TimeOffset());
61899827 1041
296dd262 1042 Printf("\nTotal available events:");
528640ed 1043
296dd262 1044 triggerAnalysis->PrintTriggerClasses();
1045 }
528640ed 1046
0c6c629b 1047 if (fHistStatistics[kStatIdxAll])
b43e01ed 1048 {
0c6c629b 1049 for (Int_t i=0; i<fCollTrigClasses.GetEntries(); i++)
1050 {
1051 Printf("\nSelection statistics for collision trigger %s:", ((TObjString*) fCollTrigClasses.At(i))->String().Data());
e9247450 1052 msg += Form("\nSelection statistics for collision trigger %s:\n", ((TObjString*) fCollTrigClasses.At(i))->String().Data());
0c6c629b 1053
1054 Printf("Total events with correct trigger class: %d", (Int_t) fHistStatistics[kStatIdxAll]->GetBinContent(1, i+1));
e9247450 1055 msg += Form("Total events with correct trigger class: %d\n", (Int_t) fHistStatistics[kStatIdxAll]->GetBinContent(1, i+1));
0c6c629b 1056 Printf("Selected collision candidates: %d", (Int_t) fHistStatistics[kStatIdxAll]->GetBinContent(fHistStatistics[kStatIdxAll]->GetXaxis()->FindBin("Accepted"), i+1));
e9247450 1057 msg += Form("Selected collision candidates: %d\n", (Int_t) fHistStatistics[kStatIdxAll]->GetBinContent(fHistStatistics[kStatIdxAll]->GetXaxis()->FindBin("Accepted"), i+1));
0c6c629b 1058 }
17ba346c 1059 }
1060
1061 if (fHistBunchCrossing)
1062 {
1063 Printf("\nBunch crossing statistics:");
e9247450 1064 msg += "\nBunch crossing statistics:\n";
17ba346c 1065
1066 for (Int_t i=1; i<=fHistBunchCrossing->GetNbinsY(); i++)
1067 {
1068 TString str;
1069 str.Form("Trigger %s has accepted events in the bunch crossings: ", fHistBunchCrossing->GetYaxis()->GetBinLabel(i));
1070
1071 for (Int_t j=1; j<=fHistBunchCrossing->GetNbinsX(); j++)
1072 if (fHistBunchCrossing->GetBinContent(j, i) > 0)
1073 str += Form("%d, ", (Int_t) fHistBunchCrossing->GetXaxis()->GetBinCenter(j));
e9247450 1074
17ba346c 1075 Printf("%s", str.Data());
e9247450 1076 msg += str;
1077 msg += "\n";
17ba346c 1078 }
1079
1080 for (Int_t j=1; j<=fHistBunchCrossing->GetNbinsX(); j++)
1081 {
0c6c629b 1082 Int_t countColl = 0;
1083 Int_t countBG = 0;
17ba346c 1084 for (Int_t i=1; i<=fHistBunchCrossing->GetNbinsY(); i++)
1085 {
1086 if (fHistBunchCrossing->GetBinContent(j, i) > 0)
0c6c629b 1087 {
1088 if (fCollTrigClasses.FindObject(fHistBunchCrossing->GetYaxis()->GetBinLabel(i)))
1089 countColl++;
1090 if (fBGTrigClasses.FindObject(fHistBunchCrossing->GetYaxis()->GetBinLabel(i)))
1091 countBG++;
1092 }
17ba346c 1093 }
0c6c629b 1094 if (countColl > 0 && countBG > 0)
1095 Printf("WARNING: Bunch crossing %d has collision and BG trigger classes active. Check BPTX functioning for this run!", (Int_t) fHistBunchCrossing->GetXaxis()->GetBinCenter(j));
17ba346c 1096 }
b43e01ed 1097 }
91bea6e7 1098
1099 if (fUsingCustomClasses)
1100 Printf("WARNING: Using custom trigger classes!");
1101 if (fSkipTriggerClassSelection)
1102 Printf("WARNING: Skipping trigger class selection!");
85c71ba7 1103 if (fSkipV0)
1104 Printf("WARNING: Ignoring V0 information in selection");
1105 if(!fBin0CallBack)
1106 Printf("WARNING: Callback not set: will not fill the statistics for the bin 0");
e9247450 1107 TString opt(option);
1108 opt.ToUpper();
1109 if (opt == "STAT") {
1110 AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
1111 if (mgr) mgr->AddStatisticsMsg(msg);
1112 }
61899827 1113}
1114
1115Long64_t AliPhysicsSelection::Merge(TCollection* list)
1116{
1117 // Merge a list of AliMultiplicityCorrection objects with this (needed for
1118 // PROOF).
1119 // Returns the number of merged objects (including this).
1120
1121 if (!list)
1122 return 0;
1123
1124 if (list->IsEmpty())
1125 return 1;
1126
1127 TIterator* iter = list->MakeIterator();
1128 TObject* obj;
17ba346c 1129
61899827 1130 // collections of all histograms
1131 const Int_t nHists = 9;
1132 TList collections[nHists];
1133
1134 Int_t count = 0;
1135 while ((obj = iter->Next())) {
1136
1137 AliPhysicsSelection* entry = dynamic_cast<AliPhysicsSelection*> (obj);
1138 if (entry == 0)
1139 continue;
518b6d59 1140 // Update run number. If this one is not initialized (-1) take the one from
1141 // the next physics selection to be merged with. In case of 2 different run
1142 // numbers issue a warning (should physics selections from different runs be
1143 // merged together) A.G.
1144 Int_t currentRun = entry->GetCurrentRun();
1145 // Nothing to merge with since run number was not initialized.
1146 if (currentRun < 0) continue;
c2ba5a61 1147 if (fCurrentRun < 0)
1148 {
1149 fCurrentRun = currentRun;
1150 fMC = entry->fMC;
1151 for (Int_t i=0; i<entry->fCollTrigClasses.GetEntries(); i++)
1152 fCollTrigClasses.Add(entry->fCollTrigClasses.At(i)->Clone());
1153 for (Int_t i=0; i<entry->fBGTrigClasses.GetEntries(); i++)
1154 fBGTrigClasses.Add(entry->fBGTrigClasses.At(i)->Clone());
1155 }
518b6d59 1156 if (fCurrentRun != currentRun)
1157 AliWarning(Form("Current run %d not matching the one to be merged with %d", fCurrentRun, currentRun));
1158
c2ba5a61 1159 if (entry->fTriggerAnalysis.GetEntries() > 0)
1160 collections[0].Add(&(entry->fTriggerAnalysis));
85c71ba7 1161 if (entry->fHistStatistics[0])
1162 collections[1].Add(entry->fHistStatistics[0]);
1163 if (entry->fHistStatistics[1])
1164 collections[2].Add(entry->fHistStatistics[1]);
17ba346c 1165 if (entry->fHistBunchCrossing)
85c71ba7 1166 collections[3].Add(entry->fHistBunchCrossing);
73cc8654 1167 if (entry->fHistTriggerPattern)
1168 collections[4].Add(entry->fHistTriggerPattern);
296dd262 1169 if (entry->fBackgroundIdentification)
73cc8654 1170 collections[5].Add(entry->fBackgroundIdentification);
61899827 1171
1172 count++;
1173 }
1174
c2ba5a61 1175 if (fTriggerAnalysis.GetEntries() == 0 && collections[0].GetEntries() > 0)
1176 {
1177 TList* firstList = (TList*) collections[0].First();
1178 for (Int_t i=0; i<firstList->GetEntries(); i++)
1179 fTriggerAnalysis.Add(firstList->At(i)->Clone());
1180
1181 collections[0].RemoveAt(0);
1182 }
17ba346c 1183 fTriggerAnalysis.Merge(&collections[0]);
c2ba5a61 1184
1185 // if this instance is empty (not initialized) nothing would be merged here --> copy first entry
1186 if (!fHistStatistics[0] && collections[1].GetEntries() > 0)
1187 {
1188 fHistStatistics[0] = (TH2F*) collections[1].First()->Clone();
1189 collections[1].RemoveAt(0);
1190 }
85c71ba7 1191 if (fHistStatistics[0])
1192 fHistStatistics[0]->Merge(&collections[1]);
c2ba5a61 1193
1194 if (!fHistStatistics[1] && collections[2].GetEntries() > 0)
1195 {
1196 fHistStatistics[1] = (TH2F*) collections[2].First()->Clone();
1197 collections[2].RemoveAt(0);
1198 }
85c71ba7 1199 if (fHistStatistics[1])
1200 fHistStatistics[1]->Merge(&collections[2]);
c2ba5a61 1201
1202 if (!fHistBunchCrossing && collections[3].GetEntries() > 0)
1203 {
1204 fHistBunchCrossing = (TH2F*) collections[3].First()->Clone();
1205 collections[3].RemoveAt(0);
1206 }
17ba346c 1207 if (fHistBunchCrossing)
85c71ba7 1208 fHistBunchCrossing->Merge(&collections[3]);
c2ba5a61 1209
1210 if (!fHistTriggerPattern && collections[4].GetEntries() > 0)
1211 {
1212 fHistTriggerPattern = (TH1F*) collections[4].First()->Clone();
1213 collections[4].RemoveAt(0);
1214 }
73cc8654 1215 if (fHistTriggerPattern)
1216 fHistTriggerPattern->Merge(&collections[4]);
c2ba5a61 1217
1218 if (!fBackgroundIdentification && collections[5].GetEntries() > 0)
1219 {
1220 fBackgroundIdentification = (AliAnalysisCuts*) collections[5].First()->Clone();
1221 collections[5].RemoveAt(0);
1222 }
296dd262 1223 if (fBackgroundIdentification)
73cc8654 1224 fBackgroundIdentification->Merge(&collections[5]);
61899827 1225
1226 delete iter;
1227
1228 return count+1;
1229}
1230
1231void AliPhysicsSelection::SaveHistograms(const char* folder) const
1232{
1233 // write histograms to current directory
1234
85c71ba7 1235 if (!fHistStatistics[0] || !fHistStatistics[1])
61899827 1236 return;
1237
1238 if (folder)
1239 {
1240 gDirectory->mkdir(folder);
1241 gDirectory->cd(folder);
1242 }
1243
85c71ba7 1244
1245 // Fill the last rows of fHistStatistics before saving
1246 if (fComputeBG) {
1247 Int_t triggerScheme = GetTriggerScheme(UInt_t(fCurrentRun));
1248 if(triggerScheme != 1){
1249 AliWarning("BG estimate only supported for trigger scheme \"1\" (CINT1 suite)");
73cc8654 1250 } else if (fUseMuonTriggers) {
1251 AliWarning("BG estimate with muon triggers to be implemented");
141265a2 1252 } else {
1253 AliInfo("BG estimate assumes that for a given run you only have A and C triggers separately or"
1254 " toghether as a AC class! Make sure this assumption holds in your case");
1255
1256 // use an anum for the different trigger classes, to make loops easier to read
1257 enum {kClassB =0 , kClassA, kClassC, kClassAC, kClassE, kNClasses};
1258 const char * classFlags[] = {"B", "A", "C", "AC", "E"}; // labels
1259
1260 UInt_t * rows[kNClasses] = {0}; // Array of matching rows
1261 Int_t nrows[kNClasses] = {0};
1262 // Get rows matching the requested trigger bits for all trigger classes
1263 for(Int_t iTrigClass = 0; iTrigClass < kNClasses; iTrigClass++){
1264 nrows[iTrigClass] = GetStatRow(classFlags[iTrigClass],fComputeBG,&rows[iTrigClass]);
1265 }
1266
5f337f3d 1267 // 0. Determine the ratios of triggers E/B, A/B, C/B from the stat histogram
1268 // Those are used to rescale the different classes to the same number of bx ids
141265a2 1269 // TODO: pass names of the rows for B, CA and E and look names of the rows. How do I handle the case in which both AC are in the same row?
1270 Int_t nBXIds[kNClasses] = {0};
1271 cout <<"Computing BG:" << endl;
1272
1273 for(Int_t iTrigClass = 0; iTrigClass < kNClasses; iTrigClass++){
1274 for(Int_t irow = 0; irow < nrows[iTrigClass]; irow++) {
1275 if(irow==0) cout << "- Class " << classFlags[iTrigClass] << endl;
1276 for (Int_t j=1; j<=fHistBunchCrossing->GetNbinsX(); j++) {
1277 if (fHistBunchCrossing->GetBinContent(j, rows[iTrigClass][irow]) > 0) {
1278 nBXIds[iTrigClass]++;
1279 }
1280 }
1281 if(nBXIds[iTrigClass]>0) cout << " Using row " << rows[iTrigClass][irow] << ": "
1282 << fHistBunchCrossing->GetYaxis()->GetBinLabel(rows[iTrigClass][irow])
1283 << " (nBXID "<< nBXIds[iTrigClass] << ")"<< endl;
1284
1285 }
1286
5f337f3d 1287 }
141265a2 1288
1289 Float_t ratioToB[kNClasses];
1290 ratioToB[kClassE] = nBXIds[kClassE] >0 ? Float_t(nBXIds[kClassB])/nBXIds[kClassE] : 0;
1291 ratioToB[kClassA] = nBXIds[kClassA] >0 ? Float_t(nBXIds[kClassB])/nBXIds[kClassA] : 0;
1292 ratioToB[kClassC] = nBXIds[kClassC] >0 ? Float_t(nBXIds[kClassB])/nBXIds[kClassC] : 0;
1293 ratioToB[kClassAC] = nBXIds[kClassAC] >0 ? Float_t(nBXIds[kClassB])/nBXIds[kClassAC] : 0;
5f337f3d 1294 Printf("Ratio between the BX ids in the different trigger classes:");
141265a2 1295 Printf(" B/E = %d/%d = %f", nBXIds[kClassB],nBXIds[kClassE], ratioToB[kClassE] );
1296 Printf(" B/A = %d/%d = %f", nBXIds[kClassB],nBXIds[kClassA], ratioToB[kClassA] );
1297 Printf(" B/C = %d/%d = %f", nBXIds[kClassB],nBXIds[kClassC], ratioToB[kClassC] );
1298 Printf(" B/AC = %d/%d = %f", nBXIds[kClassB],nBXIds[kClassAC],ratioToB[kClassAC]);
85c71ba7 1299 Int_t nHistStat = 2;
85c71ba7 1300
141265a2 1301 // 1. loop over all cols
85c71ba7 1302 for(Int_t iHistStat = 0; iHistStat < nHistStat; iHistStat++){
85c71ba7 1303 Int_t ncol = fHistStatistics[iHistStat]->GetNbinsX();
141265a2 1304 Float_t good1 = 0;
85c71ba7 1305 for(Int_t icol = 1; icol <= ncol; icol++) {
141265a2 1306 Int_t nEvents[kNClasses] = {0}; // number of events should be reset at every column
1307 // For all trigger classes, add up over row matching trigger mask (as selected before)
1308 for(Int_t iTrigClass = 0; iTrigClass < kNClasses; iTrigClass++){
1309 for(Int_t irow = 0; irow < nrows[iTrigClass]; irow++) {
1310 nEvents[iTrigClass] += (Int_t) fHistStatistics[iHistStat]->GetBinContent(icol,rows[iTrigClass][irow]);
1311 }
1312 // cout << "Events " << classFlags[iTrigClass] << " ("<<icol<<") " << nEvents[iTrigClass] << endl;
1313 }
1314 if (nEvents[kClassB]>0) {
a67fbdc0 1315 Float_t acc = ratioToB[kClassE]*nEvents[kClassE];
141265a2 1316 Double_t acc_err = TMath::Sqrt(ratioToB[kClassE]*ratioToB[kClassE]*nEvents[kClassE]);
85c71ba7 1317 // Int_t bg = cint1A + cint1C - 2*acc;
141265a2 1318
1319 // Assuming that for a given class the triggers are either recorded as A+C or AC
1320 Float_t bg = nEvents[kClassAC] > 0 ?
1321 fBIFactorAC*(ratioToB[kClassAC]*nEvents[kClassAC] - 2*acc):
1322 fBIFactorA* (ratioToB[kClassA]*nEvents[kClassA]-acc) +
1323 fBIFactorC* (ratioToB[kClassC]*nEvents[kClassC]-acc) ;
1324
8fc17e87 1325 // cout << "-----------------------" << endl;
1326 // cout << "Factors: " << fBIFactorA << " " << fBIFactorC << " " << fBIFactorAC << endl;
1327 // cout << "Ratios: " << ratioToB[kClassA] << " " << ratioToB[kClassC] << " " << ratioToB[kClassAC] << endl;
1328 // cout << "Evts: " << nEvents[kClassA] << " " << nEvents[kClassC] << " " << nEvents[kClassAC] << " " << nEvents[kClassB] << endl;
1329 // cout << "Acc: " << acc << endl;
141265a2 1330 // cout << "BG: " << bg << endl;
1331 // cout << " " << fBIFactorA* (ratioToB[kClassA]*nEvents[kClassA]-acc) <<endl;
1332 // cout << " " << fBIFactorC* (ratioToB[kClassC]*nEvents[kClassC]-acc) <<endl;
1333 // cout << " " << fBIFactorAC*(ratioToB[kClassAC]*nEvents[kClassAC] - 2*acc) << endl;
1334 // cout << "-----------------------" << endl;
1335
1336 Float_t good = Float_t(nEvents[kClassB]) - bg - acc;
85c71ba7 1337 if (icol ==1) good1 = good;
141265a2 1338 // Float_t errGood = TMath::Sqrt(2*(nEvents[kClassA]+nEvents[kClassC]+nEvents[kClassE]));// Error on the number of goods assuming only bg fluctuates
85c71ba7 1339 // DeltaG^2 = B + FA^2 A + FC^2 C + Ratio^2 (FA+FC-1)^2 E.
141265a2 1340 Float_t errGood = nEvents[kClassAC] > 0 ?
1341 TMath::Sqrt( nEvents[kClassB] +
1342 fBIFactorAC*fBIFactorAC*ratioToB[kClassAC]*ratioToB[kClassAC]*nEvents[kClassAC] +
1343 ratioToB[kClassE] * ratioToB[kClassE] *
1344 (fBIFactorAC - 1)*(fBIFactorAC - 1)*nEvents[kClassE]) :
1345 TMath::Sqrt( nEvents[kClassB] +
1346 fBIFactorA*fBIFactorA*ratioToB[kClassA]*ratioToB[kClassA]*nEvents[kClassA] +
1347 fBIFactorC*fBIFactorC*ratioToB[kClassC]*ratioToB[kClassC]*nEvents[kClassC] +
1348 ratioToB[kClassE] * ratioToB[kClassE] *
1349 (fBIFactorA + fBIFactorC - 1)*(fBIFactorA + fBIFactorC - 1)*nEvents[kClassE]);
1350
1351 Float_t errBG = nEvents[kClassAC] > 0 ?
1352 TMath::Sqrt(fBIFactorAC*fBIFactorAC*ratioToB[kClassAC]*ratioToB[kClassAC]*nEvents[kClassAC]+
1353 4*ratioToB[kClassE]*ratioToB[kClassE]*(fBIFactorAC*fBIFactorAC)*nEvents[kClassE]) :
1354 TMath::Sqrt(fBIFactorA*fBIFactorA*ratioToB[kClassA]*ratioToB[kClassA]*nEvents[kClassA]+
1355 fBIFactorC*fBIFactorC*ratioToB[kClassC]*ratioToB[kClassC]*nEvents[kClassC]+
1356 ratioToB[kClassE]*ratioToB[kClassE]*(fBIFactorA+fBIFactorC)*(fBIFactorA+fBIFactorC)*nEvents[kClassE]);
85c71ba7 1357
1358
141265a2 1359 fHistStatistics[iHistStat]->SetBinContent(icol,fBGStatOffset+kStatRowBG,bg);
1360 fHistStatistics[iHistStat]->SetBinError (icol,fBGStatOffset+kStatRowBG,errBG);
1361 fHistStatistics[iHistStat]->SetBinContent(icol,fBGStatOffset+kStatRowAcc,acc);
1362 fHistStatistics[iHistStat]->SetBinError (icol,fBGStatOffset+kStatRowAcc,acc_err);
1363 fHistStatistics[iHistStat]->SetBinContent(icol,fBGStatOffset+kStatRowGood,good);
1364 fHistStatistics[iHistStat]->SetBinError (icol,fBGStatOffset+kStatRowGood,errGood);
85c71ba7 1365
1366#ifdef VERBOSE_STAT
141265a2 1367 //kStatRowBG=0,kStatRowAcc,kStatRowBGFrac,kStatRowAccFrac,kStatRowErrGoodFrac,kStatRowGoodFrac,kStatRowGood,kStatRowErrGood
1368 Float_t accFrac = Float_t(acc) / nEvents[kClassB] *100;
1369 Float_t errAccFrac= Float_t(acc_err) / nEvents[kClassB] *100;
1370 Float_t bgFrac = Float_t(bg) / nEvents[kClassB] *100;
85c71ba7 1371 Float_t goodFrac = Float_t(good) / good1 *100;
1372 Float_t errGoodFrac = errGood/good1 * 100;
141265a2 1373 Float_t errFracBG = bg > 0 ? TMath::Sqrt((errBG/bg)*(errBG/bg) + 1/nEvents[kClassB])*bgFrac : 0;
1374 fHistStatistics[iHistStat]->SetBinContent(icol,fBGStatOffset+kStatRowBGFrac,bgFrac);
1375 fHistStatistics[iHistStat]->SetBinError (icol,fBGStatOffset+kStatRowBGFrac,errFracBG);
1376 fHistStatistics[iHistStat]->SetBinContent(icol,fBGStatOffset+kStatRowAccFrac,accFrac);
1377 fHistStatistics[iHistStat]->SetBinError (icol,fBGStatOffset+kStatRowAccFrac,errAccFrac);
1378 fHistStatistics[iHistStat]->SetBinContent(icol,fBGStatOffset+kStatRowGoodFrac,goodFrac);
1379 fHistStatistics[iHistStat]->SetBinContent(icol,fBGStatOffset+kStatRowErrGoodFrac,errGoodFrac);
1380 fHistStatistics[iHistStat]->SetBinContent(icol,fBGStatOffset+kStatRowErrGood,errGood);
85c71ba7 1381#endif
1382 }
1383 }
1384 }
ed185ba9 1385 for (Int_t iTrigClass = 0; iTrigClass < kNClasses; iTrigClass++){
1386 delete [] rows[iTrigClass];
1387 }
85c71ba7 1388 }
1389 }
1390
1391 fHistStatistics[0]->Write();
1392 fHistStatistics[1]->Write();
8e58a7b9 1393 if(fHistBunchCrossing ) fHistBunchCrossing ->Write();
1394 if(fHistTriggerPattern) fHistTriggerPattern->Write();
61899827 1395
296dd262 1396 Int_t count = fCollTrigClasses.GetEntries() + fBGTrigClasses.GetEntries();
1397 for (Int_t i=0; i < count; i++)
1398 {
1399 TString triggerClass = "trigger_histograms_";
1400 if (i < fCollTrigClasses.GetEntries())
1401 triggerClass += ((TObjString*) fCollTrigClasses.At(i))->String();
1402 else
1403 triggerClass += ((TObjString*) fBGTrigClasses.At(i - fCollTrigClasses.GetEntries()))->String();
61899827 1404
296dd262 1405 gDirectory->mkdir(triggerClass);
1406 gDirectory->cd(triggerClass);
61899827 1407
296dd262 1408 static_cast<AliTriggerAnalysis*> (fTriggerAnalysis.At(i))->SaveHistograms();
1409
1410 gDirectory->cd("..");
1411 }
1412
1413 if (fBackgroundIdentification)
1414 {
1415 gDirectory->mkdir("background_identification");
1416 gDirectory->cd("background_identification");
1417
1418 fBackgroundIdentification->GetOutput()->Write();
1419
1420 gDirectory->cd("..");
1421 }
61899827 1422
1423 if (folder)
1424 gDirectory->cd("..");
1425}
85c71ba7 1426
141265a2 1427Int_t AliPhysicsSelection::GetStatRow(const char * triggerBXClass, UInt_t offlineTriggerType, UInt_t ** rowIDs) const {
1428 // Puts inside the array rowIDs the row number for a given offline
1429 // trigger in a given bx class. Returns the total number of lines
1430 // matching the selection
1431 // triggerBXClass can be either "A", "AC", "B" or "E"
1432 // offlineTriggerType is one of the types defined in AliVEvent
1433 // User should delete rowIDs if no longer needed
1434
1435 if(!fHistStatistics[0]) {
1436 AliWarning("Not initialized, returning 0");
1437 return 0;
1438 }
1439 const Int_t nrows = fHistStatistics[0]->GetNbinsY();
1440
1441 // allocate memory for at maximum nrows
1442 Int_t nMatches = 0;
1443 (*rowIDs) = new UInt_t[nrows];
1444
1445 // Build regular expression. look for a +, followed by the beginning
1446 // of a word. Within the word, look for the class id after any
1447 // number of any char, but at most one dash ("[^-]*-?"), followed by
1448 // a - and then any char (".*") and at the class id ("\\&(\\d)")
1449 // The class id is stored.
1450 // WARNING: please check this if the trigger classes change
1451 TPRegexp re(Form("\\+\\b[^-]*-?%s-.*\\&(\\d)",triggerBXClass));
1452 // Loop over rows and find matching ones:
1453 for(Int_t irow = 1; irow <= nrows; irow++){
1454 TObjArray * matches = re.MatchS(fHistStatistics[0]->GetYaxis()->GetBinLabel(irow));
1455 if (matches->GetEntries()) {
1456 TString s = ((TObjString*)matches->At(1))->GetString();
1457 if(UInt_t(s.Atoi()) & offlineTriggerType) { // bitwise comparison with the requested mask
1458 // cout << "Marching " << s.Data() << " " << offlineTriggerType << " " << fHistStatistics[0]->GetYaxis()->GetBinLabel(irow) << endl;
1459 (*rowIDs)[nMatches] = irow;
1460 nMatches++;
1461 }
1462 }
1463 delete matches;
85c71ba7 1464 }
141265a2 1465
1466 return nMatches;
1467}
85c71ba7 1468
141265a2 1469void AliPhysicsSelection::SetBIFactors(const AliESDEvent * aESD) {
1470 // Set factors for realtive bunch intesities
1471 if(!aESD) {
1472 AliFatal("ESD not given");
1473 }
1474 Int_t run = aESD->GetRunNumber();
1475 if (run > 105268) {
1476 // intensities stored in the ESDs
1477 const AliESDRun* esdRun = aESD->GetESDRun();
1478 Double_t intAB = esdRun->GetMeanIntensityIntecting(0);
1479 Double_t intCB = esdRun->GetMeanIntensityIntecting(1);
1480 Double_t intAA = esdRun->GetMeanIntensityNonIntecting(0);
1481 Double_t intCC = esdRun->GetMeanIntensityNonIntecting(1);
1482
1483 // cout << "INT " <<intAB <<endl;
1484 // cout << "INT " <<intCB <<endl;
1485 // cout << "INT " <<intAA <<endl;
1486 // cout << "INT " <<intCC <<endl;
1487
1488 if (intAB > -1 && intAA > -1) {
1489 fBIFactorA = intAB/intAA;
1490 } else {
1491 AliWarning("Cannot set fBIFactorA, assuming 1");
1492 }
1493
1494 if (intCB > -1 && intCC > -1) {
1495 fBIFactorC = intCB/intCC;
1496 } else {
1497 AliWarning("Cannot set fBIFactorC, assuming 1");
1498 }
1499
1500 if (intAB > -1 && intAA > -1 &&
1501 intCB > -1 && intCC > -1) {
1502 fBIFactorAC = (intAB+intCB)/(intAA+intCC);
1503 } else {
1504 AliWarning("Cannot set fBIFactorAC, assuming 1");
1505 }
1506
1507 }
1508 else {
1509 // First runs. Intensities hardcoded
1510 switch(run) {
1511 case 104155:
1512 fBIFactorA = 0.961912722908;
1513 fBIFactorC = 1.04992336081;
1514 break;
1515 case 104157:
1516 fBIFactorA = 0.947312854998;
1517 fBIFactorC = 1.01599706417;
1518 break;
1519 case 104159:
1520 fBIFactorA = 0.93659320151;
1521 fBIFactorC = 0.98580804207;
1522 break;
1523 case 104160:
1524 fBIFactorA = 0.929664189926;
1525 fBIFactorC = 0.963467679851;
1526 break;
1527 case 104315:
1528 fBIFactorA = 1.08939104979;
1529 fBIFactorC = 0.931113921925;
1530 break;
1531 case 104316:
1532 fBIFactorA = 1.08351880974;
1533 fBIFactorC = 0.916068345845;
1534 break;
1535 case 104320:
1536 fBIFactorA = 1.07669281245;
1537 fBIFactorC = 0.876818744763;
1538 break;
1539 case 104321:
1540 fBIFactorA = 1.00971079602;
1541 fBIFactorC = 0.773781299076;
1542 break;
1543 case 104792:
1544 fBIFactorA = 0.787215863962;
1545 fBIFactorC = 0.778253173071;
1546 break;
1547 case 104793:
1548 fBIFactorA = 0.692211363661;
1549 fBIFactorC = 0.733152456667;
1550 break;
1551 case 104799:
1552 fBIFactorA = 1.04027825161;
1553 fBIFactorC = 1.00530825942;
1554 break;
1555 case 104800:
1556 fBIFactorA = 1.05309910671;
1557 fBIFactorC = 1.00376801855;
1558 break;
1559 case 104801:
1560 fBIFactorA = 1.0531231922;
1561 fBIFactorC = 0.992439666758;
1562 break;
1563 case 104802:
1564 fBIFactorA = 1.04191478134;
1565 fBIFactorC = 0.979368585208;
1566 break;
1567 case 104803:
1568 fBIFactorA = 1.03121314094;
1569 fBIFactorC = 0.973379962609;
1570 break;
1571 case 104824:
1572 fBIFactorA = 0.969945926722;
1573 fBIFactorC = 0.39549745806;
1574 break;
1575 case 104825:
1576 fBIFactorA = 0.968627213937;
1577 fBIFactorC = 0.310100412205;
1578 break;
1579 case 104841:
1580 fBIFactorA = 0.991601393212;
1581 fBIFactorC = 0.83762204722;
1582 break;
1583 case 104845:
1584 fBIFactorA = 0.98040863886;
1585 fBIFactorC = 0.694824205793;
1586 break;
1587 case 104867:
1588 fBIFactorA = 1.10646173412;
1589 fBIFactorC = 0.841407246916;
1590 break;
1591 case 104876:
1592 fBIFactorA = 1.12063452421;
1593 fBIFactorC = 0.78726542895;
1594 break;
1595 case 104890:
1596 fBIFactorA = 1.02346137453;
1597 fBIFactorC = 1.03355663595;
1598 break;
1599 case 104892:
1600 fBIFactorA = 1.05406025913;
1601 fBIFactorC = 1.00029166135;
1602 break;
1603 case 105143:
1604 fBIFactorA = 0.947343384349;
1605 fBIFactorC = 0.972637444408;
1606 break;
1607 case 105160:
1608 fBIFactorA = 0.908854622177;
1609 fBIFactorC = 0.958851103977;
1610 break;
1611 case 105256:
1612 fBIFactorA = 0.810076150206;
1613 fBIFactorC = 0.884663561883;
1614 break;
1615 case 105257:
1616 fBIFactorA = 0.80974912303;
1617 fBIFactorC = 0.878859123479;
1618 break;
1619 case 105268:
1620 fBIFactorA = 0.809052110679;
1621 fBIFactorC = 0.87233890989;
1622 break;
1623 default:
1624 fBIFactorA = 1;
1625 fBIFactorC = 1;
1626 }
1627 }
85c71ba7 1628
1629}