AliPhysics  vAN-20150328 (7755d27)
 All Classes Namespaces Files Functions Variables Enumerations Enumerator Macros
AliCaloTrackReader.cxx
Go to the documentation of this file.
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 
16 // --- ROOT system ---
17 #include <TFile.h>
18 #include <TGeoManager.h>
19 #include <TStreamerInfo.h>
20 
21 // ---- ANALYSIS system ----
22 #include "AliMCEvent.h"
23 #include "AliAODMCHeader.h"
24 #include "AliGenPythiaEventHeader.h"
25 #include "AliGenCocktailEventHeader.h"
26 #include "AliGenHijingEventHeader.h"
27 #include "AliESDEvent.h"
28 #include "AliAODEvent.h"
29 #include "AliVTrack.h"
30 #include "AliVParticle.h"
31 #include "AliMixedEvent.h"
32 //#include "AliTriggerAnalysis.h"
33 #include "AliESDVZERO.h"
34 #include "AliVCaloCells.h"
35 #include "AliAnalysisManager.h"
36 #include "AliInputEventHandler.h"
37 #include "AliAODMCParticle.h"
38 #include "AliStack.h"
39 #include "AliLog.h"
40 
41 // ---- Detectors ----
42 #include "AliPHOSGeoUtils.h"
43 #include "AliEMCALGeometry.h"
44 #include "AliEMCALRecoUtils.h"
45 
46 // ---- CaloTrackCorr ---
47 #include "AliCalorimeterUtils.h"
48 #include "AliCaloTrackReader.h"
49 // ---- Jets ----
50 #include "AliAODJet.h"
51 #include "AliAODJetEventBackground.h"
52 
56 
57 //________________________________________
59 //________________________________________
61 TObject(), fEventNumber(-1), //fCurrentFileName(""),
62 fDataType(0), fDebug(0),
63 fFiducialCut(0x0), fCheckFidCut(kFALSE),
64 fComparePtHardAndJetPt(0), fPtHardAndJetPtFactor(0),
65 fComparePtHardAndClusterPt(0),fPtHardAndClusterPtFactor(0),
66 fCTSPtMin(0), fEMCALPtMin(0), fPHOSPtMin(0),
67 fCTSPtMax(0), fEMCALPtMax(0), fPHOSPtMax(0),
68 fUseEMCALTimeCut(1), fUseParamTimeCut(0),
69 fUseTrackTimeCut(0), fAccessTrackTOF(0),
70 fEMCALTimeCutMin(-10000), fEMCALTimeCutMax(10000),
71 fEMCALParamTimeCutMin(), fEMCALParamTimeCutMax(),
72 fTrackTimeCutMin(-10000), fTrackTimeCutMax(10000),
73 fUseTrackDCACut(0),
74 fAODBranchList(0x0),
75 fCTSTracks(0x0), fEMCALClusters(0x0),
76 fDCALClusters(0x0), fPHOSClusters(0x0),
77 fEMCALCells(0x0), fPHOSCells(0x0),
78 fInputEvent(0x0), fOutputEvent(0x0),fMC(0x0),
79 fFillCTS(0), fFillEMCAL(0),
80 fFillDCAL(0), fFillPHOS(0),
81 fFillEMCALCells(0), fFillPHOSCells(0),
82 fRecalculateClusters(kFALSE),fCorrectELinearity(kTRUE),
83 fSelectEmbeddedClusters(kFALSE),
84 fSmearShowerShape(0), fSmearShowerShapeWidth(0), fRandom(),
85 fTrackStatus(0), fSelectSPDHitTracks(0),
86 fTrackMult(0), fTrackMultEtaCut(0.9),
87 fReadStack(kFALSE), fReadAODMCParticles(kFALSE),
88 fDeltaAODFileName(""), fFiredTriggerClassName(""),
89 
90 fEventTriggerMask(0), fMixEventTriggerMask(0), fEventTriggerAtSE(0),
91 fEventTrigMinBias(0), fEventTrigCentral(0),
92 fEventTrigSemiCentral(0), fEventTrigEMCALL0(0),
93 fEventTrigEMCALL1Gamma1(0), fEventTrigEMCALL1Gamma2(0),
94 fEventTrigEMCALL1Jet1(0), fEventTrigEMCALL1Jet2(0),
95 fBitEGA(0), fBitEJE(0),
96 
97 fAnaLED(kFALSE),
98 fTaskName(""), fCaloUtils(0x0),
99 fMixedEvent(NULL), fNMixedEvent(0), fVertex(NULL),
100 fListMixedTracksEvents(), fListMixedCaloEvents(),
101 fLastMixedTracksEvent(-1), fLastMixedCaloEvent(-1),
102 fWriteOutputDeltaAOD(kFALSE),
103 fEMCALClustersListName(""), fZvtxCut(0.),
104 fAcceptFastCluster(kFALSE), fRemoveLEDEvents(kTRUE),
105 //Trigger rejection
106 fRemoveBadTriggerEvents(0), fTriggerPatchClusterMatch(1),
107 fTriggerPatchTimeWindow(), fTriggerL0EventThreshold(0),
108 fTriggerL1EventThreshold(0), fTriggerL1EventThresholdFix(0),
109 fTriggerClusterBC(0), fTriggerClusterIndex(0), fTriggerClusterId(0),
110 fIsExoticEvent(0), fIsBadCellEvent(0), fIsBadMaxCellEvent(0),
111 fIsTriggerMatch(0), fIsTriggerMatchOpenCut(),
112 fTriggerClusterTimeRecal(kTRUE), fRemoveUnMatchedTriggers(kTRUE),
113 fDoPileUpEventRejection(kFALSE), fDoV0ANDEventSelection(kFALSE),
114 fDoVertexBCEventSelection(kFALSE),
115 fDoRejectNoTrackEvents(kFALSE),
116 fUseEventsWithPrimaryVertex(kFALSE),
117 //fTriggerAnalysis (0x0),
118 fTimeStampEventSelect(0),
119 fTimeStampEventFracMin(0), fTimeStampEventFracMax(0),
120 fTimeStampRunMin(0), fTimeStampRunMax(0),
121 fNPileUpClusters(-1), fNNonPileUpClusters(-1), fNPileUpClustersCut(3),
122 fVertexBC(-200), fRecalculateVertexBC(0),
123 fCentralityClass(""), fCentralityOpt(0),
124 fEventPlaneMethod(""),
125 fAcceptOnlyHIJINGLabels(0), fNMCProducedMin(0), fNMCProducedMax(0),
126 fFillInputNonStandardJetBranch(kFALSE),
127 fNonStandardJets(new TClonesArray("AliAODJet",100)),fInputNonStandardJetBranchName("jets"),
128 fFillInputBackgroundJetBranch(kFALSE),
129 fBackgroundJets(0x0),fInputBackgroundJetBranchName("jets"),
130 fAcceptEventsWithBit(0), fRejectEventsWithBit(0), fRejectEMCalTriggerEventsWith2Tresholds(0),
131 fMomentum()
132 {
133  InitParameters();
134 }
135 
136 //_______________________________________
138 //_______________________________________
140 {
141  DeletePointers();
142 }
143 
144 //_______________________________________
147 //_______________________________________
149 {
150  delete fFiducialCut ;
151 
152  if(fAODBranchList)
153  {
154  fAODBranchList->Delete();
155  delete fAODBranchList ;
156  }
157 
158  if(fCTSTracks)
159  {
160  if(fDataType!=kMC)fCTSTracks->Clear() ;
161  else fCTSTracks->Delete() ;
162  delete fCTSTracks ;
163  }
164 
165  if(fEMCALClusters)
166  {
167  if(fDataType!=kMC)fEMCALClusters->Clear("C") ;
168  else fEMCALClusters->Delete() ;
169  delete fEMCALClusters ;
170  }
171 
172  if(fDCALClusters)
173  {
174  if(fDataType!=kMC)fDCALClusters->Clear("C") ;
175  else fDCALClusters->Delete() ;
176  delete fDCALClusters ;
177  }
178 
179  if(fPHOSClusters)
180  {
181  if(fDataType!=kMC)fPHOSClusters->Clear("C") ;
182  else fPHOSClusters->Delete() ;
183  delete fPHOSClusters ;
184  }
185 
186  if(fVertex)
187  {
188  for (Int_t i = 0; i < fNMixedEvent; i++)
189  {
190  delete [] fVertex[i] ;
191 
192  }
193  delete [] fVertex ;
194  }
195 
196  //delete fTriggerAnalysis;
197 
198  if(fNonStandardJets)
199  {
200  if(fDataType!=kMC) fNonStandardJets->Clear("C") ;
201  else fNonStandardJets->Delete() ;
202  delete fNonStandardJets ;
203  }
204  delete fBackgroundJets ;
205 
206  fRejectEventsWithBit.Reset();
207  fAcceptEventsWithBit.Reset();
208 
209  // Pointers not owned, done by the analysis frame
210  // if(fInputEvent) delete fInputEvent ;
211  // if(fOutputEvent) delete fOutputEvent ;
212  // if(fMC) delete fMC ;
213  // Pointer not owned, deleted by maker
214  // if (fCaloUtils) delete fCaloUtils ;
215 }
216 
217 //____________________________________________________________
221 //____________________________________________________________
222 Bool_t AliCaloTrackReader::AcceptDCA(Float_t pt, Float_t dca)
223 {
224  Float_t cut = fTrackDCACut[0]+fTrackDCACut[1]/TMath::Power(pt,fTrackDCACut[2]);
225 
226  if(TMath::Abs(dca) < cut)
227  return kTRUE;
228  else
229  return kFALSE;
230 }
231 
232 //_____________________________________________________
235 //_____________________________________________________
237 {
238  Int_t nAccept = fAcceptEventsWithBit.GetSize();
239 
240  //printf("N accept %d\n", nAccept);
241 
242  if( nAccept <= 0 )
243  return kTRUE ; // accept the event
244 
245  UInt_t trigFired = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
246 
247  for(Int_t ibit = 0; ibit < nAccept; ibit++)
248  {
249  Bool_t accept = (trigFired & fAcceptEventsWithBit.At(ibit));
250 
251  //printf("accept %d, ibit %d, bit %d \n",accept, ibit,fAcceptEventsWithBit.At(ibit));
252  if(accept) return kTRUE ; // accept the event
253  }
254 
255  return kFALSE ; // reject the event
256 }
257 
258 //_____________________________________________________
261 //_____________________________________________________
263 {
264  Int_t nReject = fRejectEventsWithBit.GetSize();
265 
266  //printf("N reject %d\n", nReject);
267 
268  if( nReject <= 0 )
269  return kTRUE ; // accept the event
270 
271  UInt_t trigFired = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
272 
273  for(Int_t ibit = 0; ibit < nReject; ibit++)
274  {
275  Bool_t reject = (trigFired & fRejectEventsWithBit.At(ibit));
276 
277  //printf("reject %d, ibit %d, bit %d \n",reject, ibit,fRejectEventsWithBit.At(ibit));
278  if(reject) return kFALSE ; // reject the event
279  }
280 
281  return kTRUE ; // accept the event
282 }
283 
284 //_____________________________________________
288 //_____________________________________________
290 {
291  //-----------------------------------------------------------
292  // Reject events depending on the trigger name and event type
293  //-----------------------------------------------------------
294 
295  Int_t eventType = 0;
296  if(fInputEvent->GetHeader())
297  eventType = ((AliVHeader*)fInputEvent->GetHeader())->GetEventType();
298 
299  if( fFiredTriggerClassName !="" && !fAnaLED)
300  {
301  //printf("Event type %d\n",eventType);
302  if(eventType!=7)
303  return kFALSE; //Only physics event, do not use for simulated events!!!
304 
305  AliDebug(1,Form("FiredTriggerClass <%s>, selected class <%s>, compare name %d",
308 
309  if( !GetFiredTriggerClasses().Contains(fFiredTriggerClassName) ) return kFALSE;
310  else AliDebug(1,"Accepted triggered event");
311  }
312  else if(fAnaLED)
313  {
314  // kStartOfRun = 1, // START_OF_RUN
315  // kEndOfRun = 2, // END_OF_RUN
316  // kStartOfRunFiles = 3, // START_OF_RUN_FILES
317  // kEndOfRunFiles = 4, // END_OF_RUN_FILES
318  // kStartOfBurst = 5, // START_OF_BURST
319  // kEndOfBurst = 6, // END_OF_BURST
320  // kPhysicsEvent = 7, // PHYSICS_EVENT
321  // kCalibrationEvent = 8, // CALIBRATION_EVENT
322  // kFormatError = 9, // EVENT_FORMAT_ERROR
323  // kStartOfData = 10, // START_OF_DATA
324  // kEndOfData = 11, // END_OF_DATA
325  // kSystemSoftwareTriggerEvent = 12, // SYSTEM_SOFTWARE_TRIGGER_EVENT
326  // kDetectorSoftwareTriggerEvent = 13 // DETECTOR_SOFTWARE_TRIGGER_EVENT
327 
328  //if(eventType!=7 && fDebug > 1 )printf("AliCaloTrackReader::CheckEventTriggers() - DO LED, Event Type <%d>, 8 Calibration \n", eventType);
329  if(eventType!=8) return kFALSE;
330  }
331 
332  AliDebug(1,"Pass Trigger name rejection");
333 
334  //-----------------------------------------------------------------
335  // In case of mixing analysis, select here the trigger of the event
336  //-----------------------------------------------------------------
337  UInt_t isTrigger = kFALSE;
338  UInt_t isMB = kFALSE;
339  if(!fEventTriggerAtSE)
340  {
341  // In case of mixing analysis, accept MB events, not only Trigger
342  // Track and cluster arrays filled for MB in order to create the pool in the corresponding analysis
343  // via de method in the base class FillMixedEventPool()
344 
345  AliAnalysisManager *manager = AliAnalysisManager::GetAnalysisManager();
346  AliInputEventHandler *inputHandler = dynamic_cast<AliInputEventHandler*>(manager->GetInputEventHandler());
347 
348  if(!inputHandler) return kFALSE ; // to content coverity
349 
350  isTrigger = inputHandler->IsEventSelected() & fEventTriggerMask;
351  isMB = inputHandler->IsEventSelected() & fMixEventTriggerMask;
352 
353  if(!isTrigger && !isMB) return kFALSE;
354 
355  //printf("Selected triggered event : %s\n",GetFiredTriggerClasses().Data());
356  AliDebug(0,"Pass uninteresting triggered events rejection in case of mixing analysis");
357  }
358 
359  //-------------------------------------------------------------------------------------
360  // Reject or accept events depending on the trigger bit
361  //-------------------------------------------------------------------------------------
362 
363  Bool_t okA = AcceptEventWithTriggerBit();
364  Bool_t okR = RejectEventWithTriggerBit();
365 
366  //printf("AliCaloTrackReader::FillInputEvent() - Accept event? %d, Reject event %d? \n",okA,okR);
367 
368  if(!okA || !okR) return kFALSE;
369 
370  AliDebug(1,"Pass event bit rejection");
371 
372  //----------------------------------------------------------------------
373  // Do not count events that were likely triggered by an exotic cluster
374  // or out BC cluster
375  //----------------------------------------------------------------------
376 
377  // Set a bit with the event kind, MB, L0, L1 ...
379 
380  // In case of Mixing, avoid checking the triggers in the min bias events
381  if(!fEventTriggerAtSE && (isMB && !isTrigger)) return kTRUE;
382 
384  {
386  {
387  // Reject triggered events when there is coincidence on both EMCal trigger thresholds,
388  // but the requested trigger is the low trigger threshold
389  if(IsEventEMCALL1Jet1 () && IsEventEMCALL1Jet2 () && fFiredTriggerClassName.Contains("EJ2")) return kFALSE;
390  if(IsEventEMCALL1Gamma1() && IsEventEMCALL1Gamma2() && fFiredTriggerClassName.Contains("EG2")) return kFALSE;
391  }
392 
393  //Get Patches that triggered
395 
396  MatchTriggerCluster(patches);
397 
398  patches.Reset();
399 
400  // If requested, remove badly triggeed events, but only when the EMCal trigger bit is set
402  {
403  AliDebug(1,Form("ACCEPT triggered event? \n exotic? %d - bad cell %d - bad Max cell %d - BC %d - Matched %d\n",
405  if (fIsExoticEvent) return kFALSE;
406  else if(fIsBadCellEvent) return kFALSE;
407  else if(fRemoveUnMatchedTriggers && !fIsTriggerMatch) return kFALSE ;
408  else if(fTriggerClusterBC != 0) return kFALSE;
409  AliDebug(1,Form("\t *** YES for %s",GetFiredTriggerClasses().Data()));
410  }
411 
412  AliDebug(1,"Pass EMCal triggered event rejection \n");
413  }
414 
415  //-------------------------------------------------------------------------------------
416  //Select events only fired by a certain trigger configuration if it is provided
417  //-------------------------------------------------------------------------------------
418  if (GetFiredTriggerClasses().Contains("FAST") && !GetFiredTriggerClasses().Contains("ALL") && !fAcceptFastCluster)
419  {
420  AliDebug(1,Form("Do not count events from fast cluster, trigger name %s\n",fFiredTriggerClassName.Data()));
421  return kFALSE;
422  }
423 
424  //-------------------------------------------------------------------------------------
425  // Reject event if large clusters with large energy
426  // Use only for LHC11a data for the moment, and if input is clusterizer V1 or V1+unfolding
427  // If clusterzer NxN or V2 it does not help
428  //-------------------------------------------------------------------------------------
429  Int_t run = fInputEvent->GetRunNumber();
430  if( fRemoveLEDEvents && run > 146857 && run < 146861 )
431  {
432  Bool_t reject = RejectLEDEvents();
433 
434  if(reject) return kFALSE;
435 
436  AliDebug(1,"Pass LED event rejection");
437 
438  }// Remove LED events
439 
440  return kTRUE;
441 }
442 
443 //________________________________________________
448 //________________________________________________
450 {
451  //printf("AliCaloTrackReader::ComparePtHardAndJetPt() - GenHeaderName : %s\n",GetGenEventHeader()->ClassName());
452 
453  if(!strcmp(GetGenEventHeader()->ClassName(), "AliGenPythiaEventHeader"))
454  {
455  TParticle * jet = 0;
456  AliGenPythiaEventHeader* pygeh= (AliGenPythiaEventHeader*) GetGenEventHeader();
457  Int_t nTriggerJets = pygeh->NTriggerJets();
458  Float_t ptHard = pygeh->GetPtHard();
459 
460  AliDebug(1,Form("Njets: %d, pT Hard %f",nTriggerJets, ptHard));
461 
462  Float_t tmpjet[]={0,0,0,0};
463  for(Int_t ijet = 0; ijet< nTriggerJets; ijet++)
464  {
465  pygeh->TriggerJet(ijet, tmpjet);
466  jet = new TParticle(94, 21, -1, -1, -1, -1, tmpjet[0],tmpjet[1],tmpjet[2],tmpjet[3], 0,0,0,0);
467 
468  AliDebug(1,Form("jet %d; pycell jet pT %f",ijet, jet->Pt()));
469 
470  //Compare jet pT and pt Hard
471  if(jet->Pt() > fPtHardAndJetPtFactor * ptHard)
472  {
473  AliInfo(Form("Reject jet event with : pT Hard %2.2f, pycell jet pT %2.2f, rejection factor %1.1f\n",
474  ptHard, jet->Pt(), fPtHardAndJetPtFactor));
475  return kFALSE;
476  }
477  }
478 
479  if(jet) delete jet;
480  }
481 
482  return kTRUE ;
483 }
484 
485 //____________________________________________________
490 //____________________________________________________
492 {
493  if(!strcmp(GetGenEventHeader()->ClassName(), "AliGenPythiaEventHeader"))
494  {
495  AliGenPythiaEventHeader* pygeh= (AliGenPythiaEventHeader*) GetGenEventHeader();
496  Float_t ptHard = pygeh->GetPtHard();
497 
498  Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
499  for (Int_t iclus = 0; iclus < nclusters; iclus++)
500  {
501  AliVCluster * clus = fInputEvent->GetCaloCluster(iclus) ;
502  Float_t ecluster = clus->E();
503 
504  if(ecluster > fPtHardAndClusterPtFactor*ptHard)
505  {
506  AliInfo(Form("Reject : ecluster %2.2f, calo %d, factor %2.2f, ptHard %f",ecluster,clus->GetType(),fPtHardAndClusterPtFactor,ptHard));
507 
508  return kFALSE;
509  }
510  }
511 
512  }
513 
514  return kTRUE ;
515 }
516 
517 //____________________________________________
519 //____________________________________________
521 {
522  if(fMC)
523  return fMC->Stack();
524  else
525  {
526  AliDebug(1,"Stack is not available");
527  return 0x0 ;
528  }
529 }
530 
531 //______________________________________________
533 //______________________________________________
535 {
536  if(fMC)
537  {
538  return fMC->Header();
539  }
540  else
541  {
542  AliInfo("Header is not available");
543  return 0x0 ;
544  }
545 }
546 
547 //____________________________________________________
551 //____________________________________________________
553 {
554  fNMCProducedMin = 0;
555  fNMCProducedMax = 0;
556 
557  if (ReadStack() && fMC)
558  {
559  AliGenEventHeader * eventHeader = fMC->GenEventHeader();
560 
561  if(!fAcceptOnlyHIJINGLabels) return ;
562 
563  // TODO Check if it works from here ...
564 
565  AliGenCocktailEventHeader *cocktail = dynamic_cast<AliGenCocktailEventHeader *>(eventHeader);
566 
567  if(!cocktail) return ;
568 
569  TList *genHeaders = cocktail->GetHeaders();
570 
571  Int_t nGenerators = genHeaders->GetEntries();
572  //printf("N generators %d \n", nGenerators);
573 
574  for(Int_t igen = 0; igen < nGenerators; igen++)
575  {
576  AliGenEventHeader * eventHeader2 = (AliGenEventHeader*)genHeaders->At(igen) ;
577  TString name = eventHeader2->GetName();
578 
579  //printf("Generator %d: Class Name %s, Name %s, title %s \n",igen, eventHeader2->ClassName(), name.Data(), eventHeader2->GetTitle());
580 
582  fNMCProducedMax+= eventHeader2->NProduced();
583 
584  if(name.Contains("Hijing",TString::kIgnoreCase)) return ;
585  }
586 
587  }
588  else if(ReadAODMCParticles() && GetAODMCHeader())
589  {
590  Int_t nGenerators = GetAODMCHeader()->GetNCocktailHeaders();
591  //printf("AliCaloTrackReader::GetGenEventHeader() - N headers %d\n",nGenerators);
592 
593  if( nGenerators <= 0) return ;
594 
595  if(!fAcceptOnlyHIJINGLabels) return ;
596 
597  for(Int_t igen = 0; igen < nGenerators; igen++)
598  {
599  AliGenEventHeader * eventHeader = GetAODMCHeader()->GetCocktailHeader(igen) ;
600  TString name = eventHeader->GetName();
601 
602  //printf("Generator %d: Class Name %s, Name %s, title %s \n",igen, eventHeader->ClassName(), name.Data(), eventHeader->GetTitle());
603 
605  fNMCProducedMax+= eventHeader->NProduced();
606 
607  if(name.Contains("Hijing",TString::kIgnoreCase)) return ;
608  }
609 
610  }
611 }
612 
613 //______________________________________________________________
616 //______________________________________________________________
617 AliGenEventHeader* AliCaloTrackReader::GetGenEventHeader() const
618 {
619  if (ReadStack() && fMC)
620  {
621  AliGenEventHeader * eventHeader = fMC->GenEventHeader();
622 
623  if(!fAcceptOnlyHIJINGLabels) return eventHeader ;
624 
625  // TODO Check if it works from here ...
626 
627  AliGenCocktailEventHeader *cocktail = dynamic_cast<AliGenCocktailEventHeader *>(eventHeader);
628 
629  if(!cocktail) return 0x0 ;
630 
631  TList *genHeaders = cocktail->GetHeaders();
632 
633  Int_t nGenerators = genHeaders->GetEntries();
634  //printf("N generators %d \n", nGenerators);
635 
636  for(Int_t igen = 0; igen < nGenerators; igen++)
637  {
638  AliGenEventHeader * eventHeader2 = (AliGenEventHeader*)genHeaders->At(igen) ;
639  TString name = eventHeader2->GetName();
640 
641  //printf("Generator %d: Class Name %s, Name %s, title %s \n",igen, eventHeader2->ClassName(), name.Data(), eventHeader2->GetTitle());
642 
643  if(name.Contains("Hijing",TString::kIgnoreCase)) return eventHeader2 ;
644  }
645 
646  return 0x0;
647 
648  }
649  else if(ReadAODMCParticles() && GetAODMCHeader())
650  {
651  Int_t nGenerators = GetAODMCHeader()->GetNCocktailHeaders();
652  //printf("AliCaloTrackReader::GetGenEventHeader() - N headers %d\n",nGenerators);
653 
654  if( nGenerators <= 0) return 0x0;
655 
656  if(!fAcceptOnlyHIJINGLabels) return GetAODMCHeader()->GetCocktailHeader(0);
657 
658  for(Int_t igen = 0; igen < nGenerators; igen++)
659  {
660  AliGenEventHeader * eventHeader = GetAODMCHeader()->GetCocktailHeader(igen) ;
661  TString name = eventHeader->GetName();
662 
663  //printf("Generator %d: Class Name %s, Name %s, title %s \n",igen, eventHeader->ClassName(), name.Data(), eventHeader->GetTitle());
664 
665  if(name.Contains("Hijing",TString::kIgnoreCase)) return eventHeader ;
666  }
667 
668  return 0x0;
669 
670  }
671  else
672  {
673  //printf("AliCaloTrackReader::GetGenEventHeader() - MC header not available! \n");
674  return 0x0;
675  }
676 }
677 
678 //_________________________________________________________
681 //_________________________________________________________
683 {
684  AliInfo("Input are not AODs");
685 
686  return NULL ;
687 }
688 
689 //________________________________________________________
692 //________________________________________________________
693 AliAODMCHeader* AliCaloTrackReader::GetAODMCHeader() const
694 {
695  AliInfo("Input are not AODs");
696 
697  return NULL ;
698 }
699 
700 //___________________________________________________________
708 //___________________________________________________________
709 Int_t AliCaloTrackReader::GetVertexBC(const AliVVertex * vtx)
710 {
711  Int_t vertexBC=vtx->GetBC();
712 
713  if(!fRecalculateVertexBC) return vertexBC;
714 
715  // Value not available, recalculate it.
716 
717  Double_t bz = fInputEvent->GetMagneticField();
718  Bool_t bc0 = kFALSE;
719  Int_t ntr = GetCTSTracks()->GetEntriesFast();
720  //printf("N Tracks %d\n",ntr);
721 
722  for(Int_t i = 0 ; i < ntr ; i++)
723  {
724  AliVTrack * track = (AliVTrack*) (GetCTSTracks()->At(i));
725 
726  //Check if has TOF info, if not skip
727  ULong_t status = track->GetStatus();
728  Bool_t okTOF = (status & AliVTrack::kTOFout) == AliVTrack::kTOFout ;
729  vertexBC = track->GetTOFBunchCrossing(bz);
730  Float_t pt = track->Pt();
731 
732  if(!okTOF) continue;
733 
734  // Get DCA x, y
735  Double_t dca[2] = {1e6,1e6};
736  Double_t covar[3] = {1e6,1e6,1e6};
737  track->PropagateToDCA(vtx,bz,100.,dca,covar);
738 
739  if(AcceptDCA(pt,dca[0]))
740  {
741  if (vertexBC !=0 && fVertexBC != AliVTrack::kTOFBCNA) return vertexBC;
742  else if(vertexBC == 0) bc0 = kTRUE;
743  }
744  }
745 
746  if( bc0 ) vertexBC = 0 ;
747  else vertexBC = AliVTrack::kTOFBCNA ;
748 
749  return vertexBC;
750 }
751 
752 //_____________________________
755 //_____________________________
757 {
759  {
760  AliInfo("Cannot access stack and mcparticles at the same time, change them");
761  fReadStack = kFALSE;
762  fReadAODMCParticles = kFALSE;
763  }
764 }
765 
766 //_______________________________________
768 //_______________________________________
770 {
771  fDataType = kESD ;
772  fCTSPtMin = 0.1 ;
773  fEMCALPtMin = 0.1 ;
774  fPHOSPtMin = 0.1 ;
775  fCTSPtMax = 1000. ;
776  fEMCALPtMax = 1000. ;
777  fPHOSPtMax = 1000. ;
778 
779  //Track DCA cuts
780  // dca_xy cut = 0.0105+0.0350/TMath::Power(pt,1.1);
781  fTrackDCACut[0] = 0.0105;
782  fTrackDCACut[1] = 0.0350;
783  fTrackDCACut[2] = 1.1;
784 
785  //Do not filter the detectors input by default.
786  fFillEMCAL = kFALSE;
787  fFillDCAL = kFALSE;
788  fFillPHOS = kFALSE;
789  fFillCTS = kFALSE;
790  fFillEMCALCells = kFALSE;
791  fFillPHOSCells = kFALSE;
792 
793  fReadStack = kFALSE; // Check in the constructor of the other readers if it was set or in the configuration file
794  fReadAODMCParticles = kFALSE; // Check in the constructor of the other readers if it was set or in the configuration file
795  fDeltaAODFileName = "deltaAODPartCorr.root";
797  fEventTriggerMask = AliVEvent::kAny;
798  fMixEventTriggerMask = AliVEvent::kAnyINT;
799  fEventTriggerAtSE = kTRUE; // Use only events that pass event selection at SE base class
800 
801  fAcceptFastCluster = kTRUE;
802  fAnaLED = kFALSE;
803 
804  //We want tracks fitted in the detectors:
805  //fTrackStatus=AliESDtrack::kTPCrefit;
806  //fTrackStatus|=AliESDtrack::kITSrefit;
807  fTrackStatus = 0;
808 
809  fV0ADC[0] = 0; fV0ADC[1] = 0;
810  fV0Mul[0] = 0; fV0Mul[1] = 0;
811 
812  fZvtxCut = 10.;
813 
814  fNMixedEvent = 1;
815 
818 
819  //Centrality
820  fCentralityClass = "V0M";
821  fCentralityOpt = 100;
822  fCentralityBin[0] = fCentralityBin[1]=-1;
823 
824  fEventPlaneMethod = "V0";
825 
826  // Allocate memory (not sure this is the right place)
827  fCTSTracks = new TObjArray();
828  fEMCALClusters = new TObjArray();
829  fDCALClusters = new TObjArray();
830  fPHOSClusters = new TObjArray();
831  //fTriggerAnalysis = new AliTriggerAnalysis;
832  fAODBranchList = new TList ;
833 
834  fPileUpParamSPD[0] = 3 ; fPileUpParamSPD[1] = 0.8 ;
835  fPileUpParamSPD[2] = 3.0 ; fPileUpParamSPD[3] = 2.0 ; fPileUpParamSPD[4] = 5.0;
836 
837  // Parametrized time cut (LHC11d)
840 
841  // Parametrized time cut (LHC11c)
842  //fEMCALParamTimeCutMin[0] =-5; fEMCALParamTimeCutMin[1] =-1 ; fEMCALParamTimeCutMin[2] = 1.87; fEMCALParamTimeCutMin[3] = 0.4;
843  //fEMCALParamTimeCutMax[0] = 3.5; fEMCALParamTimeCutMax[1] = 50; fEMCALParamTimeCutMax[2] = 0.15; fEMCALParamTimeCutMax[3] = 1.6;
844 
845  fTimeStampRunMin = -1;
846  fTimeStampRunMax = 1e12;
849 
850  for(Int_t i = 0; i < 19; i++)
851  {
852  fEMCalBCEvent [i] = 0;
853  fEMCalBCEventCut[i] = 0;
854  fTrackBCEvent [i] = 0;
855  fTrackBCEventCut[i] = 0;
856  }
857 
858  // Trigger match-rejection
861 
862  fTriggerClusterBC = -10000 ;
866  fTriggerClusterId = -1;
867 
868  //Jets
871  if(!fNonStandardJets) fNonStandardJets = new TClonesArray("AliAODJet",100);
874  if(!fBackgroundJets) fBackgroundJets = new AliAODJetEventBackground();
875 
876  fSmearShowerShapeWidth = 0.005;
877 }
878 
879 //___________________________________________________
884 //___________________________________________________
886 {
887  AliGenHijingEventHeader* hijingHeader = dynamic_cast<AliGenHijingEventHeader *> (GetGenEventHeader());
888 
889  //printf("header %p, label %d\n",hijingHeader,label);
890 
891  if(!hijingHeader || label < 0 ) return kFALSE;
892 
893 
894  //printf("pass a), N produced %d\n",nproduced);
895 
896  if(label >= fNMCProducedMin && label < fNMCProducedMax)
897  {
898  //printf(" accept!, label is smaller than produced, N %d\n",nproduced);
899 
900  return kTRUE;
901  }
902 
903  if(ReadStack())
904  {
905  if(!GetStack()) return kFALSE;
906 
907  Int_t nprimaries = GetStack()->GetNtrack();
908 
909  if(label > nprimaries) return kFALSE;
910 
911  TParticle * mom = GetStack()->Particle(label);
912 
913  Int_t iMom = label;
914  Int_t iParent = mom->GetFirstMother();
915  while(iParent!=-1)
916  {
917  if(iParent >= fNMCProducedMin && iParent < fNMCProducedMax)
918  {
919  //printf("\t accept, mother is %d \n",iParent)
920  return kTRUE;
921  }
922 
923  iMom = iParent;
924  mom = GetStack()->Particle(iMom);
925  iParent = mom->GetFirstMother();
926  }
927 
928  return kFALSE ;
929 
930  } // ESD
931  else
932  {
933  TClonesArray* mcparticles = GetAODMCParticles();
934 
935  if(!mcparticles) return kFALSE;
936 
937  Int_t nprimaries = mcparticles->GetEntriesFast();
938 
939  if(label > nprimaries) return kFALSE;
940 
941  //printf("pass b) N primaries %d \n",nprimaries);
942 
943  if(label >= fNMCProducedMin && label < fNMCProducedMax)
944  {
945  return kTRUE;
946  }
947 
948  // Find grand parent, check if produced in the good range
949  AliAODMCParticle * mom = (AliAODMCParticle *) mcparticles->At(label);
950 
951  Int_t iMom = label;
952  Int_t iParent = mom->GetMother();
953  while(iParent!=-1)
954  {
955  if(iParent >= fNMCProducedMin && iParent < fNMCProducedMax)
956  {
957  //printf("\t accept, mother is %d, with nProduced %d \n",iParent, nproduced);
958  return kTRUE;
959  }
960 
961  iMom = iParent;
962  mom = (AliAODMCParticle *) mcparticles->At(iMom);
963  iParent = mom->GetMother();
964 
965  }
966 
967  //printf("pass c), no match found \n");
968 
969  return kFALSE ;
970 
971  }//AOD
972 }
973 
974 //__________________________________________________________________________
977 //__________________________________________________________________________
978 Bool_t AliCaloTrackReader::IsInTimeWindow(Double_t tof, Float_t energy) const
979 {
980  // Parametrized cut depending on E
981  if(fUseParamTimeCut)
982  {
983  Float_t minCut= fEMCALParamTimeCutMin[0]+fEMCALParamTimeCutMin[1]*TMath::Exp(-(energy-fEMCALParamTimeCutMin[2])/fEMCALParamTimeCutMin[3]);
984  Float_t maxCut= fEMCALParamTimeCutMax[0]+fEMCALParamTimeCutMax[1]*TMath::Exp(-(energy-fEMCALParamTimeCutMax[2])/fEMCALParamTimeCutMax[3]);
985  //printf("tof %f, minCut %f, maxCut %f\n",tof,minCut,maxCut);
986  if( tof < minCut || tof > maxCut ) return kFALSE ;
987  }
988 
989  //In any case, the time should to be larger than the fixed window ...
990  if( tof < fEMCALTimeCutMin || tof > fEMCALTimeCutMax ) return kFALSE ;
991 
992  return kTRUE ;
993 }
994 
995 //________________________________________________
998 //________________________________________________
1000 {
1001  return fInputEvent->IsPileupFromSPD((Int_t) fPileUpParamSPD[0] , fPileUpParamSPD[1] ,
1002  fPileUpParamSPD[2] , fPileUpParamSPD[3] , fPileUpParamSPD[4] );
1003  //printf("Param : %d, %2.2f, %2.2f, %2.2f, %2.2f\n",(Int_t) fPileUpParamSPD[0], fPileUpParamSPD[1], fPileUpParamSPD[2], fPileUpParamSPD[3], fPileUpParamSPD[4]);
1004 }
1005 
1006 //__________________________________________________
1008 //__________________________________________________
1010 {
1011  if(fNPileUpClusters > fNPileUpClustersCut) return kTRUE ;
1012  else return kFALSE;
1013 }
1014 
1015 //________________________________________________________
1017 //________________________________________________________
1019 {
1020  if( IsPileUpFromSPD() && IsPileUpFromEMCal()) return kTRUE ;
1021  else return kFALSE;
1022 }
1023 
1024 //_______________________________________________________
1026 //_______________________________________________________
1028 {
1029  if( IsPileUpFromSPD() || IsPileUpFromEMCal()) return kTRUE ;
1030  else return kFALSE;
1031 }
1032 
1033 //___________________________________________________________
1035 //___________________________________________________________
1037 {
1038  if( IsPileUpFromSPD() && !IsPileUpFromEMCal()) return kTRUE ;
1039  else return kFALSE;
1040 }
1041 
1042 //___________________________________________________________
1044 //___________________________________________________________
1046 {
1047  if( !IsPileUpFromSPD() && IsPileUpFromEMCal()) return kTRUE ;
1048  else return kFALSE;
1049 }
1050 
1051 //______________________________________________________________
1053 //______________________________________________________________
1055 {
1056  if( !IsPileUpFromSPD() && !IsPileUpFromEMCal()) return kTRUE ;
1057  else return kFALSE;
1058 }
1059 
1060 //___________________________________________________________________________________
1068 //___________________________________________________________________________________
1069 Bool_t AliCaloTrackReader::FillInputEvent(Int_t iEntry, const char * /*curFileName*/)
1070 {
1071  fEventNumber = iEntry;
1072  fTriggerClusterIndex = -1;
1073  fTriggerClusterId = -1;
1074  fIsTriggerMatch = kFALSE;
1075  fTriggerClusterBC = -10000;
1076  fIsExoticEvent = kFALSE;
1077  fIsBadCellEvent = kFALSE;
1078  fIsBadMaxCellEvent = kFALSE;
1079 
1080  fIsTriggerMatchOpenCut[0] = kFALSE ;
1081  fIsTriggerMatchOpenCut[1] = kFALSE ;
1082  fIsTriggerMatchOpenCut[2] = kFALSE ;
1083 
1084  //fCurrentFileName = TString(currentFileName);
1085  if(!fInputEvent)
1086  {
1087  AliInfo("Input event not available, skip event analysis");
1088  return kFALSE;
1089  }
1090 
1091  //-----------------------------------------------
1092  // Select the event depending on the trigger type
1093  // and other event characteristics
1094  // like the goodness of the EMCal trigger
1095  //-----------------------------------------------
1096 
1097  Bool_t accept = CheckEventTriggers();
1098  if(!accept) return kFALSE;
1099 
1100  AliDebug(1,"Pass Event trigger selection");
1101 
1102  //---------------------------------------------------------------------------
1103  // In case of analysis of events with jets, skip those with jet pt > 5 pt hard
1104  // To be used on for MC data in pT hard bins
1105  //---------------------------------------------------------------------------
1106 
1108  {
1109  if(!ComparePtHardAndJetPt()) return kFALSE ;
1110  AliDebug(1,"Pass Pt Hard - Jet rejection");
1111  }
1112 
1114  {
1115  if(!ComparePtHardAndClusterPt()) return kFALSE ;
1116  AliDebug(1,"Pass Pt Hard - Cluster rejection");
1117  }
1118 
1119  //------------------------------------------------------
1120  // Event rejection depending on time stamp
1121  //------------------------------------------------------
1122 
1124  {
1125  AliESDEvent* esd = dynamic_cast<AliESDEvent*> (fInputEvent);
1126  if(esd)
1127  {
1128  Int_t timeStamp = esd->GetTimeStamp();
1129  Float_t timeStampFrac = 1.*(timeStamp-fTimeStampRunMin) / (fTimeStampRunMax-fTimeStampRunMin);
1130 
1131  //printf("stamp0 %d, max0 %d, frac %f\n", timeStamp-fTimeStampRunMin,fTimeStampRunMax-fTimeStampRunMin, timeStampFrac);
1132 
1133  if(timeStampFrac < fTimeStampEventFracMin || timeStampFrac > fTimeStampEventFracMax) return kFALSE;
1134  }
1135 
1136  AliDebug(1,"Pass Time Stamp rejection");
1137  }
1138 
1139  //------------------------------------------------------
1140  // Event rejection depending on vertex, pileup, v0and
1141  //------------------------------------------------------
1142 
1143  FillVertexArray();
1144 
1145  //Reject events with Z vertex too large, only for SE analysis, if not, cut on the analysis code
1146  if(!GetMixedEvent() && TMath::Abs(fVertex[0][2]) > fZvtxCut) return kFALSE;
1147 
1149  {
1150  if( !CheckForPrimaryVertex() ) return kFALSE; // algorithm in ESD/AOD Readers
1151  if( TMath::Abs(fVertex[0][0] ) < 1.e-6 &&
1152  TMath::Abs(fVertex[0][1] ) < 1.e-6 &&
1153  TMath::Abs(fVertex[0][2] ) < 1.e-6 ) return kFALSE;
1154  }
1155 
1156  AliDebug(1,"Pass primary vertex rejection");
1157 
1158  //printf("Reader : IsPileUp %d, Multi %d\n",IsPileUpFromSPD(),fInputEvent->IsPileupFromSPDInMultBins());
1159 
1161  {
1162  // Do not analyze events with pileup
1163  Bool_t bPileup = IsPileUpFromSPD();
1164  //IsPileupFromSPDInMultBins() // method to try
1165  //printf("pile-up %d, %d, %2.2f, %2.2f, %2.2f, %2.2f\n",bPileup, (Int_t) fPileUpParamSPD[0], fPileUpParamSPD[1], fPileUpParamSPD[2], fPileUpParamSPD[3], fPileUpParamSPD[4]);
1166  if(bPileup) return kFALSE;
1167 
1168  AliDebug(1,"Pass Pile-Up event rejection");
1169  }
1170 
1172  {
1173  AliVVZERO* v0 = fInputEvent->GetVZEROData();
1174 
1175  Bool_t bV0AND = ((v0->GetV0ADecision()==1) && (v0->GetV0CDecision()==1));
1176  //bV0AND = fTriggerAnalysis->IsOfflineTriggerFired((AliESDEvent*)fInputEvent, AliTriggerAnalysis::kV0AND);
1177  //printf("V0AND event? %d\n",bV0AND);
1178 
1179  if(!bV0AND)
1180  {
1181  AliDebug(1,"Reject event by V0AND");
1182  return kFALSE;
1183  }
1184 
1185  AliDebug(1,"Pass V0AND event rejection");
1186  }
1187 
1188  //------------------------------------------------------
1189  // Check if there is a centrality value, PbPb analysis,
1190  // and if a centrality bin selection is requested
1191  //------------------------------------------------------
1192 
1193  //If we need a centrality bin, we select only those events in the corresponding bin.
1194  if(GetCentrality() && fCentralityBin[0]>=0 && fCentralityBin[1]>=0 && fCentralityOpt==100)
1195  {
1196  Int_t cen = GetEventCentrality();
1197  if(cen > fCentralityBin[1] || cen < fCentralityBin[0]) return kFALSE; //reject events out of bin.
1198 
1199  AliDebug(1,"Pass centrality rejection");
1200  }
1201 
1202  //---------------------------------------------------------------------------
1203  // In case of MC analysis, set the range of interest in the MC particles list
1204  // The particle selection is done inside the FillInputDetector() methods
1205  //---------------------------------------------------------------------------
1207 
1208  //printf("N min %d, N max %d\n",fNMCProducedMin,fNMCProducedMax);
1209 
1210  //-------------------------------------------------------
1211  // Get the main vertex BC, in case not available
1212  // it is calculated in FillCTS checking the BC of tracks
1213  //------------------------------------------------------
1214  fVertexBC = fInputEvent->GetPrimaryVertex()->GetBC();
1215 
1216  //-----------------------------------------------
1217  // Fill the arrays with cluster/tracks/cells data
1218  //-----------------------------------------------
1219 
1220  if(fFillCTS)
1221  {
1222  FillInputCTS();
1223  //Accept events with at least one track
1224  if(fTrackMult == 0 && fDoRejectNoTrackEvents) return kFALSE ;
1225 
1226  AliDebug(1,"Pass rejection of null track events");
1227  }
1228 
1230  {
1231  if(fVertexBC != 0 && fVertexBC != AliVTrack::kTOFBCNA) return kFALSE ;
1232 
1233  AliDebug(1,"Pass rejection of events with vertex at BC!=0");
1234  }
1235 
1236  if(fFillEMCALCells)
1238 
1239  if(fFillPHOSCells)
1241 
1242  if(fFillEMCAL || fFillDCAL)
1243  FillInputEMCAL();
1244 
1245  if(fFillPHOS)
1246  FillInputPHOS();
1247 
1248  FillInputVZERO();
1249 
1250  //one specified jet branch
1255 
1256  AliDebug(1,"Event accepted for analysis");
1257 
1258  return kTRUE ;
1259 }
1260 
1261 //__________________________________________________
1264 //__________________________________________________
1266 {
1267  if( !GetCentrality() ) return -1;
1268 
1269  if (fCentralityOpt==100) return (Int_t) GetCentrality()->GetCentralityPercentile(fCentralityClass); // 100 bins max
1270  else if(fCentralityOpt==10) return GetCentrality()->GetCentralityClass10(fCentralityClass);// 10 bins max
1271  else if(fCentralityOpt==20) return GetCentrality()->GetCentralityClass5(fCentralityClass); // 20 bins max
1272  else
1273  {
1274  AliInfo(Form("Unknown centrality option %d, use 10, 20 or 100\n",fCentralityOpt));
1275  return -1;
1276  }
1277 }
1278 
1279 //_____________________________________________________
1282 //_____________________________________________________
1284 {
1285  if( !GetEventPlane() ) return -1000;
1286 
1287  Float_t ep = GetEventPlane()->GetEventplane(GetEventPlaneMethod(), GetInputEvent());
1288 
1289  if(GetEventPlaneMethod()=="Q" && (ep < 0 || ep > TMath::Pi()))
1290  {
1291  AliDebug(1,Form("Bad EP for <Q> method : %f\n",ep));
1292  return -1000;
1293  }
1294  else if(GetEventPlaneMethod().Contains("V0") )
1295  {
1296  if((ep > TMath::Pi()/2 || ep < -TMath::Pi()/2))
1297  {
1298  AliDebug(1,Form("Bad EP for <%s> method : %f\n",GetEventPlaneMethod().Data(), ep));
1299  return -1000;
1300  }
1301 
1302  ep+=TMath::Pi()/2; // put same range as for <Q> method
1303  }
1304 
1305  AliDebug(3,Form("Event plane angle %f",ep));
1306 
1307 // if(fDebug > 0 )
1308 // {
1309 // if (ep > TMath::Pi()) printf("AliCaloTrackReader::GetEventPlaneAngle() - Too large angle = %f\n",ep);
1310 // else if(ep < 0 ) printf("AliCaloTrackReader::GetEventPlaneAngle() - Negative angle = %f\n" ,ep);
1311 // }
1312 
1313  return ep;
1314 }
1315 
1316 //__________________________________________________________
1318 //__________________________________________________________
1319 void AliCaloTrackReader::GetVertex(Double_t vertex[3]) const
1320 {
1321  vertex[0] = fVertex[0][0];
1322  vertex[1] = fVertex[0][1];
1323  vertex[2] = fVertex[0][2];
1324 }
1325 
1326 //__________________________________________________________________________
1328 //__________________________________________________________________________
1329 void AliCaloTrackReader::GetVertex(Double_t vertex[3], Int_t evtIndex) const
1330 {
1331  vertex[0] = fVertex[evtIndex][0];
1332  vertex[1] = fVertex[evtIndex][1];
1333  vertex[2] = fVertex[evtIndex][2];
1334 }
1335 
1336 //________________________________________
1339 //________________________________________
1341 {
1342  // Delete previous vertex
1343  if(fVertex)
1344  {
1345  for (Int_t i = 0; i < fNMixedEvent; i++)
1346  {
1347  delete [] fVertex[i] ;
1348  }
1349  delete [] fVertex ;
1350  }
1351 
1352  fVertex = new Double_t*[fNMixedEvent] ;
1353  for (Int_t i = 0; i < fNMixedEvent; i++)
1354  {
1355  fVertex[i] = new Double_t[3] ;
1356  fVertex[i][0] = 0.0 ;
1357  fVertex[i][1] = 0.0 ;
1358  fVertex[i][2] = 0.0 ;
1359  }
1360 
1361  if (!fMixedEvent)
1362  { // Single event analysis
1363  if(fDataType!=kMC)
1364  {
1365 
1366  if(fInputEvent->GetPrimaryVertex())
1367  {
1368  fInputEvent->GetPrimaryVertex()->GetXYZ(fVertex[0]);
1369  }
1370  else
1371  {
1372  AliWarning("NULL primary vertex");
1373  fVertex[0][0]=0.; fVertex[0][1]=0.; fVertex[0][2]=0.;
1374  }//Primary vertex pointer do not exist
1375 
1376  } else
1377  {// MC read event
1378  fVertex[0][0]=0.; fVertex[0][1]=0.; fVertex[0][2]=0.;
1379  }
1380 
1381  AliDebug(1,Form("Single Event Vertex : %f,%f,%f\n",fVertex[0][0],fVertex[0][1],fVertex[0][2]));
1382 
1383  } else
1384  { // MultiEvent analysis
1385  for (Int_t iev = 0; iev < fNMixedEvent; iev++)
1386  {
1387  if (fMixedEvent->GetVertexOfEvent(iev))
1388  fMixedEvent->GetVertexOfEvent(iev)->GetXYZ(fVertex[iev]);
1389  else
1390  AliWarning("No vertex found");
1391 
1392  AliDebug(1,Form("Multi Event %d Vertex : %f,%f,%f",iev,fVertex[iev][0],fVertex[iev][1],fVertex[iev][2]));
1393  }
1394  }
1395 }
1396 
1397 //_____________________________________
1403 //_____________________________________
1405 {
1406  AliDebug(1,"Begin");
1407 
1408  Double_t pTrack[3] = {0,0,0};
1409 
1410  Int_t nTracks = fInputEvent->GetNumberOfTracks() ;
1411  fTrackMult = 0;
1412  Int_t nstatus = 0;
1413  Double_t bz = GetInputEvent()->GetMagneticField();
1414 
1415  for(Int_t i = 0; i < 19; i++)
1416  {
1417  fTrackBCEvent [i] = 0;
1418  fTrackBCEventCut[i] = 0;
1419  }
1420 
1421  Bool_t bc0 = kFALSE;
1422  if(fRecalculateVertexBC) fVertexBC = AliVTrack::kTOFBCNA;
1423 
1424  for (Int_t itrack = 0; itrack < nTracks; itrack++)
1425  {
1426  AliVTrack * track = (AliVTrack*)fInputEvent->GetTrack(itrack) ; // retrieve track from esd
1427 
1428  if(fAcceptOnlyHIJINGLabels && !IsHIJINGLabel(track->GetLabel())) continue ;
1429 
1430  //Select tracks under certain conditions, TPCrefit, ITSrefit ... check the set bits
1431  ULong_t status = track->GetStatus();
1432 
1433  if (fTrackStatus && !((status & fTrackStatus) == fTrackStatus))
1434  continue ;
1435 
1436  nstatus++;
1437 
1438  // Select the tracks depending on cuts of AOD or ESD
1439  if(!SelectTrack(track, pTrack)) continue ;
1440 
1441  // TOF cuts
1442  Bool_t okTOF = ( (status & AliVTrack::kTOFout) == AliVTrack::kTOFout ) ;
1443  Double_t tof = -1000;
1444  Int_t trackBC = -1000 ;
1445 
1446  if(fAccessTrackTOF)
1447  {
1448  if(okTOF)
1449  {
1450  trackBC = track->GetTOFBunchCrossing(bz);
1451  SetTrackEventBC(trackBC+9);
1452 
1453  tof = track->GetTOFsignal()*1e-3;
1454 
1455  //After selecting tracks with small DCA, pointing to vertex, set vertex BC depeding on tracks BC
1457  {
1458  if (trackBC != 0 && trackBC != AliVTrack::kTOFBCNA) fVertexBC = trackBC;
1459  else if(trackBC == 0) bc0 = kTRUE;
1460  }
1461 
1462  //In any case, the time should to be larger than the fixed window ...
1463  if( fUseTrackTimeCut && (trackBC !=0 || tof < fTrackTimeCutMin || tof > fTrackTimeCutMax) )
1464  {
1465  //printf("Remove track time %f and bc = %d\n",tof,trackBC);
1466  continue ;
1467  }
1468  //else printf("Accept track time %f and bc = %d\n",tof,trackBC);
1469  }
1470  }
1471 
1472  fMomentum.SetPxPyPzE(pTrack[0],pTrack[1],pTrack[2],0);
1473 
1474  // DCA cuts
1475 
1476  if(fUseTrackDCACut)
1477  {
1478  Float_t dcaTPC =-999;
1479  //In case of AODs, TPC tracks cannot be propagated back to primary vertex,
1480  if( fDataType == kAOD ) dcaTPC = ((AliAODTrack*) track)->DCA();
1481 
1482  //normal way to get the dca, cut on dca_xy
1483  if(dcaTPC==-999)
1484  {
1485  Double_t dca[2] = {1e6,1e6};
1486  Double_t covar[3] = {1e6,1e6,1e6};
1487  Bool_t okDCA = track->PropagateToDCA(fInputEvent->GetPrimaryVertex(),bz,100.,dca,covar);
1488  if( okDCA) okDCA = AcceptDCA(fMomentum.Pt(),dca[0]);
1489  if(!okDCA)
1490  {
1491  //printf("AliCaloTrackReader::FillInputCTS() - Reject track pt %2.2f, dca_xy %2.4f\n",fMomentum.Pt(),dca[0]);
1492  continue ;
1493  }
1494  }
1495  }// DCA cuts
1496 
1497  //Count the tracks in eta < 0.9
1498  if(TMath::Abs(track->Eta())< fTrackMultEtaCut) fTrackMult++;
1499 
1500  if(fCTSPtMin > fMomentum.Pt() || fCTSPtMax < fMomentum.Pt()) continue ;
1501 
1502  // Check effect of cuts on track BC
1503  if(fAccessTrackTOF && okTOF) SetTrackEventBCcut(trackBC+9);
1504 
1505  if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(fMomentum.Eta(),fMomentum.Phi(),kCTS)) continue;
1506 
1507  AliDebug(2,Form("Selected tracks pt %3.2f, phi %3.2f, eta %3.2f",
1508  fMomentum.Pt(),fMomentum.Phi()*TMath::RadToDeg(),fMomentum.Eta()));
1509 
1510  if (fMixedEvent) track->SetID(itrack);
1511 
1512 
1513  fCTSTracks->Add(track);
1514 
1515  }// track loop
1516 
1517  if( fRecalculateVertexBC && (fVertexBC == 0 || fVertexBC == AliVTrack::kTOFBCNA))
1518  {
1519  if( bc0 ) fVertexBC = 0 ;
1520  else fVertexBC = AliVTrack::kTOFBCNA ;
1521  }
1522 
1523  AliDebug(1,Form("AOD entries %d, input tracks %d, pass status %d, multipliticy %d", fCTSTracks->GetEntriesFast(), nTracks, nstatus, fTrackMult));//fCTSTracksNormalInputEntries);
1524 }
1525 
1526 //_______________________________________________________________________________
1533 //_______________________________________________________________________________
1534 void AliCaloTrackReader::FillInputEMCALAlgorithm(AliVCluster * clus, Int_t iclus)
1535 {
1536  // Accept clusters with the proper label, TO DO, use the new more General methods!!!
1537  if(fAcceptOnlyHIJINGLabels && !IsHIJINGLabel( clus->GetLabel() )) return ;
1538 
1539  Int_t vindex = 0 ;
1540  if (fMixedEvent)
1541  vindex = fMixedEvent->EventIndexForCaloCluster(iclus);
1542 
1543 
1544  clus->GetMomentum(fMomentum, fVertex[vindex]);
1545 
1546  //if( (fDebug > 2 && fMomentum.E() > 0.1) || fDebug > 10 )
1547  AliDebug(2,Form("Input cluster E %3.2f, pt %3.2f, phi %3.2f, eta %3.2f",
1548  fMomentum.E(),fMomentum.Pt(),fMomentum.Phi()*TMath::RadToDeg(),fMomentum.Eta()));
1549 
1551  {
1552  //Recalibrate the cluster energy
1553  if(GetCaloUtils()->IsRecalibrationOn())
1554  {
1555  Float_t energy = GetCaloUtils()->RecalibrateClusterEnergy(clus, GetEMCALCells());
1556 
1557  clus->SetE(energy);
1558  //printf("Recalibrated Energy %f\n",clus->E());
1559 
1562 
1563  } // recalculate E
1564 
1565  //Recalculate distance to bad channels, if new list of bad channels provided
1567 
1568  //Recalculate cluster position
1569  if(GetCaloUtils()->IsRecalculationOfClusterPositionOn())
1570  {
1572  //clus->GetPosition(pos);
1573  //printf("After Corrections: e %f, x %f, y %f, z %f\n",clus->E(),pos[0],pos[1],pos[2]);
1574  }
1575 
1576  // Recalculate TOF
1577  if(GetCaloUtils()->GetEMCALRecoUtils()->IsTimeRecalibrationOn())
1578  {
1579  Double_t tof = clus->GetTOF();
1580  Float_t frac =-1;
1581  Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(fEMCALCells, clus,frac);
1582 
1584  {
1585  tof = fEMCALCells->GetCellTime(absIdMax);
1586  }
1587 
1588  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absIdMax,fInputEvent->GetBunchCrossNumber(),tof);
1589 
1590  clus->SetTOF(tof);
1591 
1592  }// Time recalibration
1593  }
1594 
1595  //Reject clusters with bad channels, close to borders and exotic;
1596  Bool_t goodCluster = GetCaloUtils()->GetEMCALRecoUtils()->IsGoodCluster(clus,
1597  GetCaloUtils()->GetEMCALGeometry(),
1598  GetEMCALCells(),fInputEvent->GetBunchCrossNumber());
1599 
1600  if(!goodCluster)
1601  {
1602  //if( (fDebug > 2 && fMomentum.E() > 0.1) || fDebug > 10 )
1603  AliDebug(2,Form("Bad cluster E %3.2f, pt %3.2f, phi %3.2f, eta %3.2f",
1604  fMomentum.E(),fMomentum.Pt(),fMomentum.Phi()*TMath::RadToDeg(),fMomentum.Eta()));
1605 
1606  return;
1607  }
1608 
1609  //Mask all cells in collumns facing ALICE thick material if requested
1610  if(GetCaloUtils()->GetNMaskCellColumns())
1611  {
1612  Int_t absId = -1;
1613  Int_t iSupMod = -1;
1614  Int_t iphi = -1;
1615  Int_t ieta = -1;
1616  Bool_t shared = kFALSE;
1617  GetCaloUtils()->GetEMCALRecoUtils()->GetMaxEnergyCell(GetCaloUtils()->GetEMCALGeometry(), GetEMCALCells(),clus,absId,iSupMod,ieta,iphi,shared);
1618  if(GetCaloUtils()->MaskFrameCluster(iSupMod, ieta)) return;
1619  }
1620 
1622  {
1623  if(clus->GetNLabels()==0 || clus->GetLabel() < 0) return;
1624  //else printf("Embedded cluster, %d, n label %d label %d \n",iclus,clus->GetNLabels(),clus->GetLabel());
1625  }
1626 
1627  //Float_t pos[3];
1628  //clus->GetPosition(pos);
1629  //printf("Before Corrections: e %f, x %f, y %f, z %f\n",clus->E(),pos[0],pos[1],pos[2]);
1630 
1631  //Correct non linearity or smear energy
1633  {
1635 
1636  //if( (fDebug > 5 && fMomentum.E() > 0.1) || fDebug > 10 )
1637  AliDebug(5,Form("Correct Non Lin: Old E %3.2f, New E %3.2f",
1638  fMomentum.E(),clus->E()));
1639 
1640  // In case of MC analysis, to match resolution/calibration in real data
1641  // Not needed anymore, just leave for MC studies on systematics
1642  if( GetCaloUtils()->GetEMCALRecoUtils()->IsClusterEnergySmeared() )
1643  {
1644  Float_t rdmEnergy = GetCaloUtils()->GetEMCALRecoUtils()->SmearClusterEnergy(clus);
1645 
1646  //if( (fDebug > 5 && fMomentum.E() > 0.1) || fDebug > 10 )
1647  AliDebug(5,Form("Smear energy: Old E %3.2f, New E %3.2f",clus->E(),rdmEnergy));
1648 
1649  clus->SetE(rdmEnergy);
1650  }
1651  }
1652 
1653  Double_t tof = clus->GetTOF()*1e9;
1654 
1655  Int_t bc = TMath::Nint(tof/50) + 9;
1656  //printf("tof %2.2f, bc+5=%d\n",tof,bc);
1657 
1658  SetEMCalEventBC(bc);
1659 
1660  if(fEMCALPtMin > clus->E() || fEMCALPtMax < clus->E()) return ;
1661 
1662  clus->GetMomentum(fMomentum, fVertex[vindex]);
1663 
1664  SetEMCalEventBCcut(bc);
1665 
1666  if(!IsInTimeWindow(tof,clus->E()))
1667  {
1668  fNPileUpClusters++ ;
1669  if(fUseEMCALTimeCut) return ;
1670  }
1671  else
1673 
1674  //if((fDebug > 2 && fMomentum.E() > 0.1) || fDebug > 10)
1675  AliDebug(2,Form("Selected clusters E %3.2f, pt %3.2f, phi %3.2f, eta %3.2f",
1676  fMomentum.E(),fMomentum.Pt(),fMomentum.Phi()*TMath::RadToDeg(),fMomentum.Eta()));
1677 
1678  if (fMixedEvent)
1679  clus->SetID(iclus) ;
1680 
1681  // Select cluster fiducial region
1682  Bool_t bEMCAL = kFALSE;
1683  Bool_t bDCAL = kFALSE;
1684  if(fCheckFidCut)
1685  {
1686  if(fFillEMCAL && fFiducialCut->IsInFiducialCut(fMomentum.Eta(),fMomentum.Phi(),kEMCAL)) bEMCAL = kTRUE ;
1687  if(fFillDCAL && fFiducialCut->IsInFiducialCut(fMomentum.Eta(),fMomentum.Phi(),kDCAL )) bDCAL = kTRUE ;
1688  }
1689  else
1690  {
1691  bEMCAL = kTRUE;
1692  }
1693 
1694  // Smear the SS to try to match data and simulations,
1695  // do it only for simulations.
1696  if( fSmearShowerShape && clus->GetNCells() > 2)
1697  {
1698  AliDebug(2,Form("Smear shower shape - Original: %2.4f\n", clus->GetM02()));
1699  clus->SetM02( clus->GetM02() + fRandom.Landau(0, fSmearShowerShapeWidth) );
1700  //clus->SetM02( fRandom.Landau(clus->GetM02(), fSmearShowerShapeWidth) );
1701  AliDebug(2,Form("Width %2.4f Smeared : %2.4f\n", fSmearShowerShapeWidth,clus->GetM02()));
1702  }
1703 
1704  // Fill the corresponding array. Usually just filling EMCal array with upper or lower clusters is enough, but maybe we want to do EMCal-DCal correlations.
1705  if (bEMCAL) fEMCALClusters->Add(clus);
1706  else if(bDCAL ) fDCALClusters ->Add(clus);
1707 }
1708 
1709 //_______________________________________
1716 //_______________________________________
1718 {
1719  AliDebug(1,"Begin");
1720 
1721  // First recalibrate cells, time or energy
1722  // if(GetCaloUtils()->IsRecalibrationOn())
1723  // GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCells(GetCaloUtils()->GetEMCALGeometry(),
1724  // GetEMCALCells(),
1725  // fInputEvent->GetBunchCrossNumber());
1726 
1727  fNPileUpClusters = 0; // Init counter
1728  fNNonPileUpClusters = 0; // Init counter
1729  for(Int_t i = 0; i < 19; i++)
1730  {
1731  fEMCalBCEvent [i] = 0;
1732  fEMCalBCEventCut[i] = 0;
1733  }
1734 
1735  //Loop to select clusters in fiducial cut and fill container with aodClusters
1736  if(fEMCALClustersListName=="")
1737  {
1738  Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
1739  for (Int_t iclus = 0; iclus < nclusters; iclus++)
1740  {
1741  AliVCluster * clus = 0;
1742  if ( (clus = fInputEvent->GetCaloCluster(iclus)) )
1743  {
1744  if (clus->IsEMCAL())
1745  {
1746  FillInputEMCALAlgorithm(clus, iclus);
1747  }//EMCAL cluster
1748  }// cluster exists
1749  }// cluster loop
1750 
1751  //Recalculate track matching
1753 
1754  }//Get the clusters from the input event
1755  else
1756  {
1757  TClonesArray * clusterList = 0x0;
1758 
1759  if (fInputEvent->FindListObject(fEMCALClustersListName))
1760  {
1761  clusterList = dynamic_cast<TClonesArray*> (fInputEvent->FindListObject(fEMCALClustersListName));
1762  }
1763  else if(fOutputEvent)
1764  {
1765  clusterList = dynamic_cast<TClonesArray*> (fOutputEvent->FindListObject(fEMCALClustersListName));
1766  }
1767 
1768  if(!clusterList)
1769  {
1770  AliWarning(Form("Wrong name of list with clusters? <%s>",fEMCALClustersListName.Data()));
1771  return;
1772  }
1773 
1774  Int_t nclusters = clusterList->GetEntriesFast();
1775  for (Int_t iclus = 0; iclus < nclusters; iclus++)
1776  {
1777  AliVCluster * clus = dynamic_cast<AliVCluster*> (clusterList->At(iclus));
1778  //printf("E %f\n",clus->E());
1779  if (clus) FillInputEMCALAlgorithm(clus, iclus);
1780  else AliWarning("Null cluster in list!");
1781  }// cluster loop
1782 
1783  // Recalculate the pile-up time, in case long time clusters removed during clusterization
1784  //printf("Input event INIT : Pile-up clusters %d, NO pile-up %d\n",fNPileUpClusters,fNNonPileUpClusters);
1785 
1786  fNPileUpClusters = 0; // Init counter
1787  fNNonPileUpClusters = 0; // Init counter
1788  for(Int_t i = 0; i < 19; i++)
1789  {
1790  fEMCalBCEvent [i] = 0;
1791  fEMCalBCEventCut[i] = 0;
1792  }
1793 
1794  for (Int_t iclus = 0; iclus < fInputEvent->GetNumberOfCaloClusters(); iclus++)
1795  {
1796  AliVCluster * clus = 0;
1797 
1798  if ( (clus = fInputEvent->GetCaloCluster(iclus)) )
1799  {
1800  if (clus->IsEMCAL())
1801  {
1802 
1803  Float_t frac =-1;
1804  Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(fEMCALCells, clus,frac);
1805  Double_t tof = clus->GetTOF();
1806  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absIdMax,fInputEvent->GetBunchCrossNumber(),tof);
1807  tof*=1e9;
1808 
1809  //printf("Input event cluster : AbsIdMax %d, E %2.2f, time %2.2f \n", absIdMax,clus->E(),tof);
1810 
1811  //Reject clusters with bad channels, close to borders and exotic;
1812  if(!GetCaloUtils()->GetEMCALRecoUtils()->IsGoodCluster(clus,GetCaloUtils()->GetEMCALGeometry(),GetEMCALCells(),fInputEvent->GetBunchCrossNumber())) continue;
1813 
1814  Int_t bc = TMath::Nint(tof/50) + 9;
1815  SetEMCalEventBC(bc);
1816 
1817  if(fEMCALPtMin > clus->E() || fEMCALPtMax < clus->E()) continue ;
1818 
1819  clus->GetMomentum(fMomentum, fVertex[0]);
1820 
1821  if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(fMomentum.Eta(),fMomentum.Phi(),kEMCAL)) return ;
1822 
1823  SetEMCalEventBCcut(bc);
1824 
1825  if(!IsInTimeWindow(tof,clus->E()))
1826  fNPileUpClusters++ ;
1827  else
1829 
1830  }
1831  }
1832  }
1833 
1834  //printf("Input event : Pile-up clusters %d, NO pile-up %d\n",fNPileUpClusters,fNNonPileUpClusters);
1835 
1836  // Recalculate track matching, not necessary if already done in the reclusterization task.
1837  // in case it was not done ...
1839 
1840  }
1841 
1842  AliDebug(1,Form("AOD entries %d, n pile-up clusters %d, n non pile-up %d", fEMCALClusters->GetEntriesFast(),fNPileUpClusters,fNNonPileUpClusters));
1843 }
1844 
1845 //_______________________________________
1847 //_______________________________________
1849 {
1850  AliDebug(1,"Begin");
1851 
1852  //Loop to select clusters in fiducial cut and fill container with aodClusters
1853  Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
1854  for (Int_t iclus = 0; iclus < nclusters; iclus++)
1855  {
1856  AliVCluster * clus = 0;
1857  if ( (clus = fInputEvent->GetCaloCluster(iclus)) )
1858  {
1859  if (clus->IsPHOS())
1860  {
1861  if(fAcceptOnlyHIJINGLabels && !IsHIJINGLabel(clus->GetLabel())) continue ;
1862 
1863  //Check if the cluster contains any bad channel and if close to calorimeter borders
1864  Int_t vindex = 0 ;
1865  if (fMixedEvent)
1866  vindex = fMixedEvent->EventIndexForCaloCluster(iclus);
1867  if( GetCaloUtils()->ClusterContainsBadChannel(kPHOS,clus->GetCellsAbsId(), clus->GetNCells()))
1868  continue;
1869  if(!GetCaloUtils()->CheckCellFiducialRegion(clus, fInputEvent->GetPHOSCells()))
1870  continue;
1871 
1873  {
1874  //Recalibrate the cluster energy
1876  {
1877  Float_t energy = GetCaloUtils()->RecalibrateClusterEnergy(clus, (AliAODCaloCells*)GetPHOSCells());
1878  clus->SetE(energy);
1879  }
1880  }
1881 
1882  clus->GetMomentum(fMomentum, fVertex[vindex]);
1883 
1884  if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(fMomentum.Eta(),fMomentum.Phi(),kPHOS)) continue;
1885 
1886  if(fPHOSPtMin > fMomentum.E() || fPHOSPtMax < fMomentum.E()) continue;
1887 
1888  //if(fDebug > 2 && fMomentum.E() > 0.1)
1889  AliDebug(2,Form("Selected clusters E %3.2f, pt %3.2f, phi %3.2f, eta %3.2f",
1890  fMomentum.E(),fMomentum.Pt(),fMomentum.Phi()*TMath::RadToDeg(),fMomentum.Eta()));
1891 
1892  if (fMixedEvent)
1893  {
1894  clus->SetID(iclus) ;
1895  }
1896 
1897  fPHOSClusters->Add(clus);
1898 
1899  }//PHOS cluster
1900  }//cluster exists
1901  }//esd cluster loop
1902 
1903  AliDebug(1,Form("AOD entries %d",fPHOSClusters->GetEntriesFast()));
1904 
1905 }
1906 
1907 //____________________________________________
1909 //____________________________________________
1911 {
1912  fEMCALCells = fInputEvent->GetEMCALCells();
1913 }
1914 
1915 //___________________________________________
1917 //___________________________________________
1919 {
1920  fPHOSCells = fInputEvent->GetPHOSCells();
1921 }
1922 
1923 //_______________________________________
1926 //_______________________________________
1928 {
1929  AliVVZERO* v0 = fInputEvent->GetVZEROData();
1930  //printf("Init V0: ADC (%d,%d), Multiplicity (%d,%d) \n",fV0ADC[0],fV0ADC[1],fV0Mul[0],fV0Mul[1]);
1931 
1932  if (v0)
1933  {
1934  AliESDVZERO* esdV0 = dynamic_cast<AliESDVZERO*> (v0);
1935  for (Int_t i = 0; i < 32; i++)
1936  {
1937  if(esdV0)
1938  {//Only available in ESDs
1939  fV0ADC[0] += (Int_t)esdV0->GetAdcV0C(i);
1940  fV0ADC[1] += (Int_t)esdV0->GetAdcV0A(i);
1941  }
1942 
1943  fV0Mul[0] += (Int_t)v0->GetMultiplicityV0C(i);
1944  fV0Mul[1] += (Int_t)v0->GetMultiplicityV0A(i);
1945  }
1946 
1947  AliDebug(1,Form("ADC (%d,%d), Multiplicity (%d,%d)",fV0ADC[0],fV0ADC[1],fV0Mul[0],fV0Mul[1]));
1948  }
1949  else
1950  {
1951  AliDebug(1,"Cannot retrieve V0 ESD! Run w/ null V0 charges");
1952  }
1953 }
1954 
1955 //_________________________________________________
1959 //_________________________________________________
1961 {
1962  AliDebug(2,"Begin");
1963 
1964  //
1965  //check if branch name is given
1966  if(!fInputNonStandardJetBranchName.Length())
1967  {
1968  fInputEvent->Print();
1969  AliFatal("No non-standard jet branch name specified. Specify among existing ones.");
1970  }
1971 
1972  fNonStandardJets = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(fInputNonStandardJetBranchName.Data()));
1973 
1974  if(!fNonStandardJets)
1975  {
1976  //check if jet branch exist; exit if not
1977  fInputEvent->Print();
1978 
1979  AliFatal(Form("%s:%d no reconstructed jet array with name %s in AOD", (char*)__FILE__,__LINE__,fInputNonStandardJetBranchName.Data()));
1980  }
1981  else
1982  {
1983  AliDebug(1,Form("AOD input jets %d", fNonStandardJets->GetEntriesFast()));
1984  }
1985 }
1986 
1987 //_________________________________________________
1991 //_________________________________________________
1993 {
1994  AliDebug(1,"Begin");
1995  //
1996  //check if branch name is given
1997  if(!fInputBackgroundJetBranchName.Length())
1998  {
1999  fInputEvent->Print();
2000 
2001  AliFatal("No background jet branch name specified. Specify among existing ones.");
2002  }
2003 
2004  fBackgroundJets = (AliAODJetEventBackground*)(fInputEvent->FindListObject(fInputBackgroundJetBranchName.Data()));
2005 
2006  if(!fBackgroundJets)
2007  {
2008  //check if jet branch exist; exit if not
2009  fInputEvent->Print();
2010 
2011  AliFatal(Form("%s:%d no reconstructed jet array with name %s in AOD", (char*)__FILE__,__LINE__,fInputBackgroundJetBranchName.Data()));
2012  }
2013  else
2014  {
2015  AliDebug(1,"FillInputBackgroundJets");
2016  fBackgroundJets->Print("");
2017  }
2018 }
2019 
2020 //____________________________________________________________________
2022 //____________________________________________________________________
2023 TArrayI AliCaloTrackReader::GetTriggerPatches(Int_t tmin, Int_t tmax )
2024 {
2025  // init some variables
2026  Int_t trigtimes[30], globCol, globRow,ntimes, i;
2027  Int_t absId = -1; //[100];
2028  Int_t nPatch = 0;
2029 
2030  TArrayI patches(0);
2031 
2032  // get object pointer
2033  AliVCaloTrigger *caloTrigger = GetInputEvent()->GetCaloTrigger( "EMCAL" );
2034 
2035  //printf("CaloTrigger Entries %d\n",caloTrigger->GetEntries() );
2036 
2037  // class is not empty
2038  if( caloTrigger->GetEntries() > 0 )
2039  {
2040  // must reset before usage, or the class will fail
2041  caloTrigger->Reset();
2042 
2043  // go throuth the trigger channels
2044  while( caloTrigger->Next() )
2045  {
2046  // get position in global 2x2 tower coordinates
2047  caloTrigger->GetPosition( globCol, globRow );
2048 
2049  //L0
2050  if(IsEventEMCALL0())
2051  {
2052  // get dimension of time arrays
2053  caloTrigger->GetNL0Times( ntimes );
2054 
2055  // no L0s in this channel
2056  // presence of the channel in the iterator still does not guarantee that L0 was produced!!
2057  if( ntimes < 1 )
2058  continue;
2059 
2060  // get timing array
2061  caloTrigger->GetL0Times( trigtimes );
2062  //printf("Get L0 patch : n times %d - trigger time window %d - %d\n",ntimes, tmin,tmax);
2063 
2064  // go through the array
2065  for( i = 0; i < ntimes; i++ )
2066  {
2067  // check if in cut - 8,9 shall be accepted in 2011
2068  if( trigtimes[i] >= tmin && trigtimes[i] <= tmax )
2069  {
2070  //printf("Accepted trigger time %d \n",trigtimes[i]);
2071  //if(nTrig > 99) continue;
2072  GetCaloUtils()->GetEMCALGeometry()->GetAbsFastORIndexFromPositionInEMCAL(globCol,globRow, absId);
2073  //printf("pass the time cut globCol %d, globRow %d absId %d\n",globCol,globRow, absId);
2074  patches.Set(nPatch+1);
2075  patches.AddAt(absId,nPatch++);
2076  }
2077  } // trigger time array
2078  }//L0
2079  else if(IsEventEMCALL1()) // L1
2080  {
2081  Int_t bit = 0;
2082  caloTrigger->GetTriggerBits(bit);
2083 
2084  Int_t sum = 0;
2085  caloTrigger->GetL1TimeSum(sum);
2086  //fBitEGA-=2;
2087  Bool_t isEGA1 = ((bit >> fBitEGA ) & 0x1) && IsEventEMCALL1Gamma1() ;
2088  Bool_t isEGA2 = ((bit >> (fBitEGA+1)) & 0x1) && IsEventEMCALL1Gamma2() ;
2089  Bool_t isEJE1 = ((bit >> fBitEJE ) & 0x1) && IsEventEMCALL1Jet1 () ;
2090  Bool_t isEJE2 = ((bit >> (fBitEJE+1)) & 0x1) && IsEventEMCALL1Jet2 () ;
2091 
2092  //if((bit>> fBitEGA )&0x1) printf("Trig Bit %d - bit %d - EG1 %d - EG2 %d\n",fBitEGA ,bit,IsEventEMCALL1Gamma1(),IsEventEMCALL1Gamma2());
2093  //if((bit>>(fBitEGA+1))&0x1) printf("Trig Bit %d - bit %d - EG1 %d - EG2 %d\n",fBitEGA+1,bit,IsEventEMCALL1Gamma1(),IsEventEMCALL1Gamma2());
2094 
2095  if(!isEGA1 && !isEJE1 && !isEGA2 && !isEJE2) continue;
2096 
2097  Int_t patchsize = -1;
2098  if (isEGA1 || isEGA2) patchsize = 2;
2099  else if (isEJE1 || isEJE2) patchsize = 16;
2100 
2101  //printf("**** Get L1 Patch: Bit %x, sum %d, patchsize %d, EGA1 %d, EGA2 %d, EJE1 %d, EJE2 %d, EGA bit %d, EJE bit %d, Trigger Gamma %d, Trigger Jet %d\n",
2102  // bit,sum,patchsize,isEGA1,isEGA2,isEJE1,isEJE2,fBitEGA,fBitEJE,IsEventEMCALL1Gamma(),IsEventEMCALL1Jet());
2103 
2104 
2105  // add 2x2 (EGA) or 16x16 (EJE) patches
2106  for(Int_t irow=0; irow < patchsize; irow++)
2107  {
2108  for(Int_t icol=0; icol < patchsize; icol++)
2109  {
2110  GetCaloUtils()->GetEMCALGeometry()->GetAbsFastORIndexFromPositionInEMCAL(globCol+icol,globRow+irow, absId);
2111  //printf("pass the time cut globCol %d, globRow %d absId %d\n",globCol,globRow, absId);
2112  patches.Set(nPatch+1);
2113  patches.AddAt(absId,nPatch++);
2114  }
2115  }
2116 
2117  } // L1
2118 
2119  } // trigger iterator
2120  } // go through triggers
2121 
2122  if(patches.GetSize()<=0) AliInfo(Form("No patch found! for triggers: %s and selected <%s>",
2124  //else printf(">>>>> N patches %d, test %d,first %d, last %d\n",patches.GetSize(), nPatch, patches.At(0), patches.At(patches.GetSize()-1));
2125 
2126  return patches;
2127 }
2128 
2129 //____________________________________________________________
2133 //____________________________________________________________
2135 {
2136  // Init info from previous event
2137  fTriggerClusterIndex = -1;
2138  fTriggerClusterId = -1;
2139  fTriggerClusterBC = -10000;
2140  fIsExoticEvent = kFALSE;
2141  fIsBadCellEvent = kFALSE;
2142  fIsBadMaxCellEvent = kFALSE;
2143 
2144  fIsTriggerMatch = kFALSE;
2145  fIsTriggerMatchOpenCut[0] = kFALSE;
2146  fIsTriggerMatchOpenCut[1] = kFALSE;
2147  fIsTriggerMatchOpenCut[2] = kFALSE;
2148 
2149  // Do only analysis for triggered events
2150  if(!IsEventEMCALL1() && !IsEventEMCALL0())
2151  {
2152  fTriggerClusterBC = 0;
2153  return;
2154  }
2155 
2156  //printf("***** Try to match trigger to cluster %d **** L0 %d, L1 %d\n",fTriggerPatchClusterMatch,IsEventEMCALL0(),IsEventEMCALL1());
2157 
2158  //Recover the list of clusters
2159  TClonesArray * clusterList = 0;
2160  if (fInputEvent->FindListObject(fEMCALClustersListName))
2161  {
2162  clusterList = dynamic_cast<TClonesArray*> (fInputEvent->FindListObject(fEMCALClustersListName));
2163  }
2164  else if(fOutputEvent)
2165  {
2166  clusterList = dynamic_cast<TClonesArray*> (fOutputEvent->FindListObject(fEMCALClustersListName));
2167  }
2168 
2169  // Get number of clusters and of trigger patches
2170  Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
2171  if(clusterList)
2172  nclusters = clusterList->GetEntriesFast();
2173 
2174  Int_t nPatch = patches.GetSize();
2175  Float_t exoDiffTime = GetCaloUtils()->GetEMCALRecoUtils()->GetExoticCellDiffTimeCut();
2176 
2177  //Init some variables used in the cluster loop
2178  Float_t tofPatchMax = 100000;
2179  Float_t ePatchMax =-1;
2180 
2181  Float_t tofMax = 100000;
2182  Float_t eMax =-1;
2183 
2184  Int_t clusMax =-1;
2185  Int_t idclusMax =-1;
2186  Bool_t badClMax = kFALSE;
2187  Bool_t badCeMax = kFALSE;
2188  Bool_t exoMax = kFALSE;
2189 // Int_t absIdMaxTrig= -1;
2190  Int_t absIdMaxMax = -1;
2191 
2192  Int_t nOfHighECl = 0 ;
2193 
2194  //
2195  // Check what is the trigger threshold
2196  // set minimu energym of candidate for trigger cluster
2197  //
2199 
2200  Float_t triggerThreshold = fTriggerL1EventThreshold;
2201  if(IsEventEMCALL0()) triggerThreshold = fTriggerL0EventThreshold;
2202  Float_t minE = triggerThreshold / 2.;
2203 
2204  // This method is not really suitable for JET trigger
2205  // but in case, reduce the energy cut since we do not trigger on high energy particle
2206  if(IsEventEMCALL1Jet() || minE < 1) minE = 1;
2207 
2208  AliDebug(1,Form("IsL1Trigger %d, IsL1JetTrigger? %d, IsL0Trigger %d, L1 threshold %2.1f, L0 threshold %2.1f, Min cluster E %2.2f",IsEventEMCALL1Jet(), IsEventEMCALL1(), IsEventEMCALL0(), fTriggerL1EventThreshold,fTriggerL0EventThreshold,minE));
2209 
2210  //
2211  // Loop on the clusters, check if there is any that falls into one of the patches
2212  //
2213  for (Int_t iclus = 0; iclus < nclusters; iclus++)
2214  {
2215  AliVCluster * clus = 0;
2216  if(clusterList) clus = (AliVCluster*) clusterList->At(iclus);
2217  else clus = fInputEvent->GetCaloCluster(iclus);
2218 
2219  if ( !clus ) continue ;
2220 
2221  if ( !clus->IsEMCAL() ) continue ;
2222 
2223  //Skip clusters with too low energy to be triggering
2224  if ( clus->E() < minE ) continue ;
2225 
2226  Float_t frac = -1;
2227  Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(fInputEvent->GetEMCALCells(), clus,frac);
2228 
2229  Bool_t badCluster = GetCaloUtils()->GetEMCALRecoUtils()->ClusterContainsBadChannel(GetCaloUtils()->GetEMCALGeometry(),
2230  clus->GetCellsAbsId(),clus->GetNCells());
2231  UShort_t cellMax[] = {(UShort_t) absIdMax};
2232  Bool_t badCell = GetCaloUtils()->GetEMCALRecoUtils()->ClusterContainsBadChannel(GetCaloUtils()->GetEMCALGeometry(),cellMax,1);
2233 
2234  // if cell is bad, it can happen that time calibration is not available,
2235  // when calculating if it is exotic, this can make it to be exotic by default
2236  // open it temporarily for this cluster
2237  if(badCell)
2238  GetCaloUtils()->GetEMCALRecoUtils()->SetExoticCellDiffTimeCut(10000000);
2239 
2240  Bool_t exotic = GetCaloUtils()->GetEMCALRecoUtils()->IsExoticCluster(clus, fInputEvent->GetEMCALCells());
2241 
2242  // Set back the time cut on exotics
2243  if(badCell)
2244  GetCaloUtils()->GetEMCALRecoUtils()->SetExoticCellDiffTimeCut(exoDiffTime);
2245 
2246  // Energy threshold for exotic Ecross typically at 4 GeV,
2247  // for lower energy, check that there are more than 1 cell in the cluster
2248  if(!exotic && clus->GetNCells() < 2) exotic = kTRUE;
2249 
2250  Float_t energy = clus->E();
2251  Int_t idclus = clus->GetID();
2252 
2253  Double_t tof = clus->GetTOF();
2254  if(GetCaloUtils()->GetEMCALRecoUtils()->IsTimeRecalibrationOn() && fTriggerClusterTimeRecal)
2255  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absIdMax,fInputEvent->GetBunchCrossNumber(),tof);
2256  tof *=1.e9;
2257 
2258  //printf("cluster %d, ID %d, E %2.2f, tof %2.2f, AbsId max %d, exotic %d, bad Cluster %d, bad Cell %d\n",
2259  // iclus,idclus, energy,tof,absIdMax, exotic, badCluster,badCell);
2260 
2261  // Find the highest energy cluster, avobe trigger threshold
2262  // in the event in case no match to trigger is found
2263  if( energy > eMax )
2264  {
2265  tofMax = tof;
2266  eMax = energy;
2267  badClMax = badCluster;
2268  badCeMax = badCell;
2269  exoMax = exotic;
2270  clusMax = iclus;
2271  idclusMax = idclus;
2272  absIdMaxMax = absIdMax;
2273  }
2274 
2275  // count the good clusters in the event avobe the trigger threshold
2276  // to check the exotic events
2277  if(!badCluster && !exotic)
2278  nOfHighECl++;
2279 
2280  // Find match to trigger
2281  if(fTriggerPatchClusterMatch && nPatch>0)
2282  {
2283  for(Int_t iabsId =0; iabsId < nPatch; iabsId++)
2284  {
2285  Int_t absIDCell[4];
2286  GetCaloUtils()->GetEMCALGeometry()->GetCellIndexFromFastORIndex(patches.At(iabsId), absIDCell);
2287  //if(tof > 75 ) printf("E %2.2f TOF %2.2f Trigger patch %d, cells : %d, %d, %d, %d\n",
2288  // clus->E(),tof,patches.At(iabsId), absIDCell[0],absIDCell[1],absIDCell[2],absIDCell[3]);
2289 
2290  for(Int_t ipatch = 0; ipatch < 4; ipatch++)
2291  {
2292  if(absIdMax == absIDCell[ipatch])
2293  {
2294  //printf("*** Patches : absId %d, E %2.1f, tof %f \n",absIdMax,clus->E(), tof);
2295  if(energy > ePatchMax)
2296  {
2297  tofPatchMax = tof;
2298  ePatchMax = energy;
2299  fIsBadCellEvent = badCluster;
2300  fIsBadMaxCellEvent = badCell;
2301  fIsExoticEvent = exotic;
2302  fTriggerClusterIndex = iclus;
2303  fTriggerClusterId = idclus;
2304  fIsTriggerMatch = kTRUE;
2305 // absIdMaxTrig = absIdMax;
2306  }
2307  }
2308  }// cell patch loop
2309  }// trigger patch loop
2310  } // Do trigger patch matching
2311 
2312  }// Cluster loop
2313 
2314  // If there was no match, assign as trigger
2315  // the highest energy cluster in the event
2316  if(!fIsTriggerMatch)
2317  {
2318  tofPatchMax = tofMax;
2319  ePatchMax = eMax;
2320  fIsBadCellEvent = badClMax;
2321  fIsBadMaxCellEvent = badCeMax;
2322  fIsExoticEvent = exoMax;
2323  fTriggerClusterIndex = clusMax;
2324  fTriggerClusterId = idclusMax;
2325  }
2326 
2327  Double_t tofPatchMaxUS = TMath::Abs(tofPatchMax);
2328 
2329  if (tofPatchMaxUS < 28 ) fTriggerClusterBC = 0 ;
2330  else if(tofPatchMaxUS < 75 ) fTriggerClusterBC = 1 ;
2331  else if(tofPatchMaxUS < 125) fTriggerClusterBC = 2 ;
2332  else if(tofPatchMaxUS < 175) fTriggerClusterBC = 3 ;
2333  else if(tofPatchMaxUS < 225) fTriggerClusterBC = 4 ;
2334  else if(tofPatchMaxUS < 275) fTriggerClusterBC = 5 ;
2335  else
2336  {
2337  //printf("AliCaloTrackReader::MatchTriggerCluster() - Large BC - tof %2.3f - Index %d\n",tofPatchMaxUS,fTriggerClusterIndex);
2339  else
2340  {
2341  fTriggerClusterIndex = -2;
2342  fTriggerClusterId = -2;
2343  }
2344  }
2345 
2346  if(tofPatchMax < 0) fTriggerClusterBC*=-1;
2347 
2348 
2349  // printf("AliCaloTrackReader::MatchTriggerCluster(TArrayI patches) - Trigger cluster: index %d, ID %d, E = %2.2f, tof = %2.2f (BC = %d), bad cluster? %d, bad cell? %d, exotic? %d, patch match? %d, n High E cluster %d, absId Max %d\n",
2350  // fTriggerClusterIndex, fTriggerClusterId,ePatchMax, tofPatchMax,
2351  // fTriggerClusterBC, fIsBadCellEvent,fIsBadMaxCellEvent,fIsExoticEvent, fIsTriggerMatch, nOfHighECl,absIdMaxMax);
2352  //
2353  // if(!fIsTriggerMatch) printf("\t highest energy cluster: index %d, ID %d, E = %2.2f, tof = %2.2f, bad cluster? %d, bad cell? %d, exotic? %d\n",
2354  // clusMax, idclusMax, eMax,tofMax, badClMax, badCeMax,exoMax);
2355 
2356  //Redo matching but open cuts
2357  if(!fIsTriggerMatch && fTriggerClusterId >= 0)
2358  {
2359  // Open time patch time
2360  TArrayI patchOpen = GetTriggerPatches(7,10);
2361 
2362  Int_t patchAbsIdOpenTime = -1;
2363  for(Int_t iabsId =0; iabsId < patchOpen.GetSize(); iabsId++)
2364  {
2365  Int_t absIDCell[4];
2366  patchAbsIdOpenTime = patchOpen.At(iabsId);
2367  GetCaloUtils()->GetEMCALGeometry()->GetCellIndexFromFastORIndex(patchAbsIdOpenTime, absIDCell);
2368  //if(tof > 75 ) printf("E %2.2f TOF %2.2f Trigger patch %d, cells : %d, %d, %d, %d\n",
2369  // clus->E(),tof,patches.At(iabsId), absIDCell[0],absIDCell[1],absIDCell[2],absIDCell[3]);
2370 
2371  for(Int_t ipatch = 0; ipatch < 4; ipatch++)
2372  {
2373  if(absIdMaxMax == absIDCell[ipatch])
2374  {
2375  fIsTriggerMatchOpenCut[0] = kTRUE;
2376  break;
2377  }
2378  }// cell patch loop
2379  }// trigger patch loop
2380 
2381  // Check neighbour patches
2382  Int_t patchAbsId = -1;
2383  Int_t globalCol = -1;
2384  Int_t globalRow = -1;
2385  GetCaloUtils()->GetEMCALGeometry()->GetFastORIndexFromCellIndex(absIdMaxMax, patchAbsId);
2386  GetCaloUtils()->GetEMCALGeometry()->GetPositionInEMCALFromAbsFastORIndex(patchAbsId,globalCol,globalRow);
2387 
2388  // Check patches with strict time cut
2389  Int_t patchAbsIdNeigh = -1;
2390  for(Int_t icol = globalCol-1; icol <= globalCol+1; icol++)
2391  {
2392  if(icol < 0 || icol > 47) continue;
2393 
2394  for(Int_t irow = globalRow; irow <= globalRow+1; irow++)
2395  {
2396  if(irow < 0 || irow > 63) continue;
2397 
2398  GetCaloUtils()->GetEMCALGeometry()->GetAbsFastORIndexFromPositionInEMCAL(icol, irow, patchAbsIdNeigh);
2399 
2400  if ( patchAbsIdNeigh < 0 ) continue;
2401 
2402  for(Int_t iabsId =0; iabsId < patches.GetSize(); iabsId++)
2403  {
2404  if(patchAbsIdNeigh == patches.At(iabsId))
2405  {
2406  fIsTriggerMatchOpenCut[1] = kTRUE;
2407  break;
2408  }
2409  }// trigger patch loop
2410 
2411  }// row
2412  }// col
2413 
2414  // Check patches with open time cut
2415  Int_t patchAbsIdNeighOpenTime = -1;
2416  for(Int_t icol = globalCol-1; icol <= globalCol+1; icol++)
2417  {
2418  if(icol < 0 || icol > 47) continue;
2419 
2420  for(Int_t irow = globalRow; irow <= globalRow+1; irow++)
2421  {
2422  if(irow < 0 || irow > 63) continue;
2423 
2424  GetCaloUtils()->GetEMCALGeometry()->GetAbsFastORIndexFromPositionInEMCAL(icol, irow, patchAbsIdNeighOpenTime);
2425 
2426  if ( patchAbsIdNeighOpenTime < 0 ) continue;
2427 
2428  for(Int_t iabsId =0; iabsId < patchOpen.GetSize(); iabsId++)
2429  {
2430  if(patchAbsIdNeighOpenTime == patchOpen.At(iabsId))
2431  {
2432  fIsTriggerMatchOpenCut[2] = kTRUE;
2433  break;
2434  }
2435  }// trigger patch loop
2436 
2437  }// row
2438  }// col
2439 
2440  // printf("No match, new match: Open time %d-%d, open Neigh %d-%d, both open %d-%d\n",fIsTriggerMatchOpenCut[0],patchAbsIdOpenTime,
2441  // fIsTriggerMatchOpenCut[1],patchAbsIdNeigh,
2442  // fIsTriggerMatchOpenCut[2],patchAbsIdNeighOpenTime);
2443 
2444  patchOpen.Reset();
2445 
2446  }// No trigger match found
2447  //printf("Trigger BC %d, Id %d, Index %d\n",fTriggerClusterBC,fTriggerClusterId,fTriggerClusterIndex);
2448 }
2449 
2450 //_________________________________________________________
2454 //_________________________________________________________
2456 {
2458  {
2459  // get object pointer
2460  AliVCaloTrigger *caloTrigger = GetInputEvent()->GetCaloTrigger( "EMCAL" );
2461 
2462  if ( fBitEGA == 6 )
2463  {
2464  if (IsEventEMCALL1Gamma1()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(1);
2465  else if(IsEventEMCALL1Gamma2()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(3);
2466  else if(IsEventEMCALL1Jet1 ()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(0);
2467  else if(IsEventEMCALL1Jet2 ()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(2);
2468 
2469  // printf("L1 trigger Threshold Jet1 %f, Gamma1 %f, Jet2 %f, Gamma2 %f\n",
2470  // 0.07874*caloTrigger->GetL1Threshold(0),
2471  // 0.07874*caloTrigger->GetL1Threshold(1),
2472  // 0.07874*caloTrigger->GetL1Threshold(2),
2473  // 0.07874*caloTrigger->GetL1Threshold(3));
2474  }
2475  else
2476  {
2477  // Old AOD data format, in such case, order of thresholds different!!!
2478  if (IsEventEMCALL1Gamma1()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(0);
2479  else if(IsEventEMCALL1Gamma2()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(2);
2480  else if(IsEventEMCALL1Jet1 ()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(1);
2481  else if(IsEventEMCALL1Jet2 ()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(3);
2482 
2483  // printf("L1 trigger Threshold Jet1 %f, Gamma1 %f, Jet2 %f, Gamma2 %f\n",
2484  // 0.07874*caloTrigger->GetL1Threshold(1),
2485  // 0.07874*caloTrigger->GetL1Threshold(0),
2486  // 0.07874*caloTrigger->GetL1Threshold(3),
2487  // 0.07874*caloTrigger->GetL1Threshold(2));
2488  }
2489  }
2490 
2491  // Set L0 threshold, if not set by user
2493  {
2494  // Revise for periods > LHC11d
2495  Int_t runNumber = fInputEvent->GetRunNumber();
2496  if (runNumber < 146861) fTriggerL0EventThreshold = 3. ;
2497  else if(runNumber < 154000) fTriggerL0EventThreshold = 4. ;
2498  else if(runNumber < 165000) fTriggerL0EventThreshold = 5.5;
2499  else fTriggerL0EventThreshold = 2 ;
2500  }
2501 }
2502 
2503 //________________________________________________________
2505 //________________________________________________________
2506 void AliCaloTrackReader::Print(const Option_t * opt) const
2507 {
2508  if(! opt)
2509  return;
2510 
2511  printf("***** Print: %s %s ******\n", GetName(), GetTitle() ) ;
2512  printf("Task name : %s\n", fTaskName.Data()) ;
2513  printf("Data type : %d\n", fDataType) ;
2514  printf("CTS Min pT : %2.1f GeV/c\n", fCTSPtMin) ;
2515  printf("EMCAL Min pT : %2.1f GeV/c\n", fEMCALPtMin) ;
2516  printf("PHOS Min pT : %2.1f GeV/c\n", fPHOSPtMin) ;
2517  printf("CTS Max pT : %2.1f GeV/c\n", fCTSPtMax) ;
2518  printf("EMCAL Max pT : %2.1f GeV/c\n", fEMCALPtMax) ;
2519  printf("PHOS Max pT : %2.1f GeV/c\n", fPHOSPtMax) ;
2520  printf("EMCAL Time Cut: %3.1f < TOF < %3.1f\n", fEMCALTimeCutMin, fEMCALTimeCutMax);
2521  printf("Use CTS = %d\n", fFillCTS) ;
2522  printf("Use EMCAL = %d\n", fFillEMCAL) ;
2523  printf("Use DCAL = %d\n", fFillDCAL) ;
2524  printf("Use PHOS = %d\n", fFillPHOS) ;
2525  printf("Use EMCAL Cells = %d\n", fFillEMCALCells) ;
2526  printf("Use PHOS Cells = %d\n", fFillPHOSCells) ;
2527  printf("Track status = %d\n", (Int_t) fTrackStatus) ;
2528  //printf("AODs Track filter mask = %d or hybrid %d (if filter bit comp %d), select : SPD hit %d, primary %d\n",
2529  // (Int_t) fTrackFilterMask, fSelectHybridTracks, (Int_t) fTrackFilterMaskComplementary, fSelectSPDHitTracks,fSelectPrimaryTracks) ;
2530  printf("Track Mult Eta Cut = %d\n", (Int_t) fTrackMultEtaCut) ;
2531  printf("Write delta AOD = %d\n", fWriteOutputDeltaAOD) ;
2532  printf("Recalculate Clusters = %d, E linearity = %d\n", fRecalculateClusters, fCorrectELinearity) ;
2533 
2534  printf("Use Triggers selected in SE base class %d; If not what Trigger Mask? %d; MB Trigger Mask for mixed %d \n",
2536 
2538  printf("Compare jet pt and pt hard to accept event, factor = %2.2f",fPtHardAndJetPtFactor);
2539 
2541  printf("Compare cluster pt and pt hard to accept event, factor = %2.2f",fPtHardAndClusterPtFactor);
2542 
2543  printf("Read Kine from, stack? %d, AOD ? %d \n", fReadStack, fReadAODMCParticles) ;
2544  printf("Delta AOD File Name = %s\n", fDeltaAODFileName.Data()) ;
2545  printf("Centrality: Class %s, Option %d, Bin [%d,%d] \n", fCentralityClass.Data(),fCentralityOpt,fCentralityBin[0], fCentralityBin[1]) ;
2546 
2547  printf(" \n") ;
2548 }
2549 
2550 //__________________________________________
2554 //__________________________________________
2556 {
2557  // Count number of cells with energy larger than 0.1 in SM3, cut on this number
2558  Int_t ncellsSM3 = 0;
2559  for(Int_t icell = 0; icell < fInputEvent->GetEMCALCells()->GetNumberOfCells(); icell++)
2560  {
2561  Int_t absID = fInputEvent->GetEMCALCells()->GetCellNumber(icell);
2562  Int_t sm = GetCaloUtils()->GetEMCALGeometry()->GetSuperModuleNumber(absID);
2563  if(fInputEvent->GetEMCALCells()->GetAmplitude(icell) > 0.1 && sm==3) ncellsSM3++;
2564  }
2565 
2566  Int_t ncellcut = 21;
2567  if(GetFiredTriggerClasses().Contains("EMC")) ncellcut = 35;
2568 
2569  if(ncellsSM3 >= ncellcut)
2570  {
2571  AliDebug(1,Form("Reject event with ncells in SM3 %d, cut %d, trig %s",
2572  ncellsSM3,ncellcut,GetFiredTriggerClasses().Data()));
2573  return kTRUE;
2574  }
2575 
2576  return kFALSE;
2577 }
2578 
2579 //_________________________________________________________
2583 //_________________________________________________________
2585 {
2586  if(label < 0) return ;
2587 
2588  AliAODEvent * evt = dynamic_cast<AliAODEvent*> (fInputEvent) ;
2589  if(!evt) return ;
2590 
2591  TClonesArray * arr = dynamic_cast<TClonesArray*>(evt->FindListObject("mcparticles")) ;
2592  if(!arr) return ;
2593 
2594  if(label < arr->GetEntriesFast())
2595  {
2596  AliAODMCParticle * particle = dynamic_cast<AliAODMCParticle *>(arr->At(label));
2597  if(!particle) return ;
2598 
2599  if(label == particle->Label()) return ; // label already OK
2600  //else printf("AliCaloTrackReader::RemapMCLabelForAODs() - Label %d - AOD stack %d \n",label, particle->Label());
2601  }
2602  //else printf("AliCaloTrackReader::RemapMCLabelForAODs() - Label %d > AOD labels %d \n",label, arr->GetEntriesFast());
2603 
2604  // loop on the particles list and check if there is one with the same label
2605  for(Int_t ind = 0; ind < arr->GetEntriesFast(); ind++ )
2606  {
2607  AliAODMCParticle * particle = dynamic_cast<AliAODMCParticle *>(arr->At(ind));
2608  if(!particle) continue ;
2609 
2610  if(label == particle->Label())
2611  {
2612  label = ind;
2613  //printf("AliAnalysisTaskEMCALClusterize::RemapMCLabelForAODs() - New Label Index %d \n",label);
2614  return;
2615  }
2616  }
2617 
2618  label = -1;
2619 
2620  //printf("AliCaloTrackReader::RemapMCLabelForAODs() - Label not found set to -1 \n");
2621 }
2622 
2623 //___________________________________
2625 //___________________________________
2627 {
2628  if(fCTSTracks) fCTSTracks -> Clear();
2629  if(fEMCALClusters) fEMCALClusters -> Clear("C");
2630  if(fPHOSClusters) fPHOSClusters -> Clear("C");
2631 
2632  fV0ADC[0] = 0; fV0ADC[1] = 0;
2633  fV0Mul[0] = 0; fV0Mul[1] = 0;
2634 
2635  if(fNonStandardJets) fNonStandardJets -> Clear("C");
2636  fBackgroundJets->Reset();
2637 }
2638 
2639 //___________________________________________
2643 //___________________________________________
2645 {
2646  fEventTrigMinBias = kFALSE;
2647  fEventTrigCentral = kFALSE;
2648  fEventTrigSemiCentral = kFALSE;
2649  fEventTrigEMCALL0 = kFALSE;
2650  fEventTrigEMCALL1Gamma1 = kFALSE;
2651  fEventTrigEMCALL1Gamma2 = kFALSE;
2652  fEventTrigEMCALL1Jet1 = kFALSE;
2653  fEventTrigEMCALL1Jet2 = kFALSE;
2654 
2655  AliDebug(1,Form("Select trigger mask bit %d - Trigger Event %s - Select <%s>",
2657 
2658  if(fEventTriggerMask <=0 )// in case no mask set
2659  {
2660  // EMC triggered event? Which type?
2661  if( GetFiredTriggerClasses().Contains("-B-") || GetFiredTriggerClasses().Contains("-S-") || GetFiredTriggerClasses().Contains("-I-") )
2662  {
2663  if ( GetFiredTriggerClasses().Contains("EGA" ) ||
2664  GetFiredTriggerClasses().Contains("EG1" ) )
2665  {
2666  fEventTrigEMCALL1Gamma1 = kTRUE;
2667  if( GetFiredTriggerClasses().Contains("EG1" ) && !fFiredTriggerClassName.Contains("EG1") ) fEventTrigEMCALL1Gamma1 = kFALSE;
2668  }
2669  else if( GetFiredTriggerClasses().Contains("EG2" ) )
2670  {
2671  fEventTrigEMCALL1Gamma2 = kTRUE;
2672  if( !fFiredTriggerClassName.Contains("EG2") ) fEventTrigEMCALL1Gamma2 = kFALSE;
2673  }
2674  else if( GetFiredTriggerClasses().Contains("EJE" ) ||
2675  GetFiredTriggerClasses().Contains("EJ1" ) )
2676  {
2677  fEventTrigEMCALL1Jet1 = kTRUE;
2678  if( GetFiredTriggerClasses().Contains("EJ1" ) && !fFiredTriggerClassName.Contains("EJ1") )
2679  fEventTrigEMCALL1Jet1 = kFALSE;
2680  }
2681  else if( GetFiredTriggerClasses().Contains("EJ2" ) )
2682  {
2683  fEventTrigEMCALL1Jet2 = kTRUE;
2684  if( !fFiredTriggerClassName.Contains("EJ2") ) fEventTrigEMCALL1Jet2 = kFALSE;
2685  }
2686  else if( GetFiredTriggerClasses().Contains("CEMC") &&
2687  !GetFiredTriggerClasses().Contains("EGA" ) &&
2688  !GetFiredTriggerClasses().Contains("EJE" ) &&
2689  !GetFiredTriggerClasses().Contains("EG1" ) &&
2690  !GetFiredTriggerClasses().Contains("EJ1" ) &&
2691  !GetFiredTriggerClasses().Contains("EG2" ) &&
2692  !GetFiredTriggerClasses().Contains("EJ2" ) )
2693  {
2694  fEventTrigEMCALL0 = kTRUE;
2695  }
2696 
2697  //Min bias event trigger?
2698  if (GetFiredTriggerClasses().Contains("CCENT_R2-B-NOPF-ALLNOTRD"))
2699  {
2700  fEventTrigCentral = kTRUE;
2701  }
2702  else if(GetFiredTriggerClasses().Contains("CSEMI_R1-B-NOPF-ALLNOTRD"))
2703  {
2704  fEventTrigSemiCentral = kTRUE;
2705  }
2706  else if((GetFiredTriggerClasses().Contains("CINT") || GetFiredTriggerClasses().Contains("CPBI2_B1") ) &&
2707  GetFiredTriggerClasses().Contains("-NOPF-ALLNOTRD") )
2708  {
2709  fEventTrigMinBias = kTRUE;
2710  }
2711  }
2712  }
2713  else
2714  {
2715  // EMC L1 Gamma
2716  if ( fEventTriggerMask & AliVEvent::kEMCEGA )
2717  {
2718  //printf("EGA trigger bit\n");
2719  if (GetFiredTriggerClasses().Contains("EG"))
2720  {
2721  if (GetFiredTriggerClasses().Contains("EGA")) fEventTrigEMCALL1Gamma1 = kTRUE;
2722  else
2723  {
2724  if(GetFiredTriggerClasses().Contains("EG1")) fEventTrigEMCALL1Gamma1 = kTRUE;
2725  if(GetFiredTriggerClasses().Contains("EG2")) fEventTrigEMCALL1Gamma2 = kTRUE;
2726  }
2727  }
2728  }
2729  // EMC L1 Jet
2730  else if( fEventTriggerMask & AliVEvent::kEMCEJE )
2731  {
2732  //printf("EGA trigger bit\n");
2733  if (GetFiredTriggerClasses().Contains("EJ"))
2734  {
2735  if (GetFiredTriggerClasses().Contains("EJE")) fEventTrigEMCALL1Jet1 = kTRUE;
2736  else
2737  {
2738  if(GetFiredTriggerClasses().Contains("EJ1")) fEventTrigEMCALL1Jet1 = kTRUE;
2739  if(GetFiredTriggerClasses().Contains("EJ2")) fEventTrigEMCALL1Jet2 = kTRUE;
2740  }
2741  }
2742  }
2743  // EMC L0
2744  else if((fEventTriggerMask & AliVEvent::kEMC7) ||
2745  (fEventTriggerMask & AliVEvent::kEMC1) )
2746  {
2747  //printf("L0 trigger bit\n");
2748  fEventTrigEMCALL0 = kTRUE;
2749  }
2750  // Min Bias Pb-Pb
2751  else if( fEventTriggerMask & AliVEvent::kCentral )
2752  {
2753  //printf("MB semi central trigger bit\n");
2754  fEventTrigSemiCentral = kTRUE;
2755  }
2756  // Min Bias Pb-Pb
2757  else if( fEventTriggerMask & AliVEvent::kSemiCentral )
2758  {
2759  //printf("MB central trigger bit\n");
2760  fEventTrigCentral = kTRUE;
2761  }
2762  // Min Bias pp, PbPb, pPb
2763  else if((fEventTriggerMask & AliVEvent::kMB ) ||
2764  (fEventTriggerMask & AliVEvent::kINT7) ||
2765  (fEventTriggerMask & AliVEvent::kINT8) ||
2766  (fEventTriggerMask & AliVEvent::kAnyINT) )
2767  {
2768  //printf("MB trigger bit\n");
2769  fEventTrigMinBias = kTRUE;
2770  }
2771  }
2772 
2773  AliDebug(1,Form("Event bits: \n \t MB %d, Cen %d, Sem %d, L0 %d, L1G1 %d, L1G2 %d, L1J1 %d, L1J2 %d",
2777 
2778  // L1 trigger bit
2779  if( fBitEGA == 0 && fBitEJE == 0 )
2780  {
2781  // Init the trigger bit once, correct depending on AliESD(AOD)CaloTrigger header version
2782 
2783  // Simpler way to do it ...
2784 // if( fInputEvent->GetRunNumber() < 172000 )
2785 // reader->SetEventTriggerL1Bit(4,5); // current LHC11 data
2786 // else
2787 // reader->SetEventTriggerL1Bit(6,8); // LHC12-13 data
2788 
2789  // Old values
2790  fBitEGA = 4;
2791  fBitEJE = 5;
2792 
2793  TFile* file = AliAnalysisManager::GetAnalysisManager()->GetTree()->GetCurrentFile();
2794 
2795  const TList *clist = file->GetStreamerInfoCache();
2796 
2797  if(clist)
2798  {
2799  TStreamerInfo *cinfo = (TStreamerInfo*)clist->FindObject("AliESDCaloTrigger");
2800  Int_t verid = 5; // newer ESD header version
2801  if(!cinfo)
2802  {
2803  cinfo = (TStreamerInfo*)clist->FindObject("AliAODCaloTrigger");
2804  verid = 3; // newer AOD header version
2805  }
2806 
2807  if(cinfo)
2808  {
2809  Int_t classversionid = cinfo->GetClassVersion();
2810  //printf("********* Header class version %d *********** \n",classversionid);
2811 
2812  if (classversionid >= verid)
2813  {
2814  fBitEGA = 6;
2815  fBitEJE = 8;
2816  }
2817  } else AliInfo("AliCaloTrackReader()::SetEventTriggerBit() - Streamer info for trigger class not available, bit not changed");
2818  } else AliInfo("AliCaloTrackReader::SetEventTriggerBit() - Streamer list not available!, bit not changed");
2819 
2820  } // set once the EJE, EGA trigger bit
2821 }
2822 
2823 //____________________________________________________________
2826 //____________________________________________________________
2827 void AliCaloTrackReader::SetInputEvent(AliVEvent* const input)
2828 {
2829  fInputEvent = input;
2830  fMixedEvent = dynamic_cast<AliMixedEvent*>(GetInputEvent()) ;
2831  if (fMixedEvent)
2832  fNMixedEvent = fMixedEvent->GetNumberOfEvents() ;
2833 
2834  //Delete previous vertex
2835  if(fVertex)
2836  {
2837  for (Int_t i = 0; i < fNMixedEvent; i++)
2838  {
2839  delete [] fVertex[i] ;
2840  }
2841  delete [] fVertex ;
2842  }
2843 
2844  fVertex = new Double_t*[fNMixedEvent] ;
2845  for (Int_t i = 0; i < fNMixedEvent; i++)
2846  {
2847  fVertex[i] = new Double_t[3] ;
2848  fVertex[i][0] = 0.0 ;
2849  fVertex[i][1] = 0.0 ;
2850  fVertex[i][2] = 0.0 ;
2851  }
2852 }
2853 
2854 
Bool_t IsPileUpFromSPD() const
Bool_t fUseTrackDCACut
Do DCA selection.
TArrayI GetTriggerPatches(Int_t tmin, Int_t tmax)
Recover the patches that triggered, either L0 or L1.
ClassImp(AliAnalysisTaskTriggerRates) AliAnalysisTaskTriggerRates
virtual void FillInputVZERO()
Int_t fV0ADC[2]
Integrated V0 signal.
Bool_t fComparePtHardAndClusterPt
In MonteCarlo, jet events, reject events with too large cluster energy.
AliCalorimeterUtils * GetCaloUtils() const
Float_t fTimeStampEventFracMin
Minimum value of time stamp fraction event.
Bool_t fReadAODMCParticles
Access kine information from filtered AOD MC particles.
virtual void FillInputNonStandardJets()
Double_t fEMCALTimeCutMax
Remove clusters/cells with time larger than this value, in ns.
Int_t fBitEJE
Trigger bit on VCaloTrigger for EJE.
TLorentzVector fMomentum
! Temporal TLorentzVector container, avoid declaration of TLorentzVectors per event.
Double_t fTimeStampRunMin
Minimum value of time stamp in run.
Bool_t fDoPileUpEventRejection
Select pile-up events by SPD.
Bool_t fAcceptOnlyHIJINGLabels
Select clusters or tracks that where generated by HIJING, reject other generators in case of cocktail...
TObjArray * fPHOSClusters
Temporal array with PHOS CaloClusters.
Double_t fTrackDCACut[3]
Remove tracks with DCA larger than cut, parameters of function stored here.
Bool_t fUseEventsWithPrimaryVertex
Select events with primary vertex.
virtual AliVCaloCells * GetPHOSCells() const
Bool_t fIsBadCellEvent
Bad cell triggered event flag, any cell in cluster is bad.
Bool_t IsCorrectionOfClusterEnergyOn() const
Bool_t fEventTrigEMCALL1Gamma1
Event is L1-Gamma, threshold 1 on its name, it should correspond kEMCEGA.
virtual AliHeader * GetHeader() const
Bool_t IsEventEMCALL1Jet1() const
Bool_t IsEventEMCALL1Gamma1() const
AliEMCALRecoUtils * GetEMCALRecoUtils() const
Bool_t ReadAODMCParticles() const
TString fEventPlaneMethod
Name of event plane method, by default "Q".
Int_t fTrackBCEventCut[19]
Fill one entry per event if there is a track in a given BC, depend on track pT, acceptance cut...
Bool_t fDoVertexBCEventSelection
Select events with vertex on BC=0 or -100.
AliMixedEvent * fMixedEvent
! Mixed event object. This class is not the owner.
virtual AliVEvent * GetInputEvent() const
Bool_t IsPileUpFromSPDAndNotEMCal() const
Check if event is from pile-up determined by SPD and not by EMCal.
Bool_t IsEventEMCALL1() const
Float_t fPtHardAndClusterPtFactor
Factor between ptHard and cluster pT to reject/accept event.
Bool_t fAcceptFastCluster
Accept events from fast cluster, exclude these events for LHC11a.
Bool_t fSmearShowerShape
Smear shower shape (use in MC).
AliVEvent * fInputEvent
! pointer to esd or aod input.
virtual void SetInputEvent(AliVEvent *input)
Double_t fTrackTimeCutMax
Remove tracks with time larger than this value, in ns.
Int_t fTriggerClusterBC
Event triggered by a cluster in BC -5 0 to 5.
TString fEMCALClustersListName
Alternative list of clusters produced elsewhere and not from InputEvent.
void RecalculateClusterTrackMatching(AliVEvent *event, TObjArray *clusterArray=0x0)
Recalculate track matching.
Bool_t fWriteOutputDeltaAOD
Write the created delta AOD objects into file.
Int_t fV0Mul[2]
Integrated V0 Multiplicity.
void RecalculateClusterDistanceToBadChannel(AliVCaloCells *cells, AliVCluster *clu)
Float_t fEMCALPtMin
pT Threshold on emcal clusters.
Bool_t fTriggerClusterTimeRecal
In case cluster already calibrated, do not try to recalibrate even if recalib on in AliEMCALRecoUtils...
Double_t ** fVertex
! Vertex array 3 dim for each mixed event buffer.
UInt_t fEventTriggerMask
Select this triggerered event.
virtual Int_t GetEventCentrality() const
virtual Bool_t IsHIJINGLabel(Int_t label)
virtual AliGenEventHeader * GetGenEventHeader() const
virtual AliCentrality * GetCentrality() const
AliVCaloCells * fPHOSCells
! Temporal array with PHOS AliVCaloCells.
Bool_t IsEventEMCALL1Gamma2() const
virtual AliMixedEvent * GetMixedEvent() const
Bool_t fFillInputBackgroundJetBranch
Flag to use data from background jets.
void SetTrackEventBC(Int_t bc)
virtual Bool_t ComparePtHardAndClusterPt()
TList * fAODBranchList
List with AOD branches created and needed in analysis.
void RemapMCLabelForAODs(Int_t &label)
Bool_t fEventTrigSemiCentral
Event is AliVEvent::kSemiCentral on its name, it should correspond to PbPb.
TRandom3 fRandom
! Random generator.
Bool_t AcceptDCA(Float_t pt, Float_t dca)
Bool_t IsPileUpFromNotSPDAndNotEMCal() const
Check if event not from pile-up determined neither by SPD nor by EMCal.
Bool_t fFillCTS
Use data from CTS.
virtual void InitParameters()
Initialize the parameters with default.
virtual void GetVertex(Double_t v[3]) const
Bool_t IsRecalibrationOn() const
virtual Bool_t FillInputEvent(Int_t iEntry, const char *currentFileName)
Float_t fPHOSPtMin
pT Threshold on phos clusters.
Bool_t fSelectEmbeddedClusters
Use only simulated clusters that come from embedding.
Bool_t fEventTrigMinBias
Event is min bias on its name, it should correspond to AliVEvent::kMB, AliVEvent::kAnyInt.
Float_t fEMCALParamTimeCutMin[4]
Remove clusters/cells with time smaller than parametrized value, in ns.
Int_t fEventNumber
Event number.
void RecalculateClusterPID(AliVCluster *clu)
TString fTaskName
Name of task that executes the analysis.
Bool_t fRemoveBadTriggerEvents
Remove triggered events because trigger was exotic, bad, or out of BC.
Bool_t IsPileUpFromSPDOrEMCal() const
Check if event is from pile-up determined by SPD or EMCal.
Bool_t IsPileUpFromEMCalAndNotSPD() const
Check if event is from pile-up determined by EMCal, not by SPD.
Bool_t fEventTrigEMCALL1Gamma2
Event is L1-Gamma, threshold 2 on its name, it should correspond kEMCEGA.
Bool_t fEventTrigEMCALL1Jet1
Event is L1-Gamma, threshold 1 on its name, it should correspond kEMCEGA.
Bool_t fFillEMCALCells
Use data from EMCAL.
Int_t fTriggerClusterId
Id of trigger cluster (cluster->GetID()).
Float_t fCTSPtMax
pT Threshold on charged particles.
Double_t fEMCALParamTimeCutMax[4]
Remove clusters/cells with time larger than parametrized value, in ns.
virtual void FillInputEMCALCells()
Connects the array with EMCAL cells and the pointer.
Bool_t fFillPHOSCells
Use data from PHOS.
virtual AliVCaloCells * GetEMCALCells() const
Bool_t fRemoveLEDEvents
Remove events where LED was wrongly firing - EMCAL LHC11a.
AliEMCALGeometry * GetEMCALGeometry() const
Bool_t fIsBadMaxCellEvent
Bad cell triggered event flag, only max energy cell is bad.
void RecalculateClusterShowerShapeParameters(AliVCaloCells *cells, AliVCluster *clu)
TString fFiredTriggerClassName
Name of trigger event type used to do the analysis.
virtual TClonesArray * GetAODMCParticles() const
TString GetFiredTriggerClasses() const
Bool_t ReadStack() const
TObjArray * fEMCALClusters
Temporal array with EMCAL CaloClusters.
Bool_t IsEventEMCALL1Jet() const
Float_t fPHOSPtMax
pT Threshold on phos clusters.
Int_t GetMaxEnergyCell(AliVCaloCells *cells, const AliVCluster *clu, Float_t &fraction) const
For a given CaloCluster, it gets the absId of the cell with maximum energy deposit.
Bool_t fAccessTrackTOF
Access the track TOF, in case of problems when accessing GetTOFBunchCrossing.
virtual Bool_t CheckForPrimaryVertex() const
Int_t fEMCalBCEvent[19]
Fill one entry per event if there is a cluster in a given BC.
TString fInputBackgroundJetBranchName
Name of background jet branch.
virtual AliEventplane * GetEventPlane() const
Bool_t fEventTriggerAtSE
Select triggered event at SE base task or here.
Float_t fEMCALPtMax
pT Threshold on emcal clusters.
virtual Bool_t ComparePtHardAndJetPt()
Bool_t IsInFiducialCut(Float_t eta, Float_t phi, Int_t det) const
Bool_t fEventTrigCentral
Event is AliVEvent::kCentral on its name, it should correspond to PbPb.
Float_t fTrackMultEtaCut
Track multiplicity eta cut.
Int_t fCentralityBin[2]
Minimum and maximum value of the centrality for the analysis.
virtual void FillInputPHOSCells()
Connects the array with PHOS cells and the pointer.
Base class for event, clusters and tracks filtering and preparation for the analysis.
TObjArray * fDCALClusters
Temporal array with DCAL CaloClusters, not needed in the normal case, use just EMCal array with DCal ...
Bool_t fUseTrackTimeCut
Do time cut selection.
virtual TString GetEventPlaneMethod() const
TArrayI fAcceptEventsWithBit
Accept events if trigger bit is on.
virtual void Print(const Option_t *opt) const
Print some relevant parameters set for the analysis.
void MatchTriggerCluster(TArrayI patches)
Bool_t fTimeStampEventSelect
Select events within a fraction of data taking time.
Double_t fTimeStampRunMax
Maximum value of time stamp in run.
Double_t fEMCALTimeCutMin
Remove clusters/cells with time smaller than this value, in ns.
Bool_t fEventTrigEMCALL1Jet2
Event is L1-Gamma, threshold 2 on its name, it should correspond kEMCEGA.
void SetEMCalEventBCcut(Int_t bc)
Int_t fNPileUpClustersCut
Cut to select event as pile-up.
Bool_t fCheckFidCut
Do analysis for clusters in defined region.
virtual ~AliCaloTrackReader()
Destructor.
Bool_t fEventTrigEMCALL0
Event is EMCal L0 on its name, it should correspond to AliVEvent::kEMC7, AliVEvent::kEMC1.
Int_t fTrackBCEvent[19]
Fill one entry per event if there is a track in a given BC.
virtual void FillInputEMCALAlgorithm(AliVCluster *clus, Int_t iclus)
Bool_t fComparePtHardAndJetPt
In MonteCarlo, jet events, reject fake events with wrong jet energy.
Int_t fNNonPileUpClusters
Number of clusters with time below 20 ns.
Int_t fNMixedEvent
Number of events in mixed event buffer.
Int_t fCentralityOpt
Option for the returned value of the centrality, possible options 5, 10, 100.
Float_t fZvtxCut
Cut on vertex position.
void RecalculateClusterPosition(AliVCaloCells *cells, AliVCluster *clu)
Recalculate EMCAL cluster position.
AliCaloTrackReader()
Constructor. Initialize parameters.
Bool_t fTriggerPatchClusterMatch
Search for the trigger patch and check if associated cluster was the trigger.
Bool_t fUseParamTimeCut
Use simple or parametrized time cut.
AliFiducialCut * fFiducialCut
Acceptance cuts.
Bool_t fRejectEMCalTriggerEventsWith2Tresholds
Reject events EG2 also triggered by EG1 or EJ2 also triggered by EJ1.
virtual void FillInputEMCAL()
Bool_t fFillEMCAL
Use data from EMCAL.
Bool_t fFillPHOS
Use data from PHOS.
AliVCaloCells * fEMCALCells
! Temporal array with EMCAL AliVCaloCells.
virtual AliAODMCHeader * GetAODMCHeader() const
void CorrectClusterEnergy(AliVCluster *cl)
Correct cluster energy non linearity.
Float_t fCTSPtMin
pT Threshold on charged particles.
virtual AliStack * GetStack() const
Float_t fTriggerL0EventThreshold
L0 Threshold to look for triggered events, set outside.
Float_t fTimeStampEventFracMax
Maximum value of time stamp fraction event.
virtual void ResetLists()
Reset lists, called in AliAnaCaloTrackCorrMaker.
Int_t fNPileUpClusters
Number of clusters with time avobe 20 ns.
Int_t fTrackMult
Track multiplicity.
AliAODEvent * fOutputEvent
! pointer to aod output.
Bool_t fRecalculateClusters
Correct clusters, recalculate them if recalibration parameters is given.
Bool_t fDoRejectNoTrackEvents
Reject events with no selected tracks in event.
Bool_t IsEventEMCALL1Jet2() const
Bool_t fReadStack
Access kine information from stack.
Bool_t fIsTriggerMatchOpenCut[3]
Could not match the event to a trigger patch?, retry opening cuts.
UInt_t fMixEventTriggerMask
Select this triggerered event for mixing, tipically kMB or kAnyINT.
Bool_t fAnaLED
Analyze LED data only.
virtual void FillVertexArray()
Bool_t IsInTimeWindow(Double_t tof, Float_t energy) const
Bool_t fRemoveUnMatchedTriggers
Analyze events where trigger patch and cluster where found or not.
virtual Double_t GetEventPlaneAngle() const
TObjArray * fCTSTracks
Temporal array with tracks.
Int_t fTriggerPatchTimeWindow[2]
Trigger patch selection window.
AliMCEvent * fMC
! Monte Carlo Event Handler.
Float_t fPtHardAndJetPtFactor
Factor between ptHard and jet pT to reject/accept event.
TString fDeltaAODFileName
Delta AOD file name.
Int_t GetVertexBC() const
Bool_t IsPileUpFromEMCal() const
Check if event is from pile-up determined by EMCal.
Bool_t fFillInputNonStandardJetBranch
Flag to use data from non standard jets.
TClonesArray * fNonStandardJets
! Temporal array with jets.
Int_t fEMCalBCEventCut[19]
Fill one entry per event if there is a cluster in a given BC, depend on cluster E, acceptance cut.
Int_t fTriggerClusterIndex
Index in clusters array of trigger cluster.
Bool_t fTriggerL1EventThresholdFix
L1 Threshold is fix and set outside.
Float_t fSmearShowerShapeWidth
Smear shower shape landau function "width" (use in MC).
Int_t fBitEGA
Trigger bit on VCaloTrigger for EGA.
virtual void FillInputBackgroundJets()
Bool_t fDoV0ANDEventSelection
Select events depending on V0AND.
AliAODJetEventBackground * fBackgroundJets
! Background jets.
Bool_t fCorrectELinearity
Correct cluster linearity, always on.
Int_t fVertexBC
Vertex BC.
Bool_t fIsExoticEvent
Exotic trigger event flag.
TString fInputNonStandardJetBranchName
Name of non standard jet branch.
TArrayI fRejectEventsWithBit
Reject events if trigger bit is on.
Bool_t fUseEMCALTimeCut
Do time cut selection.
Int_t fNMCProducedMax
In case of cocktail, select particles in the list with label up to this value.
Int_t fDataType
Select MC: Kinematics, Data: ESD/AOD, MCData: Both.
virtual Bool_t SelectTrack(AliVTrack *, Double_t *)
Bool_t fFillDCAL
Use data from DCAL, not needed in the normal case, use just EMCal array with DCal limits...
void SetTrackEventBCcut(Int_t bc)
Bool_t fIsTriggerMatch
Could match the event to a trigger patch?
Float_t fTriggerL1EventThreshold
L1 Threshold to look for triggered events, set in data.
TString fCentralityClass
Name of selected centrality class.
ULong_t fTrackStatus
Track selection bit, select tracks refitted in TPC, ITS ...
Bool_t IsPileUpFromSPDAndEMCal() const
Check if event is from pile-up determined by SPD and EMCal.
virtual void FillInputPHOS()
Fill the array with PHOS filtered clusters.
Int_t fNMCProducedMin
In case of cocktail, select particles in the list with label from this value.
virtual TObjArray * GetCTSTracks() const
Bool_t fRecalculateVertexBC
Recalculate vertex BC from tracks pointing to vertex.
Float_t RecalibrateClusterEnergy(AliVCluster *cluster, AliVCaloCells *cells)
Recalibrate the cluster energy, considering the recalibration map and the energy of the cells that co...
void SetEMCalEventBC(Int_t bc)
Bool_t CheckCellFiducialRegion(AliVCluster *cluster, AliVCaloCells *cells) const
Bool_t IsEventEMCALL0() const