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