]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG/muon/AliAnalysisTaskMuMu.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Added Tender to PHOSPi0Flow single-task macros.
[u/mrichter/AliRoot.git] / PWG / muon / AliAnalysisTaskMuMu.cxx
CommitLineData
2f331ac9 1#include "AliAnalysisTaskMuMu.h"
2
3#include "AliAnalysisManager.h"
4#include "AliCentrality.h"
5#include "AliCounterCollection.h"
6#include "AliHistogramCollection.h"
7#include "AliInputEventHandler.h"
8#include "AliLog.h"
9#include "AliMuonTrackCuts.h"
10#include "Riostream.h"
11#include "TCanvas.h"
12#include "TDatabasePDG.h"
13#include "TFormula.h"
14#include "TH1.h"
15#include "TH2.h"
16#include "THashList.h"
17#include "TList.h"
18#include "TMath.h"
19#include "TObjString.h"
20#include "TPaveText.h"
21#include "TROOT.h"
22#include <algorithm>
23#include "AliAODEvent.h"
24#include "AliESDEvent.h"
25#include "AliESDMuonTrack.h"
26#include "AliAODTrack.h"
27#include "AliAODTZERO.h"
28#include "AliESDTZERO.h"
29#include "AliCodeTimer.h"
30
31//
32// AliAnalysisTaskMuMu : base class for mu pairs analysis
33//
34// Mainly invariant mass (for J/psi and Upsilon) but also
35// some single control histograms.
36//
37// This base class contains common things for ESD-based
38// and AOD-based analysis
39//
40// The output contains an AliHistogramCollection and
41// an AliCounterCollection
42//
43// author: L. Aphecetche (Subatech)
44//
45
46using std::cout;
47using std::endl;
48
49ClassImp(AliAnalysisTaskMuMu)
50ClassImp(AliAnalysisTaskMuMu::PairCut)
51
52namespace
53{
54 Int_t GetNbins(Double_t xmin, Double_t xmax, Double_t xstep)
55 {
56 if ( TMath::AreEqualRel(xstep,0.0,1E-9) ) return 1;
57
58 return TMath::Nint(TMath::Abs((xmax-xmin)/xstep));
59 }
60
61 TObjArray* GetMuonTriggerList()
62 {
63 TObjArray* a = new TObjArray;
64 a->SetOwner(kTRUE);
65
66 a->Add(new TObjString("CMUL"));
67 a->Add(new TObjString("CMLL"));
68 a->Add(new TObjString("C0MUL"));
69 a->Add(new TObjString("CMSL"));
70 a->Add(new TObjString("CMSH"));
71
72 return a;
73 }
74
75 TObjArray* GetEmcalTriggerList()
76 {
77 TObjArray* a = new TObjArray;
78 a->SetOwner(kTRUE);
79
80 a->Add(new TObjString("CEMC"));
81
82 return a;
83 }
84
85 TObjArray* GetMBTriggerList()
86 {
87 TObjArray* a = new TObjArray;
88 a->SetOwner(kTRUE);
89
90 a->Add(new TObjString("CINT7-S-"));
91 a->Add(new TObjString("CINT8-S-"));
92
93 return a;
94 }
95
96}
97
98//_____________________________________________________________________________
99void AliAnalysisTaskMuMu::PairCut::Print(Option_t* /*opt*/) const
100{
101 std::cout << Form("PAIR CUT %20s SINGLE MASK %x PAIR MASK %x",GetName(),MaskForOneOrBothTrack(),MaskForTrackPair())
102 << std::endl;
103}
104
105//_____________________________________________________________________________
106AliAnalysisTaskMuMu::AliAnalysisTaskMuMu() : AliAnalysisTaskSE("AliAnalysisTaskMuMu"),
107fOutput(0),
108fHistogramCollection(0),
109fEventCounters(0),
110fMuonTrackCuts(0x0),
111fPrecomputedTrackMasks(),
112fIsFromESD(kFALSE),
113fIsDynamicTriggerClasses(kFALSE),
114fShouldSeparatePlusAndMinus(kFALSE),
115fIsHistogrammingDisabled(kFALSE),
116fBeamYear("pp"),
117fCentralityLimits(),
118fTriggerClasses(0),
119fSingleTrackCutNames(0x0),
120fPairTrackCutNames(0x0),
121fCentralityNames(0x0),
122fEventCutNames(0x0),
91dbbc20 123fUseBackgroundTriggers(kFALSE),
124fTriggerInputBitMap()
2f331ac9 125{
126 /// default ctor
127}
128
129//_____________________________________________________________________________
130AliAnalysisTaskMuMu::AliAnalysisTaskMuMu(Bool_t fromESD, const char* beamYear, TArrayF* centralities)
131: AliAnalysisTaskSE(Form("AliAnalysisTaskMuMu-from%s",fromESD ? "ESD":"AOD")),
132fOutput(0),
133fHistogramCollection(0),
134fEventCounters(0),
135fMuonTrackCuts(0x0),
136fPrecomputedTrackMasks(),
137fIsFromESD(fromESD),
138fIsDynamicTriggerClasses(kTRUE),
139fShouldSeparatePlusAndMinus(kFALSE),
140fIsHistogrammingDisabled(kFALSE),
141fBeamYear(beamYear),
142fCentralityLimits(),
143fTriggerClasses(new THashList),
144fSingleTrackCutNames(0x0),
145fPairTrackCutNames(0x0),
146fCentralityNames(new TObjArray),
147fEventCutNames(0x0),
91dbbc20 148fUseBackgroundTriggers(kFALSE),
149fTriggerInputBitMap()
2f331ac9 150{
151 /// Constructor
152 /// The list of triggers to be considered will be updated on the fly
153 /// (see method AddTriggerClasses)
154
155 fTriggerClasses->SetOwner(kTRUE);
156
157 DefineOutput(1, TList::Class());
158
159 DefineCentralityClasses(centralities);
160}
161
162//_____________________________________________________________________________
163AliAnalysisTaskMuMu::AliAnalysisTaskMuMu(Bool_t fromESD, TList* triggerClasses, const char* beamYear, TArrayF* centralities)
164: AliAnalysisTaskSE(Form("AliAnalysisTaskMuMu-from%s",fromESD ? "ESD":"AOD")),
165fOutput(0),
166fHistogramCollection(0),
167fEventCounters(0),
168fMuonTrackCuts(0x0),
169fPrecomputedTrackMasks(),
170fIsFromESD(fromESD),
171fIsDynamicTriggerClasses(kFALSE),
172fShouldSeparatePlusAndMinus(kFALSE),
173fIsHistogrammingDisabled(kFALSE),
174fBeamYear(beamYear),
175fCentralityLimits(),
176fTriggerClasses(new THashList),
177fSingleTrackCutNames(0x0),
178fPairTrackCutNames(0x0),
179fCentralityNames(new TObjArray),
180fEventCutNames(0x0),
91dbbc20 181fUseBackgroundTriggers(kFALSE),
182fTriggerInputBitMap()
2f331ac9 183{
184 /// Constructor with a predefined list of triggers to consider
185
186 fTriggerClasses->SetOwner(kTRUE);
187
188 DefineOutput(1, TList::Class());
189
190 TObjString* tname;
191 TIter next(triggerClasses);
192
193 while ( ( tname = static_cast<TObjString*>(next()) ) )
194 {
195 fTriggerClasses->Add(new TObjString(*tname));
196 if ( !beamYear )
197 {
198 if ( tname->String().BeginsWith("CMB") )
199 {
200 fBeamYear = "PbPb2010";
201 }
202 if ( tname->String().BeginsWith("CPBI") )
203 {
204 fBeamYear = "PbPb2011";
205 }
206 }
207 }
208
209 DefineCentralityClasses(centralities);
210}
211
212//_____________________________________________________________________________
213AliAnalysisTaskMuMu::~AliAnalysisTaskMuMu()
214{
215 /// dtor
216
217 if (fOutput && ! AliAnalysisManager::GetAnalysisManager()->IsProofMode())
218 {
219 delete fOutput;
220 }
221
222 delete fMuonTrackCuts;
223
224 delete fTriggerClasses;
225 delete fSingleTrackCutNames;
226 delete fPairTrackCutNames;
227 delete fCentralityNames;
228 delete fEventCutNames;
229}
230
231//_____________________________________________________________________________
91dbbc20 232void AliAnalysisTaskMuMu::AddTriggerClasses(const char* triggerlist, const char* sep)
2f331ac9 233{
91dbbc20 234 /// Given a list of trigger names (triggerlist) separated by sep
2f331ac9 235 /// add those which we don't know yet (selecting only the few ones
236 /// we're interested in, e.g CINT* CMU* CMB*
237 ///
238
239 TString slist(triggerlist);
240
91dbbc20 241 TObjArray* a = slist.Tokenize(sep);
2f331ac9 242 TObjString* s;
243 TIter next(a);
244
245 while ( ( s = static_cast<TObjString*>(next()) ) )
246 {
247 TString trigger(s->String());
248 Bool_t add(kFALSE);
249
250 if ( trigger.Contains("WU") )
251 {
252 add = kFALSE;
253 continue;
254 }
255
256 // MB triggers
257 if ( trigger.BeginsWith("CMB") && trigger.Contains("-B-") ) add = kTRUE; // MB PbPb 2010
258
259 if ( trigger.BeginsWith("CPBI") && trigger.Contains("-B-") ) add = kTRUE; // MB PbPb 2011
260 if ( trigger.BeginsWith("CVHN") && trigger.Contains("-B-") ) add = kTRUE; // Central PbPb 2011
261 if ( trigger.BeginsWith("CVLN") && trigger.Contains("-B-") ) add = kTRUE; // Semi Central PbPb 2011
262
263 if ( ( trigger.BeginsWith("CINT") /* || trigger.BeginsWith("CTRUE") */ ) && TriggerSBACECondition(trigger) ) add = kTRUE;
264
265 // main muon triggers
266 if ( ( trigger.BeginsWith("CMU") || trigger.BeginsWith("CMS") || trigger.BeginsWith("CML") ) &&
267 TriggerSBACECondition(trigger) ) add = kTRUE;
268
269 // muon triggers not in coincidence with anything else
270 if ( trigger.BeginsWith("C0MUL") ) add = kTRUE;
271
272 /*
273 if ( trigger.BeginsWith("CMUP") ||
274 trigger.BeginsWith("CCUP") )
275 {
276 // discard ultra-peripheral triggers
277 add = kFALSE;
278 }
279 */
280
281 // triggers in Monte-Carlo
282 if ( trigger.BeginsWith("MB") ||
283 trigger.BeginsWith("SC") ||
284 ( trigger.BeginsWith("CE") && !trigger.BeginsWith("CEMC") ) ||
285 trigger.BeginsWith("DMU") ||
286 trigger.BeginsWith("DML") ||
287 trigger.BeginsWith("MU") ||
288 trigger.BeginsWith("CMSNGL") ||
289 trigger.BeginsWith("CMLK") )
290 {
291 add = kTRUE; // for MC
292 }
293
294 if ( add && !fTriggerClasses->FindObject(trigger.Data()) )
295 {
296 AliDebug(1,Form("Adding %s to considered trigger classes",trigger.Data()));
297 fTriggerClasses->Add(new TObjString(trigger));
298 }
299
300 }
301
302 delete a;
303}
304
305//_____________________________________________________________________________
306void AliAnalysisTaskMuMu::AddPairCut(const char* cutName, UInt_t maskForOneOrBothTrack, UInt_t maskForTrackPair)
307{
308 /// Add a cut for a pair.
309 /// maskForOneOrBothTrack is the mask of cuts that at least one track must satisfy
310 /// maskForTrackPair is the mask of cuts that *both* tracks must satisfy.
311 /// if maskForTrackPair is 0, then no (extra) condition is applied to the pair
312
313 if ( !fPairTrackCutNames )
314 {
315 fPairTrackCutNames = new TObjArray;
316 fPairTrackCutNames->SetOwner(kTRUE);
317 }
318 fPairTrackCutNames->Add(new AliAnalysisTaskMuMu::PairCut(cutName,maskForOneOrBothTrack,maskForTrackPair));
319}
320
321//_____________________________________________________________________________
322void AliAnalysisTaskMuMu::AddSingleCut(const char* name, UInt_t mask)
323{
324 /// Add a cut for single tracks
325 if ( !fSingleTrackCutNames )
326 {
327 fSingleTrackCutNames = new TObjArray;
328 fSingleTrackCutNames->SetOwner(kTRUE);
329 }
330 TObjString* oname = new TObjString(Form("s%s",name));
331 oname->SetUniqueID(mask);
332 fSingleTrackCutNames->Add(oname);
333}
334
335//_____________________________________________________________________________
336void AliAnalysisTaskMuMu::AddEventCut(const char* name, UInt_t mask)
337{
338 /// Add a cut at event level
339 if ( !fEventCutNames )
340 {
341 fEventCutNames = new TObjArray;
342 fEventCutNames->SetOwner(kTRUE);
343 }
344 TObjString* oname = new TObjString(Form("%s",name));
345 oname->SetUniqueID(mask);
346 fEventCutNames->Add(oname);
347}
348
349
350//_____________________________________________________________________________
351void AliAnalysisTaskMuMu::AssertHistogramCollection(const char* physics, const char* triggerClassName)
352{
353 // insure that a given set of histogram is created
354 if (!fHistogramCollection->Histo(physics,triggerClassName,DefaultCentralityName(),"Zvertex"))
355 {
356 FillHistogramCollection(physics,triggerClassName);
357 }
358}
359
360//_____________________________________________________________________________
361void AliAnalysisTaskMuMu::BeautifyHistos()
362{
363 /// Put the titles, marker sizes, color, etc...
364
365 TIter next(fHistogramCollection->CreateIterator());
366 TH1* h;
367
368 while ( ( h = static_cast<TH1*>(next()) ) )
369 {
370 TString name(h->GetName());
371
372 if ( name.Contains("Plus") )
373 {
374 h->SetMarkerStyle(kFullCircle);
375 h->SetMarkerColor(kBlue);
376 h->SetLineColor(kBlue);
377 }
378 else
379 {
380 h->SetMarkerStyle(kOpenCircle);
381 h->SetMarkerColor(kRed);
382 h->SetLineColor(kRed);
383 }
384 }
385}
386
387//_____________________________________________________________________________
388const char*
389AliAnalysisTaskMuMu::CentralityName(Double_t centrality) const
390{
391 /// Get centrality name corresponding to the floating ^point value
392
393 if ( centrality > 0 && centrality <= 100.0 )
394 {
395 return Form("CENT%02d",TMath::Nint(centrality));
396 }
397 else
398 {
399 return DefaultCentralityName();
400 }
401}
402
403//_____________________________________________________________________________
404Bool_t
405AliAnalysisTaskMuMu::CheckTriggerClass(const TString& toCheck,
91dbbc20 406 const TString& firedTriggerClasses,
407 UInt_t l0Inputs) const
2f331ac9 408{
91dbbc20 409 // Check if the "toCheck" class (or logical combination of classes and L0 inputs)
2f331ac9 410 // are within the "firedTriggerClasses"
411
412 Bool_t ok(kFALSE);
413
414 TString tn(toCheck);
415 TString comp("");
416
417 if ( tn.Contains("(") || tn.Contains(")") || tn.Contains("&") || tn.Contains("|") )
418 {
419 comp=tn;
420 tn.ReplaceAll("(","");
421 tn.ReplaceAll(")","");
422 tn.ReplaceAll("&","");
423 tn.ReplaceAll("|","");
424 TObjArray* a = tn.Tokenize(" ");
425 TIter nextA(a);
426 TObjString* an;
427 while ( ( an = static_cast<TObjString*>(nextA()) ) )
428 {
429 // AliDebug(1,Form(" an %s => %d",an->String().Data(),firedTriggerClasses.Contains(an->String().Data())));
91dbbc20 430 if ( an->String().BeginsWith("0") )
431 {
432 // that's an input
433 UInt_t bit = GetTriggerInputBitMaskFromInputName(an->String().Data());
434 comp.ReplaceAll(an->String().Data(),( (l0Inputs & (bit)) == bit) ? "1" : "0");
435 }
436 else
437 {
438 comp.ReplaceAll(an->String().Data(),Form("%d",firedTriggerClasses.Contains(an->String().Data())));
439 }
2f331ac9 440 }
441 delete a;
442
443 TFormula formula("TriggerClassFormulaCheck", comp.Data());
444 if (formula.Compile() > 0)
445 {
446 AliError(Form("Could not evaluate formula %s",comp.Data()));
447 }
448 else
449 {
450 ok = formula.Eval(0);
451 }
452 }
453 else
454 {
455 ok = firedTriggerClasses.Contains(toCheck);
456 }
457
458 AliDebug(1,Form("tname %s => %d comp=%s",toCheck.Data(),ok,comp.Data()));
459
460 return ok;
461}
462
463//_____________________________________________________________________________
464void
465AliAnalysisTaskMuMu::CreateSingleHisto(const char* physics,
466 const char* triggerClassName,
467 const char* hname, const char* htitle,
468 Int_t nbinsx, Double_t xmin, Double_t xmax,
469 Int_t nbinsy, Double_t ymin, Double_t ymax,
470 Bool_t separatePlusAndMinus) const
471{
472 /// Append histograms for single track to our histogram collection
473
474 if ( separatePlusAndMinus )
475 {
476 const char* suffix[] = { "Plus", "Minus" };
477 const char* symbol[] = { "+", "-" };
478
479 for ( Int_t i = 0; i < 2; ++i )
480 {
481 TString shtitle(htitle);
482 TString shname(hname);
483
484 shtitle.ReplaceAll("#mu",Form("#mu^{%s}",symbol[i]));
485
486 shname += suffix[i];
487
488 CreateHisto(fSingleTrackCutNames,physics,triggerClassName,shname.Data(),shtitle.Data(),
489 nbinsx,xmin,xmax,nbinsy,ymin,ymax);
490 }
491 }
492 else
493 {
494 CreateHisto(fSingleTrackCutNames,physics,triggerClassName,hname,htitle,
495 nbinsx,xmin,xmax,nbinsy,ymin,ymax);
496 }
497}
498
499//_____________________________________________________________________________
500void
501AliAnalysisTaskMuMu::CreatePairHisto(const char* physics,
502 const char* triggerClassName,
503 const char* hname, const char* htitle,
504 Int_t nbinsx, Double_t xmin, Double_t xmax,
505 Int_t nbinsy, Double_t ymin, Double_t ymax) const
506{
507 /// Append histograms for track pairs to our histogram collection
508
509 CreateHisto(fPairTrackCutNames,physics,triggerClassName,hname,htitle,
510 nbinsx,xmin,xmax,nbinsy,ymin,ymax);
511}
512
513//_____________________________________________________________________________
514void
515AliAnalysisTaskMuMu::CreateEventHisto(const char* physics,
516 const char* triggerClassName,
517 const char* hname, const char* htitle,
518 Int_t nbinsx, Double_t xmin, Double_t xmax,
519 Int_t nbinsy, Double_t ymin, Double_t ymax) const
520{
521 /// Append histograms at the event level
522
523 TIter next(fCentralityNames);
524 TObjString* cent;
525
526 while ( ( cent = static_cast<TObjString*>(next()) ) )
527 {
528 TH1* h(0x0);
529
530 if ( nbinsy > 0 )
531 {
532 h = new TH2F(hname,htitle,nbinsx,xmin,xmax,nbinsy,ymin,ymax);
533 }
534 else
535 {
536 h = new TH1F(hname,htitle,nbinsx,xmin,xmax);
537 }
538
539 fHistogramCollection->Adopt(physics,triggerClassName,cent->String().Data(),h);
540 }
541}
542
543//_____________________________________________________________________________
544void
545AliAnalysisTaskMuMu::CreateHisto(TObjArray* array,
546 const char* physics,
547 const char* triggerClassName,
548 const char* hname, const char* htitle,
549 Int_t nbinsx, Double_t xmin, Double_t xmax,
550 Int_t nbinsy, Double_t ymin, Double_t ymax) const
551{
552 /// Create a bunch of histograms for all centralities
553 /// FIXME: have a way to specify the histo precision (i.e. F vs I vs S ...)
554
555 TIter next(array);
556 TObjString* tcn;
557 while ( ( tcn = static_cast<TObjString*>(next()) ) )
558 {
559 TIter nextCent(fCentralityNames);
560 TObjString* cent;
561
562 while ( ( cent = static_cast<TObjString*>(nextCent()) ) )
563 {
564 TH1* h(0x0);
565
566 if ( nbinsy > 0 )
567 {
568 h = new TH2F(hname,htitle,nbinsx,xmin,xmax,nbinsy,ymin,ymax);
569 }
570 else
571 {
572 h = new TH1F(hname,htitle,nbinsx,xmin,xmax);
573 }
574
575 fHistogramCollection->Adopt(physics,triggerClassName,cent->String().Data(),tcn->String().Data(),h);
576 }
577 }
578}
579
580//_____________________________________________________________________________
581const char*
582AliAnalysisTaskMuMu::DefaultCentralityName() const
583{
584 /// Get default centrality name
585 if ( !fBeamYear.Contains("pp") ) return "CENTX";
586 else return "PP";
587}
588
589//_____________________________________________________________________________
590void AliAnalysisTaskMuMu::DefineCentralityClasses(TArrayF* centralities)
591{
592 /// Define the default centrality classes that will be used.
593
594 if ( !fBeamYear.Contains("pp") )
595 {
596 if ( !centralities )
597 {
598 // default values
599 fCentralityLimits.push_back(10.0);
600 fCentralityLimits.push_back(30.0);
601 fCentralityLimits.push_back(50.0);
602 fCentralityLimits.push_back(80.0);
603 }
604 else
605 {
606 for ( Int_t i = 0; i < centralities->GetSize(); ++i )
607 {
608 fCentralityLimits.push_back(centralities->At(i));
609 }
610 }
611 }
612
613 for ( std::vector<double>::size_type i = 0; i < fCentralityLimits.size(); ++i )
614 {
615 Double_t limit = fCentralityLimits[i];
616 fCentralityNames->Add(new TObjString(CentralityName(limit)));
617 }
618
619 fCentralityNames->Add(new TObjString(DefaultCentralityName()));
620 fCentralityNames->SetOwner(kTRUE);
621}
622
623//_____________________________________________________________________________
624void AliAnalysisTaskMuMu::EAComputeTrackMasks(const AliVEvent& event)
625{
626 // compute the track masks for the event
627
628 fPrecomputedTrackMasks.Reset();
629 Int_t n = EAGetNumberOfMuonTracks(event);
630 fPrecomputedTrackMasks.Set(n);
631
632 if ( event.IsA() == AliESDEvent::Class() )
633 {
634 const AliESDEvent& esd = static_cast<const AliESDEvent&>(event);
635 AliESDEvent& ncesd = const_cast<AliESDEvent&>(esd);
636
637 for ( Int_t i = 0; i < n; ++i )
638 {
639 ComputeTrackMask(*(ncesd.GetMuonTrack(i)),i);
640 }
641 }
642 else if ( event.IsA() == AliAODEvent::Class() )
643 {
644 const AliAODEvent& aod = static_cast<const AliAODEvent&>(event);
645 for ( Int_t i = 0; i < n; ++i )
646 {
647 ComputeTrackMask(*(aod.GetTrack(i)),i);
648 }
649 }
650
651}
652
653//_____________________________________________________________________________
654TString AliAnalysisTaskMuMu::EAGetFiredTriggerClasses(const AliVEvent& event) const
655{
656 // hack (this method should really be part of AliVEvent itself) FIXME
657
658 if ( event.IsA() == AliESDEvent::Class() )
659 {
660 return static_cast<const AliESDEvent&>(event).GetFiredTriggerClasses();
661 }
662 else if ( event.IsA() == AliAODEvent::Class() )
663 {
664 return static_cast<const AliAODEvent&>(event).GetFiredTriggerClasses();
665 }
666 else
667 {
668 AliError(Form("Unknown class for the event = %s",event.ClassName()));
669 return "";
670 }
671}
672
673//_____________________________________________________________________________
674UInt_t AliAnalysisTaskMuMu::EAGetL0TriggerInputs(const AliVEvent& event) const
675{
676 // Get the list of level 0 trigger inputs
677
678 if ( event.IsA() == AliESDEvent::Class() )
679 {
91dbbc20 680 return static_cast<const AliESDEvent&>(event).GetHeader()->GetL0TriggerInputs();
2f331ac9 681 }
682 else if ( event.IsA() == AliAODEvent::Class() )
683 {
684 return static_cast<const AliAODEvent&>(event).GetHeader()->GetL0TriggerInputs();
685 }
686 else
687 {
688 AliError(Form("Unknown class for the event = %s",event.ClassName()));
689 return 0;
690 }
691
692 return 0;
693}
694
695//_____________________________________________________________________________
696Int_t AliAnalysisTaskMuMu::EAGetNumberOfMuonTracks(const AliVEvent& event) const
697{
698 // get the number of muon tracks from the event
699 if ( event.IsA() == AliESDEvent::Class() )
700 {
91dbbc20 701 Int_t n(0);
702 const AliESDEvent& esd = static_cast<const AliESDEvent&>(event);
703 for ( Int_t i = 0; i < esd.GetNumberOfMuonTracks(); ++i )
704 {
705 AliESDMuonTrack* muon = const_cast<AliESDEvent&>(esd).GetMuonTrack(i);
706 if ( muon->ContainTrackerData() )
707 {
708 ++n;
709 }
710 }
711 return n;
2f331ac9 712 }
713 else if ( event.IsA() == AliAODEvent::Class() )
714 {
715 return static_cast<const AliAODEvent&>(event).GetNumberOfTracks();
716 }
717 else
718 {
719 AliError(Form("Unknown class for the event = %s",event.ClassName()));
720 return 0;
721 }
722}
723
724//_____________________________________________________________________________
725AliVParticle* AliAnalysisTaskMuMu::EAGetTrack(const AliVEvent& event, Int_t i) const
726{
727 // get the i-th track from the event
728 if ( event.IsA() == AliESDEvent::Class() )
729 {
730 AliESDMuonTrack* track = static_cast<AliESDEvent&>(const_cast<AliVEvent&>(event)).GetMuonTrack(i);
731 if ( track->ContainTrackerData() )
732 {
733 return track;
734 }
735 return 0x0;
736 }
737 else if ( event.IsA() == AliAODEvent::Class() )
738 {
739 AliAODTrack* part = static_cast<const AliAODEvent&>(event).GetTrack(i);
740 if ( part->IsMuonTrack() )
741 {
742 return part;
743 }
744 return 0x0;
745 }
746 else
747 {
748 AliError(Form("Unknown class for the event = %s",event.ClassName()));
749 }
750 return 0x0;
751}
752
753//_____________________________________________________________________________
754Double_t AliAnalysisTaskMuMu::EAGetTrackDCA(const AliVParticle& track) const
755{
756 // Get track DCA
757
758 Double_t xdca(1E9);
759 Double_t ydca(xdca);
760
761 if ( track.IsA() == AliAODTrack::Class() )
762 {
763 xdca = static_cast<const AliAODTrack&>(track).XAtDCA();
764 ydca = static_cast<const AliAODTrack&>(track).YAtDCA();
765 }
766 else if ( track.IsA() == AliESDMuonTrack::Class() )
767 {
768 xdca = static_cast<const AliESDMuonTrack&>(track).GetNonBendingCoorAtDCA();
769 ydca = static_cast<const AliESDMuonTrack&>(track).GetBendingCoorAtDCA();
770 }
771 else
772 {
773 AliError(Form("Unknown track class: %s",track.ClassName()));
774 }
775
776 return TMath::Sqrt(xdca*xdca+ydca*ydca);
777}
778
779//_____________________________________________________________________________
780Double_t AliAnalysisTaskMuMu::EAGetTrackNormalizedChi2(const AliVParticle& track) const
781{
782 // get the chi2 per NDF of the track
783
784 Double_t chi2(1E9);
785
786 if ( track.IsA() == AliAODTrack::Class() )
787 {
788 chi2 = static_cast<const AliAODTrack&>(track).Chi2perNDF();
789 }
790 else if ( track.IsA() == AliESDMuonTrack::Class() )
791 {
792 chi2 = static_cast<const AliESDMuonTrack&>(track).GetNormalizedChi2();
793 }
794 else
795 {
796 AliError(Form("Unknown track class: %s",track.ClassName()));
797 }
798
799 return chi2;
800
801}
802
803//_____________________________________________________________________________
804Double_t AliAnalysisTaskMuMu::EAGetTrackChi2MatchTrigger(const AliVParticle& track) const
805{
806 /// Get the Chi2 of the matching MCH - MTR
807
808 if ( track.IsA() == AliAODTrack::Class() )
809 {
810 return static_cast<const AliAODTrack&>(track).GetChi2MatchTrigger();
811 }
812 else
813 {
814 return static_cast<const AliESDMuonTrack&>(track).GetChi2MatchTrigger();
815 }
816}
817
818//_____________________________________________________________________________
819Double_t AliAnalysisTaskMuMu::GetTrackTheta(const AliVParticle& track) const
820{
821 // Get track theta (in radian)
822
823 return TMath::ATan(EAGetTrackRabs(track)/AbsZEnd());
824}
825
826//_____________________________________________________________________________
827Double_t AliAnalysisTaskMuMu::EAGetTrackRabs(const AliVParticle& track) const
828{
829 // Get track DCA
830
831 Double_t rabs(1E9);
832
833 if ( track.IsA() == AliAODTrack::Class() )
834 {
835 rabs = static_cast<const AliAODTrack&>(track).GetRAtAbsorberEnd();
836 }
837 else if ( track.IsA() == AliESDMuonTrack::Class() )
838 {
839 rabs = static_cast<const AliESDMuonTrack&>(track).GetRAtAbsorberEnd();
840 }
841 else
842 {
843 AliError(Form("Unknown track class: %s",track.ClassName()));
844 }
845
846 return TMath::ATan(rabs)/AbsZEnd();
847}
848
849//_____________________________________________________________________________
850Bool_t AliAnalysisTaskMuMu::EAGetTZEROFlags(const AliVEvent& event, Bool_t& backgroundFlag, Bool_t& pileupFlag, Bool_t& satelliteFlag) const
851{
852 // get the TZERO decisions
853 // return false if there's no tzero information in this event
854
855 Bool_t rv(kFALSE);
856
857 if ( event.IsA() == AliESDEvent::Class() )
858 {
859 const AliESDTZERO* tzero = static_cast<AliESDEvent&>(const_cast<AliVEvent&>(event)).GetESDTZERO();
860 if ( tzero )
861 {
862 backgroundFlag = tzero->GetBackgroundFlag();
863 pileupFlag = tzero->GetPileupFlag();
864 satelliteFlag = tzero->GetSatellite();
865 rv = kTRUE;
866 }
867 }
868 else if ( event.IsA() == AliAODEvent::Class() )
869 {
870 AliAODTZERO* tzero = static_cast<const AliAODEvent&>(event).GetTZEROData();
871 if ( tzero )
872 {
873 backgroundFlag = tzero->GetBackgroundFlag();
874 pileupFlag = tzero->GetPileupFlag();
875 satelliteFlag = tzero->GetSatellite();
876 rv = kTRUE;
877 }
878 }
879 else
880 {
881 AliError(Form("Unknown class for the event = %s",event.ClassName()));
882 }
883
884 return rv;
885
886}
887
888//_____________________________________________________________________________
889Bool_t
890AliAnalysisTaskMuMu::AtLeastOneMBTrigger(const TString& firedTriggerClasses) const
891{
892 // whether or not we have a least one MB trigger in the fired trigger classes
893 static TObjArray* triggerList = GetMBTriggerList();
894 TIter next(triggerList);
895 TObjString* str;
896
897 while ( ( str = static_cast<TObjString*>(next())))
898 {
899 if ( firedTriggerClasses.Contains(str->String().Data())) return kTRUE;
900 }
901
902 return kFALSE;
903}
904
905//_____________________________________________________________________________
906Bool_t
907AliAnalysisTaskMuMu::AtLeastOneMuonTrigger(const TString& firedTriggerClasses) const
908{
909 // whether or not we have a least one muon trigger in the fired trigger classes
910 static TObjArray* triggerList = GetMuonTriggerList();
911 TIter next(triggerList);
912 TObjString* str;
913
914 while ( ( str = static_cast<TObjString*>(next())))
915 {
916 if ( firedTriggerClasses.Contains(str->String().Data())) return kTRUE;
917 }
918
919 return kFALSE;
920}
921
922//_____________________________________________________________________________
923Bool_t
924AliAnalysisTaskMuMu::AtLeastOneEmcalTrigger(const TString& firedTriggerClasses) const
925{
926 // whether or not we have a least one emcal trigger in the fired trigger classes
927 static TObjArray* triggerList = GetEmcalTriggerList();
928 TIter next(triggerList);
929 TObjString* str;
930
931 while ( ( str = static_cast<TObjString*>(next())))
932 {
933 if ( firedTriggerClasses.Contains(str->String().Data())) return kTRUE;
934 }
935
936 return kFALSE;
937}
938
939//_____________________________________________________________________________
940void AliAnalysisTaskMuMu::Fill(const char* eventtype,
941 TObjString* tname,
942 const char* centrality,
943 float fcent,
944 const AliVEvent& event)
945{
946 // Fill one set of histograms
947
948 TString seventtype(eventtype);
949 seventtype.ToLower();
950
951 fEventCounters->Count(Form("event:%s/trigger:%s/run:%d", seventtype.Data(), tname->GetName(), fCurrentRunNumber));
952
953 if ( !IsHistogrammingDisabled() )
954 {
955 AssertHistogramCollection(eventtype,tname->String().Data());
956 FillHistos(eventtype,tname->String().Data(),centrality,event);
957 fHistogramCollection->Histo(eventtype,tname->String().Data(),"Centrality")->Fill(fcent);
958 }
959}
960
961//_____________________________________________________________________________
962void AliAnalysisTaskMuMu::FillHistogramCollection(const char* physics, const char* triggerClassName)
963{
964 /// Actually create the histograms for phyics/triggerClassName
965
966 AliDebug(1,Form("(%s,%s)",physics,triggerClassName));
967
968 Double_t ptMin = 0;
969 Double_t ptMax = 12*3;
970 Int_t nbinsPt = GetNbins(ptMin,ptMax,0.5);
971 Double_t pMin = 0;
972 Double_t pMax = 100*3;
973 Int_t nbinsP = GetNbins(pMin,pMax,2.0);
974 Double_t etaMin = -5;
975 Double_t etaMax = -2;
976 Int_t nbinsEta = GetNbins(etaMin,etaMax,0.05);
977
978 Double_t rapidityMin = -5;
979 Double_t rapidityMax = -2;
980 Int_t nbinsRapidity = GetNbins(rapidityMin,rapidityMax,0.05);
981
982 CreateSingleHisto(physics,triggerClassName,"Chi2MatchTrigger","Chi2 Match Trigger",72,0,72);
983
984 CreateSingleHisto(physics,triggerClassName,"EtaRapidityMu", "Eta distribution vs Rapidity for #mu", nbinsRapidity,rapidityMin,rapidityMax,nbinsEta,etaMin,etaMax, fShouldSeparatePlusAndMinus);
985
986 CreateSingleHisto(physics,triggerClassName,"PtEtaMu", "P_{T} distribution vs Eta for #mu", nbinsEta,etaMin,etaMax, nbinsPt,ptMin,ptMax,fShouldSeparatePlusAndMinus);
987
988 CreateSingleHisto(physics,triggerClassName,"PtRapidityMu", "P_{T} distribution vs Rapidity for #mu", nbinsRapidity,rapidityMin,rapidityMax, nbinsPt,ptMin,ptMax,fShouldSeparatePlusAndMinus);
989
990 CreateSingleHisto(physics,triggerClassName,"PEtaMu", "P distribution for #mu",nbinsEta,etaMin,etaMax,nbinsP,pMin,pMax,fShouldSeparatePlusAndMinus);
991
992 Double_t chi2min = 0;
993 Double_t chi2max = 20;
994 Int_t nbinchi2 = GetNbins(chi2min,chi2max,0.05);
995
996 CreateSingleHisto(physics, triggerClassName, "Chi2Mu", "chisquare per NDF #mu", nbinchi2, chi2min, chi2max,fShouldSeparatePlusAndMinus);
997
998 Double_t minvMin = 0;
999 Double_t minvMax = 16;
1000 Int_t nMinvBins = GetNbins(minvMin,minvMax,0.025);
1001
1002 CreatePairHisto(physics,triggerClassName,"Chi12","Chi2MatchTrigger of muon 1 vs muon 2",72,0,72,72,0,72);
1003 CreatePairHisto(physics,triggerClassName,"Rabs12","Rabs of muon 1 vs muon ",100,0,100,100,0,100);
1004
1005 CreatePairHisto(physics,triggerClassName,"MinvUSPt", "#mu+#mu- inv. mass vs Pt",nbinsPt,ptMin,ptMax,nMinvBins,minvMin,minvMax);
1006 CreatePairHisto(physics,triggerClassName,"MinvUSRapidity", "#mu+#mu- inv. mass vs rapidity",nbinsRapidity,rapidityMin,rapidityMax,nMinvBins,minvMin,minvMax);
1007 CreatePairHisto(physics,triggerClassName,"USPtRapidity", "#mu+#mu- Pt vs Rapidity",nbinsRapidity,rapidityMin,rapidityMax,nbinsPt,ptMin,ptMax);
1008
1009 CreatePairHisto(physics,triggerClassName,"MinvLSPt", "unlike sign inv. mass vs Pt",nbinsPt,ptMin,ptMax,nMinvBins,minvMin,minvMax);
1010
91dbbc20 1011 Double_t xmin = -35;
1012 Double_t xmax = +35;
1013 Int_t nbins = GetNbins(xmin,xmax,0.4);
2f331ac9 1014
1015 CreateEventHisto(physics,triggerClassName,"Zvertex","z vertex",nbins,xmin,xmax);
1016
91dbbc20 1017 CreateEventHisto(physics,triggerClassName,"T0Zvertex","T0 zvertex",nbins,xmin,xmax);
1018
2f331ac9 1019 xmin = -5;
1020 xmax = 5;
1021 nbins = GetNbins(xmin,xmax,0.01);
1022
1023 CreateEventHisto(physics,triggerClassName,"Xvertex","x vertex",nbins,xmin,xmax);
1024 CreateEventHisto(physics,triggerClassName,"Yvertex","y vertex",nbins,xmin,xmax);
1025
1026// CreateEventHisto(physics,triggerClassName,"YXvertex","y vs x vertex",nbins,xmin,xmax,nbins,xmin,xmax);
1027
1028// if (!fIsFromESD)
1029// {
1030//
1031// CreateEventHisto(physics,triggerClassName,"PileUpZvertex","pileup z vertex",nbins,xmin,xmax);
1032//
1033// CreateEventHisto(physics,triggerClassName,"PileUpXvertex","pileup x vertex",nbins,xmin,xmax);
1034// CreateEventHisto(physics,triggerClassName,"PileUpYvertex","pileup y vertex",nbins,xmin,xmax);
1035//
1036// CreateEventHisto(physics,triggerClassName,"PileUpYXvertex","pileup y vs x vertex",nbins,xmin,xmax,
1037// nbins,xmin,xmax);
1038//
1039// }
1040
1041 CreateEventHisto(physics,triggerClassName,"Nevents","number of events",2,-0.5,1.5);
1042
1043 xmin = 0;
1044 xmax = 3564;
1045 nbins = GetNbins(xmin,xmax,1.0);
1046
1047 CreateEventHisto(physics,triggerClassName,"BCX","bunch-crossing ids",nbins,xmin-0.5,xmax-0.5);
1048
1049 CreateEventHisto(physics,triggerClassName,"BCXD","bunch-crossing distances",nbins,xmin-0.5,xmax-0.5);
1050
1051 xmin = -200;
1052 xmax = +200;
1053 nbins = GetNbins(xmin,xmax,1.0);
1054
1055 xmin = 0;
1056 xmax = 150;
1057 nbins = GetNbins(xmin,xmax,2.0);
1058
1059 CreateSingleHisto(physics,triggerClassName,"dcaP23Mu","#mu DCA vs P for 2-3 degrees;P (GeV);DCA (cm)",nbinsP,pMin,pMax,nbins,xmin,xmax,fShouldSeparatePlusAndMinus);
1060
1061 CreateSingleHisto(physics,triggerClassName,"dcaP310Mu","#mu DCA vs P for 3-10 degrees;P (GeV);DCA (cm)",nbinsP,pMin,pMax,nbins,xmin,xmax,fShouldSeparatePlusAndMinus);
1062
1063 CreateSingleHisto(physics,triggerClassName,"dcaPwPtCut23Mu","#mu DCA vs P for 2-3 degrees with Pt Cut;P (GeV);DCA (cm)",nbinsP,pMin,pMax,nbins,xmin,xmax,fShouldSeparatePlusAndMinus);
1064
1065 CreateSingleHisto(physics,triggerClassName,"dcaPwPtCut310Mu","#mu DCA vs P for 3-10 degrees with Pt Cut;P (GeV);DCA (cm)",nbinsP,pMin,pMax,nbins,xmin,xmax,fShouldSeparatePlusAndMinus);
1066
1067 xmin = -30;
1068 xmax = +30;
1069 nbins = GetNbins(xmin,xmax,0.1);
1070
1071 CreateEventHisto(physics,triggerClassName,"V02D","V0C+V0A versus V0A-V0C;Time V0A - V0C (ns);Time V0A+V0C (ns)",nbins,xmin,xmax,nbins,xmin,xmax);
1072
1073 CreateEventHisto(physics,triggerClassName,"V02DwT0BB","V0C+V0A versus V0A-V0C with T0 BB;Time V0A - V0C (ns);Time V0A+V0C (ns)",nbins,xmin,xmax,nbins,xmin,xmax);
1074
1075 CreateEventHisto(physics,triggerClassName,"V02DwT0BG","V0C+V0A versus V0A-V0C with T0 background flag on;Time V0A - V0C (ns);Time V0A+V0C (ns)",nbins,xmin,xmax,nbins,xmin,xmax);
1076
1077 CreateEventHisto(physics,triggerClassName,"V02DwT0PU","V0C+V0A versus V0A-V0C with T0 pile up flag on;Time V0A - V0C (ns);Time V0A+V0C (ns)",nbins,xmin,xmax,nbins,xmin,xmax);
1078
1079 CreateEventHisto(physics,triggerClassName,"V02DwT0SAT","V0C+V0A versus V0A-V0C with T0 satellite flag on;Time V0A - V0C (ns);Time V0A+V0C (ns)",nbins,xmin,xmax,nbins,xmin,xmax);
1080
1081 /*
1082 CreateEventHisto(physics,triggerClassName,"T02D","(T0C+T0A)/2 versus (T0A-T0C)/2;Time (T0A-T0C)/2 (ns);Time (T0A+T0C)/2 (ns)",nbins,xmin,xmax,nbins,xmin,xmax);
91dbbc20 1083 CreateEventHisto(physics,triggerClassName,"T0Flags","T0 flags",3,0,3);
2f331ac9 1084 */
91dbbc20 1085
2f331ac9 1086
1087
1088 TH1* h = new TH1F("Centrality","Centrality",12,-10,110);
1089
1090 fHistogramCollection->Adopt(physics,triggerClassName,h);
1091
1092 xmin = 0;
1093 xmax = 5000;
1094 nbins = GetNbins(xmin,xmax,10);
1095
1096 CreateEventHisto(physics,triggerClassName,"Tracklets","Number of tracklets",nbins,xmin,xmax);
1097
1098}
1099
1100//_____________________________________________________________________________
1101void AliAnalysisTaskMuMu::FillHistosForTrack(const char* physics,
1102 const char* triggerClassName,
1103 const char* centrality,
1104 const AliVParticle& track,
1105 Int_t trackIndex)
1106{
1107 /// Fill histograms for one track
1108
1109 TLorentzVector p(track.Px(),track.Py(),track.Pz(),
1110 TMath::Sqrt(MuonMass2()+track.P()*track.P()));
1111
1112
1113 TString charge("");
1114
1115 if ( ShouldSeparatePlusAndMinus() )
1116 {
1117 if ( track.Charge() < 0 )
1118 {
1119 charge = "Minus";
1120 }
1121 else
1122 {
1123 charge = "Plus";
1124 }
1125 }
1126
1127 UInt_t mask = GetTrackMask(trackIndex);
1128
1129 Double_t dca = EAGetTrackDCA(track);
1130
1131 Double_t theta = GetTrackTheta(track);
1132
1133 TIter next(fSingleTrackCutNames);
1134 TObjString* str;
1135
1136 while ( ( str = static_cast<TObjString*>(next()) ) )
1137 {
1138 Bool_t test = ( ( str->GetUniqueID() & mask ) == str->GetUniqueID() );
1139
1140 if ( test )
1141 {
1142 Histo(physics,triggerClassName,centrality,str->String().Data(),"Chi2MatchTrigger")->Fill(EAGetTrackChi2MatchTrigger(track));
1143
1144 Histo(physics,triggerClassName,centrality,str->String().Data(),Form("EtaRapidityMu%s",charge.Data()))->Fill(p.Rapidity(),p.Eta());
1145
1146 Histo(physics,triggerClassName,centrality,str->String().Data(),Form("PtEtaMu%s",charge.Data()))->Fill(p.Eta(),p.Pt());
1147 Histo(physics,triggerClassName,centrality,str->String().Data(),Form("PtRapidityMu%s",charge.Data()))->Fill(p.Rapidity(),p.Pt());
1148
1149 Histo(physics,triggerClassName,centrality,str->String().Data(),Form("PEtaMu%s",charge.Data()))->Fill(p.Eta(),p.P());
1150
1151
1152 // Histo(physics,triggerClassName,centrality,str->String().Data(),Form("XYdcaMu%s",charge.Data()))->Fill(ydca,xdca);
1153
1154 Histo(physics,triggerClassName,centrality,str->String().Data(),Form("Chi2Mu%s",charge.Data()))->Fill(EAGetTrackNormalizedChi2(track));
1155
1156 if ( theta >= Deg2() && theta < Deg3() )
1157 {
1158 Histo(physics,triggerClassName,centrality,str->String().Data(),Form("dcaP23Mu%s",charge.Data()))->Fill(p.P(),dca);
1159 if ( p.Pt() > 2 )
1160 {
1161 Histo(physics,triggerClassName,centrality,str->String().Data(),Form("dcaPwPtCut23Mu%s",charge.Data()))->Fill(p.P(),dca);
1162 }
1163 }
1164 else if ( theta >= Deg3() && theta < Deg10() )
1165 {
1166 Histo(physics,triggerClassName,centrality,str->String().Data(),Form("dcaP310Mu%s",charge.Data()))->Fill(p.P(),dca);
1167 if ( p.Pt() > 2 )
1168 {
1169 Histo(physics,triggerClassName,centrality,str->String().Data(),Form("dcaPwPtCut310Mu%s",charge.Data()))->Fill(p.P(),dca);
1170 }
1171 }
1172 }
1173 }
1174
1175
1176}
1177
1178//_____________________________________________________________________________
1179void AliAnalysisTaskMuMu::FillHistos(const char* physics, const char* triggerClassName,
1180 const char* centrality, const AliVEvent& event)
1181{
1182 /// Fill histograms for /physics/triggerClassName/centrality
1183
1184 Histo(physics,triggerClassName,centrality,"BCX")->Fill(1.0*event.GetBunchCrossNumber());
1185
1186 Histo(physics,triggerClassName,centrality,"Nevents")->Fill(1.0);
1187
1188 const AliVVertex* vertex = event.GetPrimaryVertex();
1189
1190 if ( vertex )
1191 {
1192 Histo(physics,triggerClassName,centrality,"Xvertex")->Fill(vertex->GetX());
1193 Histo(physics,triggerClassName,centrality,"Yvertex")->Fill(vertex->GetY());
1194 Histo(physics,triggerClassName,centrality,"Zvertex")->Fill(vertex->GetZ());
1195 }
1196
91dbbc20 1197 if ( !fIsFromESD )
1198 {
1199 const AliAODTZERO* tzero = static_cast<const AliAODEvent&>(event).GetTZEROData();
2f331ac9 1200
91dbbc20 1201 if (tzero)
1202 {
1203 Histo(physics,triggerClassName,centrality,"T0Zvertex")->Fill(tzero->GetT0VertexRaw());
1204 }
1205 }
1206 else
1207 {
1208 const AliESDTZERO* tzero = static_cast<const AliESDEvent&>(event).GetESDTZERO();
1209
1210 if (tzero)
1211 {
1212 Histo(physics,triggerClassName,centrality,"T0Zvertex")->Fill(tzero->GetT0zVertex());
1213 }
1214 }
1215
1216 AliVVZERO* vzero = event.GetVZEROData();
1217
2f331ac9 1218 if (vzero)
1219 {
1220 Float_t v0a = vzero->GetV0ATime();
1221 Float_t v0c = vzero->GetV0CTime();
1222
1223 Float_t x = v0a-v0c;
1224 Float_t y = v0a+v0c;
1225
1226 Histo(physics,triggerClassName,centrality,"V02D")->Fill(x,y);
1227
1228 Bool_t background,pileup,satellite;
1229
1230 Bool_t tzero = EAGetTZEROFlags(event,background,pileup,satellite);
1231
1232 if (tzero)
1233 {
1234 if ( background )
1235 {
1236 Histo(physics,triggerClassName,centrality,"V02DwT0BG")->Fill(x,y);
1237 }
1238
1239 if ( pileup )
1240 {
1241 Histo(physics,triggerClassName,centrality,"V02DwT0PU")->Fill(x,y);
1242 }
1243
1244 if ( satellite )
1245 {
1246 Histo(physics,triggerClassName,centrality,"V02DwT0SAT")->Fill(x,y);
1247 }
1248
1249 if ( !background && !pileup && !satellite )
1250 {
1251 Histo(physics,triggerClassName,centrality,"V02DwT0BB")->Fill(x,y);
1252 }
1253 }
1254 }
1255
1256 /* FIXME : how to properly get multiplicity from AOD and ESD consistently ?
1257 is is doable at all ?
1258 Int_t ntracklets(0);
1259 AliAODTracklets* tracklets = aod.GetTracklets();
1260 if ( tracklets )
1261 {
1262 ntracklets = tracklets->GetNumberOfTracklets();
1263 }
1264 Histo(physics,triggerClassName,centrality,"Tracklets")->Fill(ntracklets);
1265 */
1266
1267 // Track loop
1268
1269 Int_t nMuonTracks = EAGetNumberOfMuonTracks(event);
1270
1271 for (Int_t i = 0; i < nMuonTracks; ++i)
1272 {
1273 AliVParticle* tracki = EAGetTrack(event,i);
1274
1275 if (!tracki) continue;
1276
1277 FillHistosForTrack(physics,triggerClassName,centrality,*tracki,i);
1278
1279 TLorentzVector pi(tracki->Px(),tracki->Py(),tracki->Pz(),
1280 TMath::Sqrt(MuonMass2()+tracki->P()*tracki->P()));
1281
1282 for (Int_t j = i+1; j < nMuonTracks; ++j)
1283 {
1284 AliVParticle* trackj = EAGetTrack(event,j);
1285
1286 if (!trackj) continue;
1287
1288 TLorentzVector pj(trackj->Px(),trackj->Py(),trackj->Pz(),
1289 TMath::Sqrt(MuonMass2()+trackj->P()*trackj->P()));
1290
1291 pj += pi;
1292
1293 TIter next(fPairTrackCutNames);
1294 AliAnalysisTaskMuMu::PairCut* str;
1295
1296 UInt_t maski(0),maskj(0),maskij(0);
1297
1298 GetPairMask(*tracki,*trackj,i,j,maski,maskj,maskij);
1299
1300 while ( ( str = static_cast<AliAnalysisTaskMuMu::PairCut*>(next()) ) )
1301 {
1302 UInt_t singleTrackMask = str->MaskForOneOrBothTrack();
1303 UInt_t pairMask = str->MaskForTrackPair();
1304
1305 Bool_t testi = ( ( maski & singleTrackMask ) == singleTrackMask ) ;
1306 Bool_t testj = ( ( maskj & singleTrackMask ) == singleTrackMask ) ;
1307 Bool_t testij(kTRUE);
1308
1309 if (pairMask>0) testij = ( ( maskij & pairMask ) == pairMask ) ;
1310
1311 if ( ( testi || testj ) && testij )
1312 {
1313 Histo(physics,triggerClassName,centrality,str->GetName(),"Chi12")->Fill(EAGetTrackChi2MatchTrigger(*tracki),
1314 EAGetTrackChi2MatchTrigger(*trackj));
1315
1316 Histo(physics,triggerClassName,centrality,str->GetName(),"Rabs12")->Fill(EAGetTrackRabs(*tracki),
1317 EAGetTrackRabs(*trackj));
1318
1319 if ( tracki->Charge() != trackj->Charge() )
1320 {
1321 Histo(physics,triggerClassName,centrality,str->GetName(),"MinvUSPt")->Fill(pj.Pt(),pj.M());
1322 Histo(physics,triggerClassName,centrality,str->GetName(),"MinvUSRapidity")->Fill(pj.Rapidity(),pj.M());
1323 Histo(physics,triggerClassName,centrality,str->GetName(),"USPtRapidity")->Fill(pj.Rapidity(),pj.Pt());
1324 }
1325 else
1326 {
1327 Histo(physics,triggerClassName,centrality,str->GetName(),"MinvLSPt")->Fill(pj.Pt(),pj.M());
1328 }
1329 }
1330 }
1331 }
1332 } //track loop
1333}
1334
1335//_____________________________________________________________________________
1336void AliAnalysisTaskMuMu::FinishTaskOutput()
1337{
1338 /// prune empty histograms BEFORE mergin, in order to save some bytes...
1339
1340 if ( fHistogramCollection )
1341 {
1342 fHistogramCollection->PruneEmptyHistograms();
1343 }
1344}
1345
1346//_____________________________________________________________________________
1347UInt_t AliAnalysisTaskMuMu::GetEventMask(const AliVEvent& event) const
1348{
1349 /// Compute the event mask
1350
1351 /*
1352
1353 kEventAll = all events
1354 kEventPS = physics selected events
1355 kEventTVX = events with 0TVX input present
1356 kEventV0AND = events with 0VBA and 0VBC present
1357 kEventV0UP = a check of events within a narrow square range of v0a+c vs v0a-c
1358 kEventZSPD = events with a vertex computed by SPD
1359 kEventZ7 = events with | zvertex | < 7cm
1360 kEventZ10 = events with | zvertex | < 10 cm
1361 kEventSD2 = events with 0SD2 input present (was for PbPb 2010)
91dbbc20 1362 kEventMSL = events with 0MSL input present
2f331ac9 1363
1364 */
1365
1366 UInt_t m(AliAnalysisTaskMuMu::kEventAll);
1367
1368 if ( PassPhysicsSelection() ) m |= AliAnalysisTaskMuMu::kEventPS;
1369
1370 UInt_t trigger = EAGetL0TriggerInputs(event);
1371
91dbbc20 1372 UInt_t l0TVXBIT = GetTriggerInputBitMaskFromInputName("0TVX");
2f331ac9 1373
1374 if ( ( trigger & (l0TVXBIT) ) == l0TVXBIT )
1375 {
1376 m |= AliAnalysisTaskMuMu::kEventTVX;
1377 }
1378
91dbbc20 1379 UInt_t l0VBABIT = GetTriggerInputBitMaskFromInputName("0VBA");
1380 UInt_t l0VBCBIT = GetTriggerInputBitMaskFromInputName("0VBC");
2f331ac9 1381
1382 if ( ( ( trigger & (l0VBABIT ) ) == l0VBABIT ) &&
1383 ( ( trigger & (l0VBCBIT ) ) == l0VBCBIT ) )
1384 {
1385 m |= AliAnalysisTaskMuMu::kEventV0AND;
1386 }
1387
91dbbc20 1388 UInt_t l0MSLBIT = GetTriggerInputBitMaskFromInputName("0MSL");
1389
1390 if ( ( trigger & (l0MSLBIT) ) == l0MSLBIT )
1391 {
1392 m |= AliAnalysisTaskMuMu::kEventMSL;
1393 }
1394
2f331ac9 1395 if ( fBeamYear == "PbPb2010" )
1396 {
1397 // consider only events with OSM2 fired
91dbbc20 1398 UInt_t sd2 = GetTriggerInputBitMaskFromInputName("0SM2");
2f331ac9 1399 if ( ( trigger & sd2 ) == sd2 )
1400 {
1401 m |= AliAnalysisTaskMuMu::kEventSD2;
1402 }
1403 }
1404
1405 // Bool_t hasPileUp = aod->IsPileupFromSPD(3,0.8);
1406 // Bool_t hasPileUp2 = aod->IsPileupFromSPD(5,0.8);
1407 // Bool_t isIsolated = ( aod->GetBunchCrossNumber() > 1000 && aod->GetBunchCrossNumber() < 2900 );
1408
1409 // Bool_t hasPileUp2 = aod->IsPileupFromSPDInMultBins();
1410
1411 // TIter nextV(aod->GetVertices());
1412 // AliAODVertex* v;
1413 // while ( ( v = static_cast<AliAODVertex*>(nextV())) && !hasPileUp2 )
1414 // {
1415 // if ( v->GetType() == AliAODVertex::kPileupSPD ) hasPileUp2 = kTRUE;
1416 // }
1417
1418 // Bool_t isIsolated = false;//( aod->GetClosestBunchCrossingDistance() > 10 );
1419
1420 const AliVVertex* vertex = event.GetPrimaryVertex();
1421
1422 if ( vertex->IsA() == AliAODVertex::Class() )
1423 {
1424 AliAODVertex* spdVertex = static_cast<const AliAODEvent&>(event).GetPrimaryVertexSPD();
1425
1426 if ( spdVertex && spdVertex->GetNContributors() > 0 )
1427 {
1428 m |= AliAnalysisTaskMuMu::kEventZSPD;
1429 }
1430 }
1431
1432 if ( TMath::Abs(vertex->GetZ()) < 10.0 )
1433 {
1434 m |= AliAnalysisTaskMuMu::kEventZ10;
1435 }
1436
1437 if ( TMath::Abs(vertex->GetZ()) < 7.0 )
1438 {
1439 m |= AliAnalysisTaskMuMu::kEventZ7;
1440 }
1441
1442 AliVVZERO* vzero = event.GetVZEROData();
1443
1444 if (vzero)
1445 {
1446 Float_t v0a = vzero->GetV0ATime();
1447 Float_t v0c = vzero->GetV0CTime();
1448
1449 Float_t x = v0a-v0c;
1450 Float_t y = v0a+v0c;
1451
1452 if ( ( x > 6 && x < 10 ) && y > 20 )
1453 {
1454 m |= AliAnalysisTaskMuMu::kEventV0UP;
1455 }
1456 }
1457
1458 return m;
1459}
1460
91dbbc20 1461//_____________________________________________________________________________
1462UInt_t AliAnalysisTaskMuMu::GetTriggerInputBitMaskFromInputName(const char* inputName) const
1463{
1464 // Get trigger input bit from its name
1465 // FIXME : this should really come directly from the trigger configuration
1466 // object, if only this one would be available in a more convenient
1467 // way than the OCDB (e.g. in RunBasedContainer ?)
1468 //
1469
1470 if ( fTriggerInputBitMap.empty() )
1471 {
1472 // nothing given to us, use the bad bad hard-coded values !
1473
1474 TString sInputName(inputName);
1475
1476 if ( sInputName == "0SM2" ) return (1<<12);
1477
1478
1479 if ( sInputName == "0VBA" ) return (1<<0);
1480 if ( sInputName == "0VBC" ) return (1<<1);
1481 if ( sInputName == "0SMB" ) return (1<<2);
1482 if ( sInputName == "0TVX" ) return (1<<3);
1483 if ( sInputName == "0VGC" ) return (1<<4);
1484 if ( sInputName == "0VGA" ) return (1<<5);
1485 if ( sInputName == "0SH1" ) return (1<<6);
1486 if ( sInputName == "0SH2" ) return (1<<7);
1487 if ( sInputName == "0HPT" ) return (1<<8);
1488 if ( sInputName == "0AMU" ) return (1<<9);
1489 if ( sInputName == "0OB0" ) return (1<<10);
1490 if ( sInputName == "0ASL" ) return (1<<11);
1491 if ( sInputName == "0MSL" ) return (1<<12);
1492 if ( sInputName == "0MSH" ) return (1<<13);
1493 if ( sInputName == "0MUL" ) return (1<<14);
1494 if ( sInputName == "0MLL" ) return (1<<15);
1495 if ( sInputName == "0EMC" ) return (1<<16);
1496 if ( sInputName == "0PH0" ) return (1<<17);
1497 if ( sInputName == "0HWU" ) return (1<<18);
1498 if ( sInputName == "0LSR" ) return (1<<19);
1499 if ( sInputName == "0T0A" ) return (1<<20);
1500 if ( sInputName == "0BPA" ) return (1<<21);
1501 if ( sInputName == "0BPC" ) return (1<<22);
1502 if ( sInputName == "0T0C" ) return (1<<23);
1503
1504 if ( sInputName == "1EJE" ) return (1<<0);
1505 if ( sInputName == "1EGA" ) return (1<<1);
1506 if ( sInputName == "1EJ2" ) return (1<<2);
1507 if ( sInputName == "1EG2" ) return (1<<3);
1508 if ( sInputName == "1PHL" ) return (1<<4);
1509 if ( sInputName == "1PHM" ) return (1<<5);
1510 if ( sInputName == "1PHH" ) return (1<<6);
1511 if ( sInputName == "1HCO" ) return (1<<8);
1512 if ( sInputName == "1HJT" ) return (1<<9);
1513 if ( sInputName == "1HSE" ) return (1<<10);
1514 if ( sInputName == "1DUM" ) return (1<<11);
1515 if ( sInputName == "1HQU" ) return (1<<12);
1516 if ( sInputName == "1H14" ) return (1<<13);
1517 if ( sInputName == "1ZMD" ) return (1<<14);
1518 if ( sInputName == "1ZMB" ) return (1<<16);
1519 if ( sInputName == "1ZED" ) return (1<<17);
1520 if ( sInputName == "1ZAC" ) return (1<<18);
1521 if ( sInputName == "1EJE" ) return (1<<19);
1522
1523 AliError(Form("Don't know this input %s",inputName));
1524
1525 return (1<<31);
1526 }
1527 else
1528 {
1529 std::map<std::string,int>::const_iterator it = fTriggerInputBitMap.find(inputName);
1530 if ( it != fTriggerInputBitMap.end() )
1531 {
1532 return ( 1 << it->second );
1533 }
1534 else
1535 {
1536 AliError(Form("Don't know this input %s",inputName));
1537
1538 return (1<<31);
1539 }
1540 }
1541}
2f331ac9 1542
1543//_____________________________________________________________________________
1544void AliAnalysisTaskMuMu::GetPairMask(const AliVParticle& t1, const AliVParticle& t2,
1545 Int_t trackIndex1, Int_t trackIndex2,
1546 UInt_t& mask1, UInt_t& mask2,
1547 UInt_t& mask12) const
1548{
1549 /// Get the mask of the track pair
1550
1551 mask1 = GetTrackMask(trackIndex1);
1552 mask2 = GetTrackMask(trackIndex2);
1553
1554 mask12 = mask1 & mask2;
1555
1556 if ( PairRapidityCut(t1,t2) ) mask12 |= kPairRapidity;
1557}
1558
1559//_____________________________________________________________________________
1560UInt_t AliAnalysisTaskMuMu::GetTrackMask(Int_t trackIndex) const
1561{
1562 /// Get the mask of all the cuts this track pass
1563
1564 return static_cast<UInt_t>(fPrecomputedTrackMasks.At(trackIndex));
1565}
1566
1567//_____________________________________________________________________________
1568void AliAnalysisTaskMuMu::ComputeTrackMask(const AliVParticle& track, Int_t trackIndex)
1569{
1570 // Compute the track mask
1571 UInt_t m(kAll);
1572
1573 UInt_t selectionMask = fMuonTrackCuts ? fMuonTrackCuts->GetSelectionMask(&track) : 0;
1574
1575 if ( ( selectionMask & AliMuonTrackCuts::kMuThetaAbs ) == AliMuonTrackCuts::kMuThetaAbs ) m |= kRabs;
1576
1577 Double_t angle = GetTrackTheta(track);
1578
1579 if ( angle >= Deg2() && angle < Deg3() ) m |= kDeg23;
1580
1581 if ( angle >= Deg3() && angle < Deg10() ) m |= kDeg310;
1582
1583 if ( selectionMask & AliMuonTrackCuts::kMuEta ) m |= kEta;
1584
1585 Double_t pt = track.Pt();
1586
1587 if ( pt > 1.0 ) m |= kPt1;
1588 if ( pt > 1.2 ) m |= kPt1dot2;
1589 if ( pt > 1.5 ) m |= kPt1dot5;
1590 if ( pt > 2.0 ) m |= kPt2;
1591
1592 if ( track.P() > 10.0 ) m |= kP10;
1593
1594 if ( pt < 4.0 ) m |= kBelowPt;
1595
1596 if ( ( selectionMask & AliMuonTrackCuts::kMuMatchApt ) == AliMuonTrackCuts::kMuMatchApt ) m |= kMatched;
1597
1598 if ( ( selectionMask & AliMuonTrackCuts::kMuMatchLpt ) == AliMuonTrackCuts::kMuMatchLpt ) m |= kMatchedLow;
1599
1600 if ( ( selectionMask & AliMuonTrackCuts::kMuMatchHpt ) == AliMuonTrackCuts::kMuMatchHpt) m |= kMatchedHigh;
1601
1602 if ( ( selectionMask & AliMuonTrackCuts::kMuTrackChiSquare ) == AliMuonTrackCuts::kMuTrackChiSquare ) m |= kChi2;
1603
1604 if ( ( selectionMask & AliMuonTrackCuts::kMuPdca ) == AliMuonTrackCuts::kMuPdca ) m |= kDCA;
1605
1606 if ( EAGetTrackChi2MatchTrigger(track) < 16.0 ) m |= kChi2MatchTrigger;
1607
1608 fPrecomputedTrackMasks.SetAt(m,trackIndex);
1609}
1610
1611
1612//_____________________________________________________________________________
1613TH1* AliAnalysisTaskMuMu::Histo(const char* physics, const char* triggerClassName, const char* histoname)
1614{
1615 /// Get one histo back
1616 return fHistogramCollection ? fHistogramCollection->Histo(physics,triggerClassName,histoname) : 0x0;
1617}
1618
1619//_____________________________________________________________________________
1620TH1* AliAnalysisTaskMuMu::Histo(const char* physics, const char* histoname)
1621{
1622 /// Get one histo back
1623 return fHistogramCollection ? fHistogramCollection->Histo(physics,histoname) : 0x0;
1624}
1625
1626//_____________________________________________________________________________
1627TH1* AliAnalysisTaskMuMu::Histo(const char* physics,
1628 const char* triggerClassName,
1629 const char* what,
1630 const char* histoname)
1631{
1632 /// Get one histo back
1633 return fHistogramCollection ? fHistogramCollection->Histo(physics,triggerClassName,what,histoname) : 0x0;
1634}
1635
1636//_____________________________________________________________________________
1637TH1* AliAnalysisTaskMuMu::Histo(const char* physics,
1638 const char* triggerClassName,
1639 const char* cent,
1640 const char* what,
1641 const char* histoname)
1642{
1643 /// Get one histo back
1644
1645 return fHistogramCollection ? fHistogramCollection->Histo(physics,triggerClassName,cent,what,histoname) : 0x0;
1646}
1647
1648//_____________________________________________________________________________
1649Bool_t AliAnalysisTaskMuMu::IsPP() const
1650{
1651 // whether we're dealing with proton proton collisions
1652 return fBeamYear.Contains("pp");
1653}
1654
1655//_____________________________________________________________________________
1656void AliAnalysisTaskMuMu::MergeCentralities(AliHistogramCollection* histogramCollection)
1657{
1658 /// Merge CENT10 + CENT20 + ... into CENTMB
1659
1660 TList* listA = histogramCollection->CreateListOfKeysA();
1661 TList* listB = histogramCollection->CreateListOfKeysB();
1662 TList* listC = histogramCollection->CreateListOfKeysC();
1663 TList* listD = histogramCollection->CreateListOfKeysD();
1664
1665 if (!listA)
1666 {
1667 AliErrorClass("listA=0x0");
1668 return;
1669 }
1670
1671 if (!listB)
1672 {
1673 AliErrorClass("listB=0x0");
1674 return;
1675 }
1676
1677 if (!listC)
1678 {
1679 AliErrorClass("listC=0x0");
1680 return;
1681 }
1682
1683 if (!listD)
1684 {
1685 AliErrorClass("listD=0x0");
1686 return;
1687 }
1688
1689 for ( Int_t id = 0; id <= listD->GetLast(); ++id )
1690 {
1691 TString keyD = static_cast<TObjString*>(listD->At(id))->String();
1692
1693 for ( Int_t ia = 0; ia <= listA->GetLast(); ++ia )
1694 {
1695 TString keyA = static_cast<TObjString*>(listA->At(ia))->String();
1696
1697 for ( Int_t ib = 0; ib <= listB->GetLast(); ++ib )
1698 {
1699 TString keyB = static_cast<TObjString*>(listB->At(ib))->String();
1700
1701 TList* list = new TList;
1702 list->SetOwner(kTRUE);
1703
1704 AliHistogramCollection* hmerge(0x0);
1705
1706 for ( Int_t ic = 0; ic <= listC->GetLast(); ++ic )
1707 {
1708 TString keyC = static_cast<TObjString*>(listC->At(ic))->String();
1709
1710 if ( keyC != "CENTX" && keyC != "CENTMB" )
1711 {
1712 AliHistogramCollection* hc = histogramCollection->Project(keyA.Data(),keyB.Data(),keyC.Data(),keyD.Data());
1713 if (!hmerge)
1714 {
1715 hmerge = hc;
1716 }
1717 else
1718 {
1719 list->Add(hc);
1720 }
1721 }
1722 }
1723 if (hmerge)
1724 {
1725 hmerge->Merge(list);
1726 TIter next(hmerge->CreateIterator());
1727 TH1* h;
1728 while ( ( h = static_cast<TH1*>(next()) ) )
1729 {
1730 histogramCollection->Adopt(keyA.Data(),keyB.Data(),"CENTMB",keyD.Data(),static_cast<TH1*>(h->Clone()));
1731 }
1732 }
1733 delete list;
1734 delete hmerge;
1735 }
1736 }
1737 }
1738
1739 delete listA;
1740 delete listB;
1741 delete listC;
1742 delete listD;
1743}
1744
1745//_____________________________________________________________________________
1746Double_t AliAnalysisTaskMuMu::MuonMass2() const
1747{
1748 /// A usefull constant
1749 static Double_t m2 = 1.11636129640000012e-02; // using a constant here as the line below is a problem for CINT...
1750 // static Double_t m2 = TDatabasePDG::Instance()->GetParticle("mu-")->Mass()*TDatabasePDG::Instance()->GetParticle("mu-")->Mass();
1751 return m2;
1752}
1753
1754//_____________________________________________________________________________
1755void AliAnalysisTaskMuMu::NotifyRun()
1756{
1757 /// Called at each change of run
1758
1759 AliDebug(1,Form("Run %09d File %s",fCurrentRunNumber,CurrentFileName()));
1760
1761 if (!fMuonTrackCuts)
1762 {
1763 fMuonTrackCuts = new AliMuonTrackCuts;
1764
1765 fMuonTrackCuts->SetAllowDefaultParams(kTRUE );
1766
1767
1768 fMuonTrackCuts->SetFilterMask(AliMuonTrackCuts::kMuEta |
1769 AliMuonTrackCuts::kMuThetaAbs |
1770 AliMuonTrackCuts::kMuPdca |
1771 AliMuonTrackCuts::kMuMatchApt |
1772 AliMuonTrackCuts::kMuMatchLpt |
1773 AliMuonTrackCuts::kMuMatchHpt |
1774 AliMuonTrackCuts::kMuTrackChiSquare);
1775 }
1776
1777 fMuonTrackCuts->SetRun(fInputHandler);
1778}
1779
1780//_____________________________________________________________________________
1781Bool_t AliAnalysisTaskMuMu::PairRapidityCut(const AliVParticle& t1, const AliVParticle& t2) const
1782{
1783 /// Whether the pair passes the rapidity cut
1784
1785 TLorentzVector p1(t1.Px(),t1.Py(),t1.Pz(),TMath::Sqrt(MuonMass2()+t1.P()*t1.P()));
1786 TLorentzVector p2(t2.Px(),t2.Py(),t2.Pz(),TMath::Sqrt(MuonMass2()+t2.P()*t2.P()));
1787
1788 TLorentzVector total(p1+p2);
1789
1790 Double_t y = total.Rapidity();
1791
1792 Bool_t ok = ( y < -2.5 && y > -4.0 );
1793
1794 return ok;
1795}
1796
1797
1798//_____________________________________________________________________________
1799Bool_t AliAnalysisTaskMuMu::PassPhysicsSelection() const
1800{
1801 /// Whether the event pass the physics selection or not
1802
1803 AliInputEventHandler* inputEventHandler = static_cast<AliInputEventHandler*>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
1804
1805 Bool_t isPhysicsSelected = (inputEventHandler->IsEventSelected() & AliVEvent::kAny);
1806
1807 /*
1808 ( inputEventHandler->IsEventSelected() & AliVEvent::kMUSH7 )
1809 || ( inputEventHandler->IsEventSelected() & AliVEvent::kMUL7 )
1810 || ( inputEventHandler->IsEventSelected() & AliVEvent::kMUS7 )
1811 || ( inputEventHandler->IsEventSelected() & AliVEvent::kMUU7 )
1812 || ( inputEventHandler->IsEventSelected() & AliVEvent::kMuonSingleLowPt8 )
1813 || ( inputEventHandler->IsEventSelected() & AliVEvent::kMuonSingleHighPt8 )
1814 || ( inputEventHandler->IsEventSelected() & AliVEvent::kMuonLikeLowPt8 )
1815 || ( inputEventHandler->IsEventSelected() & AliVEvent::kMuonUnlikeLowPt8 );
1816
1817 if ( IsPP() )
1818 {
1819 isPhysicsSelected = isPhysicsSelected || ( inputEventHandler->IsEventSelected() & AliVEvent::kINT7 ) || ( inputEventHandler->IsEventSelected() & AliVEvent::kINT8 );
1820 }
1821 else
1822 {
1823 isPhysicsSelected = isPhysicsSelected || ( inputEventHandler->IsEventSelected() & AliVEvent::kAny );
1824 }
1825 */
1826
1827 return isPhysicsSelected;
1828}
1829
1830//_____________________________________________________________________________
1831void
1832AliAnalysisTaskMuMu::Print(Option_t* /*opt*/) const
1833{
1834 /// Print the definition of this analysis
1835
1836 cout << ClassName() << " - " << GetName() << " - " << fBeamYear.Data() << endl;
1837
1838 if ( !fSingleTrackCutNames || !fSingleTrackCutNames )
1839 {
1840 cout << "No single track cut defined yet" << endl;
1841 }
1842 else
1843 {
1844 TIter next(fSingleTrackCutNames);
1845 TObjString* str;
1846
1847 while ( ( str = static_cast<TObjString*>(next()) ) )
1848 {
1849 cout << Form("SINGLE CUT %20s MASK %x",str->String().Data(),str->GetUniqueID()) << endl;
1850 }
1851 }
1852
1853 if ( !fPairTrackCutNames || !fPairTrackCutNames )
1854 {
1855 cout << "No track pair cut defined yet" << endl;
1856 }
1857 else
1858 {
1859 TIter next2(fPairTrackCutNames);
1860 AliAnalysisTaskMuMu::PairCut* str;
1861
1862 while ( ( str = static_cast<AliAnalysisTaskMuMu::PairCut*>(next2()) ) )
1863 {
1864 str->Print();
1865 }
1866 }
1867
1868 if ( !fTriggerClasses || !fTriggerClasses->First() )
1869 {
1870 cout << "No trigger classes defined yet" << endl;
1871 }
1872 else
1873 {
1874 cout << "Trigger classes that will be considered:" << endl;
1875 TIter next(fTriggerClasses);
1876 TObjString* s;
1877
1878 while ( ( s = static_cast<TObjString*>(next()) ) )
1879 {
1880 cout << s->String().Data() << endl;
1881 }
1882 }
1883}
1884
1885//_____________________________________________________________________________
1886void
1887AliAnalysisTaskMuMu::Terminate(Option_t *)
1888{
1889 /// Called once at the end of the query
1890 /// Just a simple printout of the stat we analyse and how many histograms
1891 /// we got
1892
1893 fOutput = dynamic_cast<TList*>(GetOutputData(1));
1894
1895 if (!fOutput) return;
1896
1897 fHistogramCollection = dynamic_cast<AliHistogramCollection*>(fOutput->FindObject("mumu"));
1898
1899 if (!fHistogramCollection)
1900 {
1901 AliError("Could not find back histogram collection in output...");
1902 return;
1903 }
1904
1905 fHistogramCollection->PruneEmptyHistograms();
1906
1907 UInt_t size2 = fHistogramCollection->EstimateSize();
1908
1909 AliInfo(Form("size after prune histograms = %5.1f MB",size2/1024.0/1024.0));
1910
1911 if ( !IsPP() && fCentralityLimits.size() > 1 )
1912 {
1913 MergeCentralities(fHistogramCollection);
1914 }
1915
1916 BeautifyHistos();
1917
1918 fHistogramCollection->Print("*Minv*");
1919
1920 fEventCounters = dynamic_cast<AliCounterCollection*>(fOutput->FindObject("eventCounters"));
1921
1922 if (!fEventCounters)
1923 {
1924 AliError("Could not find back counters in output...");
1925 return;
1926 }
1927
1928 fEventCounters->Print("trigger");
1929}
1930
1931//_____________________________________________________________________________
1932Bool_t AliAnalysisTaskMuMu::TriggerSBACECondition(const TString& triggerName) const
1933{
1934 // check the beam condition in the trigger name
1935
1936 if ( triggerName.Contains("-S-") ) return kTRUE;
1937
1938 if ( triggerName.Contains("-B-") ) return kTRUE;
1939
1940 if ( fUseBackgroundTriggers )
1941 {
1942 if ( triggerName.Contains("-ACE-") ) return kTRUE;
1943 }
1944 return kFALSE;
1945}
1946
1947//_____________________________________________________________________________
1948void AliAnalysisTaskMuMu::UserExec(Option_t* /*opt*/)
1949{
1950 /// Executed at each event
1951
1952 AliVEvent* event = InputEvent();
1953
1954 if (!event) return;
1955
1956 Bool_t anyTrigger(kFALSE);
1957
1958 if ( IsDynamicTriggerClasses() )
1959 {
91dbbc20 1960 AddTriggerClasses(EAGetFiredTriggerClasses(*event)," ");
2f331ac9 1961 }
1962 else
1963 {
1964 if (fTriggerClasses->Contains("ANY") )
1965 {
1966 anyTrigger = kTRUE;
1967 }
1968 }
1969
1970 TString centrality(DefaultCentralityName());
1971
1972 Float_t fcent(-1);
1973
1974 AliCentrality* acent = event->GetCentrality();
1975
1976 if (acent) fcent = acent->GetCentralityPercentile("V0M");
1977
1978 Double_t cent(-100);
1979
1980 if ( fcent > 0 )
1981 {
1982 for ( std::vector<double>::size_type i = 0 ; i < fCentralityLimits.size() && cent < 0 ; ++i )
1983 {
1984 if ( fcent < fCentralityLimits[i] )
1985 {
1986 cent = fCentralityLimits[i];
1987 }
1988 }
1989 }
1990
1991 if ( cent > -1 )
1992 {
1993 centrality = CentralityName(cent);
1994 }
1995
1996 int nmu = EAGetNumberOfMuonTracks(*event);
1997
1998 EAComputeTrackMasks(*event);
1999
2000 TIter next(fTriggerClasses);
2001 TObjString* tname;
2002 TString firedTriggerClasses(EAGetFiredTriggerClasses(*event));
2003
2004 // first loop to count things not associated to a specific trigger
2005 TIter nextEventCut(fEventCutNames);
2006 TObjString* et;
2007
2008 UInt_t mask = GetEventMask(*event);
2009
2010 while ( ( et = static_cast<TObjString*>(nextEventCut()) ) )
2011 {
2012 Bool_t test = ( ( et->GetUniqueID() & mask ) == et->GetUniqueID() );
2013
2014 if ( test )
2015 {
2016 fEventCounters->Count(Form("event:%s/trigger:%s/run:%d", et->String().Data(), "EVERYTHING", fCurrentRunNumber));
2017
2018 if ( firedTriggerClasses == "" )
2019 {
2020 fEventCounters->Count(Form("event:%s/trigger:%s/run:%d", et->String().Data(), "EMPTY", fCurrentRunNumber));
2021 }
2022
2023 if (nmu)
2024 {
2025 fEventCounters->Count(Form("event:%s/trigger:%s/run:%d", et->String().Data(), "ATLEASTONEMUONTRACK", fCurrentRunNumber));
2026 }
2027
2028 if ( AtLeastOneMuonTrigger(firedTriggerClasses) )
2029 {
2030 fEventCounters->Count(Form("event:%s/trigger:%s/run:%d", et->String().Data(), "ATLEASTONEMUONTRIGGER", fCurrentRunNumber));
91dbbc20 2031 if ( AtLeastOneMBTrigger(firedTriggerClasses) )
2032 {
2033 fEventCounters->Count(Form("event:%s/trigger:%s/run:%d", et->String().Data(), "ATLEASTONEMUONORMBTRIGGER", fCurrentRunNumber));
2034 }
2f331ac9 2035 }
2036
2037 if ( AtLeastOneMBTrigger(firedTriggerClasses) )
2038 {
2039 fEventCounters->Count(Form("event:%s/trigger:%s/run:%d", et->String().Data(), "ATLEASTONEMBTRIGGER", fCurrentRunNumber));
2040 }
2041
91dbbc20 2042// if ( AtLeastOneEmcalTrigger(firedTriggerClasses) )
2043// {
2044// fEventCounters->Count(Form("event:%s/trigger:%s/run:%d", et->String().Data(), "ATLEASTONEEMCALTRIGGER", fCurrentRunNumber));
2045//
2046// if ( AtLeastOneMBTrigger(firedTriggerClasses) )
2047// {
2048// fEventCounters->Count(Form("event:%s/trigger:%s/run:%d", et->String().Data(), "ATLEASTONEEMCALORMBTRIGGER", fCurrentRunNumber));
2049// }
2050// }
2f331ac9 2051
2052 }
2053 }
2054
91dbbc20 2055 UInt_t l0Inputs = EAGetL0TriggerInputs(*event);
2056
2f331ac9 2057 // second loop to count only the triggers we're interested in
2058 while ( ( tname = static_cast<TObjString*>(next()) ) )
2059 {
91dbbc20 2060 if ( !CheckTriggerClass(tname->String(),firedTriggerClasses,l0Inputs) ) continue;
2f331ac9 2061
2062 nextEventCut.Reset();
2063
2064 while ( ( et = static_cast<TObjString*>(nextEventCut()) ) )
2065 {
2066 Bool_t test = ( ( et->GetUniqueID() & mask ) == et->GetUniqueID() );
2067
2068 if ( test )
2069 {
2070 Fill(et->String().Data(),tname,centrality,fcent,*event);
2071 }
2072 }
2073 }
2074
2075 // Post output data.
2076 PostData(1, fOutput);
2077}
2078
2079//_____________________________________________________________________________
2080void AliAnalysisTaskMuMu::UserCreateOutputObjects()
2081{
2082 // Create histograms
2083 // Called once
2084
2085 OpenFile(1);
2086
2087 fTriggerClasses->Print();
2088
2089 fOutput = new TList;
2090 fOutput->SetOwner(kTRUE);
2091
2092 fHistogramCollection = new AliHistogramCollection("mumu");
2093
2094 fOutput->Add(fHistogramCollection);
2095
2096 // initialize event counters
2097 fEventCounters = new AliCounterCollection("eventCounters");
2098
2099 TIter nextEventCutName(fEventCutNames);
2100 TObjString* str;
2101 TString eventRubric;
2102 while ( ( str = static_cast<TObjString*>(nextEventCutName()) ) )
2103 {
2104 if ( eventRubric.Length() > 0 ) eventRubric += "/";
2105 eventRubric += str->String();
2106 }
2107
2108 eventRubric += "/SCALER";
2109
2110 fEventCounters->AddRubric("event", eventRubric.Data());
2111
2112 fEventCounters->AddRubric("trigger", 100);
2113
2114 fEventCounters->AddRubric("run", 1000000);
2115
2116 fEventCounters->Init();
2117
2118 fOutput->Add(fEventCounters);
2119
2120 // Post output data.
2121 PostData(1, fOutput);
2122}