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