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