bug fixes
[u/mrichter/AliRoot.git] / PWG2 / FLOW / AliFlowTasks / AliAnalysisTaskFlowEvent.cxx
CommitLineData
1c1d4332 1/*************************************************************************
2* Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
3* *
4* Author: The ALICE Off-line Project. *
5* Contributors are mentioned in the code where appropriate. *
6* *
7* Permission to use, copy, modify and distribute this software and its *
8* documentation strictly for non-commercial purposes is hereby granted *
9* without fee, provided that the above copyright notice appears in all *
10* copies and that both the copyright notice and this permission notice *
11* appear in the supporting documentation. The authors make no claims *
12* about the suitability of this software for any purpose. It is *
13* provided "as is" without express or implied warranty. *
14**************************************************************************/
15
fdff44c8 16////////////////////////////////////////////////////
17// AliAnalysisTaskFlowEvent:
18//
19// analysis task for filling the flow event
20// from MCEvent, ESD, AOD ....
21// and put it in an output stream so it can
22// be used by the various flow analysis methods
23// for cuts the correction framework is used
24// which also outputs QA histograms to view
25// the effects of the cuts
26////////////////////////////////////////////////////
27
1c1d4332 28#include "Riostream.h" //needed as include
29#include "TChain.h"
30#include "TTree.h"
31#include "TFile.h" //needed as include
32#include "TList.h"
489fdf79 33#include "TRandom3.h"
34#include "TTimeStamp.h"
1c1d4332 35
7183fe85 36// ALICE Analysis Framework
1c1d4332 37class AliAnalysisTask;
38#include "AliAnalysisManager.h"
39
7183fe85 40// ESD interface
1c1d4332 41#include "AliESDEvent.h"
42#include "AliESDInputHandler.h"
43
7183fe85 44// AOD interface
1c1d4332 45#include "AliAODEvent.h"
46#include "AliAODInputHandler.h"
47
7183fe85 48// Monte Carlo Event
1c1d4332 49#include "AliMCEventHandler.h"
50#include "AliMCEvent.h"
51
7183fe85 52// ALICE Correction Framework
1c1d4332 53#include "AliCFManager.h"
54
7183fe85 55// Interface to Event generators to get Reaction Plane Angle
56#include "AliGenCocktailEventHeader.h"
57#include "AliGenHijingEventHeader.h"
48ad51a1 58#include "AliGenGeVSimEventHeader.h"
26f120fa 59#include "AliGenEposEventHeader.h"
7183fe85 60
61// Interface to make the Flow Event Simple used in the flow analysis methods
1c1d4332 62#include "AliFlowEventSimpleMaker.h"
63
7183fe85 64#include "AliAnalysisTaskFlowEvent.h"
65
1c1d4332 66
67ClassImp(AliAnalysisTaskFlowEvent)
73160af6 68
1c1d4332 69//________________________________________________________________________
1c1d4332 70AliAnalysisTaskFlowEvent::AliAnalysisTaskFlowEvent() :
65201059 71 AliAnalysisTaskSE(),
46bec39c 72 // fOutputFile(NULL),
1c1d4332 73 fEventMaker(NULL),
74 fAnalysisType("ESD"),
75 fCFManager1(NULL),
76 fCFManager2(NULL),
77 fQAInt(NULL),
78 fQADiff(NULL),
fdff44c8 79 fMinMult(0),
80 fMaxMult(10000000),
7a01f4a7 81 fMinA(-1.0),
82 fMaxA(-0.01),
83 fMinB(0.01),
84 fMaxB(1.0),
489fdf79 85 fQA(kFALSE),
86 fMCReactionPlaneAngle(0.),
87 fCount(0),
88 fNoOfLoops(1),
89 fEllipticFlowValue(0.),
90 fSigmaEllipticFlowValue(0.),
91 fMultiplicityOfEvent(1000000000),
92 fSigmaMultiplicityOfEvent(0),
65201059 93 fMyTRandom3(NULL),
94 fbAfterburnerOn(kFALSE)
1c1d4332 95{
96 // Constructor
97 cout<<"AliAnalysisTaskFlowEvent::AliAnalysisTaskFlowEvent()"<<endl;
98}
99
100//________________________________________________________________________
489fdf79 101AliAnalysisTaskFlowEvent::AliAnalysisTaskFlowEvent(const char *name, Bool_t on, UInt_t iseed) :
65201059 102 AliAnalysisTaskSE(name),
489fdf79 103// fOutputFile(NULL),
489fdf79 104 fEventMaker(NULL),
105 fAnalysisType("ESD"),
106 fCFManager1(NULL),
107 fCFManager2(NULL),
108 fQAInt(NULL),
109 fQADiff(NULL),
fdff44c8 110 fMinMult(0),
111 fMaxMult(10000000),
7a01f4a7 112 fMinA(-1.0),
113 fMaxA(-0.01),
114 fMinB(0.01),
115 fMaxB(1.0),
489fdf79 116 fQA(on),
117 fMCReactionPlaneAngle(0.),
118 fCount(0),
119 fNoOfLoops(1),
120 fEllipticFlowValue(0.),
121 fSigmaEllipticFlowValue(0.),
122 fMultiplicityOfEvent(1000000000),
123 fSigmaMultiplicityOfEvent(0),
65201059 124 fMyTRandom3(NULL),
125 fbAfterburnerOn(kFALSE)
489fdf79 126{
127 // Constructor
128 cout<<"AliAnalysisTaskFlowEvent::AliAnalysisTaskFlowEvent(const char *name, Bool_t on, UInt_t iseed)"<<endl;
489fdf79 129 fMyTRandom3 = new TRandom3(iseed);
130 gRandom->SetSeed(fMyTRandom3->Integer(65539));
65201059 131
132
133 // Define output slots here
489fdf79 134 // Define here the flow event output
65201059 135 DefineOutput(1, AliFlowEventSimple::Class());
489fdf79 136 if(on) {
65201059 137 DefineOutput(2, TList::Class());
138 DefineOutput(3, TList::Class()); }
489fdf79 139 // and for testing open an output file
140 // fOutputFile = new TFile("FlowEvents.root","RECREATE");
141
142}
143
144//________________________________________________________________________
1c1d4332 145AliAnalysisTaskFlowEvent::~AliAnalysisTaskFlowEvent()
146{
147 //
148 // Destructor
149 //
489fdf79 150 if (fMyTRandom3) delete fMyTRandom3;
1c1d4332 151 // objects in the output list are deleted
152 // by the TSelector dtor (I hope)
153
154}
155
156//________________________________________________________________________
65201059 157void AliAnalysisTaskFlowEvent::UserCreateOutputObjects()
1c1d4332 158{
159 // Called at every worker node to initialize
160 cout<<"AliAnalysisTaskFlowEvent::CreateOutputObjects()"<<endl;
161
162 if (!(fAnalysisType == "AOD" || fAnalysisType == "ESD" || fAnalysisType == "ESDMC0" || fAnalysisType == "ESDMC1" || fAnalysisType == "MC")) {
163 cout<<"WRONG ANALYSIS TYPE! only ESD, ESDMC0, ESDMC1, AOD and MC are allowed."<<endl;
164 exit(1);
165 }
166
167 // Flow Event maker
168 fEventMaker = new AliFlowEventSimpleMaker();
169}
170
171//________________________________________________________________________
65201059 172void AliAnalysisTaskFlowEvent::UserExec(Option_t *)
1c1d4332 173{
174 // Main loop
175 // Called for each event
89f41a6a 176 AliFlowEventSimple* fEvent = NULL;
7183fe85 177 Double_t fRP = 0.; // the monte carlo reaction plane angle
8c5b0db3 178 AliMCEvent* mcEvent = fMCEvent; // from TaskSE
179 AliESDEvent* myESD = dynamic_cast<AliESDEvent*>(fInputEvent); // from TaskSE
180 AliAODEvent* myAOD = dynamic_cast<AliAODEvent*>(fInputEvent); // from TaskSE
7183fe85 181
8c5b0db3 182 // if monte carlo event get reaction plane from monte carlo (depends on generator)
65201059 183 if (mcEvent) {
184
185 //COCKTAIL with HIJING
186 if (!strcmp(mcEvent-> GenEventHeader()->GetName(),"Cocktail Header")) { //returns 0 if matches
187 AliGenCocktailEventHeader *headerC = dynamic_cast<AliGenCocktailEventHeader *> (mcEvent-> GenEventHeader());
188 if (headerC) {
189 TList *lhd = headerC->GetHeaders();
190 if (lhd) {
191 AliGenHijingEventHeader *hdh = dynamic_cast<AliGenHijingEventHeader *> (lhd->At(0));
192 if (hdh) {
193 fRP = hdh->ReactionPlaneAngle();
194 //cout<<"The reactionPlane from Hijing is: "<< fRP <<endl;
7183fe85 195 }
196 }
c9f92d16 197 }
65201059 198 //else { cout<<"headerC is NULL"<<endl; }
199 }
26f120fa 200
65201059 201 //GEVSIM
202 else if (!strcmp(mcEvent-> GenEventHeader()->GetName(),"GeVSim header")) { //returns 0 if matches
203 AliGenGeVSimEventHeader* headerG = (AliGenGeVSimEventHeader*)(mcEvent->GenEventHeader());
204 if (headerG) {
205 fRP = headerG->GetEventPlane();
206 //cout<<"The reactionPlane from GeVSim is: "<< fRP <<endl;
7183fe85 207 }
65201059 208 //else { cout<<"headerG is NULL"<<endl; }
209 }
210
211 //HIJING
212 else if (!strcmp(mcEvent-> GenEventHeader()->GetName(),"Hijing")) { //returns 0 if matches
213 AliGenHijingEventHeader* headerH = (AliGenHijingEventHeader*)(mcEvent->GenEventHeader());
214 if (headerH) {
215 fRP = headerH->ReactionPlaneAngle();
216 //cout<<"The reactionPlane from Hijing is: "<< fRP <<endl;
d5acfb36 217 }
65201059 218 //else { cout<<"headerH is NULL"<<endl; }
219 }
220
221 //EPOS
222 else if (!strcmp(mcEvent->GenEventHeader()->GetName(),"EPOS")) {
223 AliGenEposEventHeader* headerE = (AliGenEposEventHeader*)(mcEvent->GenEventHeader());
224 if (headerE) {
225 fRP = headerE->ReactionPlaneAngle();
226 //cout<<"The reactionPlane from EPOS is: "<< fR <<endl;
26f120fa 227 }
65201059 228 //else { cout<<"headerE is NULL"<<endl; }
7183fe85 229 }
230 }
fdff44c8 231
8c5b0db3 232 // set the reaction plane angle for the flow event
f1f90185 233 fEventMaker->SetMCReactionPlaneAngle(fRP);
fdff44c8 234 //setting event cuts
235 fEventMaker->SetMinMult(fMinMult);
236 fEventMaker->SetMaxMult(fMaxMult);
7a01f4a7 237 //setting ranges for eta subevents
238 fEventMaker->SetSubeventEtaRange(fMinA,fMaxA,fMinB,fMaxB);
f1f90185 239
8c5b0db3 240 if (fbAfterburnerOn && fMyTRandom3) {
241 // set the new value of the values using a after burner
f1f90185 242 cout << "settings for afterburner in TaskFlowEvent.cxx:" << endl;
6d7734d2 243 cout << "fCount = " << fCount << endl;
244 cout << "fNoOfLoops = " << fNoOfLoops << endl;
245 cout << "fEllipticFlowValue = " << fEllipticFlowValue << endl;
246 cout << "fSigmaEllipticFlowValue = " << fSigmaEllipticFlowValue << endl;
247 cout << "fMultiplicityOfEvent = " << fMultiplicityOfEvent << endl;
248 cout << "fSigmaMultiplicityOfEvent = " << fSigmaMultiplicityOfEvent << endl;
8c5b0db3 249 Double_t xRPangle;
250 if (!mcEvent) { xRPangle = TMath::TwoPi()*(fMyTRandom3->Rndm()); }
251 else { xRPangle = fRP; }
252 Double_t xNewFlowValue = fMyTRandom3->Gaus(fEllipticFlowValue,fSigmaEllipticFlowValue);
253 Int_t nNewMultOfEvent = TMath::Nint(fMyTRandom3->Gaus(fMultiplicityOfEvent,fSigmaMultiplicityOfEvent));
f1f90185 254
f1f90185 255 cout << "xRPangle = " << xRPangle << endl;
256 cout << "xNewFlowValue = " << xNewFlowValue << endl;
257 cout << "nNewMultOfEvent = " << nNewMultOfEvent << endl;
258 cout << "settings for after burner" << endl;
259
260 fEventMaker->SetMCReactionPlaneAngle(xRPangle);
261 fEventMaker->SetNoOfLoops(fNoOfLoops);
262 fEventMaker->SetEllipticFlowValue(xNewFlowValue);
263 fEventMaker->SetMultiplicityOfEvent(nNewMultOfEvent);
264 //end settings afterburner
8c5b0db3 265 }
d5acfb36 266
62726ef0 267 // Fill the FlowEventSimple for MC input
1c1d4332 268 if (fAnalysisType == "MC") {
587b0a35 269 if (!fCFManager1) {cout << "ERROR: No pointer to correction framework cuts! " << endl; return; }
270 if (!fCFManager2) {cout << "ERROR: No pointer to correction framework cuts! " << endl; return; }
1c1d4332 271
272 // Process MC truth, therefore we receive the AliAnalysisManager and ask it for the AliMCEventHandler
273 // This handler can return the current MC event
587b0a35 274 if (!mcEvent) { Printf("ERROR: Could not retrieve MC event"); return;}
1c1d4332 275
6d7734d2 276 fCFManager1->SetMCEventInfo(mcEvent);
277 fCFManager2->SetMCEventInfo(mcEvent);
587b0a35 278
1c1d4332 279 // analysis
7183fe85 280 Printf("Number of MC particles: %d", mcEvent->GetNumberOfTracks());
89f41a6a 281 fEvent = fEventMaker->FillTracks(mcEvent,fCFManager1,fCFManager2);
1c1d4332 282 }
8c5b0db3 283
62726ef0 284 // Fill the FlowEventSimple for ESD input
1c1d4332 285 else if (fAnalysisType == "ESD") {
587b0a35 286 if (!fCFManager1) {cout << "ERROR: No pointer to correction framework cuts! " << endl; return; }
287 if (!fCFManager2) {cout << "ERROR: No pointer to correction framework cuts! " << endl; return; }
8c5b0db3 288
289 if (!myESD) { Printf("ERROR: ESD not available"); return;}
b0569f96 290 //check the offline trigger (check if the event has the correct trigger)
8c5b0db3 291 Printf("There are %d tracks in this event", fInputEvent->GetNumberOfTracks());
292 // analysis
293 fEvent = fEventMaker->FillTracks(myESD,fCFManager1,fCFManager2);
1c1d4332 294 }
62726ef0 295 // Fill the FlowEventSimple for ESD input combined with MC info
1c1d4332 296 else if (fAnalysisType == "ESDMC0" || fAnalysisType == "ESDMC1" ) {
587b0a35 297 if (!fCFManager1) {cout << "ERROR: No pointer to correction framework cuts! " << endl; return; }
298 if (!fCFManager2) {cout << "ERROR: No pointer to correction framework cuts! " << endl; return; }
8c5b0db3 299
300 if (!myESD) { Printf("ERROR: ESD not available"); return;}
65201059 301 Printf("There are %d tracks in this event", fInputEvent->GetNumberOfTracks());
1c1d4332 302
587b0a35 303 if (!mcEvent) {Printf("ERROR: Could not retrieve MC event"); return;}
1c1d4332 304
6d7734d2 305 fCFManager1->SetMCEventInfo(mcEvent);
306 fCFManager2->SetMCEventInfo(mcEvent);
1c1d4332 307
489fdf79 308
1c1d4332 309 if (fAnalysisType == "ESDMC0") {
8c5b0db3 310 fEvent = fEventMaker->FillTracks(myESD, mcEvent, fCFManager1, fCFManager2, 0); //0 = kine from ESD, 1 = kine from MC
1c1d4332 311 } else if (fAnalysisType == "ESDMC1") {
8c5b0db3 312 fEvent = fEventMaker->FillTracks(myESD, mcEvent, fCFManager1, fCFManager2, 1); //0 = kine from ESD, 1 = kine from MC
1c1d4332 313 }
1c1d4332 314 }
62726ef0 315 // Fill the FlowEventSimple for AOD input
1c1d4332 316 else if (fAnalysisType == "AOD") {
8c5b0db3 317 if (!myAOD) {Printf("ERROR: AOD not available"); return;}
318 Printf("There are %d tracks in this event", myAOD->GetNumberOfTracks());
1c1d4332 319
320 // analysis
65201059 321 //For the moment don't use CF //AliFlowEventSimple* fEvent = fEventMaker->FillTracks(fOutputAOD,fCFManager1,fCFManager2);
8c5b0db3 322 fEvent = fEventMaker->FillTracks(myAOD);
1c1d4332 323 }
324
46bec39c 325 //fListHistos->Print();
326 // fOutputFile->WriteObject(fEvent,"myFlowEventSimple");
65201059 327 PostData(1,fEvent);
1c1d4332 328 if (fQA) {
65201059 329 PostData(2,fQAInt);
330 PostData(3,fQADiff); }
1c1d4332 331}
332
333//________________________________________________________________________
334void AliAnalysisTaskFlowEvent::Terminate(Option_t *)
335{
336 // Called once at the end of the query -- do not call in case of CAF
337
338}
489fdf79 339
340