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