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 SetUseCentralityPatch(Bool_t isAOD049) {fUseAOD049CentralityPatch = isAOD049;}
41 void UseMultiplicity(const char *type) {fUseCentrality = kFALSE; fCentralityType = type; fCentralityType.ToUpper();}
42 void UseContinuousMix() {fContinuousMix = kTRUE;}
43 void UseBinnedMix() {fContinuousMix = kFALSE;}
44 void SetNMix(Int_t nmix) {fNMix = nmix;}
45 void SetMaxDiffMult (Double_t val) {fMaxDiffMult = val;}
46 void SetMaxDiffVz (Double_t val) {fMaxDiffVz = val;}
47 void SetMaxDiffAngle(Double_t val) {fMaxDiffAngle = val;}
48 void SetEventCuts(AliRsnCutSet *cuts) {fEventCuts = cuts;}
49 void SetMixPrintRefresh(Int_t n) {fMixPrintRefresh = n;}
50 Int_t AddTrackCuts(AliRsnCutSet *cuts);
51 TClonesArray *Outputs() {return &fHistograms;}
52 TClonesArray *Values() {return &fValues;}
54 virtual void UserCreateOutputObjects();
55 virtual void UserExec(Option_t *);
56 virtual void Terminate(Option_t *);
57 virtual void FinishTaskOutput();
59 Int_t ValueID(AliRsnMiniValue::EType type, Bool_t useMC = kFALSE);
60 Int_t CreateValue(AliRsnMiniValue::EType type, Bool_t useMC = kFALSE);
61 AliRsnMiniOutput *CreateOutput(const char *name, AliRsnMiniOutput::EOutputType type, AliRsnMiniOutput::EComputation src);
62 AliRsnMiniOutput *CreateOutput(const char *name, const char *outType, const char *compType);
66 Char_t CheckCurrentEvent();
67 void FillMiniEvent(Char_t evType);
68 Double_t ComputeAngle();
69 Double_t ComputeCentrality(Bool_t isESD);
70 Double_t ApplyCentralityPatchAOD049();
71 void FillTrueMotherESD(AliRsnMiniEvent *event);
72 void FillTrueMotherAOD(AliRsnMiniEvent *event);
73 void StoreTrueMother(AliRsnMiniPair *pair, AliRsnMiniEvent *event);
74 Bool_t EventsMatch(AliRsnMiniEvent *event1, AliRsnMiniEvent *event2);
76 Bool_t fUseMC; // use or not MC info
77 Int_t fEvNum; //! absolute event counter
78 Bool_t fUseCentrality; // if true, use centrality for event, otherwise use multiplicity
79 TString fCentralityType; // definition used to choose what centrality or multiplicity to use
80 Bool_t fUseAOD049CentralityPatch; //flag to enable AOD049 centrality patch
82 Bool_t fContinuousMix; // mixing --> technique chosen (continuous or binned)
83 Int_t fNMix; // mixing --> required number of mixes
84 Double_t fMaxDiffMult; // mixing --> max difference in multiplicity
85 Double_t fMaxDiffVz; // mixing --> max difference in Vz of prim vert
86 Double_t fMaxDiffAngle; // mixing --> max difference in reaction plane angle
88 TList *fOutput; // output list
89 TClonesArray fHistograms; // list of histogram definitions
90 TClonesArray fValues; // list of values to be computed
91 TH1F *fHEventStat; // histogram of event statistics
93 AliRsnCutSet *fEventCuts; // cuts on events
94 TObjArray fTrackCuts; // list of single track cuts
95 AliRsnEvent fRsnEvent; //! interface object to the event
96 TTree *fEvBuffer; //! mini-event buffer
97 AliTriggerAnalysis *fTriggerAna; //! trigger analysis
98 AliESDtrackCuts *fESDtrackCuts; //! quality cut for ESD tracks
99 AliRsnMiniEvent *fMiniEvent; //! mini-event cursor
100 Bool_t fBigOutput; // flag if open file for output list
101 Int_t fMixPrintRefresh; // how often info in mixing part is printed
103 ClassDef(AliRsnMiniAnalysisTask, 4); // AliRsnMiniAnalysisTask
106 inline Int_t AliRsnMiniAnalysisTask::CreateValue(AliRsnMiniValue::EType type, Bool_t useMC)
109 // Create a new value in the task,
110 // and returns its ID, which is needed for setting up histograms.
111 // If that value was already initialized, returns its ID and does not recreate it.
114 Int_t valID = ValueID(type, useMC);
115 if (valID >= 0 && valID < fValues.GetEntries()) {
116 AliInfo(Form("Value '%s' is already created in slot #%d", AliRsnMiniValue::ValueName(type, useMC), valID));
118 valID = fValues.GetEntries();
119 AliInfo(Form("Creating value '%s' in slot #%d", AliRsnMiniValue::ValueName(type, useMC), valID));
120 new (fValues[valID]) AliRsnMiniValue(type, useMC);
126 inline Int_t AliRsnMiniAnalysisTask::ValueID(AliRsnMiniValue::EType type, Bool_t useMC)
129 // Searches if a value computation is initialized
132 const char *name = AliRsnMiniValue::ValueName(type, useMC);
133 TObject *obj = fValues.FindObject(name);
135 return fValues.IndexOf(obj);
140 inline AliRsnMiniOutput *AliRsnMiniAnalysisTask::CreateOutput
141 (const char *name, AliRsnMiniOutput::EOutputType type, AliRsnMiniOutput::EComputation src)
144 // Create a new histogram definition in the task,
145 // which is then returned to the user for its configuration
148 Int_t n = fHistograms.GetEntries();
149 AliRsnMiniOutput *newDef = new (fHistograms[n]) AliRsnMiniOutput(name, type, src);
154 inline AliRsnMiniOutput *AliRsnMiniAnalysisTask::CreateOutput
155 (const char *name, const char *outType, const char *compType)
158 // Create a new histogram definition in the task,
159 // which is then returned to the user for its configuration
162 Int_t n = fHistograms.GetEntries();
163 AliRsnMiniOutput *newDef = new (fHistograms[n]) AliRsnMiniOutput(name, outType, compType);