1 #ifndef ALIRSNMINIANALYSISTASK_H
2 #define ALIRSNMINIANALYSISTASK_H
5 // Analysis task for 'mini' sub-package
6 // Contains all definitions needed for running an analysis:
8 // -- a list of track cuts (any number)
9 // -- definitions of output histograms
10 // -- values to be computed.
14 #include <TClonesArray.h>
16 #include "AliAnalysisTaskSE.h"
18 #include "AliRsnEvent.h"
19 #include "AliRsnMiniValue.h"
20 #include "AliRsnMiniOutput.h"
24 class AliTriggerAnalysis;
25 class AliRsnMiniEvent;
28 class AliRsnMiniAnalysisTask : public AliAnalysisTaskSE {
32 AliRsnMiniAnalysisTask();
33 AliRsnMiniAnalysisTask(const char *name, Bool_t isMC = kFALSE);
34 AliRsnMiniAnalysisTask(const AliRsnMiniAnalysisTask ©);
35 AliRsnMiniAnalysisTask& operator=(const AliRsnMiniAnalysisTask ©);
36 virtual ~AliRsnMiniAnalysisTask();
38 void UseMC(Bool_t yn = kTRUE) {fUseMC = yn;}
39 void UseCentrality(const char *type) {fUseCentrality = kTRUE; fCentralityType = type; fCentralityType.ToUpper();}
40 void UseMultiplicity(const char *type) {fUseCentrality = kFALSE; fCentralityType = type; fCentralityType.ToUpper();}
41 void SetNMix(Int_t nmix) {fNMix = nmix;}
42 void SetMaxDiffMult (Double_t val) {fMaxDiffMult = val;}
43 void SetMaxDiffVz (Double_t val) {fMaxDiffVz = val;}
44 void SetMaxDiffAngle(Double_t val) {fMaxDiffAngle = val;}
45 void SetEventCuts(AliRsnCutSet *cuts) {fEventCuts = cuts;}
46 Int_t AddTrackCuts(AliRsnCutSet *cuts);
47 TClonesArray *Outputs() {return &fHistograms;}
48 TClonesArray *Values() {return &fValues;}
50 virtual void UserCreateOutputObjects();
51 virtual void UserExec(Option_t*);
52 virtual void Terminate(Option_t*);
53 virtual void FinishTaskOutput();
55 Int_t ValueID(AliRsnMiniValue::EType type, Bool_t useMC = kFALSE);
56 Int_t CreateValue(AliRsnMiniValue::EType type, Bool_t useMC = kFALSE);
57 AliRsnMiniOutput *CreateOutput(const char *name, AliRsnMiniOutput::EOutputType type, AliRsnMiniOutput::EComputation src);
58 AliRsnMiniOutput *CreateOutput(const char *name, const char *outType, const char *compType);
62 Char_t CheckCurrentEvent();
63 Int_t FillMiniEvent(Char_t evType);
64 Double_t ComputeAngle();
65 Double_t ComputeCentrality(Bool_t isESD);
66 void FillTrueMotherESD(AliRsnMiniEvent *event);
67 void FillTrueMotherAOD(AliRsnMiniEvent *event);
68 void StoreTrueMother(AliRsnMiniPair *pair, AliRsnMiniEvent *event);
69 void ProcessEvents(AliRsnMiniEvent *evMain, AliRsnMiniEvent *evMix = 0x0);
71 Bool_t fUseMC; // use or not MC info
72 Int_t fEvNum; //! absolute event counter
73 Bool_t fUseCentrality; // if true, use centrality for event, otherwise use multiplicity
74 TString fCentralityType; // definition used to choose what centrality or multiplicity to use
76 Int_t fNMix; // mixing --> required number of mixes
77 Double_t fMaxDiffMult; // mixing --> max difference in multiplicity
78 Double_t fMaxDiffVz; // mixing --> max difference in Vz of prim vert
79 Double_t fMaxDiffAngle; // mixing --> max difference in reaction plane angle
81 TList *fOutput; // output list
82 TClonesArray fHistograms; // list of histogram definitions
83 TClonesArray fValues; // list of values to be computed
84 TH1F *fHEventStat; // histogram of event statistics
86 AliRsnCutSet *fEventCuts; // cuts on events
87 TObjArray fTrackCuts; // list of single track cuts
88 AliRsnEvent fRsnEvent; //! interface object to the event
89 TTree *fEvBuffer; //! mini-event buffer
90 TArrayI fNMixed; //! array to keep trace of how many times an event was mixed
91 AliTriggerAnalysis *fTriggerAna; //! trigger analysis
92 AliESDtrackCuts *fESDtrackCuts; //! quality cut for ESD tracks
93 AliRsnMiniEvent *fMiniEvent; //! mini-event cursor
95 ClassDef(AliRsnMiniAnalysisTask, 1); // AliRsnMiniAnalysisTask
98 inline Int_t AliRsnMiniAnalysisTask::CreateValue(AliRsnMiniValue::EType type, Bool_t useMC)
101 // Create a new value in the task,
102 // and returns its ID, which is needed for setting up histograms.
103 // If that value was already initialized, returns its ID and does not recreate it.
106 Int_t valID = ValueID(type, useMC);
107 if (valID >= 0 && valID < fValues.GetEntries()) {
108 AliInfo(Form("Value '%s' is already created in slot #%d", AliRsnMiniValue::ValueName(type, useMC), valID));
110 valID = fValues.GetEntries();
111 AliInfo(Form("Creating value '%s' in slot #%d", AliRsnMiniValue::ValueName(type, useMC), valID));
112 new (fValues[valID]) AliRsnMiniValue(type, useMC);
118 inline Int_t AliRsnMiniAnalysisTask::ValueID(AliRsnMiniValue::EType type, Bool_t useMC)
121 // Searches if a value computation is initialized
124 const char *name = AliRsnMiniValue::ValueName(type, useMC);
125 TObject *obj = fValues.FindObject(name);
127 return fValues.IndexOf(obj);
132 inline AliRsnMiniOutput* AliRsnMiniAnalysisTask::CreateOutput
133 (const char *name, AliRsnMiniOutput::EOutputType type, AliRsnMiniOutput::EComputation src)
136 // Create a new histogram definition in the task,
137 // which is then returned to the user for its configuration
140 Int_t n = fHistograms.GetEntries();
141 AliRsnMiniOutput *newDef = new (fHistograms[n]) AliRsnMiniOutput(name, type, src);
146 inline AliRsnMiniOutput* AliRsnMiniAnalysisTask::CreateOutput
147 (const char *name, const char *outType, const char *compType)
150 // Create a new histogram definition in the task,
151 // which is then returned to the user for its configuration
154 Int_t n = fHistograms.GetEntries();
155 AliRsnMiniOutput *newDef = new (fHistograms[n]) AliRsnMiniOutput(name, outType, compType);