AliPhysics  d497afb (d497afb)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
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 "AliGenEventHeader.h"
26 #include "AliESDEvent.h"
27 #include "AliAODEvent.h"
28 #include "AliVTrack.h"
29 #include "AliVParticle.h"
30 #include "AliMixedEvent.h"
31 //#include "AliTriggerAnalysis.h"
32 #include "AliESDVZERO.h"
33 #include "AliVCaloCells.h"
34 #include "AliAnalysisManager.h"
35 #include "AliInputEventHandler.h"
36 #include "AliAODMCParticle.h"
37 #include "AliMultSelection.h"
38 
39 // ---- Detectors ----
40 #include "AliPHOSGeoUtils.h"
41 #include "AliEMCALGeometry.h"
42 #include "AliEMCALRecoUtils.h"
43 
44 // ---- CaloTrackCorr ---
45 #include "AliCalorimeterUtils.h"
46 #include "AliCaloTrackReader.h"
47 
48 // ---- Jets ----
49 #include "AliAODJet.h"
50 #include "AliAODJetEventBackground.h"
51 
53 ClassImp(AliCaloTrackReader) ;
55 
56 //________________________________________
58 //________________________________________
60 TObject(), fEventNumber(-1), //fCurrentFileName(""),
61 fDataType(0), fDebug(0),
62 fFiducialCut(0x0), fCheckFidCut(kFALSE),
63 fComparePtHardAndJetPt(0), fPtHardAndJetPtFactor(0),
64 fComparePtHardAndClusterPt(0),fPtHardAndClusterPtFactor(0),
65 fCTSPtMin(0), fEMCALPtMin(0), fPHOSPtMin(0),
66 fCTSPtMax(0), fEMCALPtMax(0), fPHOSPtMax(0),
67 fEMCALBadChMinDist(0), fPHOSBadChMinDist (0),
68 fEMCALNCellsCut(0), fPHOSNCellsCut(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 fSmearingFunction(0), fSmearNLMMin(0), fSmearNLMMax(0),
87 fTrackStatus(0), fSelectSPDHitTracks(0),
88 fTrackMultNPtCut(0), fTrackMultEtaCut(0.9),
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(0),
107 //Trigger rejection
108 fRemoveBadTriggerEvents(0), fTriggerPatchClusterMatch(0),
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 fTimeStampEventCTPBCCorrExclude(0),
124 fTimeStampEventCTPBCCorrMin(0), fTimeStampEventCTPBCCorrMax(0),
125 fNPileUpClusters(-1), fNNonPileUpClusters(-1), fNPileUpClustersCut(3),
126 fVertexBC(-200), fRecalculateVertexBC(0),
127 fUseAliCentrality(0), fCentralityClass(""), fCentralityOpt(0),
128 fEventPlaneMethod(""),
129 fFillInputNonStandardJetBranch(kFALSE),
130 fNonStandardJets(new TClonesArray("AliAODJet",100)), fInputNonStandardJetBranchName("jets"),
131 fFillInputBackgroundJetBranch(kFALSE),
132 fBackgroundJets(0x0),fInputBackgroundJetBranchName("jets"),
133 fAcceptEventsWithBit(0), fRejectEventsWithBit(0), fRejectEMCalTriggerEventsWith2Tresholds(0),
134 fMomentum(), fOutputContainer(0x0),
135 fhEMCALClusterEtaPhi(0), fhEMCALClusterTimeE(0),
136 fEnergyHistogramNbins(0),
137 fhNEventsAfterCut(0), fNMCGenerToAccept(0), fMCGenerEventHeaderToAccept("")
138 {
139  for(Int_t i = 0; i < 8; i++) fhEMCALClusterCutsE [i]= 0x0 ;
140  for(Int_t i = 0; i < 7; i++) fhPHOSClusterCutsE [i]= 0x0 ;
141  for(Int_t i = 0; i < 6; i++) fhCTSTrackCutsPt [i]= 0x0 ;
142  for(Int_t j = 0; j < 5; j++) { fMCGenerToAccept [j] = ""; fMCGenerIndexToAccept[j] = -1; }
143 
144  InitParameters();
145 }
146 
147 //_______________________________________
149 //_______________________________________
151 {
152  DeletePointers();
153 }
154 
155 //_______________________________________
158 //_______________________________________
160 {
161  delete fFiducialCut ;
162 
163  if(fAODBranchList)
164  {
165  fAODBranchList->Delete();
166  delete fAODBranchList ;
167  }
168 
169  if(fCTSTracks)
170  {
171  if(fDataType!=kMC)fCTSTracks->Clear() ;
172  else fCTSTracks->Delete() ;
173  delete fCTSTracks ;
174  }
175 
176  if(fEMCALClusters)
177  {
178  if(fDataType!=kMC)fEMCALClusters->Clear("C") ;
179  else fEMCALClusters->Delete() ;
180  delete fEMCALClusters ;
181  }
182 
183  if(fDCALClusters)
184  {
185  if(fDataType!=kMC)fDCALClusters->Clear("C") ;
186  else fDCALClusters->Delete() ;
187  delete fDCALClusters ;
188  }
189 
190  if(fPHOSClusters)
191  {
192  if(fDataType!=kMC)fPHOSClusters->Clear("C") ;
193  else fPHOSClusters->Delete() ;
194  delete fPHOSClusters ;
195  }
196 
197  if(fVertex)
198  {
199  for (Int_t i = 0; i < fNMixedEvent; i++)
200  {
201  delete [] fVertex[i] ;
202 
203  }
204  delete [] fVertex ;
205  }
206 
207  //delete fTriggerAnalysis;
208 
209  if(fNonStandardJets)
210  {
211  if(fDataType!=kMC) fNonStandardJets->Clear("C") ;
212  else fNonStandardJets->Delete() ;
213  delete fNonStandardJets ;
214  }
215  delete fBackgroundJets ;
216 
217  fRejectEventsWithBit.Reset();
218  fAcceptEventsWithBit.Reset();
219 
220  if ( fWeightUtils ) delete fWeightUtils ;
221 
222  // Pointers not owned, done by the analysis frame
223  // if(fInputEvent) delete fInputEvent ;
224  // if(fOutputEvent) delete fOutputEvent ;
225  // if(fMC) delete fMC ;
226  // Pointer not owned, deleted by maker
227  // if (fCaloUtils) delete fCaloUtils ;
228 }
229 
230 //____________________________________________________________
234 //____________________________________________________________
236 {
237  Float_t cut = fTrackDCACut[0]+fTrackDCACut[1]/TMath::Power(pt,fTrackDCACut[2]);
238 
239  if(TMath::Abs(dca) < cut)
240  return kTRUE;
241  else
242  return kFALSE;
243 }
244 
245 //_____________________________________________________
248 //_____________________________________________________
250 {
251  Int_t nAccept = fAcceptEventsWithBit.GetSize();
252 
253  //printf("N accept %d\n", nAccept);
254 
255  if( nAccept <= 0 )
256  return kTRUE ; // accept the event
257 
258  UInt_t trigFired = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
259 
260  for(Int_t ibit = 0; ibit < nAccept; ibit++)
261  {
262  Bool_t accept = (trigFired & fAcceptEventsWithBit.At(ibit));
263 
264  //printf("accept %d, ibit %d, bit %d \n",accept, ibit,fAcceptEventsWithBit.At(ibit));
265  if(accept) return kTRUE ; // accept the event
266  }
267 
268  return kFALSE ; // reject the event
269 }
270 
271 //_____________________________________________________
276 //_____________________________________________________
278 {
279  if( !fMC || fNMCGenerToAccept <= 0 ) return kTRUE;
280 
281  TString genName;
282  Int_t genIndex;
283  genIndex = GetCocktailGeneratorAndIndex(mcLabel, genName);
284  //fMC->GetCocktailGenerator(mcLabel,genName);
285 
286  Bool_t generOK = kFALSE;
287  for(Int_t ig = 0; ig < fNMCGenerToAccept; ig++)
288  {
289  if ( fMCGenerToAccept[ig].Contains(genName) ) generOK = kTRUE;
290 
291  if ( generOK && fMCGenerIndexToAccept[ig] >= 0 && fMCGenerToAccept[ig] != genIndex) generOK = kFALSE;
292  }
293 
294  if ( !generOK ) AliDebug(1, Form("skip label %d, gen %s",mcLabel,genName.Data()) );
295 
296  return generOK;
297 }
298 
299 //_____________________________________________________________________
307 //_____________________________________________________________________
309 {
310  //method that gives the generator for a given particle with label index (or that of the corresponding primary)
311  AliVParticle* mcpart0 = (AliVParticle*) GetMC()->GetTrack(index);
312  Int_t genIndex = -1;
313 
314  if(!mcpart0)
315  {
316  printf("AliMCEvent-BREAK: No valid AliMCParticle at label %i\n",index);
317  return -1;
318  }
319 
320  nameGen = GetGeneratorNameAndIndex(index,genIndex);
321 
322  if(nameGen.Contains("nococktailheader") ) return -1;
323 
324  Int_t lab=index;
325 
326  while(nameGen.IsWhitespace())
327  {
328  AliVParticle* mcpart = (AliVParticle*) GetMC()->GetTrack(lab);
329 
330  if(!mcpart)
331  {
332  printf("AliMCEvent-BREAK: No valid AliMCParticle at label %i\n",lab);
333  break;
334  }
335 
336  Int_t mother=0;
337  mother = mcpart->GetMother();
338 
339  if(mother<0)
340  {
341  printf("AliMCEvent - BREAK: Reached primary particle without valid mother\n");
342  break;
343  }
344 
345  AliVParticle* mcmom = (AliVParticle*) GetMC()->GetTrack(mother);
346  if(!mcmom)
347  {
348  printf("AliMCEvent-BREAK: No valid AliMCParticle mother at label %i\n",mother);
349  break;
350  }
351 
352  lab=mother;
353 
354  nameGen = GetGeneratorNameAndIndex(mother,genIndex);
355  }
356 
357  return genIndex;
358 }
359 //_____________________________________________________________________
367 //_____________________________________________________________________
369 {
370  Int_t nsumpart = GetMC()->GetNumberOfPrimaries();
371 
372  genIndex = -1;
373 
374  TList* lh = GetMC()->GetCocktailList();
375  if(!lh)
376  {
377  TString noheader="nococktailheader";
378  return noheader;
379  }
380 
381  Int_t nh = lh->GetEntries();
382 
383  for (Int_t i = nh-1; i >= 0; i--)
384  {
385  AliGenEventHeader* gh = (AliGenEventHeader*)lh->At(i);
386 
387  TString genname = gh->GetName();
388 
389  Int_t npart=gh->NProduced();
390 
391  if (i == 0) npart = nsumpart;
392 
393  if(index < nsumpart && index >= (nsumpart-npart))
394  {
395  genIndex = i ;
396  return genname;
397  }
398 
399  nsumpart-=npart;
400  }
401 
402  TString empty="";
403  return empty;
404 }
405 
406 
407 //_____________________________________________________
410 //_____________________________________________________
412 {
413  Int_t nReject = fRejectEventsWithBit.GetSize();
414 
415  //printf("N reject %d\n", nReject);
416 
417  if( nReject <= 0 )
418  return kTRUE ; // accept the event
419 
420  UInt_t trigFired = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
421 
422  for(Int_t ibit = 0; ibit < nReject; ibit++)
423  {
424  Bool_t reject = (trigFired & fRejectEventsWithBit.At(ibit));
425 
426  //printf("reject %d, ibit %d, bit %d \n",reject, ibit,fRejectEventsWithBit.At(ibit));
427  if(reject) return kFALSE ; // reject the event
428  }
429 
430  return kTRUE ; // accept the event
431 }
432 
433 //_____________________________________________
437 //_____________________________________________
439 {
440  //-----------------------------------------------------------
441  // Reject events depending on the event species type
442  //-----------------------------------------------------------
443 
444  // Event types:
445  // kStartOfRun = 1, // START_OF_RUN
446  // kEndOfRun = 2, // END_OF_RUN
447  // kStartOfRunFiles = 3, // START_OF_RUN_FILES
448  // kEndOfRunFiles = 4, // END_OF_RUN_FILES
449  // kStartOfBurst = 5, // START_OF_BURST
450  // kEndOfBurst = 6, // END_OF_BURST
451  // kPhysicsEvent = 7, // PHYSICS_EVENT
452  // kCalibrationEvent = 8, // CALIBRATION_EVENT
453  // kFormatError = 9, // EVENT_FORMAT_ERROR
454  // kStartOfData = 10, // START_OF_DATA
455  // kEndOfData = 11, // END_OF_DATA
456  // kSystemSoftwareTriggerEvent = 12, // SYSTEM_SOFTWARE_TRIGGER_EVENT
457  // kDetectorSoftwareTriggerEvent = 13 // DETECTOR_SOFTWARE_TRIGGER_EVENT
458 
459  Int_t eventType = 0;
460  if(fInputEvent->GetHeader())
461  eventType = ((AliVHeader*)fInputEvent->GetHeader())->GetEventType();
462 
463  AliDebug(3,Form("Event type %d",eventType));
464 
465  // Select only Physics events in data, eventType = 7,
466  // usually done by physics selection
467  // MC not set, eventType =0, I think, so do not apply a value to fEventType by default
468  // LED events have eventType = 8, implemented a selection cut in the past
469  // but not really useful for EMCal data since LED are not reconstructed (unless wrongly assumed as physics)
470 
471  if ( fEventType >= 0 && eventType != fEventType ) return kFALSE ;
472 
473  AliDebug(1,"Pass event species selection");
474 
475  fhNEventsAfterCut->Fill(1.5);
476 
477  //-----------------------------------------------------------------
478  // In case of mixing analysis, select here the trigger of the event
479  //-----------------------------------------------------------------
480 
481  UInt_t isTrigger = kFALSE;
482  UInt_t isMB = kFALSE;
483 
484  if(!fEventTriggerAtSE)
485  {
486  // In case of mixing analysis, accept MB events, not only Trigger
487  // Track and cluster arrays filled for MB in order to create the pool in the corresponding analysis
488  // via de method in the base class FillMixedEventPool()
489 
490  AliAnalysisManager *manager = AliAnalysisManager::GetAnalysisManager();
491  AliInputEventHandler *inputHandler = dynamic_cast<AliInputEventHandler*>(manager->GetInputEventHandler());
492 
493  if(!inputHandler) return kFALSE ; // to content coverity
494 
495  isTrigger = inputHandler->IsEventSelected() & fEventTriggerMask;
496  isMB = inputHandler->IsEventSelected() & fMixEventTriggerMask;
497 
498  if(!isTrigger && !isMB) return kFALSE;
499 
500  //printf("Selected triggered event : %s\n",GetFiredTriggerClasses().Data());
501  AliDebug(1,"Pass uninteresting triggered events rejection in case of mixing analysis");
502 
503  fhNEventsAfterCut->Fill(2.5);
504  }
505 
506  //-----------------------------------------------------------
507  // Reject events depending on the trigger name
508  // Careful!, if a special MB event string is selected but the option
509  // to select events via the mask in the reader is done, it will not
510  // be taken into account.
511  //-----------------------------------------------------------
512 
513  AliDebug(1,Form("FiredTriggerClass <%s>, selected class <%s>, compare name %d",
516 
517  if ( fFiredTriggerClassName != "" && !isMB )
518  {
520  return kFALSE;
521 
522  AliDebug(1,"Accepted triggered event");
523 
524  fhNEventsAfterCut->Fill(3.5);
525  }
526 
527  //-------------------------------------------------------------------------------------
528  // Reject or accept events depending on the trigger bit
529  //-------------------------------------------------------------------------------------
530 
533 
534  //printf("AliCaloTrackReader::FillInputEvent() - Accept event? %d, Reject event %d? \n",okA,okR);
535 
536  if(!okA || !okR) return kFALSE;
537 
538  AliDebug(1,"Pass event bit rejection");
539 
540  fhNEventsAfterCut->Fill(4.5);
541 
542  //----------------------------------------------------------------------
543  // Do not count events that were likely triggered by an exotic cluster
544  // or out BC cluster
545  //----------------------------------------------------------------------
546 
547  // Set a bit with the event kind, MB, L0, L1 ...
549 
550  // In case of Mixing, avoid checking the triggers in the min bias events
551  if(!fEventTriggerAtSE && (isMB && !isTrigger)) return kTRUE;
552 
554  {
556  {
557  // Reject triggered events when there is coincidence on both EMCal trigger thresholds,
558  // but the requested trigger is the low trigger threshold
559  if(IsEventEMCALL1Jet1 () && IsEventEMCALL1Jet2 () && fFiredTriggerClassName.Contains("EJ2")) return kFALSE;
560  if(IsEventEMCALL1Gamma1() && IsEventEMCALL1Gamma2() && fFiredTriggerClassName.Contains("EG2")) return kFALSE;
561  }
562 
563  //Get Patches that triggered
565 
566  MatchTriggerCluster(patches);
567 
568  patches.Reset();
569 
570  // If requested, remove badly triggeed events, but only when the EMCal trigger bit is set
572  {
573  AliDebug(1,Form("ACCEPT triggered event? \n exotic? %d - bad cell %d - bad Max cell %d - BC %d - Matched %d\n",
575 
576  if (fIsExoticEvent) return kFALSE;
577  else if(fIsBadCellEvent) return kFALSE;
578  else if(fRemoveUnMatchedTriggers && !fIsTriggerMatch) return kFALSE ;
579  else if(fTriggerClusterBC != 0) return kFALSE;
580 
581  AliDebug(1,Form("\t *** YES for %s",GetFiredTriggerClasses().Data()));
582  }
583 
584  AliDebug(1,"Pass EMCal triggered event rejection \n");
585 
586  fhNEventsAfterCut->Fill(5.5);
587  }
588 
589  //-------------------------------------------------------------------------------------
590  // Select events only fired by a certain trigger configuration if it is provided
591  //-------------------------------------------------------------------------------------
592 
593  if (GetFiredTriggerClasses().Contains("FAST") && !GetFiredTriggerClasses().Contains("ALL") && !fAcceptFastCluster)
594  {
595  AliDebug(1,Form("Do not count events from fast cluster, trigger name %s\n",fFiredTriggerClassName.Data()));
596  return kFALSE;
597 
598  fhNEventsAfterCut->Fill(6.5);
599  }
600 
601  //-------------------------------------------------------------------------------------
602  // Reject event if large clusters with large energy
603  // Use only for LHC11a data for the moment, and if input is clusterizer V1 or V1+unfolding
604  // If clusterzer NxN or V2 it does not help
605  //-------------------------------------------------------------------------------------
606 
607  //Int_t run = fInputEvent->GetRunNumber();
608 
609  if ( fRemoveLEDEvents )
610  {
612 
613  if(reject) return kFALSE;
614 
615  AliDebug(1,"Pass LED event rejection");
616 
617  fhNEventsAfterCut->Fill(7.5);
618  } // Remove LED events
619 
620  // All selection criteria passed, accept the event
621  return kTRUE ;
622 }
623 
624 //________________________________________________
629 //________________________________________________
631 {
632  //printf("AliCaloTrackReader::ComparePtHardAndJetPt() - GenHeaderName : %s\n",GetGenEventHeader()->ClassName());
633 
634  if ( !GetGenEventHeader() )
635  {
636  AliError("Skip event, event header is not available!");
637  return kFALSE;
638  }
639 
640  if ( !strcmp(GetGenEventHeader()->ClassName(), "AliGenPythiaEventHeader") )
641  {
642  AliGenPythiaEventHeader* pygeh= (AliGenPythiaEventHeader*) GetGenEventHeader();
643 
644  // Do this check only for gamma-jet or jet-jet productions
645  Int_t process = pygeh->ProcessType();
646  //printf("Process %d \n", process);
647  if ( process != 14 && process != 18 && process != 29 && process != 114 && process != 115 && // prompt gamma
648  process != 11 && process != 12 && process != 13 && process != 28 && process != 53 && // jet-jet
649  process != 68 && process != 96 // jet-jet
650  ) return kTRUE ;
651  //printf("\t yes \n");
652 
653  Int_t nTriggerJets = pygeh->NTriggerJets();
654  Float_t ptHard = pygeh->GetPtHard();
655  TParticle * jet = 0;
656 
657  AliDebug(1,Form("Njets: %d, pT Hard %f",nTriggerJets, ptHard));
658 
659  Float_t tmpjet[]={0,0,0,0};
660  for(Int_t ijet = 0; ijet< nTriggerJets; ijet++)
661  {
662  pygeh->TriggerJet(ijet, tmpjet);
663  jet = new TParticle(94, 21, -1, -1, -1, -1, tmpjet[0],tmpjet[1],tmpjet[2],tmpjet[3], 0,0,0,0);
664 
665  AliDebug(1,Form("jet %d; pycell jet pT %f",ijet, jet->Pt()));
666 
667  //Compare jet pT and pt Hard
668  if(jet->Pt() > fPtHardAndJetPtFactor * ptHard)
669  {
670  AliInfo(Form("Reject jet event with : process %d, pT Hard %2.2f, pycell jet pT %2.2f, rejection factor %1.1f\n",
671  process, ptHard, jet->Pt(), fPtHardAndJetPtFactor));
672  return kFALSE;
673  }
674  }
675 
676  if(jet) delete jet;
677  }
678 
679  return kTRUE ;
680 }
681 
682 //____________________________________________________
687 //____________________________________________________
689 {
690  if ( !GetGenEventHeader() )
691  {
692  AliError("Skip event, event header is not available!");
693  return kFALSE;
694  }
695 
696  if ( !strcmp(GetGenEventHeader()->ClassName(), "AliGenPythiaEventHeader") )
697  {
698  AliGenPythiaEventHeader* pygeh= (AliGenPythiaEventHeader*) GetGenEventHeader();
699 
700  // Do this check only for gamma-jet productions
701  Int_t process = pygeh->ProcessType();
702  //printf("Process %d \n", process);
703  if ( process != 14 && process != 18 && process != 29 && process != 114 && process != 115 // && // prompt gamma
704  //process != 11 && process != 12 && process != 13 && process != 28 && process != 53 && // jet-jet
705  //process != 68 && process != 96 // jet-jet
706  ) return kTRUE ;
707  //printf("\t yes \n");
708 
709  Float_t ptHard = pygeh->GetPtHard();
710 
711  Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
712  for (Int_t iclus = 0; iclus < nclusters; iclus++)
713  {
714  AliVCluster * clus = fInputEvent->GetCaloCluster(iclus) ;
715  Float_t ecluster = clus->E();
716 
717  if(ecluster > fPtHardAndClusterPtFactor*ptHard)
718  {
719  AliInfo(Form("Reject : process %d, ecluster %2.2f, calo %d, factor %2.2f, ptHard %f",
720  process, ecluster,clus->GetType(),fPtHardAndClusterPtFactor,ptHard));
721 
722  return kFALSE;
723  }
724  }
725 
726  }
727 
728  return kTRUE ;
729 }
730 
731 //___________________________________________________
734 //___________________________________________________
736 {
737  fhNEventsAfterCut = new TH1I("hNEventsAfterCut", "Number of analyzed events", 20, 0, 20) ;
738  //fhNEventsAfterCut->SetXTitle("Selection");
739  fhNEventsAfterCut->SetYTitle("# events");
740  fhNEventsAfterCut->GetXaxis()->SetBinLabel(1 ,"1=Input");
741  fhNEventsAfterCut->GetXaxis()->SetBinLabel(2 ,"2=Event Type");
742  fhNEventsAfterCut->GetXaxis()->SetBinLabel(3 ,"3=Mixing Event");
743  fhNEventsAfterCut->GetXaxis()->SetBinLabel(4 ,"4=Trigger string");
744  fhNEventsAfterCut->GetXaxis()->SetBinLabel(5 ,"5=Trigger Bit");
745  fhNEventsAfterCut->GetXaxis()->SetBinLabel(6 ,"6=Good EMC Trigger");
746  fhNEventsAfterCut->GetXaxis()->SetBinLabel(7 ,"7=!Fast Cluster");
747  fhNEventsAfterCut->GetXaxis()->SetBinLabel(8 ,"8=!LED");
748  fhNEventsAfterCut->GetXaxis()->SetBinLabel(9 ,"9=Time stamp");
749  fhNEventsAfterCut->GetXaxis()->SetBinLabel(10,"10=Primary vertex");
750  fhNEventsAfterCut->GetXaxis()->SetBinLabel(11,"11=Null 3 vertex");
751  fhNEventsAfterCut->GetXaxis()->SetBinLabel(12,"12=Z vertex window");
752  fhNEventsAfterCut->GetXaxis()->SetBinLabel(13,"13=Pile-up");
753  fhNEventsAfterCut->GetXaxis()->SetBinLabel(14,"14=V0AND");
754  fhNEventsAfterCut->GetXaxis()->SetBinLabel(15,"15=Centrality");
755  fhNEventsAfterCut->GetXaxis()->SetBinLabel(16,"16=GenHeader");
756  fhNEventsAfterCut->GetXaxis()->SetBinLabel(17,"17=PtHard-Jet");
757  fhNEventsAfterCut->GetXaxis()->SetBinLabel(18,"18=PtHard-Cluster");
758  fhNEventsAfterCut->GetXaxis()->SetBinLabel(19,"19=N Track>0");
759  fhNEventsAfterCut->GetXaxis()->SetBinLabel(20,"20=TOF BC");
761 
762  if(fFillEMCAL)
763  {
764  for(Int_t i = 0; i < 8; i++)
765  {
766  TString names[] =
767  { "NoCut", "Corrected", "GoodCluster", "NonLinearity",
768  "EnergyAndFidutial", "NCells", "BadDist", "Time" } ;
769 
770  fhEMCALClusterCutsE[i] = new TH1F(Form("hEMCALReaderClusterCuts_%d_%s",i,names[i].Data()),
771  Form("EMCal %d, %s",i,names[i].Data()),
773  fhEMCALClusterCutsE[i]->SetYTitle("# clusters");
774  fhEMCALClusterCutsE[i]->SetXTitle("#it{E} (GeV)");
776  }
777 
779  ("hEMCALReaderTimeE","time vs #it{E} after cuts (if no calib, shifted -615 ns)", 100,0,100,400,-400,400);
780  fhEMCALClusterTimeE->SetXTitle("#it{E} (GeV)");
781  fhEMCALClusterTimeE->SetYTitle("#it{time} (ns)");
783 
785  ("hEMCALReaderEtaPhi","#eta vs #varphi",40,-2, 2,50, 0,10);
786  // Very open limits to check problems
787  fhEMCALClusterEtaPhi->SetXTitle("#eta");
788  fhEMCALClusterEtaPhi->SetYTitle("#varphi (rad)");
790  }
791 
792  if(fFillPHOS)
793  {
794  for(Int_t i = 0; i < 7; i++)
795  {
796  TString names[] = {"NoCut", "ExcludeCPV", "BorderCut", "FiducialCut", "EnergyCut", "NCells", "BadDist"};
797 
798  fhPHOSClusterCutsE[i] = new TH1F(Form("hPHOSReaderClusterCuts_%d_%s",i,names[i].Data()),
799  Form("PHOS Cut %d, %s",i,names[i].Data()),
801  fhPHOSClusterCutsE[i]->SetYTitle("# clusters");
802  fhPHOSClusterCutsE[i]->SetXTitle("#it{E} (GeV)");
804  }
805  }
806 
807  if(fFillCTS)
808  {
809  for(Int_t i = 0; i < 6; i++)
810  {
811  TString names[] = {"NoCut", "Status", "ESD_AOD", "TOF", "DCA","PtAcceptanceMult"};
812 
813  fhCTSTrackCutsPt[i] = new TH1F(Form("hCTSReaderClusterCuts_%d_%s",i,names[i].Data()),
814  Form("CTS Cut %d, %s",i,names[i].Data()),
816  fhCTSTrackCutsPt[i]->SetYTitle("# tracks");
817  fhCTSTrackCutsPt[i]->SetXTitle("#it{p}_{T} (GeV)");
819  }
820  }
821 
822  return fOutputContainer ;
823 }
824 
825 //_____________________________________________________
827 //_____________________________________________________
829 {
830  TString parList ; //this will be list of parameters used for this analysis.
831 
832  const Int_t buffersize = 255;
833  char onePar[buffersize] ;
834 
835  snprintf(onePar,buffersize,"--- Reader ---:") ;
836  parList+=onePar ;
837  snprintf(onePar,buffersize,"Data type: %d; zvertex cut %2.2f; EMC cluster name: <%s> ",fDataType, fZvtxCut, fEMCALClustersListName.Data()) ;
838  parList+=onePar ;
839  snprintf(onePar,buffersize,"Use detector: EMC %d, DCA %d, PHOS %d, CTS %d, EMCcells %d, PHOScells %d ; ",
841  snprintf(onePar,buffersize,"E-pT window: EMC (%2.1f,%2.1f), PHOS (%2.1f,%2.1f), CTS (%2.1f,%2.1f); ",
843  parList+=onePar ;
844  snprintf(onePar,buffersize,"Dist to bad channel: EMC > %2.1f, PHOS > %2.1f; ",fEMCALBadChMinDist,fPHOSBadChMinDist) ;
845  parList+=onePar ;
846  snprintf(onePar,buffersize,"N cells: EMC > %d, PHOS > %d; ",fEMCALNCellsCut,fPHOSNCellsCut) ;
847  parList+=onePar ;
848  snprintf(onePar,buffersize,"EMC time cut single window (%2.2f,%2.2f); ",fEMCALTimeCutMin,fEMCALTimeCutMax) ;
849  parList+=onePar ;
850  snprintf(onePar,buffersize,"Check: calo fid cut %d; ",fCheckFidCut) ;
851  parList+=onePar ;
852  snprintf(onePar,buffersize,"Track: status %d, SPD hit %d; ",(Int_t) fTrackStatus, fSelectSPDHitTracks) ;
853  parList+=onePar ;
854  snprintf(onePar,buffersize,"multip. eta cut %1.1f; npt cuts %d;",fTrackMultEtaCut, fTrackMultNPtCut) ;
855  parList+=onePar ;
856 
857  if(fUseTrackDCACut)
858  {
859  snprintf(onePar,buffersize,"DCA cut ON, param (%2.4f,%2.4f,%2.4f); ",fTrackDCACut[0],fTrackDCACut[1],fTrackDCACut[2]) ;
860  parList+=onePar ;
861  }
862 
863  snprintf(onePar,buffersize,"Recalculate Clusters = %d, E linearity = %d; ",fRecalculateClusters, fCorrectELinearity) ;
864  parList+=onePar ;
865 
866  snprintf(onePar,buffersize,"SE trigger sel. %d, not? trigger Mask? %d, MB Trigger Mask for mixed %d; ",
868  parList+=onePar ;
869 
870  snprintf(onePar,buffersize,"Select fired trigger %s; Remove Bad trigger event %d, unmatched %d; Accept fastcluster %d, Reject LED %d ",
872  parList+=onePar ;
873 
875  {
876  snprintf(onePar,buffersize,"Accept only labels from: ");
877  parList+=onePar ;
878  for(Int_t i = 0; i< fNMCGenerToAccept; i++)
879  parList+=(fMCGenerToAccept[i]+" ") ;
880 
881  snprintf(onePar,buffersize,"; ");
882  parList+=onePar ;
883  }
884 
886  {
887  snprintf(onePar,buffersize,"EMC M02 smear ON, function %d, param %2.4f, %d<=NLM<=%d; ",
889  parList+=onePar ;
890  }
891 
893  {
894  snprintf(onePar,buffersize,"jet pt / pt hard < %2.1f; ",fPtHardAndJetPtFactor);
895  parList+=onePar ;
896  }
897 
899  {
900  snprintf(onePar,buffersize,"cluster pt / pt hard < %2.2f",fPtHardAndClusterPtFactor);
901  parList+=onePar ;
902  }
903 
904  snprintf(onePar,buffersize,"Centrality: Class %s, Option %d, Bin [%d,%d]; New centrality %d; Event plane method %s; ",
906  parList+=onePar ;
907 
908  return new TObjString(parList) ;
909 }
910 
911 
912 //______________________________________________
914 //______________________________________________
916 {
917  if(fMC)
918  {
919  return fMC->Header();
920  }
921  else
922  {
923  AliInfo("Header is not available");
924  return 0x0 ;
925  }
926 }
927 
928 //_________________________________________________________
931 //_________________________________________________________
933 {
934  AliInfo("Input are not AODs");
935 
936  return NULL ;
937 }
938 
939 //________________________________________________________
942 //________________________________________________________
943 AliAODMCHeader* AliCaloTrackReader::GetAODMCHeader() const
944 {
945  AliInfo("Input are not AODs");
946 
947  return NULL ;
948 }
949 
950 //___________________________________________________________
958 //___________________________________________________________
959 Int_t AliCaloTrackReader::GetVertexBC(const AliVVertex * vtx)
960 {
961  Int_t vertexBC=vtx->GetBC();
962 
963  if(!fRecalculateVertexBC) return vertexBC;
964 
965  // Value not available, recalculate it.
966 
967  Double_t bz = fInputEvent->GetMagneticField();
968  Bool_t bc0 = kFALSE;
969  Int_t ntr = GetCTSTracks()->GetEntriesFast();
970  //printf("N Tracks %d\n",ntr);
971 
972  for(Int_t i = 0 ; i < ntr ; i++)
973  {
974  AliVTrack * track = (AliVTrack*) (GetCTSTracks()->At(i));
975 
976  //Check if has TOF info, if not skip
977  ULong_t status = track->GetStatus();
978  Bool_t okTOF = (status & AliVTrack::kTOFout) == AliVTrack::kTOFout ;
979  vertexBC = track->GetTOFBunchCrossing(bz);
980  Float_t pt = track->Pt();
981 
982  if(!okTOF) continue;
983 
984  // Get DCA x, y
985  Double_t dca[2] = {1e6,1e6};
986  Double_t covar[3] = {1e6,1e6,1e6};
987  track->PropagateToDCA(vtx,bz,100.,dca,covar);
988 
989  if(AcceptDCA(pt,dca[0]))
990  {
991  if (vertexBC !=0 && fVertexBC != AliVTrack::kTOFBCNA) return vertexBC;
992  else if(vertexBC == 0) bc0 = kTRUE;
993  }
994  }
995 
996  if( bc0 ) vertexBC = 0 ;
997  else vertexBC = AliVTrack::kTOFBCNA ;
998 
999  return vertexBC;
1000 }
1001 
1002 //_____________________________
1009 //_____________________________
1011 {
1012  return track->GetID();
1013 }
1014 
1015 //_____________________________
1018 //_____________________________
1020 {
1021  // Activate debug level in AliAnaWeights
1022  if( fWeightUtils->GetDebug() >= 0 )
1023  (AliAnalysisManager::GetAnalysisManager())->AddClassDebug(fWeightUtils->ClassName(), fWeightUtils->GetDebug());
1024 }
1025 
1026 //_______________________________________
1028 //_______________________________________
1030 {
1031  fDataType = kESD ;
1032  fCTSPtMin = 0.1 ;
1033  fEMCALPtMin = 0.1 ;
1034  fPHOSPtMin = 0.1 ;
1035  fCTSPtMax = 1000. ;
1036  fEMCALPtMax = 1000. ;
1037  fPHOSPtMax = 1000. ;
1038 
1039  fEMCALBadChMinDist = 0; // open, 2; // standard
1040  fPHOSBadChMinDist = 0; // open, 2; // standard
1041 
1042  fEMCALNCellsCut = 0; // open, 1; // standard
1043  fPHOSNCellsCut = 0; // open, 2; // standard
1044 
1045  //Track DCA cuts
1046  // dca_xy cut = 0.0105+0.0350/TMath::Power(pt,1.1);
1047  fTrackDCACut[0] = 0.0105;
1048  fTrackDCACut[1] = 0.0350;
1049  fTrackDCACut[2] = 1.1;
1050 
1051  //Do not filter the detectors input by default.
1052  fFillEMCAL = kFALSE;
1053  fFillDCAL = kFALSE;
1054  fFillPHOS = kFALSE;
1055  fFillCTS = kFALSE;
1056  fFillEMCALCells = kFALSE;
1057  fFillPHOSCells = kFALSE;
1058 
1059  fDeltaAODFileName = "deltaAODPartCorr.root";
1061  fEventTriggerMask = AliVEvent::kAny;
1062  fMixEventTriggerMask = AliVEvent::kAnyINT;
1063  fEventTriggerAtSE = kTRUE; // Use only events that pass event selection at SE base class
1064 
1065  fAcceptFastCluster = kTRUE;
1066  fEventType = -1;
1067 
1068  //We want tracks fitted in the detectors:
1069  //fTrackStatus=AliESDtrack::kTPCrefit;
1070  //fTrackStatus|=AliESDtrack::kITSrefit;
1071  fTrackStatus = 0;
1072 
1073  fV0ADC[0] = 0; fV0ADC[1] = 0;
1074  fV0Mul[0] = 0; fV0Mul[1] = 0;
1075 
1076  fZvtxCut = 10.;
1077 
1078  fNMixedEvent = 1;
1079 
1080  fPtHardAndJetPtFactor = 4. ;
1082 
1083  //Centrality
1084  fUseAliCentrality = kFALSE;
1085  fCentralityClass = "V0M";
1086  fCentralityOpt = 100;
1087  fCentralityBin[0] = fCentralityBin[1]=-1;
1088 
1089  fEventPlaneMethod = "V0";
1090 
1091  // Allocate memory (not sure this is the right place)
1092  fCTSTracks = new TObjArray();
1093  fEMCALClusters = new TObjArray();
1094  fDCALClusters = new TObjArray();
1095  fPHOSClusters = new TObjArray();
1096  //fTriggerAnalysis = new AliTriggerAnalysis;
1097  fAODBranchList = new TList ;
1098  fOutputContainer = new TList ;
1099 
1100  fEnergyHistogramNbins = 200;
1101  fEnergyHistogramLimit[0] = 0 ;
1102  fEnergyHistogramLimit[1] = 100;
1103 
1104  fPileUpParamSPD[0] = 3 ; fPileUpParamSPD[1] = 0.8 ;
1105  fPileUpParamSPD[2] = 3.0 ; fPileUpParamSPD[3] = 2.0 ; fPileUpParamSPD[4] = 5.0;
1106 
1107  // Parametrized time cut (LHC11d)
1110 
1111  // Parametrized time cut (LHC11c)
1112  //fEMCALParamTimeCutMin[0] =-5; fEMCALParamTimeCutMin[1] =-1 ; fEMCALParamTimeCutMin[2] = 1.87; fEMCALParamTimeCutMin[3] = 0.4;
1113  //fEMCALParamTimeCutMax[0] = 3.5; fEMCALParamTimeCutMax[1] = 50; fEMCALParamTimeCutMax[2] = 0.15; fEMCALParamTimeCutMax[3] = 1.6;
1114 
1115  fTimeStampRunMin = -1;
1116  fTimeStampRunMax = 1e12;
1119 
1122 
1123  for(Int_t i = 0; i < 19; i++)
1124  {
1125  fEMCalBCEvent [i] = 0;
1126  fEMCalBCEventCut[i] = 0;
1127  fTrackBCEvent [i] = 0;
1128  fTrackBCEventCut[i] = 0;
1129  }
1130 
1131  // Trigger match-rejection
1132  fTriggerPatchTimeWindow[0] = 8;
1133  fTriggerPatchTimeWindow[1] = 9;
1134 
1135  fTriggerClusterBC = -10000 ;
1138  fTriggerClusterIndex = -1;
1139  fTriggerClusterId = -1;
1140 
1141  //Jets
1144  if(!fNonStandardJets) fNonStandardJets = new TClonesArray("AliAODJet",100);
1147  if(!fBackgroundJets) fBackgroundJets = new AliAODJetEventBackground();
1148 
1149  fSmearShowerShapeWidth = 0.005;
1150  fSmearNLMMin = 1;
1151  fSmearNLMMax = 1;
1152 
1153  fWeightUtils = new AliAnaWeights() ;
1154  fEventWeight = 1 ;
1155 
1156  fTrackMultNPtCut = 8;
1157  fTrackMultPtCut[0] = 0.15; fTrackMultPtCut[1] = 0.5; fTrackMultPtCut[2] = 1.0;
1158  fTrackMultPtCut[3] = 2.0 ; fTrackMultPtCut[4] = 4.0; fTrackMultPtCut[5] = 6.0;
1159  fTrackMultPtCut[6] = 8.0 ; fTrackMultPtCut[7] = 10.;
1160  fTrackMultPtCut[8] = 15.0; fTrackMultPtCut[9] = 20.;
1161 }
1162 
1163 //__________________________________________________________________________
1166 //__________________________________________________________________________
1168 {
1169  // Parametrized cut depending on E
1170  if(fUseParamTimeCut)
1171  {
1174  //printf("tof %f, minCut %f, maxCut %f\n",tof,minCut,maxCut);
1175  if( tof < minCut || tof > maxCut ) return kFALSE ;
1176  }
1177 
1178  //In any case, the time should to be larger than the fixed window ...
1179  if( tof < fEMCALTimeCutMin || tof > fEMCALTimeCutMax ) return kFALSE ;
1180 
1181  return kTRUE ;
1182 }
1183 
1184 //________________________________________________
1187 //________________________________________________
1189 {
1190  return fInputEvent->IsPileupFromSPD((Int_t) fPileUpParamSPD[0] , fPileUpParamSPD[1] ,
1191  fPileUpParamSPD[2] , fPileUpParamSPD[3] , fPileUpParamSPD[4] );
1192  //printf("Param : %d, %2.2f, %2.2f, %2.2f, %2.2f\n",(Int_t) fPileUpParamSPD[0], fPileUpParamSPD[1], fPileUpParamSPD[2], fPileUpParamSPD[3], fPileUpParamSPD[4]);
1193 }
1194 
1195 //__________________________________________________
1197 //__________________________________________________
1199 {
1200  if(fNPileUpClusters > fNPileUpClustersCut) return kTRUE ;
1201  else return kFALSE;
1202 }
1203 
1204 //________________________________________________________
1206 //________________________________________________________
1208 {
1209  if( IsPileUpFromSPD() && IsPileUpFromEMCal()) return kTRUE ;
1210  else return kFALSE;
1211 }
1212 
1213 //_______________________________________________________
1215 //_______________________________________________________
1217 {
1218  if( IsPileUpFromSPD() || IsPileUpFromEMCal()) return kTRUE ;
1219  else return kFALSE;
1220 }
1221 
1222 //___________________________________________________________
1224 //___________________________________________________________
1226 {
1227  if( IsPileUpFromSPD() && !IsPileUpFromEMCal()) return kTRUE ;
1228  else return kFALSE;
1229 }
1230 
1231 //___________________________________________________________
1233 //___________________________________________________________
1235 {
1236  if( !IsPileUpFromSPD() && IsPileUpFromEMCal()) return kTRUE ;
1237  else return kFALSE;
1238 }
1239 
1240 //______________________________________________________________
1242 //______________________________________________________________
1244 {
1245  if( !IsPileUpFromSPD() && !IsPileUpFromEMCal()) return kTRUE ;
1246  else return kFALSE;
1247 }
1248 
1249 //___________________________________________________________________________________
1257 //___________________________________________________________________________________
1258 Bool_t AliCaloTrackReader::FillInputEvent(Int_t iEntry, const char * /*curFileName*/)
1259 {
1260  fEventNumber = iEntry;
1261  fTriggerClusterIndex = -1;
1262  fTriggerClusterId = -1;
1263  fIsTriggerMatch = kFALSE;
1264  fTriggerClusterBC = -10000;
1265  fIsExoticEvent = kFALSE;
1266  fIsBadCellEvent = kFALSE;
1267  fIsBadMaxCellEvent = kFALSE;
1268 
1269  fIsTriggerMatchOpenCut[0] = kFALSE ;
1270  fIsTriggerMatchOpenCut[1] = kFALSE ;
1271  fIsTriggerMatchOpenCut[2] = kFALSE ;
1272 
1273  //fCurrentFileName = TString(currentFileName);
1274  if(!fInputEvent)
1275  {
1276  AliInfo("Input event not available, skip event analysis");
1277  return kFALSE;
1278  }
1279 
1280  fhNEventsAfterCut->Fill(0.5);
1281 
1282  //-----------------------------------------------
1283  // Select the event depending on the trigger type
1284  // and other event characteristics
1285  // like the goodness of the EMCal trigger
1286  //-----------------------------------------------
1287 
1288  Bool_t accept = CheckEventTriggers();
1289  if(!accept) return kFALSE;
1290 
1291  AliDebug(1,"Pass Event trigger selection");
1292 
1293 
1294  //------------------------------------------------------
1295  // Event rejection depending on time stamp
1296  //------------------------------------------------------
1297 
1298  if ( fTimeStampEventSelect )
1299  {
1300  Int_t timeStamp = fInputEvent->GetTimeStamp();
1301  Float_t timeStampFrac = 1.*(timeStamp-fTimeStampRunMin) / (fTimeStampRunMax-fTimeStampRunMin);
1302 
1303  //printf("stamp0 %d, max0 %d, frac %f\n", timeStamp-fTimeStampRunMin,fTimeStampRunMax-fTimeStampRunMin, timeStampFrac);
1304 
1305  if(timeStampFrac < fTimeStampEventFracMin || timeStampFrac > fTimeStampEventFracMax) return kFALSE;
1306 
1307  AliDebug(1,"Pass Time Stamp rejection");
1308 
1309  fhNEventsAfterCut->Fill(8.5);
1310  }
1311 
1313  {
1314  AliESDEvent* esd = dynamic_cast<AliESDEvent*> (fInputEvent);
1315  if(esd)
1316  {
1317  Int_t timeStamp = esd->GetTimeStampCTPBCCorr();
1318 
1319  if(timeStamp > fTimeStampEventCTPBCCorrMin && timeStamp <= fTimeStampEventCTPBCCorrMax) return kFALSE;
1320  }
1321 
1322  AliDebug(1,"Pass Time Stamp CTPBCCorr rejection");
1323 
1324  fhNEventsAfterCut->Fill(8.5);
1325  }
1326 
1327 
1328  //------------------------------------------------------
1329  // Event rejection depending on vertex, pileup, v0and
1330  //------------------------------------------------------
1331 
1332  FillVertexArray();
1333 
1335  {
1336  if( !CheckForPrimaryVertex() ) return kFALSE; // algorithm in ESD/AOD Readers
1337 
1338  fhNEventsAfterCut->Fill(9.5);
1339 
1340  if( TMath::Abs(fVertex[0][0] ) < 1.e-6 &&
1341  TMath::Abs(fVertex[0][1] ) < 1.e-6 &&
1342  TMath::Abs(fVertex[0][2] ) < 1.e-6 ) return kFALSE;
1343 
1344  AliDebug(1,"Pass primary vertex/null rejection");
1345 
1346  fhNEventsAfterCut->Fill(10.5);
1347  }
1348 
1349  //Reject events with Z vertex too large, only for SE analysis, if not, cut on the analysis code
1350  if(!GetMixedEvent() && TMath::Abs(fVertex[0][2]) > fZvtxCut) return kFALSE;
1351 
1352  fhNEventsAfterCut->Fill(11.5);
1353 
1354  AliDebug(1,"Pass z vertex rejection");
1355 
1356 
1357  //printf("Reader : IsPileUp %d, Multi %d\n",IsPileUpFromSPD(),fInputEvent->IsPileupFromSPDInMultBins());
1358 
1360  {
1361  // Do not analyze events with pileup
1362  Bool_t bPileup = IsPileUpFromSPD();
1363  //IsPileupFromSPDInMultBins() // method to try
1364  //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]);
1365  if(bPileup) return kFALSE;
1366 
1367  AliDebug(1,"Pass Pile-Up event rejection");
1368 
1369  fhNEventsAfterCut->Fill(12.5);
1370  }
1371 
1373  {
1374  AliVVZERO* v0 = fInputEvent->GetVZEROData();
1375 
1376  Bool_t bV0AND = ((v0->GetV0ADecision()==1) && (v0->GetV0CDecision()==1));
1377  //bV0AND = fTriggerAnalysis->IsOfflineTriggerFired((AliESDEvent*)fInputEvent, AliTriggerAnalysis::kV0AND);
1378  //printf("V0AND event? %d\n",bV0AND);
1379 
1380  if(!bV0AND)
1381  {
1382  AliDebug(1,"Reject event by V0AND");
1383  return kFALSE;
1384  }
1385 
1386  AliDebug(1,"Pass V0AND event rejection");
1387 
1388  fhNEventsAfterCut->Fill(13.5);
1389  }
1390 
1391  //------------------------------------------------------
1392  // Check if there is a centrality value, PbPb analysis,
1393  // and if a centrality bin selection is requested
1394  //------------------------------------------------------
1395 
1396  //If we need a centrality bin, we select only those events in the corresponding bin.
1397  Int_t cen = -1;
1398  if ( fCentralityBin[0] >= 0 && fCentralityBin[1] >= 0 )
1399  {
1400  cen = GetEventCentrality();
1401 
1402  if(cen > fCentralityBin[1] || cen <= fCentralityBin[0]) return kFALSE; //reject events out of bin.
1403 
1404  AliDebug(1,"Pass centrality rejection");
1405 
1406  fhNEventsAfterCut->Fill(14.5);
1407  }
1408 
1409  //----------------------------------------------------------------
1410  // Reject the event if the event header name is not
1411  // the one requested among the possible generators.
1412  // Needed in case of cocktail MC generation with multiple options.
1413  //----------------------------------------------------------------
1415  {
1416  if(!GetGenEventHeader()) return kFALSE;
1417 
1418  AliDebug(1,"Pass Event header selection");
1419 
1420  fhNEventsAfterCut->Fill(15.5);
1421  }
1422 
1423  //---------------------------------------------------------------------------
1424  // In case of analysis of events with jets, skip those with jet pt > 5 pt hard
1425  // To be used on for MC data in pT hard bins
1426  //---------------------------------------------------------------------------
1427 
1429  {
1430  if(!ComparePtHardAndJetPt()) return kFALSE ;
1431 
1432  AliDebug(1,"Pass Pt Hard - Jet rejection");
1433 
1434  fhNEventsAfterCut->Fill(16.5);
1435  }
1436 
1438  {
1439  if(!ComparePtHardAndClusterPt()) return kFALSE ;
1440 
1441  AliDebug(1,"Pass Pt Hard - Cluster rejection");
1442 
1443  fhNEventsAfterCut->Fill(17.5);
1444  }
1445 
1446  //------------------------------------------------------------------
1447  // Recover the weight assigned to the event, if provided
1448  // right now only for pT-hard bins and centrality depedent weights
1449  //------------------------------------------------------------------
1451  {
1453 
1454  fWeightUtils->SetPythiaEventHeader(((AliGenPythiaEventHeader*)GetGenEventHeader()));
1455 
1457  }
1458 
1459  //-------------------------------------------------------
1460  // Get the main vertex BC, in case not available
1461  // it is calculated in FillCTS checking the BC of tracks
1462  //------------------------------------------------------
1463  fVertexBC = fInputEvent->GetPrimaryVertex()->GetBC();
1464 
1465  //-----------------------------------------------
1466  // Fill the arrays with cluster/tracks/cells data
1467  //-----------------------------------------------
1468 
1469  if(fFillCTS)
1470  {
1471  FillInputCTS();
1472  //Accept events with at least one track
1473  if(fTrackMult[0] == 0 && fDoRejectNoTrackEvents) return kFALSE ;
1474 
1475  AliDebug(1,"Pass rejection of null track events");
1476 
1477  fhNEventsAfterCut->Fill(18.5);
1478  }
1479 
1481  {
1482  if(fVertexBC != 0 && fVertexBC != AliVTrack::kTOFBCNA) return kFALSE ;
1483 
1484  AliDebug(1,"Pass rejection of events with vertex at BC!=0");
1485 
1486  fhNEventsAfterCut->Fill(19.5);
1487  }
1488 
1489  if(fFillEMCALCells)
1491 
1492  if(fFillPHOSCells)
1494 
1495  if(fFillEMCAL || fFillDCAL)
1496  FillInputEMCAL();
1497 
1498  if(fFillPHOS)
1499  FillInputPHOS();
1500 
1501  FillInputVZERO();
1502 
1503  //one specified jet branch
1508 
1509  AliDebug(1,"Event accepted for analysis");
1510 
1511  return kTRUE ;
1512 }
1513 
1514 //__________________________________________________
1517 //__________________________________________________
1519 {
1520  if(fUseAliCentrality)
1521  {
1522  if ( !GetCentrality() ) return -1;
1523 
1524  AliDebug(1,Form("Cent. Percentile: V0M %2.2f, CL0 %2.2f, CL1 %2.2f; selected class %s",
1525  GetCentrality()->GetCentralityPercentile("V0M"),
1526  GetCentrality()->GetCentralityPercentile("CL0"),
1527  GetCentrality()->GetCentralityPercentile("CL1"),
1528  fCentralityClass.Data()));
1529 
1530  if (fCentralityOpt == 100) return (Int_t) GetCentrality()->GetCentralityPercentile(fCentralityClass); // 100 bins max
1531  else if(fCentralityOpt == 10) return GetCentrality()->GetCentralityClass10(fCentralityClass);// 10 bins max
1532  else if(fCentralityOpt == 20) return GetCentrality()->GetCentralityClass5(fCentralityClass); // 20 bins max
1533  else
1534  {
1535  AliInfo(Form("Unknown centrality option %d, use 10, 20 or 100\n",fCentralityOpt));
1536  return -1;
1537  }
1538  }
1539  else
1540  {
1541  if ( !GetMultSelCen() ) return -1;
1542 
1543  AliDebug(1,Form("Mult. Percentile: V0M %2.2f, CL0 %2.2f, CL1 %2.2f; selected class %s",
1544  GetMultSelCen()->GetMultiplicityPercentile("V0M",1),
1545  GetMultSelCen()->GetMultiplicityPercentile("CL0",1),
1546  GetMultSelCen()->GetMultiplicityPercentile("CL1",1),
1547  fCentralityClass.Data()));
1548 
1549  return GetMultSelCen()->GetMultiplicityPercentile(fCentralityClass, kTRUE); // returns centrality only for events used in calibration
1550 
1551  // equivalent to
1552  //GetMultSelCen()->GetMultiplicityPercentile("V0M", kFALSE); // returns centrality for any event
1553  //Int_t qual = GetMultSelCen()->GetEvSelCode(); if (qual ! = 0) cent = qual;
1554  }
1555 }
1556 
1557 //_____________________________________________________
1560 //_____________________________________________________
1562 {
1563  if( !GetEventPlane() ) return -1000;
1564 
1565  Float_t ep = GetEventPlane()->GetEventplane(GetEventPlaneMethod(), GetInputEvent());
1566 
1567  if(GetEventPlaneMethod()=="Q" && (ep < 0 || ep > TMath::Pi()))
1568  {
1569  AliDebug(1,Form("Bad EP for <Q> method : %f\n",ep));
1570  return -1000;
1571  }
1572  else if(GetEventPlaneMethod().Contains("V0") )
1573  {
1574  if((ep > TMath::Pi()/2 || ep < -TMath::Pi()/2))
1575  {
1576  AliDebug(1,Form("Bad EP for <%s> method : %f\n",GetEventPlaneMethod().Data(), ep));
1577  return -1000;
1578  }
1579 
1580  ep+=TMath::Pi()/2; // put same range as for <Q> method
1581  }
1582 
1583  AliDebug(3,Form("Event plane angle %f",ep));
1584 
1585 // if(fDebug > 0 )
1586 // {
1587 // if (ep > TMath::Pi()) printf("AliCaloTrackReader::GetEventPlaneAngle() - Too large angle = %f\n",ep);
1588 // else if(ep < 0 ) printf("AliCaloTrackReader::GetEventPlaneAngle() - Negative angle = %f\n" ,ep);
1589 // }
1590 
1591  return ep;
1592 }
1593 
1594 //__________________________________________________________
1596 //__________________________________________________________
1598 {
1599  vertex[0] = fVertex[0][0];
1600  vertex[1] = fVertex[0][1];
1601  vertex[2] = fVertex[0][2];
1602 }
1603 
1604 //__________________________________________________________________________
1606 //__________________________________________________________________________
1607 void AliCaloTrackReader::GetVertex(Double_t vertex[3], Int_t evtIndex) const
1608 {
1609  vertex[0] = fVertex[evtIndex][0];
1610  vertex[1] = fVertex[evtIndex][1];
1611  vertex[2] = fVertex[evtIndex][2];
1612 }
1613 
1614 //________________________________________
1617 //________________________________________
1619 {
1620  // Delete previous vertex
1621  if(fVertex)
1622  {
1623  for (Int_t i = 0; i < fNMixedEvent; i++)
1624  {
1625  delete [] fVertex[i] ;
1626  }
1627  delete [] fVertex ;
1628  }
1629 
1630  fVertex = new Double_t*[fNMixedEvent] ;
1631  for (Int_t i = 0; i < fNMixedEvent; i++)
1632  {
1633  fVertex[i] = new Double_t[3] ;
1634  fVertex[i][0] = 0.0 ;
1635  fVertex[i][1] = 0.0 ;
1636  fVertex[i][2] = 0.0 ;
1637  }
1638 
1639  if (!fMixedEvent)
1640  { // Single event analysis
1641  if(fDataType!=kMC)
1642  {
1643 
1644  if(fInputEvent->GetPrimaryVertex())
1645  {
1646  fInputEvent->GetPrimaryVertex()->GetXYZ(fVertex[0]);
1647  }
1648  else
1649  {
1650  AliWarning("NULL primary vertex");
1651  fVertex[0][0]=0.; fVertex[0][1]=0.; fVertex[0][2]=0.;
1652  }//Primary vertex pointer do not exist
1653 
1654  } else
1655  {// MC read event
1656  fVertex[0][0]=0.; fVertex[0][1]=0.; fVertex[0][2]=0.;
1657  }
1658 
1659  AliDebug(1,Form("Single Event Vertex : %f,%f,%f\n",fVertex[0][0],fVertex[0][1],fVertex[0][2]));
1660 
1661  } else
1662  { // MultiEvent analysis
1663  for (Int_t iev = 0; iev < fNMixedEvent; iev++)
1664  {
1665  if (fMixedEvent->GetVertexOfEvent(iev))
1666  fMixedEvent->GetVertexOfEvent(iev)->GetXYZ(fVertex[iev]);
1667  else
1668  AliWarning("No vertex found");
1669 
1670  AliDebug(1,Form("Multi Event %d Vertex : %f,%f,%f",iev,fVertex[iev][0],fVertex[iev][1],fVertex[iev][2]));
1671  }
1672  }
1673 }
1674 
1675 //_____________________________________
1681 //_____________________________________
1683 {
1684  AliDebug(1,"Begin");
1685 
1686  Double_t pTrack[3] = {0,0,0};
1687 
1688  Int_t nTracks = fInputEvent->GetNumberOfTracks() ;
1689  Int_t nstatus = 0;
1690  Double_t bz = GetInputEvent()->GetMagneticField();
1691 
1692  for(Int_t i = 0; i < 19; i++)
1693  {
1694  fTrackBCEvent [i] = 0;
1695  fTrackBCEventCut[i] = 0;
1696  }
1697 
1698  for(Int_t iptCut = 0; iptCut < fTrackMultNPtCut; iptCut++ )
1699  {
1700  fTrackMult [iptCut] = 0;
1701  fTrackSumPt[iptCut] = 0;
1702  }
1703 
1704  Bool_t bc0 = kFALSE;
1705  if(fRecalculateVertexBC) fVertexBC = AliVTrack::kTOFBCNA;
1706 
1707  for (Int_t itrack = 0; itrack < nTracks; itrack++)
1708  {
1709  AliVTrack * track = (AliVTrack*)fInputEvent->GetTrack(itrack) ; // retrieve track from esd
1710 
1711  if ( !AcceptParticleMCLabel( TMath::Abs(track->GetLabel()) ) ) continue ;
1712 
1713  fhCTSTrackCutsPt[0]->Fill(track->Pt());
1714 
1715  //Select tracks under certain conditions, TPCrefit, ITSrefit ... check the set bits
1716  ULong_t status = track->GetStatus();
1717 
1718  if (fTrackStatus && !((status & fTrackStatus) == fTrackStatus))
1719  continue ;
1720 
1721  fhCTSTrackCutsPt[1]->Fill(track->Pt());
1722 
1723  nstatus++;
1724 
1725  //-------------------------
1726  // Select the tracks depending on cuts of AOD or ESD
1727  if(!SelectTrack(track, pTrack)) continue ;
1728 
1729  fhCTSTrackCutsPt[2]->Fill(track->Pt());
1730 
1731  //-------------------------
1732  // TOF cuts
1733  Bool_t okTOF = ( (status & AliVTrack::kTOFout) == AliVTrack::kTOFout ) ;
1734  Double_t tof = -1000;
1735  Int_t trackBC = -1000 ;
1736 
1737  if(fAccessTrackTOF)
1738  {
1739  if(okTOF)
1740  {
1741  trackBC = track->GetTOFBunchCrossing(bz);
1742  SetTrackEventBC(trackBC+9);
1743 
1744  tof = track->GetTOFsignal()*1e-3;
1745 
1746  //After selecting tracks with small DCA, pointing to vertex, set vertex BC depeding on tracks BC
1748  {
1749  if (trackBC != 0 && trackBC != AliVTrack::kTOFBCNA) fVertexBC = trackBC;
1750  else if(trackBC == 0) bc0 = kTRUE;
1751  }
1752 
1753  //In any case, the time should to be larger than the fixed window ...
1754  if( fUseTrackTimeCut && (trackBC !=0 || tof < fTrackTimeCutMin || tof > fTrackTimeCutMax) )
1755  {
1756  //printf("Remove track time %f and bc = %d\n",tof,trackBC);
1757  continue ;
1758  }
1759  //else printf("Accept track time %f and bc = %d\n",tof,trackBC);
1760  }
1761  }
1762 
1763  fhCTSTrackCutsPt[3]->Fill(track->Pt());
1764 
1765  //---------------------
1766  // DCA cuts
1767  //
1768  fMomentum.SetPxPyPzE(pTrack[0],pTrack[1],pTrack[2],0);
1769 
1770  if(fUseTrackDCACut)
1771  {
1772  Float_t dcaTPC =-999;
1773  //In case of AODs, TPC tracks cannot be propagated back to primary vertex,
1774  if( fDataType == kAOD ) dcaTPC = ((AliAODTrack*) track)->DCA();
1775 
1776  //normal way to get the dca, cut on dca_xy
1777  if(dcaTPC==-999)
1778  {
1779  Double_t dca[2] = {1e6,1e6};
1780  Double_t covar[3] = {1e6,1e6,1e6};
1781  Bool_t okDCA = track->PropagateToDCA(fInputEvent->GetPrimaryVertex(),bz,100.,dca,covar);
1782  if( okDCA) okDCA = AcceptDCA(fMomentum.Pt(),dca[0]);
1783  if(!okDCA)
1784  {
1785  //printf("AliCaloTrackReader::FillInputCTS() - Reject track pt %2.2f, dca_xy %2.4f\n",fMomentum.Pt(),dca[0]);
1786  continue ;
1787  }
1788  }
1789  }// DCA cuts
1790 
1791  fhCTSTrackCutsPt[4]->Fill(track->Pt());
1792 
1793  //-------------------------
1794  // Kinematic/acceptance cuts
1795  //
1796  // Count the tracks in eta < 0.9 and different pT cuts
1797  Float_t ptTrack = fMomentum.Pt();
1798  if(TMath::Abs(track->Eta())< fTrackMultEtaCut)
1799  {
1800  for(Int_t iptCut = 0; iptCut < fTrackMultNPtCut; iptCut++ )
1801  {
1802  if(ptTrack > fTrackMultPtCut[iptCut])
1803  {
1804  fTrackMult [iptCut]++;
1805  fTrackSumPt[iptCut]+=ptTrack;
1806  }
1807  }
1808  }
1809 
1810  if(fCTSPtMin > ptTrack || fCTSPtMax < ptTrack) continue ;
1811 
1812  // Check effect of cuts on track BC
1813  if(fAccessTrackTOF && okTOF) SetTrackEventBCcut(trackBC+9);
1814 
1815  if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(fMomentum.Eta(),fMomentum.Phi(),kCTS)) continue;
1816 
1817  fhCTSTrackCutsPt[5]->Fill(track->Pt());
1818 
1819  // ------------------------------
1820  // Add selected tracks to array
1821  AliDebug(2,Form("Selected tracks pt %3.2f, phi %3.2f deg, eta %3.2f",
1822  fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta()));
1823 
1824  fCTSTracks->Add(track);
1825 
1826  // TODO, check if remove
1827  if (fMixedEvent) track->SetID(itrack);
1828 
1829  }// track loop
1830 
1831  if( fRecalculateVertexBC && (fVertexBC == 0 || fVertexBC == AliVTrack::kTOFBCNA))
1832  {
1833  if( bc0 ) fVertexBC = 0 ;
1834  else fVertexBC = AliVTrack::kTOFBCNA ;
1835  }
1836 
1837  AliDebug(1,Form("CTS entries %d, input tracks %d, pass status %d, multipliticy %d", fCTSTracks->GetEntriesFast(), nTracks, nstatus, fTrackMult[0]));//fCTSTracksNormalInputEntries);
1838 }
1839 
1840 //_______________________________________________________________________________
1858 //_______________________________________________________________________________
1860 {
1861  // Accept clusters with the proper label, only applicable for MC
1862  if ( clus->GetLabel() >= 0 ) // -1 corresponds to noisy MC
1863  {
1864  if ( !AcceptParticleMCLabel(clus->GetLabel()) ) return ;
1865  }
1866 
1867  // TODO, not sure if needed anymore
1868  Int_t vindex = 0 ;
1869  if (fMixedEvent)
1870  vindex = fMixedEvent->EventIndexForCaloCluster(iclus);
1871 
1872  clus->GetMomentum(fMomentum, fVertex[vindex]);
1873 
1874  // No correction/cut applied yet
1875  fhEMCALClusterCutsE[0]->Fill(clus->E());
1876 
1877  //if( (fDebug > 2 && fMomentum.E() > 0.1) || fDebug > 10 )
1878  AliDebug(2,Form("Input cluster E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f, nCells %d",
1879  fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta(), clus->GetNCells()));
1880 
1881  //---------------------------
1882  // Embedding case
1883  // TO BE REVISED
1885  {
1886  if(clus->GetNLabels()==0 || clus->GetLabel() < 0) return;
1887  //else printf("Embedded cluster, %d, n label %d label %d \n",iclus,clus->GetNLabels(),clus->GetLabel());
1888  }
1889 
1890  //--------------------------------------
1891  // Apply some corrections in the cluster
1892  //
1894  {
1895  //Recalibrate the cluster energy
1897  {
1899 
1900  clus->SetE(energy);
1901  //printf("Recalibrated Energy %f\n",clus->E());
1902 
1905 
1906  } // recalculate E
1907 
1908  //Recalculate distance to bad channels, if new list of bad channels provided
1910 
1911  //Recalculate cluster position
1912  if(GetCaloUtils()->IsRecalculationOfClusterPositionOn())
1913  {
1915  //clus->GetPosition(pos);
1916  //printf("After Corrections: e %f, x %f, y %f, z %f\n",clus->E(),pos[0],pos[1],pos[2]);
1917  }
1918 
1919  // Recalculate TOF
1920  if(GetCaloUtils()->GetEMCALRecoUtils()->IsTimeRecalibrationOn())
1921  {
1922  Double_t tof = clus->GetTOF();
1923  Float_t frac =-1;
1924  Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(fEMCALCells, clus,frac);
1925 
1926  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absIdMax,fInputEvent->GetBunchCrossNumber(),tof);
1927 
1928  //additional L1 phase shift
1929  if(GetCaloUtils()->GetEMCALRecoUtils()->IsL1PhaseInTimeRecalibrationOn())
1930  {
1931  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTimeL1Phase(GetCaloUtils()->GetEMCALGeometry()->GetSuperModuleNumber(absIdMax),
1932  fInputEvent->GetBunchCrossNumber(), tof);
1933  }
1934 
1935  clus->SetTOF(tof);
1936 
1937  }// Time recalibration
1938  }
1939 
1940  // Check effect of corrections
1941  fhEMCALClusterCutsE[1]->Fill(clus->E());
1942 
1943  //-----------------------------------------------------------------
1944  // Reject clusters with bad channels, close to borders and exotic
1945  //
1946  Bool_t goodCluster = GetCaloUtils()->GetEMCALRecoUtils()->IsGoodCluster(clus,
1947  GetCaloUtils()->GetEMCALGeometry(),
1948  GetEMCALCells(),fInputEvent->GetBunchCrossNumber());
1949 
1950  if(!goodCluster)
1951  {
1952  //if( (fDebug > 2 && fMomentum.E() > 0.1) || fDebug > 10 )
1953  AliDebug(1,Form("Bad cluster E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f",
1954  fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta()));
1955 
1956  return;
1957  }
1958 
1959  // Check effect of bad cluster removal
1960  fhEMCALClusterCutsE[2]->Fill(clus->E());
1961 
1962  //Float_t pos[3];
1963  //clus->GetPosition(pos);
1964  //printf("Before Corrections: e %f, x %f, y %f, z %f\n",clus->E(),pos[0],pos[1],pos[2]);
1965 
1966  //--------------------------------------
1967  // Correct non linearity or smear energy
1968  //
1970  {
1972 
1973  AliDebug(5,Form("Correct Non Lin: Old E %3.2f, New E %3.2f",
1974  fMomentum.E(),clus->E()));
1975 
1976  // In case of MC analysis, to match resolution/calibration in real data
1977  // Not needed anymore, just leave for MC studies on systematics
1978  if( GetCaloUtils()->GetEMCALRecoUtils()->IsClusterEnergySmeared() )
1979  {
1980  Float_t rdmEnergy = GetCaloUtils()->GetEMCALRecoUtils()->SmearClusterEnergy(clus);
1981 
1982  AliDebug(5,Form("Smear energy: Old E %3.2f, New E %3.2f",clus->E(),rdmEnergy));
1983 
1984  clus->SetE(rdmEnergy);
1985  }
1986  }
1987 
1988  clus->GetMomentum(fMomentum, fVertex[vindex]);
1989  fhEMCALClusterEtaPhi->Fill(fMomentum.Eta(),GetPhi(fMomentum.Phi()));
1990 
1991  // Check effect linearity correction, energy smearing
1992  fhEMCALClusterCutsE[3]->Fill(clus->E());
1993 
1994  // Check the event BC depending on EMCal clustr before final cuts
1995  Double_t tof = clus->GetTOF()*1e9;
1996 
1997  Int_t bc = TMath::Nint(tof/50) + 9;
1998  //printf("tof %2.2f, bc+5=%d\n",tof,bc);
1999 
2000  SetEMCalEventBC(bc);
2001 
2002  //--------------------------------------
2003  // Apply some kinematical/acceptance cuts
2004  //
2005  if(fEMCALPtMin > clus->E() || fEMCALPtMax < clus->E())
2006  {
2007  AliDebug(2,Form("Cluster E out of range, %2.2f < %2.2f < %2.2f",fEMCALPtMin,clus->E(),fEMCALPtMax));
2008  return ;
2009  }
2010 
2011  // Select cluster fiducial region
2012  //
2013  Bool_t bEMCAL = kFALSE;
2014  Bool_t bDCAL = kFALSE;
2015  if(fCheckFidCut)
2016  {
2017  if(fFillEMCAL && fFiducialCut->IsInFiducialCut(fMomentum.Eta(),fMomentum.Phi(),kEMCAL)) bEMCAL = kTRUE ;
2018  if(fFillDCAL && fFiducialCut->IsInFiducialCut(fMomentum.Eta(),fMomentum.Phi(),kDCAL )) bDCAL = kTRUE ;
2019  }
2020  else
2021  {
2022  bEMCAL = kTRUE;
2023  }
2024 
2025  //---------------------------------------------------------------------
2026  // Mask all cells in collumns facing ALICE thick material if requested
2027  //
2029  {
2030  Int_t absId = -1;
2031  Int_t iSupMod = -1;
2032  Int_t iphi = -1;
2033  Int_t ieta = -1;
2034  Bool_t shared = kFALSE;
2035  GetCaloUtils()->GetEMCALRecoUtils()->GetMaxEnergyCell(GetCaloUtils()->GetEMCALGeometry(), GetEMCALCells(),clus,absId,iSupMod,ieta,iphi,shared);
2036 
2037  AliDebug(2,Form("Masked collumn: cluster E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f",
2038  fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta()));
2039 
2040 
2041  if(GetCaloUtils()->MaskFrameCluster(iSupMod, ieta))
2042  {
2043  AliDebug(2,"Mask cluster");
2044  return;
2045  }
2046  }
2047 
2048  // Check effect of energy and fiducial cuts
2049  fhEMCALClusterCutsE[4]->Fill(clus->E());
2050 
2051  //----------------------------------------------------
2052  // Apply N cells cut
2053  //
2054  if(clus->GetNCells() <= fEMCALNCellsCut && fDataType != AliCaloTrackReader::kMC)
2055  {
2056  AliDebug(2,Form("Cluster with n cells %d < %d",clus->GetNCells(), fEMCALNCellsCut));
2057  return ;
2058  }
2059 
2060  // Check effect of n cells cut
2061  fhEMCALClusterCutsE[5]->Fill(clus->E());
2062 
2063  //----------------------------------------------------
2064  // Apply distance to bad channel cut
2065  //
2066  Double_t distBad = clus->GetDistanceToBadChannel() ; //Distance to bad channel
2067 
2068  if(distBad < 0.) distBad=9999. ; //workout strange convension dist = -1. ;
2069 
2070  if(distBad < fEMCALBadChMinDist)
2071  {
2072  AliDebug(2, Form("Cluster close to bad, dist %2.2f < %2.2f",distBad,fEMCALBadChMinDist));
2073  return ;
2074  }
2075 
2076  // Check effect distance to bad channel cut
2077  fhEMCALClusterCutsE[6]->Fill(clus->E());
2078 
2079  //------------------------------------------
2080  // Apply time cut, count EMCal BC before cut
2081  //
2082  SetEMCalEventBCcut(bc);
2083 
2084  // Shift time in case of no calibration with rough factor
2085  Double_t tofShift = tof;
2086  if(tof > 400) tofShift-=615;
2087  fhEMCALClusterTimeE->Fill(clus->E(),tofShift);
2088 
2089  if(!IsInTimeWindow(tof,clus->E()))
2090  {
2091  fNPileUpClusters++ ;
2092  if(fUseEMCALTimeCut)
2093  {
2094  AliDebug(2,Form("Out of time window E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f, time %e",
2095  fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta(),tof));
2096 
2097  return ;
2098  }
2099  }
2100  else
2102 
2103  // Check effect of time cut
2104  fhEMCALClusterCutsE[7]->Fill(clus->E());
2105 
2106 
2107  //----------------------------------------------------
2108  // Smear the SS to try to match data and simulations,
2109  // do it only for simulations.
2110  //
2111  if ( fSmearShowerShape && clus->GetNCells() > 2 )
2112  {
2113  Int_t nMaxima = GetCaloUtils()->GetNumberOfLocalMaxima(clus, GetEMCALCells());
2114 // Int_t nMaxima = clus->GetNExMax(); // For Run2
2115 
2116  if ( nMaxima >= fSmearNLMMin && nMaxima <= fSmearNLMMax )
2117  {
2118  AliDebug(2,Form("Smear shower shape - Original: %2.4f\n", clus->GetM02()));
2120  {
2121  clus->SetM02( clus->GetM02() + fRandom.Landau(0, fSmearShowerShapeWidth) );
2122  }
2124  {
2125  if(iclus%3 == 0 && clus->GetM02() > 0.1) clus->SetM02( clus->GetM02() + fRandom.Landau(0.05, fSmearShowerShapeWidth) ); //fSmearShowerShapeWidth = 0.035
2126  }
2127  else if (fSmearingFunction == kNoSmearing)
2128  {
2129  clus->SetM02( clus->GetM02() );
2130  }
2131  //clus->SetM02( fRandom.Landau(clus->GetM02(), fSmearShowerShapeWidth) );
2132  AliDebug(2,Form("Width %2.4f Smeared : %2.4f\n", fSmearShowerShapeWidth,clus->GetM02()));
2133  }
2134  }
2135 
2136  //--------------------------------------------------------
2137  // Fill the corresponding array with the selected clusters
2138  // Usually just filling EMCal array with upper or lower clusters is enough,
2139  // but maybe we want to do EMCal-DCal correlations.
2140 
2141  //if((fDebug > 2 && fMomentum.E() > 0.1) || fDebug > 10)
2142  AliDebug(2,Form("Selected clusters (EMCAL%d, DCAL%d), E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f",
2143  bEMCAL,bDCAL,fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta()));
2144 
2145 
2146  if (bEMCAL) fEMCALClusters->Add(clus);
2147  else if(bDCAL ) fDCALClusters ->Add(clus);
2148 
2149  // TODO, not sure if needed anymore
2150  if (fMixedEvent)
2151  clus->SetID(iclus) ;
2152 }
2153 
2154 //_______________________________________
2161 //_______________________________________
2163 {
2164  AliDebug(1,"Begin");
2165 
2166  // First recalibrate cells, time or energy
2167  // if(GetCaloUtils()->IsRecalibrationOn())
2168  // GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCells(GetCaloUtils()->GetEMCALGeometry(),
2169  // GetEMCALCells(),
2170  // fInputEvent->GetBunchCrossNumber());
2171 
2172  fNPileUpClusters = 0; // Init counter
2173  fNNonPileUpClusters = 0; // Init counter
2174  for(Int_t i = 0; i < 19; i++)
2175  {
2176  fEMCalBCEvent [i] = 0;
2177  fEMCalBCEventCut[i] = 0;
2178  }
2179 
2180  //Loop to select clusters in fiducial cut and fill container with aodClusters
2181  if(fEMCALClustersListName=="")
2182  {
2183  Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
2184  for (Int_t iclus = 0; iclus < nclusters; iclus++)
2185  {
2186  AliVCluster * clus = 0;
2187  if ( (clus = fInputEvent->GetCaloCluster(iclus)) )
2188  {
2189  if (clus->IsEMCAL())
2190  {
2191  FillInputEMCALAlgorithm(clus, iclus);
2192  }//EMCAL cluster
2193  }// cluster exists
2194  }// cluster loop
2195 
2196  //Recalculate track matching
2198 
2199  }//Get the clusters from the input event
2200  else
2201  {
2202  TClonesArray * clusterList = 0x0;
2203 
2204  if (fInputEvent->FindListObject(fEMCALClustersListName))
2205  {
2206  clusterList = dynamic_cast<TClonesArray*> (fInputEvent->FindListObject(fEMCALClustersListName));
2207  }
2208  else if(fOutputEvent)
2209  {
2210  clusterList = dynamic_cast<TClonesArray*> (fOutputEvent->FindListObject(fEMCALClustersListName));
2211  }
2212 
2213  if(!clusterList)
2214  {
2215  AliWarning(Form("Wrong name of list with clusters? <%s>",fEMCALClustersListName.Data()));
2216  return;
2217  }
2218 
2219  Int_t nclusters = clusterList->GetEntriesFast();
2220  for (Int_t iclus = 0; iclus < nclusters; iclus++)
2221  {
2222  AliVCluster * clus = dynamic_cast<AliVCluster*> (clusterList->At(iclus));
2223  //printf("E %f\n",clus->E());
2224  if (clus) FillInputEMCALAlgorithm(clus, iclus);
2225  else AliWarning("Null cluster in list!");
2226  }// cluster loop
2227 
2228  // Recalculate the pile-up time, in case long time clusters removed during clusterization
2229  //printf("Input event INIT : Pile-up clusters %d, NO pile-up %d\n",fNPileUpClusters,fNNonPileUpClusters);
2230 
2231  fNPileUpClusters = 0; // Init counter
2232  fNNonPileUpClusters = 0; // Init counter
2233  for(Int_t i = 0; i < 19; i++)
2234  {
2235  fEMCalBCEvent [i] = 0;
2236  fEMCalBCEventCut[i] = 0;
2237  }
2238 
2239  for (Int_t iclus = 0; iclus < fInputEvent->GetNumberOfCaloClusters(); iclus++)
2240  {
2241  AliVCluster * clus = 0;
2242 
2243  if ( (clus = fInputEvent->GetCaloCluster(iclus)) )
2244  {
2245  if (clus->IsEMCAL())
2246  {
2247 
2248  Float_t frac =-1;
2249  Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(fEMCALCells, clus,frac);
2250  Double_t tof = clus->GetTOF();
2251  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absIdMax,fInputEvent->GetBunchCrossNumber(),tof);
2252  //additional L1 phase shift
2253  if(GetCaloUtils()->GetEMCALRecoUtils()->IsL1PhaseInTimeRecalibrationOn())
2254  {
2255  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTimeL1Phase(GetCaloUtils()->GetEMCALGeometry()->GetSuperModuleNumber(absIdMax), fInputEvent->GetBunchCrossNumber(), tof);
2256  }
2257 
2258  tof*=1e9;
2259 
2260  //printf("Input event cluster : AbsIdMax %d, E %2.2f, time %2.2f \n", absIdMax,clus->E(),tof);
2261 
2262  //Reject clusters with bad channels, close to borders and exotic;
2263  if(!GetCaloUtils()->GetEMCALRecoUtils()->IsGoodCluster(clus,GetCaloUtils()->GetEMCALGeometry(),GetEMCALCells(),fInputEvent->GetBunchCrossNumber())) continue;
2264 
2265  Int_t bc = TMath::Nint(tof/50) + 9;
2266  SetEMCalEventBC(bc);
2267 
2268  if(fEMCALPtMin > clus->E() || fEMCALPtMax < clus->E()) continue ;
2269 
2270  clus->GetMomentum(fMomentum, fVertex[0]);
2271 
2272  if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(fMomentum.Eta(),fMomentum.Phi(),kEMCAL)) return ;
2273 
2274  SetEMCalEventBCcut(bc);
2275 
2276  if(!IsInTimeWindow(tof,clus->E()))
2277  fNPileUpClusters++ ;
2278  else
2280 
2281  }
2282  }
2283  }
2284 
2285  //printf("Input event : Pile-up clusters %d, NO pile-up %d\n",fNPileUpClusters,fNNonPileUpClusters);
2286 
2287  // Recalculate track matching, not necessary if already done in the reclusterization task.
2288  // in case it was not done ...
2290 
2291  }
2292 
2293  AliDebug(1,Form("EMCal selected clusters %d",
2294  fEMCALClusters->GetEntriesFast()));
2295  AliDebug(2,Form("\t n pile-up clusters %d, n non pile-up %d",
2297 }
2298 
2299 //_______________________________________
2301 //_______________________________________
2303 {
2304  AliDebug(1,"Begin");
2305 
2306  // Loop to select clusters in fiducial cut and fill container with aodClusters
2307  Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
2308  TString genName;
2309  for (Int_t iclus = 0; iclus < nclusters; iclus++)
2310  {
2311  AliVCluster * clus = fInputEvent->GetCaloCluster(iclus) ;
2312  if ( !clus ) continue ;
2313 
2314  if ( !clus->IsPHOS() ) continue ;
2315 
2316  if(clus->GetLabel() >=0 ) // -1 corresponds to noisy MC
2317  {
2318  if ( !AcceptParticleMCLabel(clus->GetLabel()) ) continue ;
2319  }
2320 
2321  fhPHOSClusterCutsE[0]->Fill(clus->E());
2322 
2323  // Skip CPV input
2324  if( clus->GetType() == AliVCluster::kPHOSCharged ) continue ;
2325 
2326  fhPHOSClusterCutsE[1]->Fill(clus->E());
2327 
2328  //---------------------------------------------
2329  // Old feature, try to rely on PHOS tender
2330  //
2332  {
2333  // Recalibrate the cluster energy
2335  {
2336  Float_t energy = GetCaloUtils()->RecalibrateClusterEnergy(clus, (AliAODCaloCells*)GetPHOSCells());
2337  clus->SetE(energy);
2338  }
2339  }
2340 
2341  //----------------------------------------------------------------------------------
2342  // Check if the cluster contains any bad channel and if close to calorimeter borders
2343  //
2344  // Old feature, try to rely on PHOS tender
2345  if( GetCaloUtils()->ClusterContainsBadChannel(kPHOS,clus->GetCellsAbsId(), clus->GetNCells()))
2346  continue;
2347 
2348  if(!GetCaloUtils()->CheckCellFiducialRegion(clus, fInputEvent->GetPHOSCells()))
2349  continue;
2350 
2351  // TODO, add exotic cut???
2352 
2353  fhPHOSClusterCutsE[2]->Fill(clus->E());
2354 
2355  // TODO Dead code? remove?
2356  Int_t vindex = 0 ;
2357  if (fMixedEvent)
2358  vindex = fMixedEvent->EventIndexForCaloCluster(iclus);
2359 
2360  clus->GetMomentum(fMomentum, fVertex[vindex]);
2361 
2362  //----------------------------------------------------------------------------------
2363  // Remove clusters close to borders
2364  //
2366  continue ;
2367 
2368  fhPHOSClusterCutsE[3]->Fill(clus->E());
2369 
2370  //----------------------------------------------------------------------------------
2371  // Remove clusters with too low energy
2372  //
2373  if (fPHOSPtMin > fMomentum.E() || fPHOSPtMax < fMomentum.E() )
2374  continue ;
2375 
2376  fhPHOSClusterCutsE[4]->Fill(clus->E());
2377 
2378  //----------------------------------------------------
2379  // Apply N cells cut
2380  //
2381  if(clus->GetNCells() <= fPHOSNCellsCut && fDataType != AliCaloTrackReader::kMC) return ;
2382 
2383  // Check effect of n cells cut
2384  fhPHOSClusterCutsE[5]->Fill(clus->E());
2385 
2386  //----------------------------------------------------
2387  // Apply distance to bad channel cut
2388  //
2389  Double_t distBad = clus->GetDistanceToBadChannel() ; //Distance to bad channel
2390 
2391  if(distBad < 0.) distBad=9999. ; //workout strange convension dist = -1. ;
2392 
2393  if(distBad < fPHOSBadChMinDist) return ;
2394 
2395  // Check effect distance to bad channel cut
2396  fhPHOSClusterCutsE[6]->Fill(clus->E());
2397 
2398  // TODO, add time cut
2399 
2400  //----------------------------------------------------------------------------------
2401  // Add selected clusters to array
2402  //
2403  //if(fDebug > 2 && fMomentum.E() > 0.1)
2404  AliDebug(2,Form("Selected clusters E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f",
2405  fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta()));
2406 
2407  fPHOSClusters->Add(clus);
2408 
2409  // TODO Dead code? remove?
2410  if (fMixedEvent)
2411  clus->SetID(iclus) ;
2412 
2413  } // esd/aod cluster loop
2414 
2415  AliDebug(1,Form("PHOS selected clusters %d",fPHOSClusters->GetEntriesFast())) ;
2416 }
2417 
2418 //____________________________________________
2420 //____________________________________________
2422 {
2423  fEMCALCells = fInputEvent->GetEMCALCells();
2424 }
2425 
2426 //___________________________________________
2428 //___________________________________________
2430 {
2431  fPHOSCells = fInputEvent->GetPHOSCells();
2432 }
2433 
2434 //_______________________________________
2437 //_______________________________________
2439 {
2440  AliVVZERO* v0 = fInputEvent->GetVZEROData();
2441  //printf("Init V0: ADC (%d,%d), Multiplicity (%d,%d) \n",fV0ADC[0],fV0ADC[1],fV0Mul[0],fV0Mul[1]);
2442 
2443  if (v0)
2444  {
2445  AliESDVZERO* esdV0 = dynamic_cast<AliESDVZERO*> (v0);
2446  for (Int_t i = 0; i < 32; i++)
2447  {
2448  if(esdV0)
2449  {//Only available in ESDs
2450  fV0ADC[0] += (Int_t)esdV0->GetAdcV0C(i);
2451  fV0ADC[1] += (Int_t)esdV0->GetAdcV0A(i);
2452  }
2453 
2454  fV0Mul[0] += (Int_t)v0->GetMultiplicityV0C(i);
2455  fV0Mul[1] += (Int_t)v0->GetMultiplicityV0A(i);
2456  }
2457 
2458  AliDebug(1,Form("ADC (%d,%d), Multiplicity (%d,%d)",fV0ADC[0],fV0ADC[1],fV0Mul[0],fV0Mul[1]));
2459  }
2460  else
2461  {
2462  AliDebug(1,"Cannot retrieve V0 ESD! Run w/ null V0 charges");
2463  }
2464 }
2465 
2466 //_________________________________________________
2470 //_________________________________________________
2472 {
2473  AliDebug(2,"Begin");
2474 
2475  //
2476  //check if branch name is given
2477  if(!fInputNonStandardJetBranchName.Length())
2478  {
2479  fInputEvent->Print();
2480  AliFatal("No non-standard jet branch name specified. Specify among existing ones.");
2481  }
2482 
2483  fNonStandardJets = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(fInputNonStandardJetBranchName.Data()));
2484 
2485  if(!fNonStandardJets)
2486  {
2487  //check if jet branch exist; exit if not
2488  fInputEvent->Print();
2489 
2490  AliFatal(Form("%s:%d no reconstructed jet array with name %s in AOD", (char*)__FILE__,__LINE__,fInputNonStandardJetBranchName.Data()));
2491  }
2492  else
2493  {
2494  AliDebug(1,Form("AOD input jets %d", fNonStandardJets->GetEntriesFast()));
2495  }
2496 }
2497 
2498 //_________________________________________________
2502 //_________________________________________________
2504 {
2505  AliDebug(1,"Begin");
2506  //
2507  //check if branch name is given
2508  if(!fInputBackgroundJetBranchName.Length())
2509  {
2510  fInputEvent->Print();
2511 
2512  AliFatal("No background jet branch name specified. Specify among existing ones.");
2513  }
2514 
2515  fBackgroundJets = (AliAODJetEventBackground*)(fInputEvent->FindListObject(fInputBackgroundJetBranchName.Data()));
2516 
2517  if(!fBackgroundJets)
2518  {
2519  //check if jet branch exist; exit if not
2520  fInputEvent->Print();
2521 
2522  AliFatal(Form("%s:%d no reconstructed jet array with name %s in AOD", (char*)__FILE__,__LINE__,fInputBackgroundJetBranchName.Data()));
2523  }
2524  else
2525  {
2526  AliDebug(1,"FillInputBackgroundJets");
2527  fBackgroundJets->Print("");
2528  }
2529 }
2530 
2531 //____________________________________________________________________
2537 //____________________________________________________________________
2539 {
2540  // init some variables
2541  Int_t trigtimes[30], globCol, globRow,ntimes, i;
2542  Int_t absId = -1; //[100];
2543  Int_t nPatch = 0;
2544 
2545  TArrayI patches(0);
2546 
2547  // get object pointer
2548  AliVCaloTrigger *caloTrigger = GetInputEvent()->GetCaloTrigger( "EMCAL" );
2549 
2550  if(!caloTrigger)
2551  {
2552  AliError("Trigger patches input (AliVCaloTrigger) not available in data!");
2553  return patches;
2554  }
2555 
2556  //printf("CaloTrigger Entries %d\n",caloTrigger->GetEntries() );
2557 
2558  // class is not empty
2559  if( caloTrigger->GetEntries() > 0 )
2560  {
2561  // must reset before usage, or the class will fail
2562  caloTrigger->Reset();
2563 
2564  // go throuth the trigger channels
2565  while( caloTrigger->Next() )
2566  {
2567  // get position in global 2x2 tower coordinates
2568  caloTrigger->GetPosition( globCol, globRow );
2569 
2570  //L0
2571  if(IsEventEMCALL0())
2572  {
2573  // get dimension of time arrays
2574  caloTrigger->GetNL0Times( ntimes );
2575 
2576  // no L0s in this channel
2577  // presence of the channel in the iterator still does not guarantee that L0 was produced!!
2578  if( ntimes < 1 )
2579  continue;
2580 
2581  // get timing array
2582  caloTrigger->GetL0Times( trigtimes );
2583  //printf("Get L0 patch : n times %d - trigger time window %d - %d\n",ntimes, tmin,tmax);
2584 
2585  // go through the array
2586  for( i = 0; i < ntimes; i++ )
2587  {
2588  // check if in cut - 8,9 shall be accepted in 2011
2589  if( trigtimes[i] >= tmin && trigtimes[i] <= tmax )
2590  {
2591  //printf("Accepted trigger time %d \n",trigtimes[i]);
2592  //if(nTrig > 99) continue;
2593  GetCaloUtils()->GetEMCALGeometry()->GetAbsFastORIndexFromPositionInEMCAL(globCol,globRow, absId);
2594  //printf("pass the time cut globCol %d, globRow %d absId %d\n",globCol,globRow, absId);
2595  patches.Set(nPatch+1);
2596  patches.AddAt(absId,nPatch++);
2597  }
2598  } // trigger time array
2599  }//L0
2600  else if(IsEventEMCALL1()) // L1
2601  {
2602  Int_t bit = 0;
2603  caloTrigger->GetTriggerBits(bit);
2604 
2605  Int_t sum = 0;
2606  caloTrigger->GetL1TimeSum(sum);
2607  //fBitEGA-=2;
2608  Bool_t isEGA1 = ((bit >> fBitEGA ) & 0x1) && IsEventEMCALL1Gamma1() ;
2609  Bool_t isEGA2 = ((bit >> (fBitEGA+1)) & 0x1) && IsEventEMCALL1Gamma2() ;
2610  Bool_t isEJE1 = ((bit >> fBitEJE ) & 0x1) && IsEventEMCALL1Jet1 () ;
2611  Bool_t isEJE2 = ((bit >> (fBitEJE+1)) & 0x1) && IsEventEMCALL1Jet2 () ;
2612 
2613  //if((bit>> fBitEGA )&0x1) printf("Trig Bit %d - bit %d - EG1 %d - EG2 %d\n",fBitEGA ,bit,IsEventEMCALL1Gamma1(),IsEventEMCALL1Gamma2());
2614  //if((bit>>(fBitEGA+1))&0x1) printf("Trig Bit %d - bit %d - EG1 %d - EG2 %d\n",fBitEGA+1,bit,IsEventEMCALL1Gamma1(),IsEventEMCALL1Gamma2());
2615 
2616  if(!isEGA1 && !isEJE1 && !isEGA2 && !isEJE2) continue;
2617 
2618  Int_t patchsize = -1;
2619  if (isEGA1 || isEGA2) patchsize = 2;
2620  else if (isEJE1 || isEJE2) patchsize = 16;
2621 
2622  //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",
2623  // bit,sum,patchsize,isEGA1,isEGA2,isEJE1,isEJE2,fBitEGA,fBitEJE,IsEventEMCALL1Gamma(),IsEventEMCALL1Jet());
2624 
2625 
2626  // add 2x2 (EGA) or 16x16 (EJE) patches
2627  for(Int_t irow=0; irow < patchsize; irow++)
2628  {
2629  for(Int_t icol=0; icol < patchsize; icol++)
2630  {
2631  GetCaloUtils()->GetEMCALGeometry()->GetAbsFastORIndexFromPositionInEMCAL(globCol+icol,globRow+irow, absId);
2632  //printf("pass the time cut globCol %d, globRow %d absId %d\n",globCol,globRow, absId);
2633  patches.Set(nPatch+1);
2634  patches.AddAt(absId,nPatch++);
2635  }
2636  }
2637 
2638  } // L1
2639 
2640  } // trigger iterator
2641  } // go through triggers
2642 
2643  if(patches.GetSize()<=0) AliInfo(Form("No patch found! for triggers: %s and selected <%s>",
2645  //else printf(">>>>> N patches %d, test %d,first %d, last %d\n",patches.GetSize(), nPatch, patches.At(0), patches.At(patches.GetSize()-1));
2646 
2647  return patches;
2648 }
2649 
2650 //____________________________________________________________
2654 //____________________________________________________________
2656 {
2657  // Init info from previous event
2658  fTriggerClusterIndex = -1;
2659  fTriggerClusterId = -1;
2660  fTriggerClusterBC = -10000;
2661  fIsExoticEvent = kFALSE;
2662  fIsBadCellEvent = kFALSE;
2663  fIsBadMaxCellEvent = kFALSE;
2664 
2665  fIsTriggerMatch = kFALSE;
2666  fIsTriggerMatchOpenCut[0] = kFALSE;
2667  fIsTriggerMatchOpenCut[1] = kFALSE;
2668  fIsTriggerMatchOpenCut[2] = kFALSE;
2669 
2670  // Do only analysis for triggered events
2671  if(!IsEventEMCALL1() && !IsEventEMCALL0())
2672  {
2673  fTriggerClusterBC = 0;
2674  return;
2675  }
2676 
2677  //printf("***** Try to match trigger to cluster %d **** L0 %d, L1 %d\n",fTriggerPatchClusterMatch,IsEventEMCALL0(),IsEventEMCALL1());
2678 
2679  //Recover the list of clusters
2680  TClonesArray * clusterList = 0;
2681  if (fInputEvent->FindListObject(fEMCALClustersListName))
2682  {
2683  clusterList = dynamic_cast<TClonesArray*> (fInputEvent->FindListObject(fEMCALClustersListName));
2684  }
2685  else if(fOutputEvent)
2686  {
2687  clusterList = dynamic_cast<TClonesArray*> (fOutputEvent->FindListObject(fEMCALClustersListName));
2688  }
2689 
2690  // Get number of clusters and of trigger patches
2691  Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
2692  if(clusterList)
2693  nclusters = clusterList->GetEntriesFast();
2694 
2695  Int_t nPatch = patches.GetSize();
2696  Float_t exoDiffTime = GetCaloUtils()->GetEMCALRecoUtils()->GetExoticCellDiffTimeCut();
2697 
2698  //Init some variables used in the cluster loop
2699  Float_t tofPatchMax = 100000;
2700  Float_t ePatchMax =-1;
2701 
2702  Float_t tofMax = 100000;
2703  Float_t eMax =-1;
2704 
2705  Int_t clusMax =-1;
2706  Int_t idclusMax =-1;
2707  Bool_t badClMax = kFALSE;
2708  Bool_t badCeMax = kFALSE;
2709  Bool_t exoMax = kFALSE;
2710 // Int_t absIdMaxTrig= -1;
2711  Int_t absIdMaxMax = -1;
2712 
2713  Int_t nOfHighECl = 0 ;
2714 
2715  //
2716  // Check what is the trigger threshold
2717  // set minimu energym of candidate for trigger cluster
2718  //
2720 
2721  Float_t triggerThreshold = fTriggerL1EventThreshold;
2722  if(IsEventEMCALL0()) triggerThreshold = fTriggerL0EventThreshold;
2723  Float_t minE = triggerThreshold / 2.;
2724 
2725  // This method is not really suitable for JET trigger
2726  // but in case, reduce the energy cut since we do not trigger on high energy particle
2727  if(IsEventEMCALL1Jet() || minE < 1) minE = 1;
2728 
2729  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));
2730 
2731  //
2732  // Loop on the clusters, check if there is any that falls into one of the patches
2733  //
2734  for (Int_t iclus = 0; iclus < nclusters; iclus++)
2735  {
2736  AliVCluster * clus = 0;
2737  if(clusterList) clus = (AliVCluster*) clusterList->At(iclus);
2738  else clus = fInputEvent->GetCaloCluster(iclus);
2739 
2740  if ( !clus ) continue ;
2741 
2742  if ( !clus->IsEMCAL() ) continue ;
2743 
2744  //Skip clusters with too low energy to be triggering
2745  if ( clus->E() < minE ) continue ;
2746 
2747  Float_t frac = -1;
2748  Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(fInputEvent->GetEMCALCells(), clus,frac);
2749 
2750  Bool_t badCluster = GetCaloUtils()->GetEMCALRecoUtils()->ClusterContainsBadChannel(GetCaloUtils()->GetEMCALGeometry(),
2751  clus->GetCellsAbsId(),clus->GetNCells());
2752  UShort_t cellMax[] = {(UShort_t) absIdMax};
2753  Bool_t badCell = GetCaloUtils()->GetEMCALRecoUtils()->ClusterContainsBadChannel(GetCaloUtils()->GetEMCALGeometry(),cellMax,1);
2754 
2755  // if cell is bad, it can happen that time calibration is not available,
2756  // when calculating if it is exotic, this can make it to be exotic by default
2757  // open it temporarily for this cluster
2758  if(badCell)
2759  GetCaloUtils()->GetEMCALRecoUtils()->SetExoticCellDiffTimeCut(10000000);
2760 
2761  Bool_t exotic = GetCaloUtils()->GetEMCALRecoUtils()->IsExoticCluster(clus, fInputEvent->GetEMCALCells());
2762 
2763  // Set back the time cut on exotics
2764  if(badCell)
2765  GetCaloUtils()->GetEMCALRecoUtils()->SetExoticCellDiffTimeCut(exoDiffTime);
2766 
2767  // Energy threshold for exotic Ecross typically at 4 GeV,
2768  // for lower energy, check that there are more than 1 cell in the cluster
2769  if(!exotic && clus->GetNCells() < 2) exotic = kTRUE;
2770 
2771  Float_t energy = clus->E();
2772  Int_t idclus = clus->GetID();
2773 
2774  Double_t tof = clus->GetTOF();
2775  if(GetCaloUtils()->GetEMCALRecoUtils()->IsTimeRecalibrationOn() && fTriggerClusterTimeRecal){
2776  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absIdMax,fInputEvent->GetBunchCrossNumber(),tof);
2777  //additional L1 phase shift
2778  if(GetCaloUtils()->GetEMCALRecoUtils()->IsL1PhaseInTimeRecalibrationOn()){
2779  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTimeL1Phase(GetCaloUtils()->GetEMCALGeometry()->GetSuperModuleNumber(absIdMax), fInputEvent->GetBunchCrossNumber(), tof);
2780  }
2781  }
2782  tof *=1.e9;
2783 
2784  //printf("cluster %d, ID %d, E %2.2f, tof %2.2f, AbsId max %d, exotic %d, bad Cluster %d, bad Cell %d\n",
2785  // iclus,idclus, energy,tof,absIdMax, exotic, badCluster,badCell);
2786 
2787  // Find the highest energy cluster, avobe trigger threshold
2788  // in the event in case no match to trigger is found
2789  if( energy > eMax )
2790  {
2791  tofMax = tof;
2792  eMax = energy;
2793  badClMax = badCluster;
2794  badCeMax = badCell;
2795  exoMax = exotic;
2796  clusMax = iclus;
2797  idclusMax = idclus;
2798  absIdMaxMax = absIdMax;
2799  }
2800 
2801  // count the good clusters in the event avobe the trigger threshold
2802  // to check the exotic events
2803  if(!badCluster && !exotic)
2804  nOfHighECl++;
2805 
2806  // Find match to trigger
2807  if(fTriggerPatchClusterMatch && nPatch>0)
2808  {
2809  for(Int_t iabsId =0; iabsId < nPatch; iabsId++)
2810  {
2811  Int_t absIDCell[4];
2812  GetCaloUtils()->GetEMCALGeometry()->GetCellIndexFromFastORIndex(patches.At(iabsId), absIDCell);
2813  //if(tof > 75 ) printf("E %2.2f TOF %2.2f Trigger patch %d, cells : %d, %d, %d, %d\n",
2814  // clus->E(),tof,patches.At(iabsId), absIDCell[0],absIDCell[1],absIDCell[2],absIDCell[3]);
2815 
2816  for(Int_t ipatch = 0; ipatch < 4; ipatch++)
2817  {
2818  if(absIdMax == absIDCell[ipatch])
2819  {
2820  //printf("*** Patches : absId %d, E %2.1f, tof %f \n",absIdMax,clus->E(), tof);
2821  if(energy > ePatchMax)
2822  {
2823  tofPatchMax = tof;
2824  ePatchMax = energy;
2825  fIsBadCellEvent = badCluster;
2826  fIsBadMaxCellEvent = badCell;
2827  fIsExoticEvent = exotic;
2828  fTriggerClusterIndex = iclus;
2829  fTriggerClusterId = idclus;
2830  fIsTriggerMatch = kTRUE;
2831 // absIdMaxTrig = absIdMax;
2832  }
2833  }
2834  }// cell patch loop
2835  }// trigger patch loop
2836  } // Do trigger patch matching
2837 
2838  }// Cluster loop
2839 
2840  // If there was no match, assign as trigger
2841  // the highest energy cluster in the event
2842  if(!fIsTriggerMatch)
2843  {
2844  tofPatchMax = tofMax;
2845  ePatchMax = eMax;
2846  fIsBadCellEvent = badClMax;
2847  fIsBadMaxCellEvent = badCeMax;
2848  fIsExoticEvent = exoMax;
2849  fTriggerClusterIndex = clusMax;
2850  fTriggerClusterId = idclusMax;
2851  }
2852 
2853  Double_t tofPatchMaxUS = TMath::Abs(tofPatchMax);
2854 
2855  if (tofPatchMaxUS < 28 ) fTriggerClusterBC = 0 ;
2856  else if(tofPatchMaxUS < 75 ) fTriggerClusterBC = 1 ;
2857  else if(tofPatchMaxUS < 125) fTriggerClusterBC = 2 ;
2858  else if(tofPatchMaxUS < 175) fTriggerClusterBC = 3 ;
2859  else if(tofPatchMaxUS < 225) fTriggerClusterBC = 4 ;
2860  else if(tofPatchMaxUS < 275) fTriggerClusterBC = 5 ;
2861  else
2862  {
2863  //printf("AliCaloTrackReader::MatchTriggerCluster() - Large BC - tof %2.3f - Index %d\n",tofPatchMaxUS,fTriggerClusterIndex);
2865  else
2866  {
2867  fTriggerClusterIndex = -2;
2868  fTriggerClusterId = -2;
2869  }
2870  }
2871 
2872  if(tofPatchMax < 0) fTriggerClusterBC*=-1;
2873 
2874 
2875  // 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",
2876  // fTriggerClusterIndex, fTriggerClusterId,ePatchMax, tofPatchMax,
2877  // fTriggerClusterBC, fIsBadCellEvent,fIsBadMaxCellEvent,fIsExoticEvent, fIsTriggerMatch, nOfHighECl,absIdMaxMax);
2878  //
2879  // 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",
2880  // clusMax, idclusMax, eMax,tofMax, badClMax, badCeMax,exoMax);
2881 
2882  //Redo matching but open cuts
2883  if(!fIsTriggerMatch && fTriggerClusterId >= 0)
2884  {
2885  // Open time patch time
2886  TArrayI patchOpen = GetTriggerPatches(7,10);
2887 
2888  Int_t patchAbsIdOpenTime = -1;
2889  for(Int_t iabsId =0; iabsId < patchOpen.GetSize(); iabsId++)
2890  {
2891  Int_t absIDCell[4];
2892  patchAbsIdOpenTime = patchOpen.At(iabsId);
2893  GetCaloUtils()->GetEMCALGeometry()->GetCellIndexFromFastORIndex(patchAbsIdOpenTime, absIDCell);
2894  //if(tof > 75 ) printf("E %2.2f TOF %2.2f Trigger patch %d, cells : %d, %d, %d, %d\n",
2895  // clus->E(),tof,patches.At(iabsId), absIDCell[0],absIDCell[1],absIDCell[2],absIDCell[3]);
2896 
2897  for(Int_t ipatch = 0; ipatch < 4; ipatch++)
2898  {
2899  if(absIdMaxMax == absIDCell[ipatch])
2900  {
2901  fIsTriggerMatchOpenCut[0] = kTRUE;
2902  break;
2903  }
2904  }// cell patch loop
2905  }// trigger patch loop
2906 
2907  // Check neighbour patches
2908  Int_t patchAbsId = -1;
2909  Int_t globalCol = -1;
2910  Int_t globalRow = -1;
2911  GetCaloUtils()->GetEMCALGeometry()->GetFastORIndexFromCellIndex(absIdMaxMax, patchAbsId);
2912  GetCaloUtils()->GetEMCALGeometry()->GetPositionInEMCALFromAbsFastORIndex(patchAbsId,globalCol,globalRow);
2913 
2914  // Check patches with strict time cut
2915  Int_t patchAbsIdNeigh = -1;
2916  for(Int_t icol = globalCol-1; icol <= globalCol+1; icol++)
2917  {
2918  if(icol < 0 || icol > 47) continue;
2919 
2920  for(Int_t irow = globalRow; irow <= globalRow+1; irow++)
2921  {
2922  if(irow < 0 || irow > 63) continue;
2923 
2924  GetCaloUtils()->GetEMCALGeometry()->GetAbsFastORIndexFromPositionInEMCAL(icol, irow, patchAbsIdNeigh);
2925 
2926  if ( patchAbsIdNeigh < 0 ) continue;
2927 
2928  for(Int_t iabsId =0; iabsId < patches.GetSize(); iabsId++)
2929  {
2930  if(patchAbsIdNeigh == patches.At(iabsId))
2931  {
2932  fIsTriggerMatchOpenCut[1] = kTRUE;
2933  break;
2934  }
2935  }// trigger patch loop
2936 
2937  }// row
2938  }// col
2939 
2940  // Check patches with open time cut
2941  Int_t patchAbsIdNeighOpenTime = -1;
2942  for(Int_t icol = globalCol-1; icol <= globalCol+1; icol++)
2943  {
2944  if(icol < 0 || icol > 47) continue;
2945 
2946  for(Int_t irow = globalRow; irow <= globalRow+1; irow++)
2947  {
2948  if(irow < 0 || irow > 63) continue;
2949 
2950  GetCaloUtils()->GetEMCALGeometry()->GetAbsFastORIndexFromPositionInEMCAL(icol, irow, patchAbsIdNeighOpenTime);
2951 
2952  if ( patchAbsIdNeighOpenTime < 0 ) continue;
2953 
2954  for(Int_t iabsId =0; iabsId < patchOpen.GetSize(); iabsId++)
2955  {
2956  if(patchAbsIdNeighOpenTime == patchOpen.At(iabsId))
2957  {
2958  fIsTriggerMatchOpenCut[2] = kTRUE;
2959  break;
2960  }
2961  }// trigger patch loop
2962 
2963  }// row
2964  }// col
2965 
2966  // printf("No match, new match: Open time %d-%d, open Neigh %d-%d, both open %d-%d\n",fIsTriggerMatchOpenCut[0],patchAbsIdOpenTime,
2967  // fIsTriggerMatchOpenCut[1],patchAbsIdNeigh,
2968  // fIsTriggerMatchOpenCut[2],patchAbsIdNeighOpenTime);
2969 
2970  patchOpen.Reset();
2971 
2972  }// No trigger match found
2973  //printf("Trigger BC %d, Id %d, Index %d\n",fTriggerClusterBC,fTriggerClusterId,fTriggerClusterIndex);
2974 }
2975 
2976 //_________________________________________________________
2980 //_________________________________________________________
2982 {
2984  {
2985  // get object pointer
2986  AliVCaloTrigger *caloTrigger = GetInputEvent()->GetCaloTrigger( "EMCAL" );
2987 
2988  if ( fBitEGA == 6 )
2989  {
2990  if (IsEventEMCALL1Gamma1()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(1);
2991  else if(IsEventEMCALL1Gamma2()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(3);
2992  else if(IsEventEMCALL1Jet1 ()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(0);
2993  else if(IsEventEMCALL1Jet2 ()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(2);
2994 
2995  // printf("L1 trigger Threshold Jet1 %f, Gamma1 %f, Jet2 %f, Gamma2 %f\n",
2996  // 0.07874*caloTrigger->GetL1Threshold(0),
2997  // 0.07874*caloTrigger->GetL1Threshold(1),
2998  // 0.07874*caloTrigger->GetL1Threshold(2),
2999  // 0.07874*caloTrigger->GetL1Threshold(3));
3000  }
3001  else
3002  {
3003  // Old AOD data format, in such case, order of thresholds different!!!
3004  if (IsEventEMCALL1Gamma1()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(0);
3005  else if(IsEventEMCALL1Gamma2()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(2);
3006  else if(IsEventEMCALL1Jet1 ()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(1);
3007  else if(IsEventEMCALL1Jet2 ()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(3);
3008 
3009  // printf("L1 trigger Threshold Jet1 %f, Gamma1 %f, Jet2 %f, Gamma2 %f\n",
3010  // 0.07874*caloTrigger->GetL1Threshold(1),
3011  // 0.07874*caloTrigger->GetL1Threshold(0),
3012  // 0.07874*caloTrigger->GetL1Threshold(3),
3013  // 0.07874*caloTrigger->GetL1Threshold(2));
3014  }
3015  }
3016 
3017  // Set L0 threshold, if not set by user
3019  {
3020  // Revise for periods > LHC11d
3021  Int_t runNumber = fInputEvent->GetRunNumber();
3022  if (runNumber < 146861) fTriggerL0EventThreshold = 3. ; // LHC11a
3023  else if(runNumber < 154000) fTriggerL0EventThreshold = 4. ; // LHC11b,c
3024  else if(runNumber < 165000) fTriggerL0EventThreshold = 5.5; // LHC11c,d,e
3025  else if(runNumber < 194000) fTriggerL0EventThreshold = 2 ; // LHC12
3026  else if(runNumber < 197400) fTriggerL0EventThreshold = 3 ; // LHC13def
3027  else if(runNumber < 197400) fTriggerL0EventThreshold = 2 ; // LHC13g
3028  else if(runNumber < 244300) fTriggerL0EventThreshold = 5 ; // LHC15 in, phys 1, 5 in phys2
3029  else if(runNumber < 266400) fTriggerL0EventThreshold = 2.5; // LHC16ir
3030  else fTriggerL0EventThreshold = 3.5; // LHC16s
3031  }
3032 }
3033 
3034 //________________________________________________________
3036 //________________________________________________________
3037 void AliCaloTrackReader::Print(const Option_t * opt) const
3038 {
3039  if(! opt)
3040  return;
3041 
3042  printf("***** Print: %s %s ******\n", GetName(), GetTitle() ) ;
3043  printf("Task name : %s\n", fTaskName.Data()) ;
3044  printf("Data type : %d\n", fDataType) ;
3045  printf("CTS Min pT : %2.1f GeV/c\n", fCTSPtMin) ;
3046  printf("EMCAL Min pT : %2.1f GeV/c\n", fEMCALPtMin) ;
3047  printf("PHOS Min pT : %2.1f GeV/c\n", fPHOSPtMin) ;
3048  printf("CTS Max pT : %2.1f GeV/c\n", fCTSPtMax) ;
3049  printf("EMCAL Max pT : %2.1f GeV/c\n", fEMCALPtMax) ;
3050  printf("PHOS Max pT : %2.1f GeV/c\n", fPHOSPtMax) ;
3051  printf("EMCAL Bad Dist > %2.1f \n" , fEMCALBadChMinDist) ;
3052  printf("PHOS Bad Dist > %2.1f \n" , fPHOSBadChMinDist) ;
3053  printf("EMCAL N cells > %d \n" , fEMCALNCellsCut) ;
3054  printf("PHOS N cells > %d \n" , fPHOSNCellsCut) ;
3055  printf("EMCAL Time Cut: %3.1f < TOF < %3.1f\n", fEMCALTimeCutMin, fEMCALTimeCutMax);
3056  printf("Use CTS = %d\n", fFillCTS) ;
3057  printf("Use EMCAL = %d\n", fFillEMCAL) ;
3058  printf("Use DCAL = %d\n", fFillDCAL) ;
3059  printf("Use PHOS = %d\n", fFillPHOS) ;
3060  printf("Use EMCAL Cells = %d\n", fFillEMCALCells) ;
3061  printf("Use PHOS Cells = %d\n", fFillPHOSCells) ;
3062  printf("Track status = %d\n", (Int_t) fTrackStatus) ;
3063 
3064  printf("Track Mult Eta Cut = %2.2f\n", fTrackMultEtaCut) ;
3065 
3066  printf("Track Mult Pt Cuts:") ;
3067  for(Int_t icut = 0; icut < fTrackMultNPtCut; icut++) printf(" %2.2f GeV;",fTrackMultPtCut[icut]);
3068  printf(" \n") ;
3069 
3070  printf("Write delta AOD = %d\n", fWriteOutputDeltaAOD) ;
3071  printf("Recalculate Clusters = %d, E linearity = %d\n", fRecalculateClusters, fCorrectELinearity) ;
3072 
3073  printf("Use Triggers selected in SE base class %d; If not what Trigger Mask? %d; MB Trigger Mask for mixed %d \n",
3075 
3077  printf("Compare jet pt and pt hard to accept event, factor = %2.2f",fPtHardAndJetPtFactor);
3078 
3080  printf("Compare cluster pt and pt hard to accept event, factor = %2.2f",fPtHardAndClusterPtFactor);
3081 
3082  printf("Delta AOD File Name = %s\n", fDeltaAODFileName.Data()) ;
3083  printf("Centrality: Class %s, Option %d, Bin [%d,%d] \n", fCentralityClass.Data(),fCentralityOpt,fCentralityBin[0], fCentralityBin[1]) ;
3084 
3085  printf(" \n") ;
3086 }
3087 
3088 //__________________________________________
3094 //__________________________________________
3096 {
3097  // For LHC11a
3098  // Count number of cells with energy larger than 0.1 in SM3, cut on this number
3099  if ( fRemoveLEDEvents == 1 )
3100  {
3101  Int_t ncellsSM3 = 0;
3102  for(Int_t icell = 0; icell < fInputEvent->GetEMCALCells()->GetNumberOfCells(); icell++)
3103  {
3104  Int_t absID = fInputEvent->GetEMCALCells()->GetCellNumber(icell);
3105  Int_t sm = GetCaloUtils()->GetEMCALGeometry()->GetSuperModuleNumber(absID);
3106  if(fInputEvent->GetEMCALCells()->GetAmplitude(icell) > 0.1 && sm==3) ncellsSM3++;
3107  }
3108 
3109  Int_t ncellcut = 21;
3110  if(GetFiredTriggerClasses().Contains("EMC")) ncellcut = 35;
3111 
3112  if(ncellsSM3 >= ncellcut)
3113  {
3114  AliDebug(1,Form("Reject event with ncells in SM3 %d, cut %d, trig %s",
3115  ncellsSM3,ncellcut,GetFiredTriggerClasses().Data()));
3116  return kTRUE;
3117  }
3118  }
3119  // For testing
3120  // Count number of cells with energy larger than 0.1 any SM, cut on this number
3121  else if ( fRemoveLEDEvents > 1 )
3122  {
3123  Int_t ncellsSM[] = {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0};
3124 
3125  for(Int_t icell = 0; icell < fInputEvent->GetEMCALCells()->GetNumberOfCells(); icell++)
3126  {
3127  Int_t absID = fInputEvent->GetEMCALCells()->GetCellNumber(icell);
3128  Int_t sm = GetCaloUtils()->GetEMCALGeometry()->GetSuperModuleNumber(absID);
3129  if(fInputEvent->GetEMCALCells()->GetAmplitude(icell) > 0.1) ncellsSM[sm]++;
3130  }
3131 
3132  Int_t ncellcut = 21;
3133  if(GetFiredTriggerClasses().Contains("EMC")) ncellcut = 35;
3134 
3135  for(Int_t ism = 0; ism < GetCaloUtils()->GetEMCALGeometry()->GetNumberOfSuperModules(); ism++)
3136  {
3137  if(ncellsSM[ism] >= ncellcut)
3138  {
3139  AliDebug(1,Form("Reject event with ncells in SM%d %d, cut %d, trig %s",
3140  ism,ncellsSM[ism],ncellcut,GetFiredTriggerClasses().Data()));
3141 
3142  printf("Reject event with ncells in SM%d %d, cut %d, trig %s\n",
3143  ism,ncellsSM[ism],ncellcut,GetFiredTriggerClasses().Data());
3144  return kTRUE;
3145  }
3146  }
3147  }
3148 
3149  return kFALSE;
3150 }
3151 
3152 //_________________________________________________________
3156 //_________________________________________________________
3158 {
3159  if(label < 0) return ;
3160 
3161  AliAODEvent * evt = dynamic_cast<AliAODEvent*> (fInputEvent) ;
3162  if(!evt) return ;
3163 
3164  TClonesArray * arr = dynamic_cast<TClonesArray*>(evt->FindListObject("mcparticles")) ;
3165  if(!arr) return ;
3166 
3167  if(label < arr->GetEntriesFast())
3168  {
3169  AliAODMCParticle * particle = dynamic_cast<AliAODMCParticle *>(arr->At(label));
3170  if(!particle) return ;
3171 
3172  if(label == particle->Label()) return ; // label already OK
3173  //else printf("AliCaloTrackReader::RemapMCLabelForAODs() - Label %d - AOD stack %d \n",label, particle->Label());
3174  }
3175  //else printf("AliCaloTrackReader::RemapMCLabelForAODs() - Label %d > AOD labels %d \n",label, arr->GetEntriesFast());
3176 
3177  // loop on the particles list and check if there is one with the same label
3178  for(Int_t ind = 0; ind < arr->GetEntriesFast(); ind++ )
3179  {
3180  AliAODMCParticle * particle = dynamic_cast<AliAODMCParticle *>(arr->At(ind));
3181  if(!particle) continue ;
3182 
3183  if(label == particle->Label())
3184  {
3185  label = ind;
3186  //printf("AliAnalysisTaskEMCALClusterize::RemapMCLabelForAODs() - New Label Index %d \n",label);
3187  return;
3188  }
3189  }
3190 
3191  label = -1;
3192 
3193  //printf("AliCaloTrackReader::RemapMCLabelForAODs() - Label not found set to -1 \n");
3194 }
3195 
3196 //___________________________________
3198 //___________________________________
3200 {
3201  if(fCTSTracks) fCTSTracks -> Clear();
3202  if(fEMCALClusters) fEMCALClusters -> Clear("C");
3203  if(fPHOSClusters) fPHOSClusters -> Clear("C");
3204 
3205  fV0ADC[0] = 0; fV0ADC[1] = 0;
3206  fV0Mul[0] = 0; fV0Mul[1] = 0;
3207 
3208  if(fNonStandardJets) fNonStandardJets -> Clear("C");
3209  fBackgroundJets->Reset();
3210 }
3211 
3212 //___________________________________________
3216 //___________________________________________
3218 {
3219  fEventTrigMinBias = kFALSE;
3220  fEventTrigCentral = kFALSE;
3221  fEventTrigSemiCentral = kFALSE;
3222  fEventTrigEMCALL0 = kFALSE;
3223  fEventTrigEMCALL1Gamma1 = kFALSE;
3224  fEventTrigEMCALL1Gamma2 = kFALSE;
3225  fEventTrigEMCALL1Jet1 = kFALSE;
3226  fEventTrigEMCALL1Jet2 = kFALSE;
3227 
3228  AliDebug(1,Form("Select trigger mask bit %d - Trigger Event %s - Select <%s>",
3230 
3231  if(fEventTriggerMask <=0 )// in case no mask set
3232  {
3233  // EMC triggered event? Which type?
3234  if( GetFiredTriggerClasses().Contains("-B-") || GetFiredTriggerClasses().Contains("-S-") || GetFiredTriggerClasses().Contains("-I-") )
3235  {
3236  if ( GetFiredTriggerClasses().Contains("EGA" ) ||
3237  GetFiredTriggerClasses().Contains("EG1" ) )
3238  {
3239  fEventTrigEMCALL1Gamma1 = kTRUE;
3240  if( GetFiredTriggerClasses().Contains("EG1" ) && !fFiredTriggerClassName.Contains("EG1") ) fEventTrigEMCALL1Gamma1 = kFALSE;
3241  }
3242  else if( GetFiredTriggerClasses().Contains("EG2" ) )
3243  {
3244  fEventTrigEMCALL1Gamma2 = kTRUE;
3245  if( !fFiredTriggerClassName.Contains("EG2") ) fEventTrigEMCALL1Gamma2 = kFALSE;
3246  }
3247  else if( GetFiredTriggerClasses().Contains("EJE" ) ||
3248  GetFiredTriggerClasses().Contains("EJ1" ) )
3249  {
3250  fEventTrigEMCALL1Jet1 = kTRUE;
3251  if( GetFiredTriggerClasses().Contains("EJ1" ) && !fFiredTriggerClassName.Contains("EJ1") )
3252  fEventTrigEMCALL1Jet1 = kFALSE;
3253  }
3254  else if( GetFiredTriggerClasses().Contains("EJ2" ) )
3255  {
3256  fEventTrigEMCALL1Jet2 = kTRUE;
3257  if( !fFiredTriggerClassName.Contains("EJ2") ) fEventTrigEMCALL1Jet2 = kFALSE;
3258  }
3259  else if( GetFiredTriggerClasses().Contains("CEMC") &&
3260  !GetFiredTriggerClasses().Contains("EGA" ) &&
3261  !GetFiredTriggerClasses().Contains("EJE" ) &&
3262  !GetFiredTriggerClasses().Contains("EG1" ) &&
3263  !GetFiredTriggerClasses().Contains("EJ1" ) &&
3264  !GetFiredTriggerClasses().Contains("EG2" ) &&
3265  !GetFiredTriggerClasses().Contains("EJ2" ) )
3266  {
3267  fEventTrigEMCALL0 = kTRUE;
3268  }
3269 
3270  //Min bias event trigger?
3271  if (GetFiredTriggerClasses().Contains("CCENT_R2-B-NOPF-ALLNOTRD"))
3272  {
3273  fEventTrigCentral = kTRUE;
3274  }
3275  else if(GetFiredTriggerClasses().Contains("CSEMI_R1-B-NOPF-ALLNOTRD"))
3276  {
3277  fEventTrigSemiCentral = kTRUE;
3278  }
3279  else if((GetFiredTriggerClasses().Contains("CINT") || GetFiredTriggerClasses().Contains("CPBI2_B1") ) &&
3280  GetFiredTriggerClasses().Contains("-NOPF-ALLNOTRD") )
3281  {
3282  fEventTrigMinBias = kTRUE;
3283  }
3284  }
3285  }
3286  else
3287  {
3288  // EMC L1 Gamma
3289  if ( fEventTriggerMask & AliVEvent::kEMCEGA )
3290  {
3291  //printf("EGA trigger bit\n");
3292  if (GetFiredTriggerClasses().Contains("EG"))
3293  {
3294  if (GetFiredTriggerClasses().Contains("EGA")) fEventTrigEMCALL1Gamma1 = kTRUE;
3295  else
3296  {
3297  if(GetFiredTriggerClasses().Contains("EG1")) fEventTrigEMCALL1Gamma1 = kTRUE;
3298  if(GetFiredTriggerClasses().Contains("EG2")) fEventTrigEMCALL1Gamma2 = kTRUE;
3299  }
3300  }
3301  }
3302  // EMC L1 Jet
3303  else if( fEventTriggerMask & AliVEvent::kEMCEJE )
3304  {
3305  //printf("EGA trigger bit\n");
3306  if (GetFiredTriggerClasses().Contains("EJ"))
3307  {
3308  if (GetFiredTriggerClasses().Contains("EJE")) fEventTrigEMCALL1Jet1 = kTRUE;
3309  else
3310  {
3311  if(GetFiredTriggerClasses().Contains("EJ1")) fEventTrigEMCALL1Jet1 = kTRUE;
3312  if(GetFiredTriggerClasses().Contains("EJ2")) fEventTrigEMCALL1Jet2 = kTRUE;
3313  }
3314  }
3315  }
3316  // EMC L0
3317  else if((fEventTriggerMask & AliVEvent::kEMC7) ||
3318  (fEventTriggerMask & AliVEvent::kEMC1) )
3319  {
3320  //printf("L0 trigger bit\n");
3321  fEventTrigEMCALL0 = kTRUE;
3322  }
3323  // Min Bias Pb-Pb
3324  else if( fEventTriggerMask & AliVEvent::kCentral )
3325  {
3326  //printf("MB semi central trigger bit\n");
3327  fEventTrigSemiCentral = kTRUE;
3328  }
3329  // Min Bias Pb-Pb
3330  else if( fEventTriggerMask & AliVEvent::kSemiCentral )
3331  {
3332  //printf("MB central trigger bit\n");
3333  fEventTrigCentral = kTRUE;
3334  }
3335  // Min Bias pp, PbPb, pPb
3336  else if((fEventTriggerMask & AliVEvent::kMB ) ||
3337  (fEventTriggerMask & AliVEvent::kINT7) ||
3338  (fEventTriggerMask & AliVEvent::kINT8) ||
3339  (fEventTriggerMask & AliVEvent::kAnyINT) )
3340  {
3341  //printf("MB trigger bit\n");
3342  fEventTrigMinBias = kTRUE;
3343  }
3344  }
3345 
3346  AliDebug(1,Form("Event bits: \n \t MB %d, Cen %d, Sem %d, L0 %d, L1G1 %d, L1G2 %d, L1J1 %d, L1J2 %d",
3350 
3351  // L1 trigger bit
3352  if( fBitEGA == 0 && fBitEJE == 0 )
3353  {
3354  // Init the trigger bit once, correct depending on AliESD(AOD)CaloTrigger header version
3355 
3356  // Simpler way to do it ...
3357 // if( fInputEvent->GetRunNumber() < 172000 )
3358 // reader->SetEventTriggerL1Bit(4,5); // current LHC11 data
3359 // else
3360 // reader->SetEventTriggerL1Bit(6,8); // LHC12-13 data
3361 
3362  // Old values
3363  fBitEGA = 4;
3364  fBitEJE = 5;
3365 
3366  TFile* file = AliAnalysisManager::GetAnalysisManager()->GetTree()->GetCurrentFile();
3367 
3368  const TList *clist = file->GetStreamerInfoCache();
3369 
3370  if(clist)
3371  {
3372  TStreamerInfo *cinfo = (TStreamerInfo*)clist->FindObject("AliESDCaloTrigger");
3373  Int_t verid = 5; // newer ESD header version
3374  if(!cinfo)
3375  {
3376  cinfo = (TStreamerInfo*)clist->FindObject("AliAODCaloTrigger");
3377  verid = 3; // newer AOD header version
3378  }
3379 
3380  if(cinfo)
3381  {
3382  Int_t classversionid = cinfo->GetClassVersion();
3383  //printf("********* Header class version %d *********** \n",classversionid);
3384 
3385  if (classversionid >= verid)
3386  {
3387  fBitEGA = 6;
3388  fBitEJE = 8;
3389  }
3390  } else AliInfo("AliCaloTrackReader()::SetEventTriggerBit() - Streamer info for trigger class not available, bit not changed");
3391  } else AliInfo("AliCaloTrackReader::SetEventTriggerBit() - Streamer list not available!, bit not changed");
3392 
3393  } // set once the EJE, EGA trigger bit
3394 }
3395 
3396 //____________________________________________________________
3399 //____________________________________________________________
3400 void AliCaloTrackReader::SetInputEvent(AliVEvent* const input)
3401 {
3402  fInputEvent = input;
3403  fMixedEvent = dynamic_cast<AliMixedEvent*>(GetInputEvent()) ;
3404  if (fMixedEvent)
3405  fNMixedEvent = fMixedEvent->GetNumberOfEvents() ;
3406 
3407  //Delete previous vertex
3408  if(fVertex)
3409  {
3410  for (Int_t i = 0; i < fNMixedEvent; i++)
3411  {
3412  delete [] fVertex[i] ;
3413  }
3414  delete [] fVertex ;
3415  }
3416 
3417  fVertex = new Double_t*[fNMixedEvent] ;
3418  for (Int_t i = 0; i < fNMixedEvent; i++)
3419  {
3420  fVertex[i] = new Double_t[3] ;
3421  fVertex[i][0] = 0.0 ;
3422  fVertex[i][1] = 0.0 ;
3423  fVertex[i][2] = 0.0 ;
3424  }
3425 }
3426 
3427 
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)
virtual void FillInputVZERO()
virtual AliMCEvent * GetMC() const
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.
Double_t fTimeStampEventCTPBCCorrMax
Maximum value of time stamp corrected by CTP in run.
AliCalorimeterUtils * GetCaloUtils() const
Float_t fTimeStampEventFracMin
Minimum value of time stamp fraction event.
double Double_t
Definition: External.C:58
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:61
TObjArray * fPHOSClusters
Temporal array with PHOS CaloClusters.
TString GetGeneratorNameAndIndex(Int_t index, Int_t &genIndex) const
Double_t fTrackDCACut[3]
Remove tracks with DCA larger than cut, parameters of function stored here.
Definition: External.C:236
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
Int_t fRemoveLEDEvents
Remove events where LED was wrongly firing - only EMCAL LHC11a for this equal to 1, generalized to any SM for larger.
AliEMCALRecoUtils * GetEMCALRecoUtils() const
Int_t fSmearNLMMin
Do smearing for clusters with at least this value.
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.
TH2F * fhEMCALClusterTimeE
! Control histogram on EMCAL timing
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.
energy
Definition: HFPtSpectrum.C:44
Double_t bz
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:34
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.
virtual AliMultSelection * GetMultSelCen() const
TString fEMCALClustersListName
Alternative list of clusters produced elsewhere and not from InputEvent.
Float_t fPHOSBadChMinDist
Minimal distance to bad channel to accept cluster in PHOS, cm.
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.
TH1F * fhCTSTrackCutsPt[6]
! Control histogram on the different CTS tracks selection cuts, pT
Int_t GetDebug() const
Definition: AliAnaWeights.h:48
UInt_t fEventTriggerMask
Select this triggerered event.
virtual Int_t GetEventCentrality() const
virtual AliCentrality * GetCentrality() const
AliVCaloCells * fPHOSCells
! Temporal array with PHOS AliVCaloCells.
TList * fOutputContainer
! Output container with cut control histograms.
Bool_t IsEventEMCALL1Gamma2() const
virtual AliMixedEvent * GetMixedEvent() const
Int_t fEnergyHistogramNbins
Binning of the control histograms, min and max window.
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 AliGenEventHeader * GetGenEventHeader() const
virtual void InitParameters()
Initialize the parameters with default.
virtual void GetVertex(Double_t v[3]) const
Bool_t IsRecalibrationOn() const
Int_t fTrackMult[10]
Track multiplicity, count for different pT cuts.
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.
Float_t fTrackMultPtCut[10]
Track multiplicity and sum pt cuts list.
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.
Double_t fTimeStampEventCTPBCCorrMin
Minimum value of time stamp corrected by CTP in run.
Bool_t fFillEMCALCells
Use data from EMCAL.
Float_t fEnergyHistogramLimit[2]
Binning of the control histograms, number of bins.
Int_t fTriggerClusterId
Id of trigger cluster (cluster->GetID()).
Int_t fPHOSNCellsCut
Accept for the analysis PHOS clusters with more than fNCellsCut cells.
Float_t fCTSPtMax
pT Threshold on charged particles.
Int_t GetCocktailGeneratorAndIndex(Int_t index, TString &nameGen) const
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
TH2F * fhEMCALClusterEtaPhi
! Control histogram on EMCAL clusters acceptance, before fiducial cuts
AliEMCALGeometry * GetEMCALGeometry() const
Bool_t fIsBadMaxCellEvent
Bad cell triggered event flag, only max energy cell is bad.
int Int_t
Definition: External.C:63
void RecalculateClusterShowerShapeParameters(AliVCaloCells *cells, AliVCluster *clu)
Bool_t fSelectSPDHitTracks
Ensure that track hits SPD layers.
TString fFiredTriggerClassName
Name of trigger event type used to do the analysis.
virtual TClonesArray * GetAODMCParticles() const
Definition: External.C:204
TString GetFiredTriggerClasses() const
unsigned int UInt_t
Definition: External.C:33
TObjArray * fEMCALClusters
Temporal array with EMCAL CaloClusters.
Bool_t IsEventEMCALL1Jet() const
Float_t fPHOSPtMax
pT Threshold on phos clusters.
Bool_t fAccessTrackTOF
Access the track TOF, in case of problems when accessing GetTOFBunchCrossing.
float Float_t
Definition: External.C:68
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.
Int_t fNMCGenerToAccept
Number of MC generators that should not be included in analysis.
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()
void RecalculateClusterTrackMatching(AliVEvent *event, TObjArray *clusterArray=0x0, AliMCEvent *mc=0x0)
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.
virtual TObjString * GetListOfParameters()
Save parameters used for analysis in a string.
Float_t fTrackMultEtaCut
Track multiplicity eta cut.
TString fMCGenerToAccept[5]
List with name of generators that should not be included.
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
virtual Int_t GetTrackID(AliVTrack *track)
TArrayI fAcceptEventsWithBit
Accept events if trigger bit is on.
unsigned long ULong_t
Definition: External.C:38
TH1F * fhEMCALClusterCutsE[8]
! Control histogram on the different EMCal cluster selection cuts, E
virtual void Print(const Option_t *opt) const
Print some relevant parameters set for the analysis.
void MatchTriggerCluster(TArrayI patches)
Bool_t fTimeStampEventSelect
Select events within a fraction of data taking time.
Double_t fTimeStampRunMax
Maximum value of time stamp in run.
Double_t fEMCALTimeCutMin
Remove clusters/cells with time smaller than this value, in ns.
Bool_t fEventTrigEMCALL1Jet2
Event is L1-Gamma, threshold 2 on its name, it should correspond kEMCEGA.
void SetEMCalEventBCcut(Int_t bc)
Int_t fNPileUpClustersCut
Cut to select event as pile-up.
Bool_t fCheckFidCut
Do analysis for clusters in defined region.
virtual ~AliCaloTrackReader()
Destructor.
Bool_t fEventTrigEMCALL0
Event is EMCal L0 on its name, it should correspond to AliVEvent::kEMC7, AliVEvent::kEMC1.
Int_t GetNumberOfLocalMaxima(AliVCluster *cluster, AliVCaloCells *cells)
Find the number of local maxima in cluster.
Bool_t IsWeightSettingOn() const
Definition: AliAnaWeights.h:78
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.
Bool_t Data(TH1F *h, Double_t *rangefit, Bool_t writefit, Double_t &sgn, Double_t &errsgn, Double_t &bkg, Double_t &errbkg, Double_t &sgnf, Double_t &errsgnf, Double_t &sigmafit, Int_t &status)
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)
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.
Int_t GetNMaskCellColumns() const
AliFiducialCut * fFiducialCut
Acceptance cuts.
Bool_t fRejectEMCalTriggerEventsWith2Tresholds
Reject events EG2 also triggered by EG1 or EJ2 also triggered by EJ1.
TH1F * fhPHOSClusterCutsE[7]
! Control histogram on the different PHOS cluster selection cuts, E
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
Int_t fEMCALNCellsCut
Accept for the analysis EMCAL clusters with more than fNCellsCut cells.
void CorrectClusterEnergy(AliVCluster *cl)
Correct cluster energy non linearity.
Float_t fCTSPtMin
pT Threshold on charged particles.
Int_t fSmearNLMMax
Do smearing for clusters with at maximum this value.
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.
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
Float_t fEMCALBadChMinDist
Minimal distance to bad channel to accept cluster in EMCal, cell units.
Bool_t fTimeStampEventCTPBCCorrExclude
Activate event selection within a range of data taking time CTP corrected. ESD only.
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.
Float_t RadToDeg(Float_t rad) const
TFile * file
TList with histograms for a given trigger.
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
void SetPythiaEventHeader(AliGenPythiaEventHeader *py)
Definition: AliAnaWeights.h:83
TObjArray * fCTSTracks
Temporal array with tracks.
Int_t fTriggerPatchTimeWindow[2]
Trigger patch selection window.
unsigned short UShort_t
Definition: External.C:28
Bool_t AcceptParticleMCLabel(Int_t mcLabel) const
AliMCEvent * fMC
! Monte Carlo Event Handler.
TH1I * fhNEventsAfterCut
! Each bin represents number of events resulting after a given selection cut: vertex, trigger, ...
Float_t fPtHardAndJetPtFactor
Factor between ptHard and jet pT to reject/accept event.
TString fMCGenerEventHeaderToAccept
Accept events that contain at least this event header name.
const char Option_t
Definition: External.C:48
TString fDeltaAODFileName
Delta AOD file name.
Int_t GetVertexBC() const
Float_t fTrackSumPt[10]
Track sum pT, count for different pT cuts.
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.
bool Bool_t
Definition: External.C:53
Float_t fSmearShowerShapeWidth
Smear shower shape landau function "width" (use in MC).
Int_t fBitEGA
Trigger bit on VCaloTrigger for EGA.
Int_t fTrackMultNPtCut
Track multiplicty, number of pt cuts.
virtual void FillInputBackgroundJets()
Bool_t fDoV0ANDEventSelection
Select events depending on V0AND.
AliAODJetEventBackground * fBackgroundJets
! Background jets.
Bool_t fCorrectELinearity
Correct cluster linearity, always on.
virtual TList * GetCreateControlHistograms()
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.
Float_t GetPhi(Float_t phi) const
Shift phi angle in case of negative value 360 degrees. Example TLorenzVector::Phi defined in -pi to p...
Int_t fDataType
Select MC: Kinematics, Data: ESD/AOD, MCData: Both.
virtual Bool_t SelectTrack(AliVTrack *, Double_t *)
Int_t GetMaxEnergyCell(AliVCaloCells *cells, AliVCluster *clu, Float_t &fraction) const
For a given CaloCluster, it gets the absId of the cell with maximum energy deposit.
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.
Int_t fMCGenerIndexToAccept[5]
List with index of generators that should not be included.
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.
virtual TObjArray * GetCTSTracks() const
Bool_t fUseAliCentrality
Select as centrality estimator AliCentrality (Run1) or AliMultSelection (Run1 and Run2) ...
Int_t fSmearingFunction
Choice of smearing function. 0 no smearing. 1 smearing from Gustavo (Landau center at 0)...
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 reject
Bool_t CheckCellFiducialRegion(AliVCluster *cluster, AliVCaloCells *cells) const
Bool_t IsEventEMCALL0() const