]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG4/PartCorrBase/AliCaloTrackReader.cxx
Update AliEMCALRecoUtils with new different method to recalculate position, posibilit...
[u/mrichter/AliRoot.git] / PWG4 / PartCorrBase / AliCaloTrackReader.cxx
CommitLineData
1c5acb87 1/**************************************************************************
2 * Copyright(c) 1998-1999, 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/* $Id: $ */
16
17//_________________________________________________________________________
18// Base class for reading data: MonteCarlo, ESD or AOD, of PHOS EMCAL and
19// Central Barrel Tracking detectors (CTS).
ff45398a 20// Not all MC particles/tracks/clusters are kept, some kinematical/fiducial restrictions are done.
591cc579 21// Mother class of : AliCaloTrackESDReader: Fills ESD data in 3 TObjArrays (PHOS, EMCAL, CTS)
22// : AliCaloTrackMCReader: Fills Kinematics data in 3 TObjArrays (PHOS, EMCAL, CTS)
f37fa8d2 23// : AliCaloTrackReader: Fills AOD data in 3 TObjArrays (PHOS, EMCAL, CTS)
1e68a3f4 24//
25// This part is commented: Mixing analysis can be done, input AOD with events
26// is opened in the AliCaloTrackReader::Init()
27
1c5acb87 28//-- Author: Gustavo Conesa (LNF-INFN)
29//////////////////////////////////////////////////////////////////////////////
30
31
32// --- ROOT system ---
591cc579 33#include "TFile.h"
1c5acb87 34
35//---- ANALYSIS system ----
36#include "AliCaloTrackReader.h"
477d6cee 37#include "AliMCEvent.h"
591cc579 38#include "AliAODMCHeader.h"
39#include "AliGenPythiaEventHeader.h"
c8fe2783 40#include "AliVEvent.h"
8dacfd76 41#include "AliAODEvent.h"
c8fe2783 42#include "AliVTrack.h"
43#include "AliVParticle.h"
44#include "AliMixedEvent.h"
0ae57829 45#include "AliESDtrack.h"
9584c261 46#include "AliEMCALRecoUtils.h"
1c5acb87 47
48ClassImp(AliCaloTrackReader)
49
50
51//____________________________________________________________________________
52 AliCaloTrackReader::AliCaloTrackReader() :
a79a2424 53 TObject(), fEventNumber(-1), fCurrentFileName(""),fDataType(0), fDebug(0),
decca433 54 fFiducialCut(0x0), fCheckFidCut(kFALSE), fComparePtHardAndJetPt(kFALSE), fPtHardAndJetPtFactor(7),
f37fa8d2 55 fCTSPtMin(0), fEMCALPtMin(0),fPHOSPtMin(0), fAODBranchList(new TList ),
591cc579 56 fAODCTS(new TObjArray()), fAODEMCAL(new TObjArray()), fAODPHOS(new TObjArray()),
1c5acb87 57 fEMCALCells(0x0), fPHOSCells(0x0),
477d6cee 58 fInputEvent(0x0), fOutputEvent(0x0),fMC(0x0),
1c5acb87 59 fFillCTS(0),fFillEMCAL(0),fFillPHOS(0),
591cc579 60 fFillEMCALCells(0),fFillPHOSCells(0),
1e68a3f4 61// fSecondInputAODTree(0x0), fSecondInputAODEvent(0x0),
62// fSecondInputFileName(""),fSecondInputFirstEvent(0),
63// fAODCTSNormalInputEntries(0), fAODEMCALNormalInputEntries(0),
64// fAODPHOSNormalInputEntries(0),
65 fTrackStatus(0),
afabc52f 66 fReadStack(kFALSE), fReadAODMCParticles(kFALSE),
1e68a3f4 67 fDeltaAODFileName("deltaAODPartCorr.root"),fFiredTriggerClassName(""),
c8fe2783 68 fAnaLED(kFALSE),fTaskName(""),fCaloUtils(0x0),
f37fa8d2 69 fMixedEvent(NULL), fNMixedEvent(1), fVertex(NULL),
9584c261 70 fWriteOutputDeltaAOD(kFALSE),fOldAOD(kFALSE){
1c5acb87 71 //Ctor
72
73 //Initialize parameters
74 InitParameters();
75}
1c5acb87 76
77//_________________________________
78AliCaloTrackReader::~AliCaloTrackReader() {
79 //Dtor
80
ff45398a 81 if(fFiducialCut) delete fFiducialCut ;
29b2ceec 82
f37fa8d2 83 if(fAODBranchList){
84 fAODBranchList->Delete();
85 delete fAODBranchList ;
86 }
87
1c5acb87 88 if(fAODCTS){
89 fAODCTS->Clear() ;
90 delete fAODCTS ;
91 }
92
93 if(fAODEMCAL){
94 fAODEMCAL->Clear() ;
95 delete fAODEMCAL ;
96 }
97
98 if(fAODPHOS){
99 fAODPHOS->Clear() ;
100 delete fAODPHOS ;
101 }
102
103 if(fEMCALCells){
104 fEMCALCells->Clear() ;
105 delete fEMCALCells ;
106 }
107
108 if(fPHOSCells){
109 fPHOSCells->Clear() ;
110 delete fPHOSCells ;
111 }
eb310db0 112
7e25653f 113 if(fVertex){
114 for (Int_t i = 0; i < fNMixedEvent; i++) {
115 delete [] fVertex[i] ;
116
117 }
118 delete [] fVertex ;
119 }
120
7787a778 121// Pointers not owned, done by the analysis frame
122// if(fInputEvent) delete fInputEvent ;
123// if(fOutputEvent) delete fOutputEvent ;
124// if(fMC) delete fMC ;
591cc579 125
7787a778 126// if(fSecondInputAODTree){
127// fSecondInputAODTree->Clear();
128// delete fSecondInputAODTree;
129// }
130//
131// if(fSecondInputAODEvent) delete fSecondInputAODEvent ;
c1ac3823 132
7787a778 133 // Pointer not owned, deleted by maker
134 //if (fCaloUtils) delete fCaloUtils ;
c1ac3823 135
136}
477d6cee 137
591cc579 138//_________________________________________________________________________
139Bool_t AliCaloTrackReader::ComparePtHardAndJetPt(){
140 // Check the event, if the requested ptHard is much larger than the jet pT, then there is a problem.
141 // Only for PYTHIA.
142 if(!fReadStack) return kTRUE; //Information not filtered to AOD
143
144 if(!strcmp(GetGenEventHeader()->ClassName(), "AliGenPythiaEventHeader")){
898c9d44 145 TParticle * jet = 0;
591cc579 146 AliGenPythiaEventHeader* pygeh= (AliGenPythiaEventHeader*) GetGenEventHeader();
147 Int_t nTriggerJets = pygeh->NTriggerJets();
148 Float_t ptHard = pygeh->GetPtHard();
149
150 //if(fDebug > 1) printf("AliMCAnalysisUtils::PythiaEventHeader: Njets: %d, pT Hard %f\n",nTriggerJets, ptHard);
898c9d44 151 Float_t tmpjet[]={0,0,0,0};
591cc579 152 for(Int_t ijet = 0; ijet< nTriggerJets; ijet++){
153 pygeh->TriggerJet(ijet, tmpjet);
154 jet = new TParticle(94, 21, -1, -1, -1, -1, tmpjet[0],tmpjet[1],tmpjet[2],tmpjet[3], 0,0,0,0);
155 //Compare jet pT and pt Hard
156 //if(fDebug > 1) printf("AliMCAnalysisUtils:: %d pycell jet pT %f\n",ijet, jet->Pt());
157 if(jet->Pt() > fPtHardAndJetPtFactor * ptHard) {
158 printf("AliMCAnalysisUtils::PythiaEventHeader: Njets: %d, pT Hard %2.2f, pycell jet pT %2.2f, rejection factor %1.1f\n",
159 nTriggerJets, ptHard, jet->Pt(), fPtHardAndJetPtFactor);
160 return kFALSE;
161 }
162 }
898c9d44 163 if(jet) delete jet;
591cc579 164 }
898c9d44 165
591cc579 166 return kTRUE ;
167
168}
169
1c5acb87 170//____________________________________________________________________________
171AliStack* AliCaloTrackReader::GetStack() const {
172 //Return pointer to stack
173 if(fMC)
174 return fMC->Stack();
175 else{
477d6cee 176 if(fDebug > 1) printf("AliCaloTrackReader::GetStack() - Stack is not available\n");
1c5acb87 177 return 0x0 ;
178 }
179}
180
181//____________________________________________________________________________
182AliHeader* AliCaloTrackReader::GetHeader() const {
183 //Return pointer to header
184 if(fMC)
185 return fMC->Header();
186 else{
477d6cee 187 printf("AliCaloTrackReader::Header is not available\n");
1c5acb87 188 return 0x0 ;
189 }
190}
191//____________________________________________________________________________
192AliGenEventHeader* AliCaloTrackReader::GetGenEventHeader() const {
193 //Return pointer to Generated event header
194 if(fMC)
195 return fMC->GenEventHeader();
196 else{
477d6cee 197 printf("AliCaloTrackReader::GenEventHeader is not available\n");
1c5acb87 198 return 0x0 ;
199 }
200}
201
591cc579 202//____________________________________________________________________________
203TClonesArray* AliCaloTrackReader::GetAODMCParticles(Int_t input) const {
1e68a3f4 204 //Return list of particles in AOD. Do it for the corresponding input event.
205
c8fe2783 206 TClonesArray * rv = NULL ;
1e68a3f4 207 if(fDataType == kAOD){
208
c8fe2783 209 if(input == 0){
1e68a3f4 210 //Normal input AOD
211 AliAODEvent * evt = dynamic_cast<AliAODEvent*> (fInputEvent) ;
212 if(evt)
213 rv = (TClonesArray*)evt->FindListObject("mcparticles");
214 else
215 printf("AliCaloTrackReader::GetAODMCParticles() - wrong AOD input index? %d, or non existing tree? \n",input);
216
217 } //else if(input == 1 && fSecondInputAODEvent){ //Second input AOD
218// rv = (TClonesArray*) fSecondInputAODEvent->FindListObject("mcparticles");
219// }
220
c8fe2783 221 } else {
1e68a3f4 222 printf("AliCaloTrackReader::GetAODMCParticles() - Input are not AODs\n");
c8fe2783 223 }
1e68a3f4 224
c8fe2783 225 return rv ;
591cc579 226}
227
228//____________________________________________________________________________
229AliAODMCHeader* AliCaloTrackReader::GetAODMCHeader(Int_t input) const {
230 //Return MC header in AOD. Do it for the corresponding input event.
1e68a3f4 231 AliAODMCHeader *mch = NULL;
591cc579 232 if(fDataType == kAOD){
233 //Normal input AOD
1e68a3f4 234 if(input == 0) {
235 mch = (AliAODMCHeader*)((AliAODEvent*)fInputEvent)->FindListObject("mcheader");
236 }
237// //Second input AOD
238// else if(input == 1){
239// mch = (AliAODMCHeader*) fSecondInputAODEvent->FindListObject("mcheader");
240// }
591cc579 241 else {
242 printf("AliCaloTrackReader::GetAODMCHeader() - wrong AOD input index, %d\n",input);
591cc579 243 }
244 }
245 else {
246 printf("AliCaloTrackReader::GetAODMCHeader() - Input are not AODs\n");
591cc579 247 }
1e68a3f4 248
249 return mch;
591cc579 250}
251
252//_______________________________________________________________
253void AliCaloTrackReader::Init()
254{
255 //Init reader. Method to be called in AliAnaPartCorrMaker
256
257 //Get the file with second input events if the filename is given
258 //Get the tree and connect the AODEvent. Only with AODs
f1f0bd84 259
591cc579 260 if(fReadStack && fReadAODMCParticles){
261 printf("AliCaloTrackReader::Init() - Cannot access stack and mcparticles at the same time, change them \n");
262 fReadStack = kFALSE;
263 fReadAODMCParticles = kFALSE;
264 }
265
1e68a3f4 266// if(fSecondInputFileName!=""){
267// if(fDataType == kAOD){
268// TFile * input2 = new TFile(fSecondInputFileName,"read");
269// printf("AliCaloTrackReader::Init() - Second input file opened: %s, size %d \n", input2->GetName(), (Int_t) input2->GetSize());
270// fSecondInputAODTree = (TTree*) input2->Get("aodTree");
271// if(fSecondInputAODTree) printf("AliCaloTrackReader::Init() - Second input tree opened: %s, entries %d \n",
272// fSecondInputAODTree->GetName(), (Int_t) fSecondInputAODTree->GetEntries());
273// else{
274// printf("AliCaloTrackReader::Init() - Second input tree not available, STOP \n");
275// abort();
276// }
277// fSecondInputAODEvent = new AliAODEvent;
278// fSecondInputAODEvent->ReadFromTree(fSecondInputAODTree);
279// if(fSecondInputFirstEvent >= fSecondInputAODTree->GetEntriesFast()){
280// printf("AliCaloTrackReader::Init() - Requested first event of second input %d, is larger than number of events %d, STOP\n",
281// fSecondInputFirstEvent, (Int_t) fSecondInputAODTree->GetEntriesFast());
282// abort();
283// }
284// }
285// else printf("AliCaloTrackReader::Init() - Second input not added, reader is not AOD\n");
286// }
591cc579 287}
765d44e7 288
1c5acb87 289//_______________________________________________________________
290void AliCaloTrackReader::InitParameters()
291{
1c5acb87 292 //Initialize the parameters of the analysis.
591cc579 293 fDataType = kESD ;
1c5acb87 294 fCTSPtMin = 0.2 ;
29b2ceec 295 fEMCALPtMin = 0.2 ;
296 fPHOSPtMin = 0.2 ;
1c5acb87 297
902aa95c 298 //Do not filter the detectors input by default.
299 fFillEMCAL = kFALSE;
300 fFillPHOS = kFALSE;
301 fFillCTS = kFALSE;
1c5acb87 302 fFillEMCALCells = kFALSE;
591cc579 303 fFillPHOSCells = kFALSE;
1c5acb87 304
1e68a3f4 305 //fSecondInputFileName = "" ;
306 //fSecondInputFirstEvent = 0 ;
591cc579 307 fReadStack = kFALSE; // Check in the constructor of the other readers if it was set or in the configuration file
308 fReadAODMCParticles = kFALSE; // Check in the constructor of the other readers if it was set or in the configuration file
42dc8e7d 309 fDeltaAODFileName = "deltaAODPartCorr.root";
72d2488e 310 fFiredTriggerClassName = "";
765d44e7 311
c1ac3823 312 fAnaLED = kFALSE;
0ae57829 313
314 //We want tracks fitted in the detectors:
315 fTrackStatus=AliESDtrack::kTPCrefit;
316 fTrackStatus|=AliESDtrack::kITSrefit;
317
4e2b43d8 318}
1c5acb87 319
320//________________________________________________________________
321void AliCaloTrackReader::Print(const Option_t * opt) const
322{
323
324 //Print some relevant parameters set for the analysis
325 if(! opt)
326 return;
327
328 printf("***** Print: %s %s ******\n", GetName(), GetTitle() ) ;
486258c9 329 printf("Task name : %s\n", fTaskName.Data()) ;
1c5acb87 330 printf("Data type : %d\n", fDataType) ;
331 printf("CTS Min pT : %2.1f GeV/c\n", fCTSPtMin) ;
332 printf("EMCAL Min pT : %2.1f GeV/c\n", fEMCALPtMin) ;
333 printf("PHOS Min pT : %2.1f GeV/c\n", fPHOSPtMin) ;
334 printf("Use CTS = %d\n", fFillCTS) ;
335 printf("Use EMCAL = %d\n", fFillEMCAL) ;
336 printf("Use PHOS = %d\n", fFillPHOS) ;
337 printf("Use EMCAL Cells = %d\n", fFillEMCALCells) ;
338 printf("Use PHOS Cells = %d\n", fFillPHOSCells) ;
591cc579 339 printf("Track status = %d\n", (Int_t) fTrackStatus) ;
f37fa8d2 340 printf("Write delta AOD = %d\n", fWriteOutputDeltaAOD) ;
1e68a3f4 341
591cc579 342 if(fComparePtHardAndJetPt)
343 printf("Compare jet pt and pt hard to accept event, factor = %2.2f",fPtHardAndJetPtFactor);
344
1e68a3f4 345// if(fSecondInputFileName!="") {
346// printf("Second Input File Name = %s\n", fSecondInputFileName.Data()) ;
347// printf("Second Input First Event = %d\n", fSecondInputFirstEvent) ;
348// }
591cc579 349
42dc8e7d 350 printf("Read Kine from, stack? %d, AOD ? %d \n", fReadStack, fReadAODMCParticles) ;
42dc8e7d 351 printf("Delta AOD File Name = %s\n", fDeltaAODFileName.Data()) ;
1c5acb87 352 printf(" \n") ;
353}
354
355//___________________________________________________
29b2ceec 356Bool_t AliCaloTrackReader::FillInputEvent(const Int_t iEntry, const char * currentFileName) {
6639984f 357 //Fill the event counter and input lists that are needed, called by the analysis maker.
902aa95c 358
6639984f 359 fEventNumber = iEntry;
a79a2424 360 fCurrentFileName = TString(currentFileName);
be1f5fa4 361 if(!fInputEvent) {
f7c2338a 362 if(fDebug >= 0) printf("AliCaloTrackReader::FillInputEvent() - Input event not available, skip event analysis\n");
be1f5fa4 363 return kFALSE;
364 }
72d2488e 365 //Select events only fired by a certain trigger configuration if it is provided
be1f5fa4 366 Int_t eventType = 0;
367 if(fInputEvent->GetHeader())
368 eventType = ((AliVHeader*)fInputEvent->GetHeader())->GetEventType();
c1ac3823 369 if( fFiredTriggerClassName !="" && !fAnaLED){
c8fe2783 370 if(eventType!=7)
371 return kFALSE; //Only physics event, do not use for simulated events!!!
7ec23b5a 372 if(fDebug > 0)
373 printf("AliCaloTrackReader::FillInputEvent() - FiredTriggerClass <%s>, selected class <%s>, compare name %d\n",
374 GetFiredTriggerClasses().Data(),fFiredTriggerClassName.Data(), GetFiredTriggerClasses().Contains(fFiredTriggerClassName));
375 if( !GetFiredTriggerClasses().Contains(fFiredTriggerClassName) ) return kFALSE;
72d2488e 376 }
c1ac3823 377 else if(fAnaLED){
378// kStartOfRun = 1, // START_OF_RUN
379// kEndOfRun = 2, // END_OF_RUN
380// kStartOfRunFiles = 3, // START_OF_RUN_FILES
381// kEndOfRunFiles = 4, // END_OF_RUN_FILES
382// kStartOfBurst = 5, // START_OF_BURST
383// kEndOfBurst = 6, // END_OF_BURST
384// kPhysicsEvent = 7, // PHYSICS_EVENT
385// kCalibrationEvent = 8, // CALIBRATION_EVENT
386// kFormatError = 9, // EVENT_FORMAT_ERROR
387// kStartOfData = 10, // START_OF_DATA
388// kEndOfData = 11, // END_OF_DATA
389// kSystemSoftwareTriggerEvent = 12, // SYSTEM_SOFTWARE_TRIGGER_EVENT
390// kDetectorSoftwareTriggerEvent = 13 // DETECTOR_SOFTWARE_TRIGGER_EVENT
391
392 if(eventType!=7 && fDebug > 1 )printf("AliCaloTrackReader::FillInputEvent() - DO LED, Event Type <%d>, 8 Calibration \n", eventType);
393 if(eventType!=8)return kFALSE;
394 }
395
29b2ceec 396 //In case of analysis of events with jets, skip those with jet pt > 5 pt hard
397 if(fComparePtHardAndJetPt && GetStack()) {
7ec23b5a 398 if(!ComparePtHardAndJetPt()) return kFALSE ;
29b2ceec 399 }
591cc579 400
401 //In case of mixing events with other AOD file
1e68a3f4 402 // if(fDataType == kAOD && fSecondInputAODTree){
403//
404// if(fDebug > 1)
405// printf("AliCaloTrackReader::FillInputEvent() - Get event %d from second input AOD file \n", iEntry+fSecondInputFirstEvent);
406// if(fSecondInputAODTree->GetEntriesFast() <= iEntry+fSecondInputFirstEvent) {
407// if(fSecondInputAODTree->GetEntriesFast() == iEntry+fSecondInputFirstEvent)
408// printf("AliCaloTrackReader::FillInputEvent() - Skip events from event %d, no more events in second AOD file \n", iEntry);
409// return kFALSE;
410// }
411//
412// //Get the Event
413// Int_t nbytes = fSecondInputAODTree->GetEvent(iEntry+fSecondInputFirstEvent);
414// if ( nbytes == 0 ) {//If nothing in AOD
415// printf("AliCaloTrackReader::FillInputEvent() - Nothing in Second AOD input, STOP\n");
416// abort() ;
417// }
418//
419// }
09e819c9 420
f8006433 421 //Fill Vertex array
422
423 FillVertexArray();
c8fe2783 424
f8006433 425 //Fill the arrays with cluster/tracks/cells data
426 if(fFillEMCALCells)
c8fe2783 427 FillInputEMCALCells();
428 if(fFillPHOSCells)
429 FillInputPHOSCells();
09e819c9 430
c8fe2783 431 if(fFillCTS)
432 FillInputCTS();
433 if(fFillEMCAL)
434 FillInputEMCAL();
435 if(fFillPHOS)
436 FillInputPHOS();
09e819c9 437
7ec23b5a 438
29b2ceec 439 return kTRUE ;
1c5acb87 440}
441
442//__________________________________________________
443void AliCaloTrackReader::ResetLists() {
444 // Reset lists, called by the analysis maker
445
7e25653f 446 if(fAODCTS) fAODCTS -> Clear();
447 if(fAODEMCAL) fAODEMCAL -> Clear();
448 if(fAODPHOS) fAODPHOS -> Clear();
1c5acb87 449 if(fEMCALCells) fEMCALCells -> Clear();
7e25653f 450 if(fPHOSCells) fPHOSCells -> Clear();
1c5acb87 451}
c8fe2783 452
453//____________________________________________________________________________
454void AliCaloTrackReader::SetInputEvent(AliVEvent* const input)
455{
456 fInputEvent = input;
457 fMixedEvent = dynamic_cast<AliMixedEvent*>(GetInputEvent()) ;
458 if (fMixedEvent) {
459 fNMixedEvent = fMixedEvent->GetNumberOfEvents() ;
460 }
eb310db0 461
462 //Delete previous vertex
463 if(fVertex){
464 for (Int_t i = 0; i < fNMixedEvent; i++) {
465 delete [] fVertex[i] ;
466 }
467 delete [] fVertex ;
468 }
469
c8fe2783 470 fVertex = new Double_t*[fNMixedEvent] ;
471 for (Int_t i = 0; i < fNMixedEvent; i++) {
472 fVertex[i] = new Double_t[3] ;
473 fVertex[i][0] = 0.0 ;
474 fVertex[i][1] = 0.0 ;
475 fVertex[i][2] = 0.0 ;
476 }
477}
478
479//____________________________________________________________________________
f8006433 480void AliCaloTrackReader::GetVertex(Double_t vertex[3]) const {
481 //Return vertex position to be used for single event analysis
482 vertex[0]=fVertex[0][0];
483 vertex[1]=fVertex[0][1];
484 vertex[2]=fVertex[0][2];
485}
486
487//____________________________________________________________________________
488void AliCaloTrackReader::GetVertex(Double_t vertex[3], const Int_t evtIndex) const {
489 //Return vertex position for mixed event, recover the vertex in a particular event.
490
491 //Int_t evtIndex = 0; // for single events only one vertex stored in position 0, default value
492 //if (fMixedEvent && clusterID >=0) {
493 // evtIndex=GetMixedEvent()->EventIndexForCaloCluster(clusterID) ;
494 //}
495
496 vertex[0]=fVertex[evtIndex][0]; vertex[1]=fVertex[evtIndex][1]; vertex[2]=fVertex[evtIndex][2];
497
498}
499//
500
501
502//____________________________________________________________________________
503void AliCaloTrackReader::FillVertexArray() {
504
505 //Fill data member with vertex
506 //In case of Mixed event, multiple vertices
507
508 //Delete previous vertex
509 if(fVertex){
510 for (Int_t i = 0; i < fNMixedEvent; i++) {
511 delete [] fVertex[i] ;
512 }
513 delete [] fVertex ;
514 }
515
516 fVertex = new Double_t*[fNMixedEvent] ;
517 for (Int_t i = 0; i < fNMixedEvent; i++) {
518 fVertex[i] = new Double_t[3] ;
519 fVertex[i][0] = 0.0 ;
520 fVertex[i][1] = 0.0 ;
521 fVertex[i][2] = 0.0 ;
522 }
523
524 if (!fMixedEvent) { //Single event analysis
edffc439 525
526 if(fDataType!=kMC)fInputEvent->GetPrimaryVertex()->GetXYZ(fVertex[0]);
527 else {
528 fVertex[0][0]=0.; fVertex[0][1]=0.; fVertex[0][2]=0.;
529 }
530
f8006433 531 if(fDebug > 1)
532 printf("AliCaloTrackReader::FillVertexArray() - Single Event Vertex : %f,%f,%f\n",fVertex[0][0],fVertex[0][1],fVertex[0][2]);
533
534 } else { // MultiEvent analysis
535 for (Int_t iev = 0; iev < fNMixedEvent; iev++) {
536 fMixedEvent->GetVertexOfEvent(iev)->GetXYZ(fVertex[iev]);
537 if(fDebug > 1)
538 printf("AliCaloTrackReader::FillVertexArray() - Multi Event %d Vertex : %f,%f,%f\n",iev,fVertex[iev][0],fVertex[iev][1],fVertex[iev][2]);
539
540 }
c8fe2783 541 }
f8006433 542
c8fe2783 543}
544
f8006433 545
f37fa8d2 546//____________________________________________________________________________
c8fe2783 547void AliCaloTrackReader::FillInputCTS() {
f37fa8d2 548 //Return array with Central Tracking System (CTS) tracks
549
550 if(fDebug > 2 ) printf("AliCaloTrackReader::FillInputCTS()\n");
551
c8fe2783 552 Int_t nTracks = fInputEvent->GetNumberOfTracks() ;
c8fe2783 553 Double_t p[3];
554 for (Int_t itrack = 0; itrack < nTracks; itrack++) {////////////// track loop
555 AliVTrack * track = (AliVTrack*)fInputEvent->GetTrack(itrack) ; // retrieve track from esd
556
557 //Select tracks under certain conditions, TPCrefit, ITSrefit ... check the set bits
558 if (fTrackStatus && !((track->GetStatus() & fTrackStatus) == fTrackStatus))
559 continue ;
560
561 track->GetPxPyPz(p) ;
562 TLorentzVector momentum(p[0],p[1],p[2],0);
563
564 if(fCTSPtMin < momentum.Pt()){
565
566 if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(momentum,"CTS"))
567 continue;
568
569 if(fDebug > 2 && momentum.Pt() > 0.1)
f37fa8d2 570 printf("AliCaloTrackReader::FillInputCTS() - Selected tracks E %3.2f, pt %3.2f, phi %3.2f, eta %3.2f\n",
c8fe2783 571 momentum.E(),momentum.Pt(),momentum.Phi()*TMath::RadToDeg(),momentum.Eta());
572
573 if (fMixedEvent) {
574 track->SetID(itrack);
575 }
f37fa8d2 576
577 fAODCTS->Add(track);
578
c8fe2783 579 }//Pt and Fiducial cut passed.
580 }// track loop
581
1e68a3f4 582 //fAODCTSNormalInputEntries = fAODCTS->GetEntriesFast();
f37fa8d2 583 if(fDebug > 1) printf("AliCaloTrackReader::FillInputCTS() - aod entries %d\n", fAODCTS->GetEntriesFast());//fAODCTSNormalInputEntries);
c8fe2783 584
585 // //If second input event available, add the clusters.
586 // if(fSecondInputAODTree && fSecondInputAODEvent){
587 // nTracks = fSecondInputAODEvent->GetNumberOfTracks() ;
f37fa8d2 588 // if(fDebug > 1) printf("AliCaloTrackReader::FillInputCTS() - Add second input tracks, entries %d\n", nTracks) ;
c8fe2783 589 // for (Int_t itrack = 0; itrack < nTracks; itrack++) {////////////// track loop
590 // AliAODTrack * track = ((AliAODEvent*)fSecondInputAODEvent)->GetTrack(itrack) ; // retrieve track from esd
591 //
592 // //Select tracks under certain conditions, TPCrefit, ITSrefit ... check the set bits
593 // if (fTrackStatus && !((track->GetStatus() & fTrackStatus) == fTrackStatus)) continue ;
594 //
595 // track->GetPxPyPz(p) ;
596 // TLorentzVector momentum(p[0],p[1],p[2],0);
597 //
598 // if(fCTSPtMin < momentum.Pt()){
599 //
600 // if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(momentum,"CTS")) continue;
601 //
f37fa8d2 602 // if(fDebug > 2 && momentum.Pt() > 0.1) printf("AliCaloTrackReader::FillInputCTS() - Selected tracks E %3.2f, pt %3.2f, phi %3.2f, eta %3.2f\n",
c8fe2783 603 // momentum.E(),momentum.Pt(),momentum.Phi()*TMath::RadToDeg(),momentum.Eta());
604 //
f37fa8d2 605 // fAODCTS->Add(track);
c8fe2783 606 //
607 // }//Pt and Fiducial cut passed.
608 // }// track loop
609 //
f37fa8d2 610 // if(fDebug > 1) printf("AliCaloTrackReader::FillInputCTS() - aod normal entries %d, after second input %d\n", fAODCTSNormalInputEntries, fAODCTS->GetEntriesFast());
c8fe2783 611 // } //second input loop
612 //
613}
614
615 //____________________________________________________________________________
616void AliCaloTrackReader::FillInputEMCAL() {
f37fa8d2 617 //Return array with EMCAL clusters in aod format
c8fe2783 618
f37fa8d2 619 if(fDebug > 2 ) printf("AliCaloTrackReader::FillInputEMCAL()\n");
620
621 //Loop to select clusters in fiducial cut and fill container with aodClusters
c8fe2783 622 Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
623 for (Int_t iclus = 0; iclus < nclusters; iclus++) {
624 AliVCluster * clus = 0;
625 if ( (clus = fInputEvent->GetCaloCluster(iclus)) ) {
f37fa8d2 626 if (IsEMCALCluster(clus)){
c8fe2783 627
1e68a3f4 628 //Check if the cluster contains any bad channel and if close to calorimeter borders
c8fe2783 629 Int_t vindex = 0 ;
630 if (fMixedEvent)
631 vindex = fMixedEvent->EventIndexForCaloCluster(iclus);
632
633 if(GetCaloUtils()->ClusterContainsBadChannel("EMCAL",clus->GetCellsAbsId(), clus->GetNCells()))
634 continue;
635 if(!GetCaloUtils()->CheckCellFiducialRegion(clus, (AliVCaloCells*)fInputEvent->GetEMCALCells(), fInputEvent, vindex))
636 continue;
637
638 TLorentzVector momentum ;
639
640 clus->GetMomentum(momentum, fVertex[vindex]);
641
642 if(fEMCALPtMin < momentum.Pt()){
643
644 if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(momentum,"EMCAL"))
645 continue;
646
647 if(fDebug > 2 && momentum.E() > 0.1)
f37fa8d2 648 printf("AliCaloTrackReader::FillInputEMCAL() - Selected clusters E %3.2f, pt %3.2f, phi %3.2f, eta %3.2f\n",
c8fe2783 649 momentum.E(),momentum.Pt(),momentum.Phi()*TMath::RadToDeg(),momentum.Eta());
650
9584c261 651 Float_t pos[3];
652 clus->GetPosition(pos);
653 //printf("Before Corrections: e %f, x %f, y %f, z %f\n",clus->E(),pos[0],pos[1],pos[2]);
654
1e68a3f4 655 //Recalibrate the cluster energy
c8fe2783 656 if(GetCaloUtils()->IsRecalibrationOn()) {
9584c261 657 Float_t energy = GetCaloUtils()->RecalibrateClusterEnergy(clus, GetEMCALCells());
c8fe2783 658 clus->SetE(energy);
9584c261 659 //printf("Recalibrated Energy %f\n",clus->E());
660 }
661
662 //Correct non linearity
663 if(GetCaloUtils()->IsCorrectionOfClusterEnergyOn()){
664 GetCaloUtils()->CorrectClusterEnergy(clus) ;
665 //printf("Linearity Corrected Energy %f\n",clus->E());
666 }
667
668 //Recalculate cluster position
669 if(GetCaloUtils()->IsRecalculationOfClusterPositionOn()){
670 GetCaloUtils()->RecalculateClusterPosition(GetEMCALCells(),clus);
671 clus->GetPosition(pos);
672 //printf("After Corrections: e %f, x %f, y %f, z %f\n",clus->E(),pos[0],pos[1],pos[2]);
c8fe2783 673 }
674
675 if (fMixedEvent) {
676 clus->SetID(iclus) ;
677 }
f37fa8d2 678
679 fAODEMCAL->Add(clus);
680
c8fe2783 681 }//Pt and Fiducial cut passed.
682 }//EMCAL cluster
683 }// cluster exists
684 }// cluster loop
1e68a3f4 685 //fAODEMCALNormalInputEntries = fAODEMCAL->GetEntriesFast();
f37fa8d2 686 if(fDebug > 1) printf("AliCaloTrackReader::FillInputEMCAL() - aod entries %d\n", fAODEMCAL->GetEntriesFast());//fAODEMCALNormalInputEntries);
c8fe2783 687
688 //If second input event available, add the clusters.
689 // if(fSecondInputAODTree && fSecondInputAODEvent){
690 // GetSecondInputAODVertex(v);
691 // nclusters = ((AliAODEvent*)fSecondInputAODEvent)->GetNumberOfCaloClusters();
f37fa8d2 692 // if(fDebug > 1) printf("AliCaloTrackReader::FillInputEMCAL() - Add second input clusters, entries %d\n", nclusters) ;
c8fe2783 693 // for (Int_t iclus = 0; iclus < nclusters; iclus++) {
694 // AliAODCaloCluster * clus = 0;
695 // if ( (clus = ((AliAODEvent*)fSecondInputAODEvent)->GetCaloCluster(iclus)) ) {
696 // if (clus->IsEMCAL()){
697 // TLorentzVector momentum ;
698 // clus->GetMomentum(momentum, v);
699 //
700 // if(fEMCALPtMin < momentum.Pt()){
701 //
702 // if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(momentum,"EMCAL")) continue;
703 //
f37fa8d2 704 // if(fDebug > 2 && momentum.E() > 0.1) printf("AliCaloTrackReader::FillInputEMCAL() - Selected clusters E %3.2f, pt %3.2f, phi %3.2f, eta %3.2f\n",
c8fe2783 705 // momentum.E(),momentum.Pt(),momentum.Phi()*TMath::RadToDeg(),momentum.Eta());
f37fa8d2 706 // fAODEMCAL->Add(clus);
c8fe2783 707 // }//Pt and Fiducial cut passed.
708 // }//EMCAL cluster
709 // }// cluster exists
710 // }// cluster loop
711 //
f37fa8d2 712 // if(fDebug > 1) printf("AliCaloTrackReader::FillInputEMCAL() - aod normal entries %d, after second input %d\n", fAODEMCALNormalInputEntries, fAODEMCAL->GetEntriesFast());
c8fe2783 713 //
714 // } //second input loop
715}
716
717 //____________________________________________________________________________
718void AliCaloTrackReader::FillInputPHOS() {
f37fa8d2 719 //Return array with PHOS clusters in aod format
c8fe2783 720
f37fa8d2 721 if(fDebug > 2 ) printf("AliCaloTrackReader::FillInputPHOS()\n");
722
723 //Loop to select clusters in fiducial cut and fill container with aodClusters
c8fe2783 724 Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
725 for (Int_t iclus = 0; iclus < nclusters; iclus++) {
726 AliVCluster * clus = 0;
727 if ( (clus = fInputEvent->GetCaloCluster(iclus)) ) {
f37fa8d2 728 if (IsPHOSCluster(clus)){
729 //Check if the cluster contains any bad channel and if close to calorimeter borders
c8fe2783 730 Int_t vindex = 0 ;
731 if (fMixedEvent)
732 vindex = fMixedEvent->EventIndexForCaloCluster(iclus);
733 if( GetCaloUtils()->ClusterContainsBadChannel("PHOS",clus->GetCellsAbsId(), clus->GetNCells()))
734 continue;
735 if(!GetCaloUtils()->CheckCellFiducialRegion(clus, fInputEvent->GetPHOSCells(), fInputEvent, vindex))
736 continue;
737
738 TLorentzVector momentum ;
739
740 clus->GetMomentum(momentum, fVertex[vindex]);
741
742 if(fPHOSPtMin < momentum.Pt()){
743
744 if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(momentum,"PHOS"))
745 continue;
746
747 if(fDebug > 2 && momentum.E() > 0.1)
f37fa8d2 748 printf("AliCaloTrackReader::FillInputPHOS() - Selected clusters E %3.2f, pt %3.2f, phi %3.2f, eta %3.2f\n",
c8fe2783 749 momentum.E(),momentum.Pt(),momentum.Phi()*TMath::RadToDeg(),momentum.Eta());
750
751 //Recalibrate the cluster energy
752 if(GetCaloUtils()->IsRecalibrationOn()) {
753 Float_t energy = GetCaloUtils()->RecalibrateClusterEnergy(clus, (AliAODCaloCells*)GetPHOSCells());
754 clus->SetE(energy);
755 }
756
757 if (fMixedEvent) {
758 clus->SetID(iclus) ;
759 }
760
f37fa8d2 761 fAODPHOS->Add(clus);
762
c8fe2783 763 }//Pt and Fiducial cut passed.
764 }//PHOS cluster
765 }//cluster exists
766 }//esd cluster loop
767
1e68a3f4 768 //fAODPHOSNormalInputEntries = fAODPHOS->GetEntriesFast() ;
f37fa8d2 769 if(fDebug > 1) printf("AliCaloTrackReader::FillInputPHOS() - aod entries %d\n", fAODPHOS->GetEntriesFast());//fAODPHOSNormalInputEntries);
c8fe2783 770
771 //If second input event available, add the clusters.
772 // if(fSecondInputAODTree && fSecondInputAODEvent){
773 // GetSecondInputAODVertex(v);
774 // nclusters = ((AliAODEvent*)fSecondInputAODEvent)->GetNumberOfCaloClusters();
f37fa8d2 775 // if(fDebug > 1) printf("AliCaloTrackReader::FillInputPHOS() - Add second input clusters, entries %d\n", nclusters);
c8fe2783 776 // for (Int_t iclus = 0; iclus < nclusters; iclus++) {
777 // AliAODCaloCluster * clus = 0;
778 // if ( (clus = ((AliAODEvent*)fSecondInputAODEvent)->GetCaloCluster(iclus)) ) {
779 // if (clus->IsPHOS()){
780 // TLorentzVector momentum ;
781 // clus->GetMomentum(momentum, v);
782 //
783 // if(fPHOSPtMin < momentum.Pt()){
784 //
785 // if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(momentum,"PHOS")) continue;
786 //
f37fa8d2 787 // if(fDebug > 2 && momentum.E() > 0.1) printf("AliCaloTrackReader::FillInputPHOS() - Selected clusters E %3.2f, pt %3.2f, phi %3.2f, eta %3.2f\n",
c8fe2783 788 // momentum.E(),momentum.Pt(),momentum.Phi()*TMath::RadToDeg(),momentum.Eta());
f37fa8d2 789 // fAODPHOS->Add(clus);
c8fe2783 790 // }//Pt and Fiducial cut passed.
791 // }//PHOS cluster
792 // }// cluster exists
793 // }// cluster loop
f37fa8d2 794 // if(fDebug > 1) printf("AliCaloTrackReader::FillInputPHOS() - aod normal entries %d, after second input %d\n", fAODPHOSNormalInputEntries, fAODPHOS->GetEntriesFast());
c8fe2783 795 // } //second input loop
796
797}
798
f37fa8d2 799//____________________________________________________________________________
c8fe2783 800void AliCaloTrackReader::FillInputEMCALCells() {
801 //Return array with EMCAL cells in aod format
802
803 fEMCALCells = fInputEvent->GetEMCALCells();
804
805}
806
f37fa8d2 807//____________________________________________________________________________
c8fe2783 808void AliCaloTrackReader::FillInputPHOSCells() {
809 //Return array with PHOS cells in aod format
810
811 fPHOSCells = fInputEvent->GetPHOSCells();
812
813}
814
f37fa8d2 815//____________________________________________________________________________
816Bool_t AliCaloTrackReader::IsEMCALCluster(AliVCluster* cluster) const {
817 // Check if it is a cluster from EMCAL. For old AODs cluster type has
818 // different number and need to patch here
819
820 if(fDataType==kAOD && fOldAOD)
821 {
822 if (cluster->GetType() == 2) return kTRUE;
823 else return kFALSE;
824 }
825 else
826 {
827 return cluster->IsEMCAL();
828 }
829
830}
831
832//____________________________________________________________________________
833Bool_t AliCaloTrackReader::IsPHOSCluster(AliVCluster * cluster) const {
834 //Check if it is a cluster from PHOS.For old AODs cluster type has
835 // different number and need to patch here
836
837 if(fDataType==kAOD && fOldAOD)
838 {
839 Int_t type = cluster->GetType();
840 if (type == 0 || type == 1) return kTRUE;
841 else return kFALSE;
842 }
843 else
844 {
845 return cluster->IsPHOS();
846 }
847
848}
849