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