1 /*************************************************************************
2 * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 #include "Riostream.h" //needed as include
23 class AliAnalysisTask;
24 #include "AliAnalysisManager.h"
26 #include "AliESDEvent.h"
27 #include "AliESDInputHandler.h"
29 #include "AliAODEvent.h"
30 #include "AliAODInputHandler.h"
32 #include "AliMCEventHandler.h"
33 #include "AliMCEvent.h"
35 #include "AliCFManager.h"
37 #include "AliAnalysisTaskLYZEventPlane.h"
38 #include "AliFlowEventSimpleMaker.h"
39 #include "AliFlowCommonHist.h"
40 #include "AliFlowCommonHistResults.h"
41 #include "AliFlowLYZEventPlane.h"
42 #include "AliFlowAnalysisWithLYZEventPlane.h"
44 // AliAnalysisTaskLYZEventPlane:
46 // analysis task for Lee Yang Zeros Event Plane
48 // Author: Naomi van der Kolk (kolk@nikhef.nl)
50 ClassImp(AliAnalysisTaskLYZEventPlane)
52 //________________________________________________________________________
53 AliAnalysisTaskLYZEventPlane::AliAnalysisTaskLYZEventPlane(const char *name, Bool_t on) :
54 AliAnalysisTask(name, ""),
70 cout<<"AliAnalysisTaskLYZEventPlane::AliAnalysisTaskLYZEventPlane(const char *name)"<<endl;
72 // Define input and output slots here
73 // Input slot #0 works with a TChain
74 DefineInput(0, TChain::Class());
75 DefineInput(1, TList::Class());
76 // Output slot #0 writes into a TList container
77 DefineOutput(0, TList::Class());
79 DefineOutput(1, TList::Class());
80 DefineOutput(2, TList::Class()); }
83 //________________________________________________________________________
84 AliAnalysisTaskLYZEventPlane::AliAnalysisTaskLYZEventPlane() :
100 cout<<"AliAnalysisTaskLYZEventPlane::AliAnalysisTaskLYZEventPlane()"<<endl;
104 //________________________________________________________________________
105 AliAnalysisTaskLYZEventPlane::~AliAnalysisTaskLYZEventPlane()
111 //________________________________________________________________________
112 void AliAnalysisTaskLYZEventPlane::ConnectInputData(Option_t *)
114 // Connect ESD or AOD here
116 cout<<"AliAnalysisTaskLYZEventPlane::ConnectInputData(Option_t *)"<<endl;
118 TTree* tree = dynamic_cast<TTree*> (GetInputData(0));
120 Printf("ERROR: Could not read chain from input slot 0");
122 // Disable all branches and enable only the needed ones
123 if (fAnalysisType == "MC") {
124 // we want to process only MC
125 tree->SetBranchStatus("*", kFALSE);
127 AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
130 Printf("ERROR: Could not get ESDInputHandler");
132 fESD = esdH->GetEvent();
135 else if (fAnalysisType == "ESD" || fAnalysisType == "ESDMC0" || fAnalysisType == "ESDMC1" ) {
136 tree->SetBranchStatus("*", kFALSE);
137 tree->SetBranchStatus("Tracks.*", kTRUE);
139 AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
142 Printf("ERROR: Could not get ESDInputHandler");
144 fESD = esdH->GetEvent();
146 else if (fAnalysisType == "AOD") {
147 AliAODInputHandler *aodH = dynamic_cast<AliAODInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
150 Printf("ERROR: Could not get AODInputHandler");
153 fAOD = aodH->GetEvent();
157 Printf("Wrong analysis type: Only ESD, ESDMC0, ESDMC1, AOD and MC types are allowed!");
164 //________________________________________________________________________
165 void AliAnalysisTaskLYZEventPlane::CreateOutputObjects()
168 cout<<"AliAnalysisTaskLYZEventPlane::CreateOutputObjects()"<<endl;
170 if (!(fAnalysisType == "AOD" || fAnalysisType == "ESD" || fAnalysisType == "ESDMC0" || fAnalysisType == "ESDMC1" || fAnalysisType == "MC")) {
171 cout<<"WRONG ANALYSIS TYPE! only ESD, , ESDMC0, ESDMC1, AOD and MC are allowed."<<endl;
176 fEventMaker = new AliFlowEventSimpleMaker();
177 //lee yang zeros event plane
178 fLyzEp = new AliFlowLYZEventPlane() ;
180 fLyz = new AliFlowAnalysisWithLYZEventPlane() ;
182 // Get data from input slot
183 TList* pSecondRunList = (TList*)GetInputData(1);
184 if (pSecondRunList) {
185 fLyzEp -> SetSecondRunList(pSecondRunList);
186 fLyz -> SetSecondRunList(pSecondRunList);
187 } else { cout<<"No Second run List!"<<endl; exit(0); }
192 if (fLyz->GetHistList()) {
193 fListHistos = fLyz->GetHistList();
194 fListHistos->Print();
196 else { cout<<"ERROR: Could not retrieve histogram list"<<endl;}
201 //________________________________________________________________________
202 void AliAnalysisTaskLYZEventPlane::Exec(Option_t *)
205 // Called for each event
207 if (fAnalysisType == "MC") {
208 // Process MC truth, therefore we receive the AliAnalysisManager and ask it for the AliMCEventHandler
209 // This handler can return the current MC event
211 AliMCEventHandler* eventHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
213 Printf("ERROR: Could not retrieve MC event handler");
217 AliMCEvent* mcEvent = eventHandler->MCEvent();
219 Printf("ERROR: Could not retrieve MC event");
223 fCFManager1->SetEventInfo(mcEvent);
224 fCFManager2->SetEventInfo(mcEvent);
226 Printf("MC particles: %d", mcEvent->GetNumberOfTracks());
228 //lee yang zeros analysis
229 AliFlowEventSimple* fEvent = fEventMaker->FillTracks(mcEvent,fCFManager1,fCFManager2);
230 fLyz->Make(fEvent,fLyzEp);
234 else if (fAnalysisType == "ESD") {
236 Printf("ERROR: fESD not available");
239 Printf("There are %d tracks in this event", fESD->GetNumberOfTracks());
241 //lee yang zeros analysis
242 AliFlowEventSimple* fEvent = fEventMaker->FillTracks(fESD,fCFManager1,fCFManager2);
243 fLyz->Make(fEvent,fLyzEp);
246 else if (fAnalysisType == "ESDMC0" || fAnalysisType == "ESDMC1") {
248 Printf("ERROR: fESD not available");
251 Printf("There are %d tracks in this event", fESD->GetNumberOfTracks());
253 AliMCEventHandler* eventHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
255 Printf("ERROR: Could not retrieve MC event handler");
259 AliMCEvent* mcEvent = eventHandler->MCEvent();
261 Printf("ERROR: Could not retrieve MC event");
265 fCFManager1->SetEventInfo(mcEvent);
266 fCFManager2->SetEventInfo(mcEvent);
268 //lee yang zeros analysis
269 AliFlowEventSimple* fEvent = NULL;
270 if (fAnalysisType == "ESDMC0") {
271 fEvent = fEventMaker->FillTracks(fESD, mcEvent,fCFManager1, fCFManager2,0); //0 = kine from ESD, 1 = kine from MC
272 } else if (fAnalysisType == "ESDMC1") {
273 fEvent = fEventMaker->FillTracks(fESD, mcEvent,fCFManager1, fCFManager2,0); //0 = kine from ESD, 1 = kine from MC
276 fLyz->Make(fEvent,fLyzEp);
281 else if (fAnalysisType == "AOD") {
283 Printf("ERROR: fAOD not available");
286 Printf("There are %d tracks in this event", fAOD->GetNumberOfTracks());
288 //lee yang zeros analysis
289 //for the moment don't use CF
290 //AliFlowEventSimple* fEvent = fEventMaker->FillTracks(fAOD,fCFManager1,fCFManager2);
291 AliFlowEventSimple* fEvent = fEventMaker->FillTracks(fAOD);
292 fLyz->Make(fEvent,fLyzEp);
296 PostData(0,fListHistos);
299 PostData(2,fQADiff); }
302 //________________________________________________________________________
303 void AliAnalysisTaskLYZEventPlane::Terminate(Option_t *)
305 // Called once at the end of the query
306 AliFlowAnalysisWithLYZEventPlane* fLyzTerm = new AliFlowAnalysisWithLYZEventPlane() ;
307 fListHistos = (TList*)GetOutputData(0);
308 //cout << "histogram list in Terminate" << endl;
310 //Get the common histograms from the output list
311 AliFlowCommonHist *pCommonHist = dynamic_cast<AliFlowCommonHist*>
312 (fListHistos->FindObject("AliFlowCommonHistLYZEP"));
313 AliFlowCommonHistResults *pCommonHistResults = dynamic_cast<AliFlowCommonHistResults*>
314 (fListHistos->FindObject("AliFlowCommonHistResultsLYZEP"));
316 TProfile* pHistProR0theta = dynamic_cast<TProfile*>
317 (fListHistos->FindObject("First_FlowPro_r0theta_LYZ"));
319 TProfile* pHistProVetaRP = dynamic_cast<TProfile*>
320 (fListHistos->FindObject("FlowPro_VetaRP_LYZEP"));
321 TProfile* pHistProVetaPOI = dynamic_cast<TProfile*>
322 (fListHistos->FindObject("FlowPro_VetaPOI_LYZEP"));
323 TProfile* pHistProVPtRP = dynamic_cast<TProfile*>
324 (fListHistos->FindObject("FlowPro_VPtRP_LYZEP"));
325 TProfile* pHistProVPtPOI = dynamic_cast<TProfile*>
326 (fListHistos->FindObject("FlowPro_VPtPOI_LYZEP"));
328 TH1F* pHistQsumforChi = dynamic_cast<TH1F*>
329 (fListHistos->FindObject("Flow_QsumforChi_LYZEP"));
331 if (pCommonHist && pCommonHistResults && pHistProR0theta &&
332 pHistProVetaRP && pHistProVetaPOI && pHistProVPtRP &&
333 pHistProVPtPOI && pHistQsumforChi ) {
334 fLyzTerm -> SetCommonHists(pCommonHist);
335 fLyzTerm -> SetCommonHistsRes(pCommonHistResults);
336 fLyzTerm -> SetFirstr0theta(pHistProR0theta);
337 fLyzTerm -> SetHistProVetaRP(pHistProVetaRP);
338 fLyzTerm -> SetHistProVetaPOI(pHistProVetaPOI);
339 fLyzTerm -> SetHistProVPtRP(pHistProVPtRP);
340 fLyzTerm -> SetHistProVPtPOI(pHistProVPtPOI);
341 fLyzTerm -> SetHistQsumforChi(pHistQsumforChi);
342 fLyzTerm -> Finish();
343 PostData(0,fListHistos);
345 cout<<"WARNING: Histograms needed to run Finish() are not accessable!"<<endl;
348 //fListHistos->Print();
349 } else { cout << "histogram list pointer is empty" << endl;}
351 //cout<<".....finished LYZ EventPlane"<<endl;