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