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