AliAODEvent::GetHeader() returns AliVHeader
[u/mrichter/AliRoot.git] / PWGCF / FEMTOSCOPY / AliFemto / AliAnalysisTaskFemto.cxx
1 //------------------------------------------------------
2 // AliAnalysisTaskFemto - A task for the analysis framework
3 // from the FEMTOSCOPY analysis of PWG2. Creates the necessary
4 // connection between the ESD or AOD input and the femtoscopic
5 // code.
6 // Author: Adam Kisiel, OSU; Adam.Kisiel@cern.ch
7 //------------------------------------------------------
8 #include "TROOT.h"
9 #include "TChain.h"
10 #include "TH1.h"
11 #include "TCanvas.h"
12 #include "TSystem.h"
13 #include "TFile.h"
14 #include "TInterpreter.h"
15
16 //#include "AliAnalysisTask.h"
17 #include "AliAnalysisTaskSE.h"
18
19 #include "AliESDEvent.h"
20
21 #include "AliFemtoAnalysis.h"
22 #include "AliAnalysisTaskFemto.h"
23 #include "AliVHeader.h"
24 #include "AliGenEventHeader.h"
25 #include "AliGenHijingEventHeader.h"
26 #include "AliGenCocktailEventHeader.h"
27
28 ClassImp(AliAnalysisTaskFemto)
29
30 // Default name for the setup macro of femto analysis  
31 // This function MUST be defined in the separate file !!!
32 // extern AliFemtoManager *ConfigFemtoAnalysis();
33
34 //________________________________________________________________________
35 AliAnalysisTaskFemto::AliAnalysisTaskFemto(TString name, TString aConfigMacro, TString aConfigParams, Bool_t aVerbose):
36 AliAnalysisTaskSE(name), //AliAnalysisTask(name,""), 
37   fESD(0), 
38   fESDpid(0),
39   fAOD(0),
40   fAODpidUtil(0),
41   fAODheader(0),
42   fStack(0),
43   fOutputList(0), 
44   fReader(0x0),
45   fManager(0x0),
46   fAnalysisType(0),
47     fConfigMacro(aConfigMacro),
48     fConfigParams(aConfigParams),
49     fVerbose(aVerbose)
50 {
51   // Constructor.
52   // Input slot #0 works with an Ntuple
53   //DefineInput(0, TChain::Class());
54   // Output slot #0 writes into a TH1 container
55   DefineOutput(0, TList::Class());
56
57 }
58 //________________________________________________________________________
59 AliAnalysisTaskFemto::AliAnalysisTaskFemto(TString name, TString aConfigMacro="ConfigFemtoAnalysis.C", Bool_t aVerbose): 
60     AliAnalysisTaskSE(name), //AliAnalysisTask(name,""), 
61     fESD(0), 
62     fESDpid(0),
63     fAOD(0),
64     fAODpidUtil(0),
65     fAODheader(0),
66     fStack(0),
67     fOutputList(0), 
68     fReader(0x0),
69     fManager(0x0),
70     fAnalysisType(0),
71     fConfigMacro(aConfigMacro),
72     fConfigParams(""),
73     fVerbose(aVerbose)
74 {
75   // Constructor.
76   // Input slot #0 works with an Ntuple
77   //DefineInput(0, TChain::Class());
78   // Output slot #0 writes into a TH1 container
79   DefineOutput(0, TList::Class());
80
81 }
82
83 AliAnalysisTaskFemto::AliAnalysisTaskFemto(const AliAnalysisTaskFemto& aFemtoTask):
84     AliAnalysisTaskSE(aFemtoTask), //AliAnalysisTask(aFemtoTask), 
85     fESD(0), 
86     fESDpid(0),
87     fAOD(0),
88     fAODpidUtil(0),
89     fAODheader(0),
90     fStack(0),
91     fOutputList(0), 
92     fReader(0x0),
93     fManager(0x0),
94     fAnalysisType(0),
95     fConfigMacro(0),
96     fConfigParams(0),
97     fVerbose(kFALSE)
98 {
99   // copy constructor
100   fESD = aFemtoTask.fESD; 
101   fESDpid = aFemtoTask.fESDpid; 
102   fAOD = aFemtoTask.fAOD; 
103   fAODpidUtil = aFemtoTask.fAODpidUtil;
104   fAODheader = aFemtoTask.fAODheader;
105   fStack = aFemtoTask.fStack;
106   fOutputList = aFemtoTask.fOutputList;   
107   fReader = aFemtoTask.fReader;       
108   fManager = aFemtoTask.fManager;      
109   fAnalysisType = aFemtoTask.fAnalysisType; 
110
111   fConfigMacro = aFemtoTask.fConfigMacro;
112   fConfigParams = aFemtoTask.fConfigParams;
113
114     fVerbose = aFemtoTask.fVerbose;
115 }
116
117
118 AliAnalysisTaskFemto& AliAnalysisTaskFemto::operator=(const AliAnalysisTaskFemto& aFemtoTask){
119   // assignment operator
120   if (this == &aFemtoTask)
121     return *this;
122
123   fESD = aFemtoTask.fESD; 
124   fESDpid = aFemtoTask.fESDpid;
125   fAOD = aFemtoTask.fAOD; 
126   fAODpidUtil = aFemtoTask.fAODpidUtil;
127   fAODheader = aFemtoTask.fAODheader;
128   fStack = aFemtoTask.fStack;
129   fOutputList = aFemtoTask.fOutputList;   
130   fReader = aFemtoTask.fReader;       
131   fManager = aFemtoTask.fManager;      
132   fAnalysisType = aFemtoTask.fAnalysisType; 
133
134   fConfigMacro = aFemtoTask.fConfigMacro;
135   fConfigParams = aFemtoTask.fConfigParams;
136     fVerbose = aFemtoTask.fVerbose;
137
138
139   return *this;
140 }
141
142 AliAnalysisTaskFemto::~AliAnalysisTaskFemto() 
143 {
144 }
145
146
147 //________________________________________________________________________
148 void AliAnalysisTaskFemto::ConnectInputData(Option_t *) {
149   AliInfo(Form("   ConnectInputData %s\n", GetName()));
150
151   fESD = 0;
152   fESDpid = 0;
153   fAOD = 0;
154   fAODpidUtil = 0;
155   fAODheader=0;
156   fAnalysisType = 0;
157
158   TTree* tree = dynamic_cast<TTree*> (GetInputData(0));
159   if (!tree) {
160     AliWarning("Could not read chain from input slot 0");
161     return;
162   } 
163
164   AliFemtoEventReaderESDChain *femtoReader = dynamic_cast<AliFemtoEventReaderESDChain *> (fReader);
165   if ((dynamic_cast<AliFemtoEventReaderESDChain *> (fReader))) {
166     AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
167     
168     if(esdH) {
169       if (fVerbose)
170         AliInfo("Selected ESD analysis");
171       fAnalysisType = 1;
172       
173       fESD = esdH->GetEvent();
174       fESDpid = esdH->GetESDpid();
175       femtoReader->SetESDPid(fESDpid);
176     }
177   }
178   else if ((dynamic_cast<AliFemtoEventReaderKinematicsChain *> (fReader))) {
179     AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
180     
181     if(esdH) {
182       if (fVerbose)
183         AliInfo("Selected ESD analysis");
184       fAnalysisType = 1;
185       fESD = esdH->GetEvent();
186     }
187   }
188  else if ((dynamic_cast<AliFemtoEventReaderKinematicsChainESD *> (fReader))) {
189     AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
190     if(esdH) {
191       if (fVerbose)
192         AliInfo("Selected ESD analysis");
193       fAnalysisType = 1;
194       fESD = esdH->GetEvent();
195     }
196   }
197
198   AliFemtoEventReaderESDChainKine *femtoReaderESDKine = dynamic_cast<AliFemtoEventReaderESDChainKine *> (fReader);
199   if ((dynamic_cast<AliFemtoEventReaderESDChainKine *> (fReader))) {
200     AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
201     
202     if(esdH) {
203       if (fVerbose)           
204         AliInfo("Selected ESD analysis");
205       fAnalysisType = 1;
206       fESD = esdH->GetEvent();
207       fESDpid = esdH->GetESDpid();
208       femtoReaderESDKine->SetESDPid(fESDpid);
209
210     }
211   }
212
213   //  AliFemtoEventReaderKinematicsChain *femtoReaderKine = dynamic_cast<AliFemtoEventReaderKinematicsChain *> (fReader);
214   if ((dynamic_cast<AliFemtoEventReaderKinematicsChain *> (fReader))) {
215     AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
216     
217     if(esdH) {
218             if (fVerbose)
219       AliInfo("Selected ESD analysis");
220       fAnalysisType = 1;
221       
222       //       if (!esdH) {
223       //        AliWarning("Could not get ESDInputHandler");
224       //       } 
225       //       else {
226       fESD = esdH->GetEvent();
227       //fESDpid = esdH->GetESDpid();
228       //femtoReader->SetESDPid(fESDpid);
229       //       }
230     }
231   }
232   
233   
234   AliFemtoEventReaderAODChain *femtoReaderAOD = dynamic_cast<AliFemtoEventReaderAODChain *> (fReader);
235   if (dynamic_cast<AliFemtoEventReaderAODChain *> (fReader)) {
236     AliAODInputHandler *aodH = dynamic_cast<AliAODInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
237     
238     if (!aodH) {
239       TObject *handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
240             if (fVerbose)
241       AliInfo("Has output handler ");
242       if( handler && handler->InheritsFrom("AliAODHandler") ) {
243                 if (fVerbose)
244         AliInfo("Selected AOD analysis");
245         
246         fAOD = ((AliAODHandler*)handler)->GetAOD();
247         fAnalysisType = 2;
248       }
249       else {
250                 if (fVerbose)
251         AliWarning("Selected AOD reader but no AOD handler found");
252       }
253     } 
254     else {
255             if (fVerbose)
256       AliInfo("Selected AOD analysis");
257       fAnalysisType = 2;
258       
259       fAOD = aodH->GetEvent();
260
261       fAODpidUtil = aodH->GetAODpidUtil(); //correct way
262       //fAODpidUtil = new AliAODpidUtil(); //not correct way
263       //      printf("aodH->GetAODpidUtil(): %x",aodH->GetAODpidUtil());
264             if (fVerbose)
265       cout<<"AliAnalysisTaskFemto::AodpidUtil:"<<fAODpidUtil<<endl;
266       femtoReaderAOD->SetAODpidUtil(fAODpidUtil);
267
268       fAODheader = dynamic_cast<AliAODHeader*>(fAOD->GetHeader());
269       if(!fAODheader) AliFatal("Not a standard AOD");
270       femtoReaderAOD->SetAODheader(fAODheader);
271    
272     }
273   }
274
275   if ((!fAOD) && (!fESD)) {
276         if (fVerbose)
277     AliWarning("Wrong analysis type: Only ESD and AOD types are allowed!");
278   }
279 }
280
281 //________________________________________________________________________
282 void AliAnalysisTaskFemto::CreateOutputObjects() {
283     if (fVerbose)
284   AliInfo("Creating Femto Analysis objects\n");
285
286   gSystem->SetIncludePath("-I$ROOTSYS/include -I./STEERBase/ -I./ESD/ -I./AOD/ -I./ANALYSIS/ -I./ANALYSISalice/ -I./PWG2AOD/AOD -I./PWG2femtoscopy/FEMTOSCOPY/AliFemto -I./PWG2femtoscopyUser/FEMTOSCOPY/AliFemtoUser");
287   //  char fcm[2000];
288 //   sprintf(fcm, "%s++", fConfigMacro);
289 //   gROOT->LoadMacro(fcm);
290   gROOT->LoadMacro(fConfigMacro);
291   //  fJetFinder = (AliJetFinder*) gInterpreter->ProcessLine("ConfigJetAnalysis()");
292   if (!fConfigParams)
293     SetFemtoManager((AliFemtoManager *) gInterpreter->ProcessLine("ConfigFemtoAnalysis()"));
294   else
295       SetFemtoManager((AliFemtoManager *) gInterpreter->ProcessLine(Form("ConfigFemtoAnalysis(%s)", fConfigParams.Data())));
296
297   TList *tOL;
298   fOutputList = fManager->Analysis(0)->GetOutputList();
299
300   for (unsigned int ian = 1; ian<fManager->AnalysisCollection()->size(); ian++) {
301     tOL = fManager->Analysis(ian)->GetOutputList();
302
303     TIter nextListCf(tOL);
304     while (TObject *obj = nextListCf()) {
305       fOutputList->Add(obj);
306     }
307
308     delete tOL;
309   }
310
311   PostData(0, fOutputList);
312 }
313
314 //________________________________________________________________________
315 void AliAnalysisTaskFemto::Exec(Option_t *) {
316   // Task making a femtoscopic analysis.
317   if(fOfflineTriggerMask){
318     Bool_t isSelected = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & fOfflineTriggerMask);
319     if(!isSelected) {
320       if (fVerbose)
321         cout << "AliAnalysisTaskFemto: is not selected" << endl; 
322       return;
323     }
324   }
325   if (fAnalysisType==1) {
326     if (!fESD) {
327       if (fVerbose)
328         AliWarning("fESD not available");
329       return;
330     }
331     //Get MC data
332     AliMCEventHandler*    mctruth = (AliMCEventHandler*) 
333       ((AliAnalysisManager::GetAnalysisManager())->GetMCtruthEventHandler());
334     
335     AliGenHijingEventHeader *hdh = 0;
336     if(mctruth) {
337       fStack = mctruth->MCEvent()->Stack();
338
339       AliGenCocktailEventHeader *hd = dynamic_cast<AliGenCocktailEventHeader *> (mctruth->MCEvent()->GenEventHeader());
340       
341       if (hd) {
342         
343         //      AliInfo ("Got MC cocktail event header %p\n", (void *) hd);
344         TList *lhd = hd->GetHeaders();
345         //      AliInfo ("Got list of headers %d\n", lhd->GetEntries());
346         
347         for (int iterh=0; iterh<lhd->GetEntries(); iterh++) 
348           {
349             hdh = dynamic_cast<AliGenHijingEventHeader *> (lhd->At(iterh));
350             //      AliInfo ("HIJING header at %i is %p\n", iterh, (void *) hdh);
351           }
352       }    
353     }
354
355     // Get ESD
356     AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
357     
358     if (!esdH) {
359       if (fVerbose)
360         AliWarning("Could not get ESDInputHandler");
361       return;
362     } 
363     else {
364       fESD = esdH->GetEvent();
365       fESDpid = esdH->GetESDpid();   
366      }
367
368     if (fVerbose)
369       AliInfo(Form("Tracks in ESD: %d \n",fESD->GetNumberOfTracks()));
370
371     if (fESD->GetNumberOfTracks() >= 0) {
372     
373       if (!fReader) {
374         if (fVerbose)
375           AliWarning("No ESD reader for ESD analysis !\n");
376       }
377       
378       AliFemtoEventReaderESDChain* fesdc = dynamic_cast<AliFemtoEventReaderESDChain *> (fReader);
379       if (fesdc)
380         {
381           // Process the event with no Kine information
382           fesdc->SetESDSource(fESD);
383           fManager->ProcessEvent();
384         }
385     }
386     AliFemtoEventReaderKinematicsChain* fkinec = dynamic_cast<AliFemtoEventReaderKinematicsChain *> (fReader);
387     if (fkinec)
388       {
389         // Process the event with Kine information only
390         fkinec->SetStackSource(fStack);
391         fManager->ProcessEvent();
392       }
393
394
395     AliFemtoEventReaderESDChainKine* fesdck = dynamic_cast<AliFemtoEventReaderESDChainKine *> (fReader);
396     if (fesdck) 
397       {
398         // Process the event with Kine information
399         fesdck->SetESDSource(fESD);
400         fesdck->SetStackSource(fStack);
401         cout<<"Set Stack:"<<fStack<<endl;
402         fesdck->SetGenEventHeader(hdh);
403         fManager->ProcessEvent();
404       }
405
406      
407     AliFemtoEventReaderKinematicsChainESD* fkcesd = dynamic_cast<AliFemtoEventReaderKinematicsChainESD *> (fReader);
408     if (fkcesd)
409       {
410         // Process the event with Kine information
411         fkcesd->SetESDSource(fESD);
412         fkcesd->SetStackSource(fStack);
413         fkcesd->SetGenEventHeader(hdh);
414         fManager->ProcessEvent();
415       }
416
417     AliFemtoEventReaderStandard* fstd = dynamic_cast<AliFemtoEventReaderStandard *> (fReader);
418     if (fstd) 
419       {
420         // Process the event with Kine information
421         fstd->SetESDSource(fESD);
422         if (mctruth) {
423           fstd->SetStackSource(fStack);
424           fstd->SetGenEventHeader(hdh);
425           fstd->SetInputType(AliFemtoEventReaderStandard::kESDKine);
426         }
427         else
428           fstd->SetInputType(AliFemtoEventReaderStandard::kESD);
429         fManager->ProcessEvent();
430       }
431   
432
433     // Post the output histogram list
434     PostData(0, fOutputList);
435   }
436
437   if (fAnalysisType==2) {    
438     if (!fAOD) {
439             if (fVerbose)
440       AliWarning("fAOD not available");
441       return;
442     }
443
444     // Get AOD
445 //     AliAODInputHandler *aodH = dynamic_cast<AliAODInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
446       
447 //     if (!aodH) {
448 //       AliWarning("Could not get AODInputHandler");
449 //       return;
450 //     } 
451 //     else {
452
453 //       fAOD = aodH->GetEvent();
454 //     }
455
456
457
458
459     if (fVerbose)
460       AliInfo(Form("Tracks in AOD: %d \n",fAOD->GetNumberOfTracks()));
461     
462     if (fAOD->GetNumberOfTracks() > 0) {
463       if (!fReader) {
464         if (fVerbose)
465           AliWarning("No AOD reader for AOD analysis! \n");
466       }
467       else {
468         AliFemtoEventReaderAODChain* faodc = dynamic_cast<AliFemtoEventReaderAODChain *> (fReader);
469
470         if (faodc) {
471           // Process the event
472           faodc->SetAODSource(fAOD);
473           fManager->ProcessEvent();
474         }
475         AliFemtoEventReaderStandard* fstd = dynamic_cast<AliFemtoEventReaderStandard *> (fReader);
476
477         if (fstd) {
478           // Process the event
479           fstd->SetAODSource(fAOD);
480           fstd->SetInputType(AliFemtoEventReaderStandard::kAOD);
481           fManager->ProcessEvent();
482         }
483       }
484     } 
485
486     // Post the output histogram list
487     PostData(0, fOutputList);
488   }
489 }      
490
491 //________________________________________________________________________
492 void AliAnalysisTaskFemto::Terminate(Option_t *) {
493   // Do the final processing
494   if (fManager) {
495     fManager->Finish();
496   }
497 }
498 //________________________________________________________________________
499 void AliAnalysisTaskFemto:: FinishTaskOutput() {
500   // Do the final processing
501   if (fManager) {
502     fManager->Finish();
503   }
504 }
505 //________________________________________________________________________
506 void AliAnalysisTaskFemto::SetFemtoReaderESD(AliFemtoEventReaderESDChain *aReader)
507 {
508     if (fVerbose)
509   AliInfo("Selecting Femto reader for ESD\n");
510   fReader = aReader;
511 }
512 //________________________________________________________________________
513 void AliAnalysisTaskFemto::SetFemtoReaderESDKine(AliFemtoEventReaderESDChainKine *aReader)
514 {
515     if (fVerbose)
516   AliInfo("Selecting Femto reader for ESD with Kinematics information\n");
517   fReader = aReader;
518 }
519 //________________________________________________________________________
520 void AliAnalysisTaskFemto::SetFemtoReaderAOD(AliFemtoEventReaderAODChain *aReader)
521 {
522     if (fVerbose)
523   AliInfo("Selecting Femto reader for AOD\n");
524   fReader = aReader;
525 }
526 void AliAnalysisTaskFemto::SetFemtoReaderStandard(AliFemtoEventReaderStandard *aReader)
527 {
528     if (fVerbose)
529   AliInfo("Selecting Standard all-purpose Femto reader\n");
530   fReader = aReader;
531 }
532 void AliAnalysisTaskFemto::SetFemtoReaderKinematics(AliFemtoEventReaderKinematicsChain *aReader)
533 {
534     if (fVerbose)
535   printf("Selecting Femto reader for Kinematics (Monte Carlo) information\n");
536   fReader = aReader;
537 }
538 void AliAnalysisTaskFemto::SetFemtoReaderKinematicsESD(AliFemtoEventReaderKinematicsChainESD *aReader)
539 {
540   if (fVerbose)
541     printf("Selecting Femto reader for Kinematics (Monte Carlo) information + ESD\n");
542   fReader = aReader;
543 }
544 //________________________________________________________________________
545 void AliAnalysisTaskFemto::SetFemtoManager(AliFemtoManager *aManager)
546 {
547   fManager = aManager;
548     if (fVerbose)
549   AliInfo(Form("Got reader %p\n", (void *) aManager->EventReader()));
550   AliFemtoEventReaderESDChain     *tReaderESDChain     = dynamic_cast<AliFemtoEventReaderESDChain *> (aManager->EventReader());
551   AliFemtoEventReaderESDChainKine *tReaderESDChainKine = dynamic_cast<AliFemtoEventReaderESDChainKine *> (aManager->EventReader());
552   AliFemtoEventReaderAODChain     *tReaderAODChain     = dynamic_cast<AliFemtoEventReaderAODChain *> (aManager->EventReader());
553   AliFemtoEventReaderStandard     *tReaderStandard     = dynamic_cast<AliFemtoEventReaderStandard *> (aManager->EventReader());
554   AliFemtoEventReaderKinematicsChain *tReaderKineChain = dynamic_cast<AliFemtoEventReaderKinematicsChain *> (aManager->EventReader());
555   AliFemtoEventReaderKinematicsChainESD *tReaderKineChainESD = dynamic_cast<AliFemtoEventReaderKinematicsChainESD *> (aManager->EventReader());
556
557  if ((!tReaderESDChain) && (!tReaderESDChainKine) && (!tReaderAODChain) && (!tReaderStandard) && (!tReaderKineChain)  && (!tReaderKineChainESD)) {
558         if (fVerbose)
559     AliWarning("No AliFemto event reader created. Will not run femto analysis.\n");
560     return;
561   }
562   if (tReaderESDChain) SetFemtoReaderESD(tReaderESDChain);
563   if (tReaderESDChainKine) SetFemtoReaderESDKine(tReaderESDChainKine);
564   if (tReaderAODChain) SetFemtoReaderAOD(tReaderAODChain);
565   if (tReaderStandard) SetFemtoReaderStandard(tReaderStandard);
566   if (tReaderKineChain) SetFemtoReaderKinematics(tReaderKineChain);
567   if (tReaderKineChainESD) SetFemtoReaderKinematicsESD(tReaderKineChainESD);
568 }
569