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