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