AliPhysics  c69cd47 (c69cd47)
 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 //_____________________________
955 //_____________________________
957 {
959  {
960  AliInfo("Cannot access stack and mcparticles at the same time, change them");
961  fReadStack = kFALSE;
962  fReadAODMCParticles = kFALSE;
963  }
964 
965  // Activate debug level in AliAnaWeights
966  if( fWeightUtils->GetDebug() >= 0 )
967  (AliAnalysisManager::GetAnalysisManager())->AddClassDebug(fWeightUtils->ClassName(), fWeightUtils->GetDebug());
968 }
969 
970 //_______________________________________
972 //_______________________________________
974 {
975  fDataType = kESD ;
976  fCTSPtMin = 0.1 ;
977  fEMCALPtMin = 0.1 ;
978  fPHOSPtMin = 0.1 ;
979  fCTSPtMax = 1000. ;
980  fEMCALPtMax = 1000. ;
981  fPHOSPtMax = 1000. ;
982 
983  fEMCALBadChMinDist = 0; // open, 2; // standard
984  fPHOSBadChMinDist = 0; // open, 2; // standard
985 
986  fEMCALNCellsCut = 0; // open, 1; // standard
987  fPHOSNCellsCut = 0; // open, 2; // standard
988 
989  //Track DCA cuts
990  // dca_xy cut = 0.0105+0.0350/TMath::Power(pt,1.1);
991  fTrackDCACut[0] = 0.0105;
992  fTrackDCACut[1] = 0.0350;
993  fTrackDCACut[2] = 1.1;
994 
995  //Do not filter the detectors input by default.
996  fFillEMCAL = kFALSE;
997  fFillDCAL = kFALSE;
998  fFillPHOS = kFALSE;
999  fFillCTS = kFALSE;
1000  fFillEMCALCells = kFALSE;
1001  fFillPHOSCells = kFALSE;
1002 
1003  fReadStack = kFALSE; // Check in the constructor of the other readers if it was set or in the configuration file
1004  fReadAODMCParticles = kFALSE; // Check in the constructor of the other readers if it was set or in the configuration file
1005  fDeltaAODFileName = "deltaAODPartCorr.root";
1007  fEventTriggerMask = AliVEvent::kAny;
1008  fMixEventTriggerMask = AliVEvent::kAnyINT;
1009  fEventTriggerAtSE = kTRUE; // Use only events that pass event selection at SE base class
1010 
1011  fAcceptFastCluster = kTRUE;
1012  fEventType = -1;
1013 
1014  //We want tracks fitted in the detectors:
1015  //fTrackStatus=AliESDtrack::kTPCrefit;
1016  //fTrackStatus|=AliESDtrack::kITSrefit;
1017  fTrackStatus = 0;
1018 
1019  fV0ADC[0] = 0; fV0ADC[1] = 0;
1020  fV0Mul[0] = 0; fV0Mul[1] = 0;
1021 
1022  fZvtxCut = 10.;
1023 
1024  fNMixedEvent = 1;
1025 
1026  fPtHardAndJetPtFactor = 7.;
1028 
1029  //Centrality
1030  fUseAliCentrality = kFALSE;
1031  fCentralityClass = "V0M";
1032  fCentralityOpt = 100;
1033  fCentralityBin[0] = fCentralityBin[1]=-1;
1034 
1035  fEventPlaneMethod = "V0";
1036 
1037  // Allocate memory (not sure this is the right place)
1038  fCTSTracks = new TObjArray();
1039  fEMCALClusters = new TObjArray();
1040  fDCALClusters = new TObjArray();
1041  fPHOSClusters = new TObjArray();
1042  //fTriggerAnalysis = new AliTriggerAnalysis;
1043  fAODBranchList = new TList ;
1044  fOutputContainer = new TList ;
1045 
1046  fEnergyHistogramNbins = 200;
1047  fEnergyHistogramLimit[0] = 0 ;
1048  fEnergyHistogramLimit[1] = 100;
1049 
1050  fPileUpParamSPD[0] = 3 ; fPileUpParamSPD[1] = 0.8 ;
1051  fPileUpParamSPD[2] = 3.0 ; fPileUpParamSPD[3] = 2.0 ; fPileUpParamSPD[4] = 5.0;
1052 
1053  // Parametrized time cut (LHC11d)
1056 
1057  // Parametrized time cut (LHC11c)
1058  //fEMCALParamTimeCutMin[0] =-5; fEMCALParamTimeCutMin[1] =-1 ; fEMCALParamTimeCutMin[2] = 1.87; fEMCALParamTimeCutMin[3] = 0.4;
1059  //fEMCALParamTimeCutMax[0] = 3.5; fEMCALParamTimeCutMax[1] = 50; fEMCALParamTimeCutMax[2] = 0.15; fEMCALParamTimeCutMax[3] = 1.6;
1060 
1061  fTimeStampRunMin = -1;
1062  fTimeStampRunMax = 1e12;
1065 
1066  for(Int_t i = 0; i < 19; i++)
1067  {
1068  fEMCalBCEvent [i] = 0;
1069  fEMCalBCEventCut[i] = 0;
1070  fTrackBCEvent [i] = 0;
1071  fTrackBCEventCut[i] = 0;
1072  }
1073 
1074  // Trigger match-rejection
1075  fTriggerPatchTimeWindow[0] = 8;
1076  fTriggerPatchTimeWindow[1] = 9;
1077 
1078  fTriggerClusterBC = -10000 ;
1081  fTriggerClusterIndex = -1;
1082  fTriggerClusterId = -1;
1083 
1084  //Jets
1087  if(!fNonStandardJets) fNonStandardJets = new TClonesArray("AliAODJet",100);
1090  if(!fBackgroundJets) fBackgroundJets = new AliAODJetEventBackground();
1091 
1092  fSmearShowerShapeWidth = 0.005;
1093 
1094  fWeightUtils = new AliAnaWeights() ;
1095  fEventWeight = 1 ;
1096 
1097 }
1098 
1099 //___________________________________________________
1104 //___________________________________________________
1106 {
1107  AliGenHijingEventHeader* hijingHeader = dynamic_cast<AliGenHijingEventHeader *> (GetGenEventHeader());
1108 
1109  //printf("header %p, label %d\n",hijingHeader,label);
1110 
1111  if(!hijingHeader || label < 0 ) return kFALSE;
1112 
1113 
1114  //printf("pass a), N produced %d\n",nproduced);
1115 
1116  if(label >= fNMCProducedMin && label < fNMCProducedMax)
1117  {
1118  //printf(" accept!, label is smaller than produced, N %d\n",nproduced);
1119 
1120  return kTRUE;
1121  }
1122 
1123  if(ReadStack())
1124  {
1125  if(!GetStack()) return kFALSE;
1126 
1127  Int_t nprimaries = GetStack()->GetNtrack();
1128 
1129  if(label > nprimaries) return kFALSE;
1130 
1131  TParticle * mom = GetStack()->Particle(label);
1132 
1133  Int_t iMom = label;
1134  Int_t iParent = mom->GetFirstMother();
1135  while(iParent!=-1)
1136  {
1137  if(iParent >= fNMCProducedMin && iParent < fNMCProducedMax)
1138  {
1139  //printf("\t accept, mother is %d \n",iParent)
1140  return kTRUE;
1141  }
1142 
1143  iMom = iParent;
1144  mom = GetStack()->Particle(iMom);
1145  iParent = mom->GetFirstMother();
1146  }
1147 
1148  return kFALSE ;
1149 
1150  } // ESD
1151  else
1152  {
1153  TClonesArray* mcparticles = GetAODMCParticles();
1154 
1155  if(!mcparticles) return kFALSE;
1156 
1157  Int_t nprimaries = mcparticles->GetEntriesFast();
1158 
1159  if(label > nprimaries) return kFALSE;
1160 
1161  //printf("pass b) N primaries %d \n",nprimaries);
1162 
1163  if(label >= fNMCProducedMin && label < fNMCProducedMax)
1164  {
1165  return kTRUE;
1166  }
1167 
1168  // Find grand parent, check if produced in the good range
1169  AliAODMCParticle * mom = (AliAODMCParticle *) mcparticles->At(label);
1170 
1171  Int_t iMom = label;
1172  Int_t iParent = mom->GetMother();
1173  while(iParent!=-1)
1174  {
1175  if(iParent >= fNMCProducedMin && iParent < fNMCProducedMax)
1176  {
1177  //printf("\t accept, mother is %d, with nProduced %d \n",iParent, nproduced);
1178  return kTRUE;
1179  }
1180 
1181  iMom = iParent;
1182  mom = (AliAODMCParticle *) mcparticles->At(iMom);
1183  iParent = mom->GetMother();
1184 
1185  }
1186 
1187  //printf("pass c), no match found \n");
1188 
1189  return kFALSE ;
1190 
1191  }//AOD
1192 }
1193 
1194 //__________________________________________________________________________
1197 //__________________________________________________________________________
1198 Bool_t AliCaloTrackReader::IsInTimeWindow(Double_t tof, Float_t energy) const
1199 {
1200  // Parametrized cut depending on E
1201  if(fUseParamTimeCut)
1202  {
1203  Float_t minCut= fEMCALParamTimeCutMin[0]+fEMCALParamTimeCutMin[1]*TMath::Exp(-(energy-fEMCALParamTimeCutMin[2])/fEMCALParamTimeCutMin[3]);
1204  Float_t maxCut= fEMCALParamTimeCutMax[0]+fEMCALParamTimeCutMax[1]*TMath::Exp(-(energy-fEMCALParamTimeCutMax[2])/fEMCALParamTimeCutMax[3]);
1205  //printf("tof %f, minCut %f, maxCut %f\n",tof,minCut,maxCut);
1206  if( tof < minCut || tof > maxCut ) return kFALSE ;
1207  }
1208 
1209  //In any case, the time should to be larger than the fixed window ...
1210  if( tof < fEMCALTimeCutMin || tof > fEMCALTimeCutMax ) return kFALSE ;
1211 
1212  return kTRUE ;
1213 }
1214 
1215 //________________________________________________
1218 //________________________________________________
1220 {
1221  return fInputEvent->IsPileupFromSPD((Int_t) fPileUpParamSPD[0] , fPileUpParamSPD[1] ,
1222  fPileUpParamSPD[2] , fPileUpParamSPD[3] , fPileUpParamSPD[4] );
1223  //printf("Param : %d, %2.2f, %2.2f, %2.2f, %2.2f\n",(Int_t) fPileUpParamSPD[0], fPileUpParamSPD[1], fPileUpParamSPD[2], fPileUpParamSPD[3], fPileUpParamSPD[4]);
1224 }
1225 
1226 //__________________________________________________
1228 //__________________________________________________
1230 {
1231  if(fNPileUpClusters > fNPileUpClustersCut) return kTRUE ;
1232  else return kFALSE;
1233 }
1234 
1235 //________________________________________________________
1237 //________________________________________________________
1239 {
1240  if( IsPileUpFromSPD() && IsPileUpFromEMCal()) return kTRUE ;
1241  else return kFALSE;
1242 }
1243 
1244 //_______________________________________________________
1246 //_______________________________________________________
1248 {
1249  if( IsPileUpFromSPD() || IsPileUpFromEMCal()) return kTRUE ;
1250  else return kFALSE;
1251 }
1252 
1253 //___________________________________________________________
1255 //___________________________________________________________
1257 {
1258  if( IsPileUpFromSPD() && !IsPileUpFromEMCal()) return kTRUE ;
1259  else return kFALSE;
1260 }
1261 
1262 //___________________________________________________________
1264 //___________________________________________________________
1266 {
1267  if( !IsPileUpFromSPD() && IsPileUpFromEMCal()) return kTRUE ;
1268  else return kFALSE;
1269 }
1270 
1271 //______________________________________________________________
1273 //______________________________________________________________
1275 {
1276  if( !IsPileUpFromSPD() && !IsPileUpFromEMCal()) return kTRUE ;
1277  else return kFALSE;
1278 }
1279 
1280 //___________________________________________________________________________________
1288 //___________________________________________________________________________________
1289 Bool_t AliCaloTrackReader::FillInputEvent(Int_t iEntry, const char * /*curFileName*/)
1290 {
1291  fEventNumber = iEntry;
1292  fTriggerClusterIndex = -1;
1293  fTriggerClusterId = -1;
1294  fIsTriggerMatch = kFALSE;
1295  fTriggerClusterBC = -10000;
1296  fIsExoticEvent = kFALSE;
1297  fIsBadCellEvent = kFALSE;
1298  fIsBadMaxCellEvent = kFALSE;
1299 
1300  fIsTriggerMatchOpenCut[0] = kFALSE ;
1301  fIsTriggerMatchOpenCut[1] = kFALSE ;
1302  fIsTriggerMatchOpenCut[2] = kFALSE ;
1303 
1304  //fCurrentFileName = TString(currentFileName);
1305  if(!fInputEvent)
1306  {
1307  AliInfo("Input event not available, skip event analysis");
1308  return kFALSE;
1309  }
1310 
1311  fhNEventsAfterCut->Fill(0.5);
1312 
1313  //-----------------------------------------------
1314  // Select the event depending on the trigger type
1315  // and other event characteristics
1316  // like the goodness of the EMCal trigger
1317  //-----------------------------------------------
1318 
1319  Bool_t accept = CheckEventTriggers();
1320  if(!accept) return kFALSE;
1321 
1322  AliDebug(1,"Pass Event trigger selection");
1323 
1324  //---------------------------------------------------------------------------
1325  // In case of analysis of events with jets, skip those with jet pt > 5 pt hard
1326  // To be used on for MC data in pT hard bins
1327  //---------------------------------------------------------------------------
1328 
1330  {
1331  if(!ComparePtHardAndJetPt()) return kFALSE ;
1332  AliDebug(1,"Pass Pt Hard - Jet rejection");
1333  fhNEventsAfterCut->Fill(8.5);
1334  }
1335 
1337  {
1338  if(!ComparePtHardAndClusterPt()) return kFALSE ;
1339  AliDebug(1,"Pass Pt Hard - Cluster rejection");
1340  fhNEventsAfterCut->Fill(9.5);
1341  }
1342 
1343  //------------------------------------------------------
1344  // Event rejection depending on time stamp
1345  //------------------------------------------------------
1346 
1348  {
1349  AliESDEvent* esd = dynamic_cast<AliESDEvent*> (fInputEvent);
1350  if(esd)
1351  {
1352  Int_t timeStamp = esd->GetTimeStamp();
1353  Float_t timeStampFrac = 1.*(timeStamp-fTimeStampRunMin) / (fTimeStampRunMax-fTimeStampRunMin);
1354 
1355  //printf("stamp0 %d, max0 %d, frac %f\n", timeStamp-fTimeStampRunMin,fTimeStampRunMax-fTimeStampRunMin, timeStampFrac);
1356 
1357  if(timeStampFrac < fTimeStampEventFracMin || timeStampFrac > fTimeStampEventFracMax) return kFALSE;
1358  }
1359 
1360  AliDebug(1,"Pass Time Stamp rejection");
1361  fhNEventsAfterCut->Fill(10.5);
1362  }
1363 
1364  //------------------------------------------------------
1365  // Event rejection depending on vertex, pileup, v0and
1366  //------------------------------------------------------
1367 
1368  FillVertexArray();
1369 
1370  //Reject events with Z vertex too large, only for SE analysis, if not, cut on the analysis code
1371  if(!GetMixedEvent() && TMath::Abs(fVertex[0][2]) > fZvtxCut) return kFALSE;
1372 
1373  fhNEventsAfterCut->Fill(11.5);
1374 
1376  {
1377  if( !CheckForPrimaryVertex() ) return kFALSE; // algorithm in ESD/AOD Readers
1378  if( TMath::Abs(fVertex[0][0] ) < 1.e-6 &&
1379  TMath::Abs(fVertex[0][1] ) < 1.e-6 &&
1380  TMath::Abs(fVertex[0][2] ) < 1.e-6 ) return kFALSE;
1381  }
1382 
1383  AliDebug(1,"Pass primary vertex rejection");
1384 
1385  fhNEventsAfterCut->Fill(12.5);
1386 
1387  //printf("Reader : IsPileUp %d, Multi %d\n",IsPileUpFromSPD(),fInputEvent->IsPileupFromSPDInMultBins());
1388 
1390  {
1391  // Do not analyze events with pileup
1392  Bool_t bPileup = IsPileUpFromSPD();
1393  //IsPileupFromSPDInMultBins() // method to try
1394  //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]);
1395  if(bPileup) return kFALSE;
1396 
1397  AliDebug(1,"Pass Pile-Up event rejection");
1398 
1399  fhNEventsAfterCut->Fill(13.5);
1400  }
1401 
1403  {
1404  AliVVZERO* v0 = fInputEvent->GetVZEROData();
1405 
1406  Bool_t bV0AND = ((v0->GetV0ADecision()==1) && (v0->GetV0CDecision()==1));
1407  //bV0AND = fTriggerAnalysis->IsOfflineTriggerFired((AliESDEvent*)fInputEvent, AliTriggerAnalysis::kV0AND);
1408  //printf("V0AND event? %d\n",bV0AND);
1409 
1410  if(!bV0AND)
1411  {
1412  AliDebug(1,"Reject event by V0AND");
1413  return kFALSE;
1414  }
1415 
1416  AliDebug(1,"Pass V0AND event rejection");
1417 
1418  fhNEventsAfterCut->Fill(14.5);
1419  }
1420 
1421  //------------------------------------------------------
1422  // Check if there is a centrality value, PbPb analysis,
1423  // and if a centrality bin selection is requested
1424  //------------------------------------------------------
1425 
1426  //If we need a centrality bin, we select only those events in the corresponding bin.
1427  Int_t cen = -1;
1428  if ( fCentralityBin[0] >= 0 && fCentralityBin[1] >= 0 )
1429  {
1430  cen = GetEventCentrality();
1431 
1432  if(cen > fCentralityBin[1] || cen < fCentralityBin[0]) return kFALSE; //reject events out of bin.
1433 
1434  AliDebug(1,"Pass centrality rejection");
1435 
1436  fhNEventsAfterCut->Fill(15.5);
1437  }
1438 
1439  // Recover the weight assigned to the event, if provided
1440  // right now only for pT-hard bins and centrality depedent weights
1442  {
1444 
1445  fWeightUtils->SetPythiaEventHeader(((AliGenPythiaEventHeader*)GetGenEventHeader()));
1446 
1448  }
1449  //---------------------------------------------------------------------------
1450  // In case of MC analysis, set the range of interest in the MC particles list
1451  // The particle selection is done inside the FillInputDetector() methods
1452  //---------------------------------------------------------------------------
1454 
1455  //printf("N min %d, N max %d\n",fNMCProducedMin,fNMCProducedMax);
1456 
1457  //-------------------------------------------------------
1458  // Get the main vertex BC, in case not available
1459  // it is calculated in FillCTS checking the BC of tracks
1460  //------------------------------------------------------
1461  fVertexBC = fInputEvent->GetPrimaryVertex()->GetBC();
1462 
1463  //-----------------------------------------------
1464  // Fill the arrays with cluster/tracks/cells data
1465  //-----------------------------------------------
1466 
1467  if(fFillCTS)
1468  {
1469  FillInputCTS();
1470  //Accept events with at least one track
1471  if(fTrackMult == 0 && fDoRejectNoTrackEvents) return kFALSE ;
1472 
1473  fhNEventsAfterCut->Fill(16.5);
1474 
1475  AliDebug(1,"Pass rejection of null track events");
1476  }
1477 
1479  {
1480  if(fVertexBC != 0 && fVertexBC != AliVTrack::kTOFBCNA) return kFALSE ;
1481 
1482  AliDebug(1,"Pass rejection of events with vertex at BC!=0");
1483 
1484  fhNEventsAfterCut->Fill(17.5);
1485  }
1486 
1487  if(fFillEMCALCells)
1489 
1490  if(fFillPHOSCells)
1492 
1493  if(fFillEMCAL || fFillDCAL)
1494  FillInputEMCAL();
1495 
1496  if(fFillPHOS)
1497  FillInputPHOS();
1498 
1499  FillInputVZERO();
1500 
1501  //one specified jet branch
1506 
1507  AliDebug(1,"Event accepted for analysis");
1508 
1509  return kTRUE ;
1510 }
1511 
1512 //__________________________________________________
1515 //__________________________________________________
1517 {
1518  if(fUseAliCentrality)
1519  {
1520  if ( !GetCentrality() ) return -1;
1521 
1522  if (fCentralityOpt == 100) return (Int_t) GetCentrality()->GetCentralityPercentile(fCentralityClass); // 100 bins max
1523  else if(fCentralityOpt == 10) return GetCentrality()->GetCentralityClass10(fCentralityClass);// 10 bins max
1524  else if(fCentralityOpt == 20) return GetCentrality()->GetCentralityClass5(fCentralityClass); // 20 bins max
1525  else
1526  {
1527  AliInfo(Form("Unknown centrality option %d, use 10, 20 or 100\n",fCentralityOpt));
1528  return -1;
1529  }
1530  }
1531  else
1532  {
1533  if ( !GetMultSelCen() ) return -1;
1534 
1535  return GetMultSelCen()->GetMultiplicityPercentile(fCentralityClass, kTRUE); // returns centrality only for events used in calibration
1536 
1537  // equivalent to
1538  //GetMultSelCen()->GetMultiplicityPercentile("V0M", kFALSE); // returns centrality for any event
1539  //Int_t qual = GetMultSelCen()->GetEvSelCode(); if (qual ! = 0) cent = qual;
1540  }
1541 }
1542 
1543 //_____________________________________________________
1546 //_____________________________________________________
1548 {
1549  if( !GetEventPlane() ) return -1000;
1550 
1551  Float_t ep = GetEventPlane()->GetEventplane(GetEventPlaneMethod(), GetInputEvent());
1552 
1553  if(GetEventPlaneMethod()=="Q" && (ep < 0 || ep > TMath::Pi()))
1554  {
1555  AliDebug(1,Form("Bad EP for <Q> method : %f\n",ep));
1556  return -1000;
1557  }
1558  else if(GetEventPlaneMethod().Contains("V0") )
1559  {
1560  if((ep > TMath::Pi()/2 || ep < -TMath::Pi()/2))
1561  {
1562  AliDebug(1,Form("Bad EP for <%s> method : %f\n",GetEventPlaneMethod().Data(), ep));
1563  return -1000;
1564  }
1565 
1566  ep+=TMath::Pi()/2; // put same range as for <Q> method
1567  }
1568 
1569  AliDebug(3,Form("Event plane angle %f",ep));
1570 
1571 // if(fDebug > 0 )
1572 // {
1573 // if (ep > TMath::Pi()) printf("AliCaloTrackReader::GetEventPlaneAngle() - Too large angle = %f\n",ep);
1574 // else if(ep < 0 ) printf("AliCaloTrackReader::GetEventPlaneAngle() - Negative angle = %f\n" ,ep);
1575 // }
1576 
1577  return ep;
1578 }
1579 
1580 //__________________________________________________________
1582 //__________________________________________________________
1583 void AliCaloTrackReader::GetVertex(Double_t vertex[3]) const
1584 {
1585  vertex[0] = fVertex[0][0];
1586  vertex[1] = fVertex[0][1];
1587  vertex[2] = fVertex[0][2];
1588 }
1589 
1590 //__________________________________________________________________________
1592 //__________________________________________________________________________
1593 void AliCaloTrackReader::GetVertex(Double_t vertex[3], Int_t evtIndex) const
1594 {
1595  vertex[0] = fVertex[evtIndex][0];
1596  vertex[1] = fVertex[evtIndex][1];
1597  vertex[2] = fVertex[evtIndex][2];
1598 }
1599 
1600 //________________________________________
1603 //________________________________________
1605 {
1606  // Delete previous vertex
1607  if(fVertex)
1608  {
1609  for (Int_t i = 0; i < fNMixedEvent; i++)
1610  {
1611  delete [] fVertex[i] ;
1612  }
1613  delete [] fVertex ;
1614  }
1615 
1616  fVertex = new Double_t*[fNMixedEvent] ;
1617  for (Int_t i = 0; i < fNMixedEvent; i++)
1618  {
1619  fVertex[i] = new Double_t[3] ;
1620  fVertex[i][0] = 0.0 ;
1621  fVertex[i][1] = 0.0 ;
1622  fVertex[i][2] = 0.0 ;
1623  }
1624 
1625  if (!fMixedEvent)
1626  { // Single event analysis
1627  if(fDataType!=kMC)
1628  {
1629 
1630  if(fInputEvent->GetPrimaryVertex())
1631  {
1632  fInputEvent->GetPrimaryVertex()->GetXYZ(fVertex[0]);
1633  }
1634  else
1635  {
1636  AliWarning("NULL primary vertex");
1637  fVertex[0][0]=0.; fVertex[0][1]=0.; fVertex[0][2]=0.;
1638  }//Primary vertex pointer do not exist
1639 
1640  } else
1641  {// MC read event
1642  fVertex[0][0]=0.; fVertex[0][1]=0.; fVertex[0][2]=0.;
1643  }
1644 
1645  AliDebug(1,Form("Single Event Vertex : %f,%f,%f\n",fVertex[0][0],fVertex[0][1],fVertex[0][2]));
1646 
1647  } else
1648  { // MultiEvent analysis
1649  for (Int_t iev = 0; iev < fNMixedEvent; iev++)
1650  {
1651  if (fMixedEvent->GetVertexOfEvent(iev))
1652  fMixedEvent->GetVertexOfEvent(iev)->GetXYZ(fVertex[iev]);
1653  else
1654  AliWarning("No vertex found");
1655 
1656  AliDebug(1,Form("Multi Event %d Vertex : %f,%f,%f",iev,fVertex[iev][0],fVertex[iev][1],fVertex[iev][2]));
1657  }
1658  }
1659 }
1660 
1661 //_____________________________________
1667 //_____________________________________
1669 {
1670  AliDebug(1,"Begin");
1671 
1672  Double_t pTrack[3] = {0,0,0};
1673 
1674  Int_t nTracks = fInputEvent->GetNumberOfTracks() ;
1675  fTrackMult = 0;
1676  Int_t nstatus = 0;
1677  Double_t bz = GetInputEvent()->GetMagneticField();
1678 
1679  for(Int_t i = 0; i < 19; i++)
1680  {
1681  fTrackBCEvent [i] = 0;
1682  fTrackBCEventCut[i] = 0;
1683  }
1684 
1685  Bool_t bc0 = kFALSE;
1686  if(fRecalculateVertexBC) fVertexBC = AliVTrack::kTOFBCNA;
1687 
1688  for (Int_t itrack = 0; itrack < nTracks; itrack++)
1689  {
1690  AliVTrack * track = (AliVTrack*)fInputEvent->GetTrack(itrack) ; // retrieve track from esd
1691 
1692  if(fAcceptOnlyHIJINGLabels && !IsHIJINGLabel(track->GetLabel())) continue ;
1693 
1694  fhCTSTrackCutsPt[0]->Fill(track->Pt());
1695 
1696  //Select tracks under certain conditions, TPCrefit, ITSrefit ... check the set bits
1697  ULong_t status = track->GetStatus();
1698 
1699  if (fTrackStatus && !((status & fTrackStatus) == fTrackStatus))
1700  continue ;
1701 
1702  fhCTSTrackCutsPt[1]->Fill(track->Pt());
1703 
1704  nstatus++;
1705 
1706  //-------------------------
1707  // Select the tracks depending on cuts of AOD or ESD
1708  if(!SelectTrack(track, pTrack)) continue ;
1709 
1710  fhCTSTrackCutsPt[2]->Fill(track->Pt());
1711 
1712  //-------------------------
1713  // TOF cuts
1714  Bool_t okTOF = ( (status & AliVTrack::kTOFout) == AliVTrack::kTOFout ) ;
1715  Double_t tof = -1000;
1716  Int_t trackBC = -1000 ;
1717 
1718  if(fAccessTrackTOF)
1719  {
1720  if(okTOF)
1721  {
1722  trackBC = track->GetTOFBunchCrossing(bz);
1723  SetTrackEventBC(trackBC+9);
1724 
1725  tof = track->GetTOFsignal()*1e-3;
1726 
1727  //After selecting tracks with small DCA, pointing to vertex, set vertex BC depeding on tracks BC
1729  {
1730  if (trackBC != 0 && trackBC != AliVTrack::kTOFBCNA) fVertexBC = trackBC;
1731  else if(trackBC == 0) bc0 = kTRUE;
1732  }
1733 
1734  //In any case, the time should to be larger than the fixed window ...
1735  if( fUseTrackTimeCut && (trackBC !=0 || tof < fTrackTimeCutMin || tof > fTrackTimeCutMax) )
1736  {
1737  //printf("Remove track time %f and bc = %d\n",tof,trackBC);
1738  continue ;
1739  }
1740  //else printf("Accept track time %f and bc = %d\n",tof,trackBC);
1741  }
1742  }
1743 
1744  fhCTSTrackCutsPt[3]->Fill(track->Pt());
1745 
1746  //---------------------
1747  // DCA cuts
1748  //
1749  fMomentum.SetPxPyPzE(pTrack[0],pTrack[1],pTrack[2],0);
1750 
1751  if(fUseTrackDCACut)
1752  {
1753  Float_t dcaTPC =-999;
1754  //In case of AODs, TPC tracks cannot be propagated back to primary vertex,
1755  if( fDataType == kAOD ) dcaTPC = ((AliAODTrack*) track)->DCA();
1756 
1757  //normal way to get the dca, cut on dca_xy
1758  if(dcaTPC==-999)
1759  {
1760  Double_t dca[2] = {1e6,1e6};
1761  Double_t covar[3] = {1e6,1e6,1e6};
1762  Bool_t okDCA = track->PropagateToDCA(fInputEvent->GetPrimaryVertex(),bz,100.,dca,covar);
1763  if( okDCA) okDCA = AcceptDCA(fMomentum.Pt(),dca[0]);
1764  if(!okDCA)
1765  {
1766  //printf("AliCaloTrackReader::FillInputCTS() - Reject track pt %2.2f, dca_xy %2.4f\n",fMomentum.Pt(),dca[0]);
1767  continue ;
1768  }
1769  }
1770  }// DCA cuts
1771 
1772  fhCTSTrackCutsPt[4]->Fill(track->Pt());
1773 
1774  //-------------------------
1775  // Kinematic/acceptance cuts
1776  //
1777  // Count the tracks in eta < 0.9
1778  if(TMath::Abs(track->Eta())< fTrackMultEtaCut) fTrackMult++;
1779 
1780  if(fCTSPtMin > fMomentum.Pt() || fCTSPtMax < fMomentum.Pt()) continue ;
1781 
1782  // Check effect of cuts on track BC
1783  if(fAccessTrackTOF && okTOF) SetTrackEventBCcut(trackBC+9);
1784 
1785  if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(fMomentum.Eta(),fMomentum.Phi(),kCTS)) continue;
1786 
1787  fhCTSTrackCutsPt[5]->Fill(track->Pt());
1788 
1789  // ------------------------------
1790  // Add selected tracks to array
1791  AliDebug(2,Form("Selected tracks pt %3.2f, phi %3.2f deg, eta %3.2f",
1792  fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta()));
1793 
1794  fCTSTracks->Add(track);
1795 
1796  // TODO, check if remove
1797  if (fMixedEvent) track->SetID(itrack);
1798 
1799  }// track loop
1800 
1801  if( fRecalculateVertexBC && (fVertexBC == 0 || fVertexBC == AliVTrack::kTOFBCNA))
1802  {
1803  if( bc0 ) fVertexBC = 0 ;
1804  else fVertexBC = AliVTrack::kTOFBCNA ;
1805  }
1806 
1807  AliDebug(1,Form("AOD entries %d, input tracks %d, pass status %d, multipliticy %d", fCTSTracks->GetEntriesFast(), nTracks, nstatus, fTrackMult));//fCTSTracksNormalInputEntries);
1808 }
1809 
1810 //_______________________________________________________________________________
1828 //_______________________________________________________________________________
1829 void AliCaloTrackReader::FillInputEMCALAlgorithm(AliVCluster * clus, Int_t iclus)
1830 {
1831  // Accept clusters with the proper label, TO DO, use the new more General methods!!!
1832  if(fAcceptOnlyHIJINGLabels && !IsHIJINGLabel( clus->GetLabel() )) return ;
1833 
1834  // TODO, not sure if needed anymore
1835  Int_t vindex = 0 ;
1836  if (fMixedEvent)
1837  vindex = fMixedEvent->EventIndexForCaloCluster(iclus);
1838 
1839  clus->GetMomentum(fMomentum, fVertex[vindex]);
1840 
1841  // No correction/cut applied yet
1842  fhEMCALClusterCutsE[0]->Fill(clus->E());
1843 
1844  //if( (fDebug > 2 && fMomentum.E() > 0.1) || fDebug > 10 )
1845  AliDebug(2,Form("Input cluster E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f",
1846  fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta()));
1847 
1848  //---------------------------
1849  // Embedding case
1850  // TO BE REVISED
1852  {
1853  if(clus->GetNLabels()==0 || clus->GetLabel() < 0) return;
1854  //else printf("Embedded cluster, %d, n label %d label %d \n",iclus,clus->GetNLabels(),clus->GetLabel());
1855  }
1856 
1857  //--------------------------------------
1858  // Apply some corrections in the cluster
1859  //
1861  {
1862  //Recalibrate the cluster energy
1864  {
1866 
1867  clus->SetE(energy);
1868  //printf("Recalibrated Energy %f\n",clus->E());
1869 
1872 
1873  } // recalculate E
1874 
1875  //Recalculate distance to bad channels, if new list of bad channels provided
1877 
1878  //Recalculate cluster position
1879  if(GetCaloUtils()->IsRecalculationOfClusterPositionOn())
1880  {
1882  //clus->GetPosition(pos);
1883  //printf("After Corrections: e %f, x %f, y %f, z %f\n",clus->E(),pos[0],pos[1],pos[2]);
1884  }
1885 
1886  // Recalculate TOF
1887  if(GetCaloUtils()->GetEMCALRecoUtils()->IsTimeRecalibrationOn())
1888  {
1889  Double_t tof = clus->GetTOF();
1890  Float_t frac =-1;
1891  Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(fEMCALCells, clus,frac);
1892 
1893  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absIdMax,fInputEvent->GetBunchCrossNumber(),tof);
1894 
1895  //additional L1 phase shift
1896  if(GetCaloUtils()->GetEMCALRecoUtils()->IsL1PhaseInTimeRecalibrationOn())
1897  {
1898  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTimeL1Phase(GetCaloUtils()->GetEMCALGeometry()->GetSuperModuleNumber(absIdMax),
1899  fInputEvent->GetBunchCrossNumber(), tof);
1900  }
1901 
1902  clus->SetTOF(tof);
1903 
1904  }// Time recalibration
1905  }
1906 
1907  // Check effect of corrections
1908  fhEMCALClusterCutsE[1]->Fill(clus->E());
1909 
1910  //-----------------------------------------------------------------
1911  // Reject clusters with bad channels, close to borders and exotic
1912  //
1913  Bool_t goodCluster = GetCaloUtils()->GetEMCALRecoUtils()->IsGoodCluster(clus,
1914  GetCaloUtils()->GetEMCALGeometry(),
1915  GetEMCALCells(),fInputEvent->GetBunchCrossNumber());
1916 
1917  if(!goodCluster)
1918  {
1919  //if( (fDebug > 2 && fMomentum.E() > 0.1) || fDebug > 10 )
1920  AliDebug(2,Form("Bad cluster E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f",
1921  fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta()));
1922 
1923  return;
1924  }
1925 
1926  // Check effect of bad cluster removal
1927  fhEMCALClusterCutsE[2]->Fill(clus->E());
1928 
1929  //Float_t pos[3];
1930  //clus->GetPosition(pos);
1931  //printf("Before Corrections: e %f, x %f, y %f, z %f\n",clus->E(),pos[0],pos[1],pos[2]);
1932 
1933  //--------------------------------------
1934  // Correct non linearity or smear energy
1935  //
1937  {
1939 
1940  //if( (fDebug > 5 && fMomentum.E() > 0.1) || fDebug > 10 )
1941  AliDebug(5,Form("Correct Non Lin: Old E %3.2f, New E %3.2f",
1942  fMomentum.E(),clus->E()));
1943 
1944  // In case of MC analysis, to match resolution/calibration in real data
1945  // Not needed anymore, just leave for MC studies on systematics
1946  if( GetCaloUtils()->GetEMCALRecoUtils()->IsClusterEnergySmeared() )
1947  {
1948  Float_t rdmEnergy = GetCaloUtils()->GetEMCALRecoUtils()->SmearClusterEnergy(clus);
1949 
1950  //if( (fDebug > 5 && fMomentum.E() > 0.1) || fDebug > 10 )
1951  AliDebug(5,Form("Smear energy: Old E %3.2f, New E %3.2f",clus->E(),rdmEnergy));
1952 
1953  clus->SetE(rdmEnergy);
1954  }
1955  }
1956 
1957  clus->GetMomentum(fMomentum, fVertex[vindex]);
1958 
1959  // Check effect linearity correction, energy smearing
1960  fhEMCALClusterCutsE[3]->Fill(clus->E());
1961 
1962  // Check the event BC depending on EMCal clustr before final cuts
1963  Double_t tof = clus->GetTOF()*1e9;
1964 
1965  Int_t bc = TMath::Nint(tof/50) + 9;
1966  //printf("tof %2.2f, bc+5=%d\n",tof,bc);
1967 
1968  SetEMCalEventBC(bc);
1969 
1970  //--------------------------------------
1971  // Apply some kinematical/acceptance cuts
1972  //
1973  if(fEMCALPtMin > clus->E() || fEMCALPtMax < clus->E()) return ;
1974 
1975  // Select cluster fiducial region
1976  //
1977  Bool_t bEMCAL = kFALSE;
1978  Bool_t bDCAL = kFALSE;
1979  if(fCheckFidCut)
1980  {
1981  if(fFillEMCAL && fFiducialCut->IsInFiducialCut(fMomentum.Eta(),fMomentum.Phi(),kEMCAL)) bEMCAL = kTRUE ;
1982  if(fFillDCAL && fFiducialCut->IsInFiducialCut(fMomentum.Eta(),fMomentum.Phi(),kDCAL )) bDCAL = kTRUE ;
1983  }
1984  else
1985  {
1986  bEMCAL = kTRUE;
1987  }
1988 
1989  //---------------------------------------------------------------------
1990  // Mask all cells in collumns facing ALICE thick material if requested
1991  //
1993  {
1994  Int_t absId = -1;
1995  Int_t iSupMod = -1;
1996  Int_t iphi = -1;
1997  Int_t ieta = -1;
1998  Bool_t shared = kFALSE;
1999  GetCaloUtils()->GetEMCALRecoUtils()->GetMaxEnergyCell(GetCaloUtils()->GetEMCALGeometry(), GetEMCALCells(),clus,absId,iSupMod,ieta,iphi,shared);
2000 
2001  AliDebug(2,Form("Masked collumn: cluster E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f",
2002  fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta()));
2003 
2004 
2005  if(GetCaloUtils()->MaskFrameCluster(iSupMod, ieta)) return;
2006  }
2007 
2008  // Check effect of energy and fiducial cuts
2009  fhEMCALClusterCutsE[4]->Fill(clus->E());
2010 
2011 
2012  //------------------------------------------
2013  // Apply time cut, count EMCal BC before cut
2014  //
2015  SetEMCalEventBCcut(bc);
2016 
2017  if(!IsInTimeWindow(tof,clus->E()))
2018  {
2019  fNPileUpClusters++ ;
2020  if(fUseEMCALTimeCut)
2021  {
2022  AliDebug(2,Form("Out of time window E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f, time %e",
2023  fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta(),tof));
2024 
2025  return ;
2026  }
2027  }
2028  else
2030 
2031  // Check effect of time cut
2032  fhEMCALClusterCutsE[5]->Fill(clus->E());
2033 
2034 
2035  //----------------------------------------------------
2036  // Apply N cells cut
2037  //
2038  if(clus->GetNCells() <= fEMCALNCellsCut && fDataType != AliCaloTrackReader::kMC) return ;
2039 
2040  // Check effect of n cells cut
2041  fhEMCALClusterCutsE[6]->Fill(clus->E());
2042 
2043  //----------------------------------------------------
2044  // Apply distance to bad channel cut
2045  //
2046  Double_t distBad = clus->GetDistanceToBadChannel() ; //Distance to bad channel
2047 
2048  if(distBad < 0.) distBad=9999. ; //workout strange convension dist = -1. ;
2049 
2050  if(distBad < fEMCALBadChMinDist) return ;
2051 
2052  // Check effect distance to bad channel cut
2053  fhEMCALClusterCutsE[7]->Fill(clus->E());
2054 
2055  //----------------------------------------------------
2056  // Smear the SS to try to match data and simulations,
2057  // do it only for simulations.
2058  //
2059  if( fSmearShowerShape && clus->GetNCells() > 2)
2060  {
2061  AliDebug(2,Form("Smear shower shape - Original: %2.4f\n", clus->GetM02()));
2063  {
2064  clus->SetM02( clus->GetM02() + fRandom.Landau(0, fSmearShowerShapeWidth) );
2065  }
2067  {
2068  if(iclus%3 == 0 && clus->GetM02() > 0.1) clus->SetM02( clus->GetM02() + fRandom.Landau(0.05, fSmearShowerShapeWidth) ); //fSmearShowerShapeWidth = 0.035
2069  }
2070  //clus->SetM02( fRandom.Landau(clus->GetM02(), fSmearShowerShapeWidth) );
2071  AliDebug(2,Form("Width %2.4f Smeared : %2.4f\n", fSmearShowerShapeWidth,clus->GetM02()));
2072  }
2073 
2074  //--------------------------------------------------------
2075  // Fill the corresponding array with the selected clusters
2076  // Usually just filling EMCal array with upper or lower clusters is enough,
2077  // but maybe we want to do EMCal-DCal correlations.
2078 
2079  //if((fDebug > 2 && fMomentum.E() > 0.1) || fDebug > 10)
2080  AliDebug(2,Form("Selected clusters (EMCAL%d, DCAL%d), E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f",
2081  bEMCAL,bDCAL,fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta()));
2082 
2083 
2084  if (bEMCAL) fEMCALClusters->Add(clus);
2085  else if(bDCAL ) fDCALClusters ->Add(clus);
2086 
2087  // TODO, not sure if needed anymore
2088  if (fMixedEvent)
2089  clus->SetID(iclus) ;
2090 }
2091 
2092 //_______________________________________
2099 //_______________________________________
2101 {
2102  AliDebug(1,"Begin");
2103 
2104  // First recalibrate cells, time or energy
2105  // if(GetCaloUtils()->IsRecalibrationOn())
2106  // GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCells(GetCaloUtils()->GetEMCALGeometry(),
2107  // GetEMCALCells(),
2108  // fInputEvent->GetBunchCrossNumber());
2109 
2110  fNPileUpClusters = 0; // Init counter
2111  fNNonPileUpClusters = 0; // Init counter
2112  for(Int_t i = 0; i < 19; i++)
2113  {
2114  fEMCalBCEvent [i] = 0;
2115  fEMCalBCEventCut[i] = 0;
2116  }
2117 
2118  //Loop to select clusters in fiducial cut and fill container with aodClusters
2119  if(fEMCALClustersListName=="")
2120  {
2121  Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
2122  for (Int_t iclus = 0; iclus < nclusters; iclus++)
2123  {
2124  AliVCluster * clus = 0;
2125  if ( (clus = fInputEvent->GetCaloCluster(iclus)) )
2126  {
2127  if (clus->IsEMCAL())
2128  {
2129  FillInputEMCALAlgorithm(clus, iclus);
2130  }//EMCAL cluster
2131  }// cluster exists
2132  }// cluster loop
2133 
2134  //Recalculate track matching
2136 
2137  }//Get the clusters from the input event
2138  else
2139  {
2140  TClonesArray * clusterList = 0x0;
2141 
2142  if (fInputEvent->FindListObject(fEMCALClustersListName))
2143  {
2144  clusterList = dynamic_cast<TClonesArray*> (fInputEvent->FindListObject(fEMCALClustersListName));
2145  }
2146  else if(fOutputEvent)
2147  {
2148  clusterList = dynamic_cast<TClonesArray*> (fOutputEvent->FindListObject(fEMCALClustersListName));
2149  }
2150 
2151  if(!clusterList)
2152  {
2153  AliWarning(Form("Wrong name of list with clusters? <%s>",fEMCALClustersListName.Data()));
2154  return;
2155  }
2156 
2157  Int_t nclusters = clusterList->GetEntriesFast();
2158  for (Int_t iclus = 0; iclus < nclusters; iclus++)
2159  {
2160  AliVCluster * clus = dynamic_cast<AliVCluster*> (clusterList->At(iclus));
2161  //printf("E %f\n",clus->E());
2162  if (clus) FillInputEMCALAlgorithm(clus, iclus);
2163  else AliWarning("Null cluster in list!");
2164  }// cluster loop
2165 
2166  // Recalculate the pile-up time, in case long time clusters removed during clusterization
2167  //printf("Input event INIT : Pile-up clusters %d, NO pile-up %d\n",fNPileUpClusters,fNNonPileUpClusters);
2168 
2169  fNPileUpClusters = 0; // Init counter
2170  fNNonPileUpClusters = 0; // Init counter
2171  for(Int_t i = 0; i < 19; i++)
2172  {
2173  fEMCalBCEvent [i] = 0;
2174  fEMCalBCEventCut[i] = 0;
2175  }
2176 
2177  for (Int_t iclus = 0; iclus < fInputEvent->GetNumberOfCaloClusters(); iclus++)
2178  {
2179  AliVCluster * clus = 0;
2180 
2181  if ( (clus = fInputEvent->GetCaloCluster(iclus)) )
2182  {
2183  if (clus->IsEMCAL())
2184  {
2185 
2186  Float_t frac =-1;
2187  Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(fEMCALCells, clus,frac);
2188  Double_t tof = clus->GetTOF();
2189  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absIdMax,fInputEvent->GetBunchCrossNumber(),tof);
2190  //additional L1 phase shift
2191  if(GetCaloUtils()->GetEMCALRecoUtils()->IsL1PhaseInTimeRecalibrationOn())
2192  {
2193  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTimeL1Phase(GetCaloUtils()->GetEMCALGeometry()->GetSuperModuleNumber(absIdMax), fInputEvent->GetBunchCrossNumber(), tof);
2194  }
2195 
2196  tof*=1e9;
2197 
2198  //printf("Input event cluster : AbsIdMax %d, E %2.2f, time %2.2f \n", absIdMax,clus->E(),tof);
2199 
2200  //Reject clusters with bad channels, close to borders and exotic;
2201  if(!GetCaloUtils()->GetEMCALRecoUtils()->IsGoodCluster(clus,GetCaloUtils()->GetEMCALGeometry(),GetEMCALCells(),fInputEvent->GetBunchCrossNumber())) continue;
2202 
2203  Int_t bc = TMath::Nint(tof/50) + 9;
2204  SetEMCalEventBC(bc);
2205 
2206  if(fEMCALPtMin > clus->E() || fEMCALPtMax < clus->E()) continue ;
2207 
2208  clus->GetMomentum(fMomentum, fVertex[0]);
2209 
2210  if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(fMomentum.Eta(),fMomentum.Phi(),kEMCAL)) return ;
2211 
2212  SetEMCalEventBCcut(bc);
2213 
2214  if(!IsInTimeWindow(tof,clus->E()))
2215  fNPileUpClusters++ ;
2216  else
2218 
2219  }
2220  }
2221  }
2222 
2223  //printf("Input event : Pile-up clusters %d, NO pile-up %d\n",fNPileUpClusters,fNNonPileUpClusters);
2224 
2225  // Recalculate track matching, not necessary if already done in the reclusterization task.
2226  // in case it was not done ...
2228 
2229  }
2230 
2231  AliDebug(1,Form("AOD entries %d, n pile-up clusters %d, n non pile-up %d", fEMCALClusters->GetEntriesFast(),fNPileUpClusters,fNNonPileUpClusters));
2232 }
2233 
2234 //_______________________________________
2236 //_______________________________________
2238 {
2239  AliDebug(1,"Begin");
2240 
2241  // Loop to select clusters in fiducial cut and fill container with aodClusters
2242  Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
2243  for (Int_t iclus = 0; iclus < nclusters; iclus++)
2244  {
2245  AliVCluster * clus = fInputEvent->GetCaloCluster(iclus) ;
2246  if ( !clus ) continue ;
2247 
2248  if ( !clus->IsPHOS() ) continue ;
2249 
2250  if(fAcceptOnlyHIJINGLabels && !IsHIJINGLabel(clus->GetLabel())) continue ;
2251 
2252  fhPHOSClusterCutsE[0]->Fill(clus->E());
2253 
2254  // Skip CPV input
2255  if( clus->GetType() == AliVCluster::kPHOSCharged ) continue ;
2256 
2257  fhPHOSClusterCutsE[1]->Fill(clus->E());
2258 
2259  //---------------------------------------------
2260  // Old feature, try to rely on PHOS tender
2261  //
2263  {
2264  // Recalibrate the cluster energy
2266  {
2267  Float_t energy = GetCaloUtils()->RecalibrateClusterEnergy(clus, (AliAODCaloCells*)GetPHOSCells());
2268  clus->SetE(energy);
2269  }
2270  }
2271 
2272  //----------------------------------------------------------------------------------
2273  // Check if the cluster contains any bad channel and if close to calorimeter borders
2274  //
2275  // Old feature, try to rely on PHOS tender
2276  if( GetCaloUtils()->ClusterContainsBadChannel(kPHOS,clus->GetCellsAbsId(), clus->GetNCells()))
2277  continue;
2278 
2279  if(!GetCaloUtils()->CheckCellFiducialRegion(clus, fInputEvent->GetPHOSCells()))
2280  continue;
2281 
2282  // TODO, add exotic cut???
2283 
2284  fhPHOSClusterCutsE[2]->Fill(clus->E());
2285 
2286  // TODO Dead code? remove?
2287  Int_t vindex = 0 ;
2288  if (fMixedEvent)
2289  vindex = fMixedEvent->EventIndexForCaloCluster(iclus);
2290 
2291  clus->GetMomentum(fMomentum, fVertex[vindex]);
2292 
2293  //----------------------------------------------------------------------------------
2294  // Remove clusters close to borders
2295  //
2297  continue ;
2298 
2299  fhPHOSClusterCutsE[3]->Fill(clus->E());
2300 
2301  //----------------------------------------------------------------------------------
2302  // Remove clusters with too low energy
2303  //
2304  if (fPHOSPtMin > fMomentum.E() || fPHOSPtMax < fMomentum.E() )
2305  continue ;
2306 
2307  fhPHOSClusterCutsE[4]->Fill(clus->E());
2308 
2309  //----------------------------------------------------
2310  // Apply N cells cut
2311  //
2312  if(clus->GetNCells() <= fPHOSNCellsCut && fDataType != AliCaloTrackReader::kMC) return ;
2313 
2314  // Check effect of n cells cut
2315  fhPHOSClusterCutsE[5]->Fill(clus->E());
2316 
2317  //----------------------------------------------------
2318  // Apply distance to bad channel cut
2319  //
2320  Double_t distBad = clus->GetDistanceToBadChannel() ; //Distance to bad channel
2321 
2322  if(distBad < 0.) distBad=9999. ; //workout strange convension dist = -1. ;
2323 
2324  if(distBad < fPHOSBadChMinDist) return ;
2325 
2326  // Check effect distance to bad channel cut
2327  fhPHOSClusterCutsE[6]->Fill(clus->E());
2328 
2329  // TODO, add time cut
2330 
2331  //----------------------------------------------------------------------------------
2332  // Add selected clusters to array
2333  //
2334  //if(fDebug > 2 && fMomentum.E() > 0.1)
2335  AliDebug(2,Form("Selected clusters E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f",
2336  fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta()));
2337 
2338  fPHOSClusters->Add(clus);
2339 
2340  // TODO Dead code? remove?
2341  if (fMixedEvent)
2342  clus->SetID(iclus) ;
2343 
2344  } // esd/aod cluster loop
2345 
2346  AliDebug(1,Form("AOD entries %d",fPHOSClusters->GetEntriesFast())) ;
2347 }
2348 
2349 //____________________________________________
2351 //____________________________________________
2353 {
2354  fEMCALCells = fInputEvent->GetEMCALCells();
2355 }
2356 
2357 //___________________________________________
2359 //___________________________________________
2361 {
2362  fPHOSCells = fInputEvent->GetPHOSCells();
2363 }
2364 
2365 //_______________________________________
2368 //_______________________________________
2370 {
2371  AliVVZERO* v0 = fInputEvent->GetVZEROData();
2372  //printf("Init V0: ADC (%d,%d), Multiplicity (%d,%d) \n",fV0ADC[0],fV0ADC[1],fV0Mul[0],fV0Mul[1]);
2373 
2374  if (v0)
2375  {
2376  AliESDVZERO* esdV0 = dynamic_cast<AliESDVZERO*> (v0);
2377  for (Int_t i = 0; i < 32; i++)
2378  {
2379  if(esdV0)
2380  {//Only available in ESDs
2381  fV0ADC[0] += (Int_t)esdV0->GetAdcV0C(i);
2382  fV0ADC[1] += (Int_t)esdV0->GetAdcV0A(i);
2383  }
2384 
2385  fV0Mul[0] += (Int_t)v0->GetMultiplicityV0C(i);
2386  fV0Mul[1] += (Int_t)v0->GetMultiplicityV0A(i);
2387  }
2388 
2389  AliDebug(1,Form("ADC (%d,%d), Multiplicity (%d,%d)",fV0ADC[0],fV0ADC[1],fV0Mul[0],fV0Mul[1]));
2390  }
2391  else
2392  {
2393  AliDebug(1,"Cannot retrieve V0 ESD! Run w/ null V0 charges");
2394  }
2395 }
2396 
2397 //_________________________________________________
2401 //_________________________________________________
2403 {
2404  AliDebug(2,"Begin");
2405 
2406  //
2407  //check if branch name is given
2408  if(!fInputNonStandardJetBranchName.Length())
2409  {
2410  fInputEvent->Print();
2411  AliFatal("No non-standard jet branch name specified. Specify among existing ones.");
2412  }
2413 
2414  fNonStandardJets = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(fInputNonStandardJetBranchName.Data()));
2415 
2416  if(!fNonStandardJets)
2417  {
2418  //check if jet branch exist; exit if not
2419  fInputEvent->Print();
2420 
2421  AliFatal(Form("%s:%d no reconstructed jet array with name %s in AOD", (char*)__FILE__,__LINE__,fInputNonStandardJetBranchName.Data()));
2422  }
2423  else
2424  {
2425  AliDebug(1,Form("AOD input jets %d", fNonStandardJets->GetEntriesFast()));
2426  }
2427 }
2428 
2429 //_________________________________________________
2433 //_________________________________________________
2435 {
2436  AliDebug(1,"Begin");
2437  //
2438  //check if branch name is given
2439  if(!fInputBackgroundJetBranchName.Length())
2440  {
2441  fInputEvent->Print();
2442 
2443  AliFatal("No background jet branch name specified. Specify among existing ones.");
2444  }
2445 
2446  fBackgroundJets = (AliAODJetEventBackground*)(fInputEvent->FindListObject(fInputBackgroundJetBranchName.Data()));
2447 
2448  if(!fBackgroundJets)
2449  {
2450  //check if jet branch exist; exit if not
2451  fInputEvent->Print();
2452 
2453  AliFatal(Form("%s:%d no reconstructed jet array with name %s in AOD", (char*)__FILE__,__LINE__,fInputBackgroundJetBranchName.Data()));
2454  }
2455  else
2456  {
2457  AliDebug(1,"FillInputBackgroundJets");
2458  fBackgroundJets->Print("");
2459  }
2460 }
2461 
2462 //____________________________________________________________________
2468 //____________________________________________________________________
2469 TArrayI AliCaloTrackReader::GetTriggerPatches(Int_t tmin, Int_t tmax )
2470 {
2471  // init some variables
2472  Int_t trigtimes[30], globCol, globRow,ntimes, i;
2473  Int_t absId = -1; //[100];
2474  Int_t nPatch = 0;
2475 
2476  TArrayI patches(0);
2477 
2478  // get object pointer
2479  AliVCaloTrigger *caloTrigger = GetInputEvent()->GetCaloTrigger( "EMCAL" );
2480 
2481  if(!caloTrigger)
2482  {
2483  AliError("Trigger patches input (AliVCaloTrigger) not available in data!");
2484  return patches;
2485  }
2486 
2487  //printf("CaloTrigger Entries %d\n",caloTrigger->GetEntries() );
2488 
2489  // class is not empty
2490  if( caloTrigger->GetEntries() > 0 )
2491  {
2492  // must reset before usage, or the class will fail
2493  caloTrigger->Reset();
2494 
2495  // go throuth the trigger channels
2496  while( caloTrigger->Next() )
2497  {
2498  // get position in global 2x2 tower coordinates
2499  caloTrigger->GetPosition( globCol, globRow );
2500 
2501  //L0
2502  if(IsEventEMCALL0())
2503  {
2504  // get dimension of time arrays
2505  caloTrigger->GetNL0Times( ntimes );
2506 
2507  // no L0s in this channel
2508  // presence of the channel in the iterator still does not guarantee that L0 was produced!!
2509  if( ntimes < 1 )
2510  continue;
2511 
2512  // get timing array
2513  caloTrigger->GetL0Times( trigtimes );
2514  //printf("Get L0 patch : n times %d - trigger time window %d - %d\n",ntimes, tmin,tmax);
2515 
2516  // go through the array
2517  for( i = 0; i < ntimes; i++ )
2518  {
2519  // check if in cut - 8,9 shall be accepted in 2011
2520  if( trigtimes[i] >= tmin && trigtimes[i] <= tmax )
2521  {
2522  //printf("Accepted trigger time %d \n",trigtimes[i]);
2523  //if(nTrig > 99) continue;
2524  GetCaloUtils()->GetEMCALGeometry()->GetAbsFastORIndexFromPositionInEMCAL(globCol,globRow, absId);
2525  //printf("pass the time cut globCol %d, globRow %d absId %d\n",globCol,globRow, absId);
2526  patches.Set(nPatch+1);
2527  patches.AddAt(absId,nPatch++);
2528  }
2529  } // trigger time array
2530  }//L0
2531  else if(IsEventEMCALL1()) // L1
2532  {
2533  Int_t bit = 0;
2534  caloTrigger->GetTriggerBits(bit);
2535 
2536  Int_t sum = 0;
2537  caloTrigger->GetL1TimeSum(sum);
2538  //fBitEGA-=2;
2539  Bool_t isEGA1 = ((bit >> fBitEGA ) & 0x1) && IsEventEMCALL1Gamma1() ;
2540  Bool_t isEGA2 = ((bit >> (fBitEGA+1)) & 0x1) && IsEventEMCALL1Gamma2() ;
2541  Bool_t isEJE1 = ((bit >> fBitEJE ) & 0x1) && IsEventEMCALL1Jet1 () ;
2542  Bool_t isEJE2 = ((bit >> (fBitEJE+1)) & 0x1) && IsEventEMCALL1Jet2 () ;
2543 
2544  //if((bit>> fBitEGA )&0x1) printf("Trig Bit %d - bit %d - EG1 %d - EG2 %d\n",fBitEGA ,bit,IsEventEMCALL1Gamma1(),IsEventEMCALL1Gamma2());
2545  //if((bit>>(fBitEGA+1))&0x1) printf("Trig Bit %d - bit %d - EG1 %d - EG2 %d\n",fBitEGA+1,bit,IsEventEMCALL1Gamma1(),IsEventEMCALL1Gamma2());
2546 
2547  if(!isEGA1 && !isEJE1 && !isEGA2 && !isEJE2) continue;
2548 
2549  Int_t patchsize = -1;
2550  if (isEGA1 || isEGA2) patchsize = 2;
2551  else if (isEJE1 || isEJE2) patchsize = 16;
2552 
2553  //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",
2554  // bit,sum,patchsize,isEGA1,isEGA2,isEJE1,isEJE2,fBitEGA,fBitEJE,IsEventEMCALL1Gamma(),IsEventEMCALL1Jet());
2555 
2556 
2557  // add 2x2 (EGA) or 16x16 (EJE) patches
2558  for(Int_t irow=0; irow < patchsize; irow++)
2559  {
2560  for(Int_t icol=0; icol < patchsize; icol++)
2561  {
2562  GetCaloUtils()->GetEMCALGeometry()->GetAbsFastORIndexFromPositionInEMCAL(globCol+icol,globRow+irow, absId);
2563  //printf("pass the time cut globCol %d, globRow %d absId %d\n",globCol,globRow, absId);
2564  patches.Set(nPatch+1);
2565  patches.AddAt(absId,nPatch++);
2566  }
2567  }
2568 
2569  } // L1
2570 
2571  } // trigger iterator
2572  } // go through triggers
2573 
2574  if(patches.GetSize()<=0) AliInfo(Form("No patch found! for triggers: %s and selected <%s>",
2576  //else printf(">>>>> N patches %d, test %d,first %d, last %d\n",patches.GetSize(), nPatch, patches.At(0), patches.At(patches.GetSize()-1));
2577 
2578  return patches;
2579 }
2580 
2581 //____________________________________________________________
2585 //____________________________________________________________
2587 {
2588  // Init info from previous event
2589  fTriggerClusterIndex = -1;
2590  fTriggerClusterId = -1;
2591  fTriggerClusterBC = -10000;
2592  fIsExoticEvent = kFALSE;
2593  fIsBadCellEvent = kFALSE;
2594  fIsBadMaxCellEvent = kFALSE;
2595 
2596  fIsTriggerMatch = kFALSE;
2597  fIsTriggerMatchOpenCut[0] = kFALSE;
2598  fIsTriggerMatchOpenCut[1] = kFALSE;
2599  fIsTriggerMatchOpenCut[2] = kFALSE;
2600 
2601  // Do only analysis for triggered events
2602  if(!IsEventEMCALL1() && !IsEventEMCALL0())
2603  {
2604  fTriggerClusterBC = 0;
2605  return;
2606  }
2607 
2608  //printf("***** Try to match trigger to cluster %d **** L0 %d, L1 %d\n",fTriggerPatchClusterMatch,IsEventEMCALL0(),IsEventEMCALL1());
2609 
2610  //Recover the list of clusters
2611  TClonesArray * clusterList = 0;
2612  if (fInputEvent->FindListObject(fEMCALClustersListName))
2613  {
2614  clusterList = dynamic_cast<TClonesArray*> (fInputEvent->FindListObject(fEMCALClustersListName));
2615  }
2616  else if(fOutputEvent)
2617  {
2618  clusterList = dynamic_cast<TClonesArray*> (fOutputEvent->FindListObject(fEMCALClustersListName));
2619  }
2620 
2621  // Get number of clusters and of trigger patches
2622  Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
2623  if(clusterList)
2624  nclusters = clusterList->GetEntriesFast();
2625 
2626  Int_t nPatch = patches.GetSize();
2627  Float_t exoDiffTime = GetCaloUtils()->GetEMCALRecoUtils()->GetExoticCellDiffTimeCut();
2628 
2629  //Init some variables used in the cluster loop
2630  Float_t tofPatchMax = 100000;
2631  Float_t ePatchMax =-1;
2632 
2633  Float_t tofMax = 100000;
2634  Float_t eMax =-1;
2635 
2636  Int_t clusMax =-1;
2637  Int_t idclusMax =-1;
2638  Bool_t badClMax = kFALSE;
2639  Bool_t badCeMax = kFALSE;
2640  Bool_t exoMax = kFALSE;
2641 // Int_t absIdMaxTrig= -1;
2642  Int_t absIdMaxMax = -1;
2643 
2644  Int_t nOfHighECl = 0 ;
2645 
2646  //
2647  // Check what is the trigger threshold
2648  // set minimu energym of candidate for trigger cluster
2649  //
2651 
2652  Float_t triggerThreshold = fTriggerL1EventThreshold;
2653  if(IsEventEMCALL0()) triggerThreshold = fTriggerL0EventThreshold;
2654  Float_t minE = triggerThreshold / 2.;
2655 
2656  // This method is not really suitable for JET trigger
2657  // but in case, reduce the energy cut since we do not trigger on high energy particle
2658  if(IsEventEMCALL1Jet() || minE < 1) minE = 1;
2659 
2660  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));
2661 
2662  //
2663  // Loop on the clusters, check if there is any that falls into one of the patches
2664  //
2665  for (Int_t iclus = 0; iclus < nclusters; iclus++)
2666  {
2667  AliVCluster * clus = 0;
2668  if(clusterList) clus = (AliVCluster*) clusterList->At(iclus);
2669  else clus = fInputEvent->GetCaloCluster(iclus);
2670 
2671  if ( !clus ) continue ;
2672 
2673  if ( !clus->IsEMCAL() ) continue ;
2674 
2675  //Skip clusters with too low energy to be triggering
2676  if ( clus->E() < minE ) continue ;
2677 
2678  Float_t frac = -1;
2679  Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(fInputEvent->GetEMCALCells(), clus,frac);
2680 
2681  Bool_t badCluster = GetCaloUtils()->GetEMCALRecoUtils()->ClusterContainsBadChannel(GetCaloUtils()->GetEMCALGeometry(),
2682  clus->GetCellsAbsId(),clus->GetNCells());
2683  UShort_t cellMax[] = {(UShort_t) absIdMax};
2684  Bool_t badCell = GetCaloUtils()->GetEMCALRecoUtils()->ClusterContainsBadChannel(GetCaloUtils()->GetEMCALGeometry(),cellMax,1);
2685 
2686  // if cell is bad, it can happen that time calibration is not available,
2687  // when calculating if it is exotic, this can make it to be exotic by default
2688  // open it temporarily for this cluster
2689  if(badCell)
2690  GetCaloUtils()->GetEMCALRecoUtils()->SetExoticCellDiffTimeCut(10000000);
2691 
2692  Bool_t exotic = GetCaloUtils()->GetEMCALRecoUtils()->IsExoticCluster(clus, fInputEvent->GetEMCALCells());
2693 
2694  // Set back the time cut on exotics
2695  if(badCell)
2696  GetCaloUtils()->GetEMCALRecoUtils()->SetExoticCellDiffTimeCut(exoDiffTime);
2697 
2698  // Energy threshold for exotic Ecross typically at 4 GeV,
2699  // for lower energy, check that there are more than 1 cell in the cluster
2700  if(!exotic && clus->GetNCells() < 2) exotic = kTRUE;
2701 
2702  Float_t energy = clus->E();
2703  Int_t idclus = clus->GetID();
2704 
2705  Double_t tof = clus->GetTOF();
2706  if(GetCaloUtils()->GetEMCALRecoUtils()->IsTimeRecalibrationOn() && fTriggerClusterTimeRecal){
2707  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absIdMax,fInputEvent->GetBunchCrossNumber(),tof);
2708  //additional L1 phase shift
2709  if(GetCaloUtils()->GetEMCALRecoUtils()->IsL1PhaseInTimeRecalibrationOn()){
2710  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTimeL1Phase(GetCaloUtils()->GetEMCALGeometry()->GetSuperModuleNumber(absIdMax), fInputEvent->GetBunchCrossNumber(), tof);
2711  }
2712  }
2713  tof *=1.e9;
2714 
2715  //printf("cluster %d, ID %d, E %2.2f, tof %2.2f, AbsId max %d, exotic %d, bad Cluster %d, bad Cell %d\n",
2716  // iclus,idclus, energy,tof,absIdMax, exotic, badCluster,badCell);
2717 
2718  // Find the highest energy cluster, avobe trigger threshold
2719  // in the event in case no match to trigger is found
2720  if( energy > eMax )
2721  {
2722  tofMax = tof;
2723  eMax = energy;
2724  badClMax = badCluster;
2725  badCeMax = badCell;
2726  exoMax = exotic;
2727  clusMax = iclus;
2728  idclusMax = idclus;
2729  absIdMaxMax = absIdMax;
2730  }
2731 
2732  // count the good clusters in the event avobe the trigger threshold
2733  // to check the exotic events
2734  if(!badCluster && !exotic)
2735  nOfHighECl++;
2736 
2737  // Find match to trigger
2738  if(fTriggerPatchClusterMatch && nPatch>0)
2739  {
2740  for(Int_t iabsId =0; iabsId < nPatch; iabsId++)
2741  {
2742  Int_t absIDCell[4];
2743  GetCaloUtils()->GetEMCALGeometry()->GetCellIndexFromFastORIndex(patches.At(iabsId), absIDCell);
2744  //if(tof > 75 ) printf("E %2.2f TOF %2.2f Trigger patch %d, cells : %d, %d, %d, %d\n",
2745  // clus->E(),tof,patches.At(iabsId), absIDCell[0],absIDCell[1],absIDCell[2],absIDCell[3]);
2746 
2747  for(Int_t ipatch = 0; ipatch < 4; ipatch++)
2748  {
2749  if(absIdMax == absIDCell[ipatch])
2750  {
2751  //printf("*** Patches : absId %d, E %2.1f, tof %f \n",absIdMax,clus->E(), tof);
2752  if(energy > ePatchMax)
2753  {
2754  tofPatchMax = tof;
2755  ePatchMax = energy;
2756  fIsBadCellEvent = badCluster;
2757  fIsBadMaxCellEvent = badCell;
2758  fIsExoticEvent = exotic;
2759  fTriggerClusterIndex = iclus;
2760  fTriggerClusterId = idclus;
2761  fIsTriggerMatch = kTRUE;
2762 // absIdMaxTrig = absIdMax;
2763  }
2764  }
2765  }// cell patch loop
2766  }// trigger patch loop
2767  } // Do trigger patch matching
2768 
2769  }// Cluster loop
2770 
2771  // If there was no match, assign as trigger
2772  // the highest energy cluster in the event
2773  if(!fIsTriggerMatch)
2774  {
2775  tofPatchMax = tofMax;
2776  ePatchMax = eMax;
2777  fIsBadCellEvent = badClMax;
2778  fIsBadMaxCellEvent = badCeMax;
2779  fIsExoticEvent = exoMax;
2780  fTriggerClusterIndex = clusMax;
2781  fTriggerClusterId = idclusMax;
2782  }
2783 
2784  Double_t tofPatchMaxUS = TMath::Abs(tofPatchMax);
2785 
2786  if (tofPatchMaxUS < 28 ) fTriggerClusterBC = 0 ;
2787  else if(tofPatchMaxUS < 75 ) fTriggerClusterBC = 1 ;
2788  else if(tofPatchMaxUS < 125) fTriggerClusterBC = 2 ;
2789  else if(tofPatchMaxUS < 175) fTriggerClusterBC = 3 ;
2790  else if(tofPatchMaxUS < 225) fTriggerClusterBC = 4 ;
2791  else if(tofPatchMaxUS < 275) fTriggerClusterBC = 5 ;
2792  else
2793  {
2794  //printf("AliCaloTrackReader::MatchTriggerCluster() - Large BC - tof %2.3f - Index %d\n",tofPatchMaxUS,fTriggerClusterIndex);
2796  else
2797  {
2798  fTriggerClusterIndex = -2;
2799  fTriggerClusterId = -2;
2800  }
2801  }
2802 
2803  if(tofPatchMax < 0) fTriggerClusterBC*=-1;
2804 
2805 
2806  // 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",
2807  // fTriggerClusterIndex, fTriggerClusterId,ePatchMax, tofPatchMax,
2808  // fTriggerClusterBC, fIsBadCellEvent,fIsBadMaxCellEvent,fIsExoticEvent, fIsTriggerMatch, nOfHighECl,absIdMaxMax);
2809  //
2810  // 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",
2811  // clusMax, idclusMax, eMax,tofMax, badClMax, badCeMax,exoMax);
2812 
2813  //Redo matching but open cuts
2814  if(!fIsTriggerMatch && fTriggerClusterId >= 0)
2815  {
2816  // Open time patch time
2817  TArrayI patchOpen = GetTriggerPatches(7,10);
2818 
2819  Int_t patchAbsIdOpenTime = -1;
2820  for(Int_t iabsId =0; iabsId < patchOpen.GetSize(); iabsId++)
2821  {
2822  Int_t absIDCell[4];
2823  patchAbsIdOpenTime = patchOpen.At(iabsId);
2824  GetCaloUtils()->GetEMCALGeometry()->GetCellIndexFromFastORIndex(patchAbsIdOpenTime, absIDCell);
2825  //if(tof > 75 ) printf("E %2.2f TOF %2.2f Trigger patch %d, cells : %d, %d, %d, %d\n",
2826  // clus->E(),tof,patches.At(iabsId), absIDCell[0],absIDCell[1],absIDCell[2],absIDCell[3]);
2827 
2828  for(Int_t ipatch = 0; ipatch < 4; ipatch++)
2829  {
2830  if(absIdMaxMax == absIDCell[ipatch])
2831  {
2832  fIsTriggerMatchOpenCut[0] = kTRUE;
2833  break;
2834  }
2835  }// cell patch loop
2836  }// trigger patch loop
2837 
2838  // Check neighbour patches
2839  Int_t patchAbsId = -1;
2840  Int_t globalCol = -1;
2841  Int_t globalRow = -1;
2842  GetCaloUtils()->GetEMCALGeometry()->GetFastORIndexFromCellIndex(absIdMaxMax, patchAbsId);
2843  GetCaloUtils()->GetEMCALGeometry()->GetPositionInEMCALFromAbsFastORIndex(patchAbsId,globalCol,globalRow);
2844 
2845  // Check patches with strict time cut
2846  Int_t patchAbsIdNeigh = -1;
2847  for(Int_t icol = globalCol-1; icol <= globalCol+1; icol++)
2848  {
2849  if(icol < 0 || icol > 47) continue;
2850 
2851  for(Int_t irow = globalRow; irow <= globalRow+1; irow++)
2852  {
2853  if(irow < 0 || irow > 63) continue;
2854 
2855  GetCaloUtils()->GetEMCALGeometry()->GetAbsFastORIndexFromPositionInEMCAL(icol, irow, patchAbsIdNeigh);
2856 
2857  if ( patchAbsIdNeigh < 0 ) continue;
2858 
2859  for(Int_t iabsId =0; iabsId < patches.GetSize(); iabsId++)
2860  {
2861  if(patchAbsIdNeigh == patches.At(iabsId))
2862  {
2863  fIsTriggerMatchOpenCut[1] = kTRUE;
2864  break;
2865  }
2866  }// trigger patch loop
2867 
2868  }// row
2869  }// col
2870 
2871  // Check patches with open time cut
2872  Int_t patchAbsIdNeighOpenTime = -1;
2873  for(Int_t icol = globalCol-1; icol <= globalCol+1; icol++)
2874  {
2875  if(icol < 0 || icol > 47) continue;
2876 
2877  for(Int_t irow = globalRow; irow <= globalRow+1; irow++)
2878  {
2879  if(irow < 0 || irow > 63) continue;
2880 
2881  GetCaloUtils()->GetEMCALGeometry()->GetAbsFastORIndexFromPositionInEMCAL(icol, irow, patchAbsIdNeighOpenTime);
2882 
2883  if ( patchAbsIdNeighOpenTime < 0 ) continue;
2884 
2885  for(Int_t iabsId =0; iabsId < patchOpen.GetSize(); iabsId++)
2886  {
2887  if(patchAbsIdNeighOpenTime == patchOpen.At(iabsId))
2888  {
2889  fIsTriggerMatchOpenCut[2] = kTRUE;
2890  break;
2891  }
2892  }// trigger patch loop
2893 
2894  }// row
2895  }// col
2896 
2897  // printf("No match, new match: Open time %d-%d, open Neigh %d-%d, both open %d-%d\n",fIsTriggerMatchOpenCut[0],patchAbsIdOpenTime,
2898  // fIsTriggerMatchOpenCut[1],patchAbsIdNeigh,
2899  // fIsTriggerMatchOpenCut[2],patchAbsIdNeighOpenTime);
2900 
2901  patchOpen.Reset();
2902 
2903  }// No trigger match found
2904  //printf("Trigger BC %d, Id %d, Index %d\n",fTriggerClusterBC,fTriggerClusterId,fTriggerClusterIndex);
2905 }
2906 
2907 //_________________________________________________________
2911 //_________________________________________________________
2913 {
2915  {
2916  // get object pointer
2917  AliVCaloTrigger *caloTrigger = GetInputEvent()->GetCaloTrigger( "EMCAL" );
2918 
2919  if ( fBitEGA == 6 )
2920  {
2921  if (IsEventEMCALL1Gamma1()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(1);
2922  else if(IsEventEMCALL1Gamma2()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(3);
2923  else if(IsEventEMCALL1Jet1 ()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(0);
2924  else if(IsEventEMCALL1Jet2 ()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(2);
2925 
2926  // printf("L1 trigger Threshold Jet1 %f, Gamma1 %f, Jet2 %f, Gamma2 %f\n",
2927  // 0.07874*caloTrigger->GetL1Threshold(0),
2928  // 0.07874*caloTrigger->GetL1Threshold(1),
2929  // 0.07874*caloTrigger->GetL1Threshold(2),
2930  // 0.07874*caloTrigger->GetL1Threshold(3));
2931  }
2932  else
2933  {
2934  // Old AOD data format, in such case, order of thresholds different!!!
2935  if (IsEventEMCALL1Gamma1()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(0);
2936  else if(IsEventEMCALL1Gamma2()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(2);
2937  else if(IsEventEMCALL1Jet1 ()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(1);
2938  else if(IsEventEMCALL1Jet2 ()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(3);
2939 
2940  // printf("L1 trigger Threshold Jet1 %f, Gamma1 %f, Jet2 %f, Gamma2 %f\n",
2941  // 0.07874*caloTrigger->GetL1Threshold(1),
2942  // 0.07874*caloTrigger->GetL1Threshold(0),
2943  // 0.07874*caloTrigger->GetL1Threshold(3),
2944  // 0.07874*caloTrigger->GetL1Threshold(2));
2945  }
2946  }
2947 
2948  // Set L0 threshold, if not set by user
2950  {
2951  // Revise for periods > LHC11d
2952  Int_t runNumber = fInputEvent->GetRunNumber();
2953  if (runNumber < 146861) fTriggerL0EventThreshold = 3. ;
2954  else if(runNumber < 154000) fTriggerL0EventThreshold = 4. ;
2955  else if(runNumber < 165000) fTriggerL0EventThreshold = 5.5;
2956  else fTriggerL0EventThreshold = 2 ;
2957  }
2958 }
2959 
2960 //________________________________________________________
2962 //________________________________________________________
2963 void AliCaloTrackReader::Print(const Option_t * opt) const
2964 {
2965  if(! opt)
2966  return;
2967 
2968  printf("***** Print: %s %s ******\n", GetName(), GetTitle() ) ;
2969  printf("Task name : %s\n", fTaskName.Data()) ;
2970  printf("Data type : %d\n", fDataType) ;
2971  printf("CTS Min pT : %2.1f GeV/c\n", fCTSPtMin) ;
2972  printf("EMCAL Min pT : %2.1f GeV/c\n", fEMCALPtMin) ;
2973  printf("PHOS Min pT : %2.1f GeV/c\n", fPHOSPtMin) ;
2974  printf("CTS Max pT : %2.1f GeV/c\n", fCTSPtMax) ;
2975  printf("EMCAL Max pT : %2.1f GeV/c\n", fEMCALPtMax) ;
2976  printf("PHOS Max pT : %2.1f GeV/c\n", fPHOSPtMax) ;
2977  printf("EMCAL Bad Dist > %2.1f \n" , fEMCALBadChMinDist) ;
2978  printf("PHOS Bad Dist > %2.1f \n" , fPHOSBadChMinDist) ;
2979  printf("EMCAL N cells > %d \n" , fEMCALNCellsCut) ;
2980  printf("PHOS N cells > %d \n" , fPHOSNCellsCut) ;
2981  printf("EMCAL Time Cut: %3.1f < TOF < %3.1f\n", fEMCALTimeCutMin, fEMCALTimeCutMax);
2982  printf("Use CTS = %d\n", fFillCTS) ;
2983  printf("Use EMCAL = %d\n", fFillEMCAL) ;
2984  printf("Use DCAL = %d\n", fFillDCAL) ;
2985  printf("Use PHOS = %d\n", fFillPHOS) ;
2986  printf("Use EMCAL Cells = %d\n", fFillEMCALCells) ;
2987  printf("Use PHOS Cells = %d\n", fFillPHOSCells) ;
2988  printf("Track status = %d\n", (Int_t) fTrackStatus) ;
2989 
2990  printf("Track Mult Eta Cut = %2.2f\n", fTrackMultEtaCut) ;
2991  printf("Write delta AOD = %d\n", fWriteOutputDeltaAOD) ;
2992  printf("Recalculate Clusters = %d, E linearity = %d\n", fRecalculateClusters, fCorrectELinearity) ;
2993 
2994  printf("Use Triggers selected in SE base class %d; If not what Trigger Mask? %d; MB Trigger Mask for mixed %d \n",
2996 
2998  printf("Compare jet pt and pt hard to accept event, factor = %2.2f",fPtHardAndJetPtFactor);
2999 
3001  printf("Compare cluster pt and pt hard to accept event, factor = %2.2f",fPtHardAndClusterPtFactor);
3002 
3003  printf("Read Kine from, stack? %d, AOD ? %d \n", fReadStack, fReadAODMCParticles) ;
3004  printf("Delta AOD File Name = %s\n", fDeltaAODFileName.Data()) ;
3005  printf("Centrality: Class %s, Option %d, Bin [%d,%d] \n", fCentralityClass.Data(),fCentralityOpt,fCentralityBin[0], fCentralityBin[1]) ;
3006 
3007  printf(" \n") ;
3008 }
3009 
3010 //__________________________________________
3014 //__________________________________________
3016 {
3017  // Count number of cells with energy larger than 0.1 in SM3, cut on this number
3018  Int_t ncellsSM3 = 0;
3019  for(Int_t icell = 0; icell < fInputEvent->GetEMCALCells()->GetNumberOfCells(); icell++)
3020  {
3021  Int_t absID = fInputEvent->GetEMCALCells()->GetCellNumber(icell);
3022  Int_t sm = GetCaloUtils()->GetEMCALGeometry()->GetSuperModuleNumber(absID);
3023  if(fInputEvent->GetEMCALCells()->GetAmplitude(icell) > 0.1 && sm==3) ncellsSM3++;
3024  }
3025 
3026  Int_t ncellcut = 21;
3027  if(GetFiredTriggerClasses().Contains("EMC")) ncellcut = 35;
3028 
3029  if(ncellsSM3 >= ncellcut)
3030  {
3031  AliDebug(1,Form("Reject event with ncells in SM3 %d, cut %d, trig %s",
3032  ncellsSM3,ncellcut,GetFiredTriggerClasses().Data()));
3033  return kTRUE;
3034  }
3035 
3036  return kFALSE;
3037 }
3038 
3039 //_________________________________________________________
3043 //_________________________________________________________
3045 {
3046  if(label < 0) return ;
3047 
3048  AliAODEvent * evt = dynamic_cast<AliAODEvent*> (fInputEvent) ;
3049  if(!evt) return ;
3050 
3051  TClonesArray * arr = dynamic_cast<TClonesArray*>(evt->FindListObject("mcparticles")) ;
3052  if(!arr) return ;
3053 
3054  if(label < arr->GetEntriesFast())
3055  {
3056  AliAODMCParticle * particle = dynamic_cast<AliAODMCParticle *>(arr->At(label));
3057  if(!particle) return ;
3058 
3059  if(label == particle->Label()) return ; // label already OK
3060  //else printf("AliCaloTrackReader::RemapMCLabelForAODs() - Label %d - AOD stack %d \n",label, particle->Label());
3061  }
3062  //else printf("AliCaloTrackReader::RemapMCLabelForAODs() - Label %d > AOD labels %d \n",label, arr->GetEntriesFast());
3063 
3064  // loop on the particles list and check if there is one with the same label
3065  for(Int_t ind = 0; ind < arr->GetEntriesFast(); ind++ )
3066  {
3067  AliAODMCParticle * particle = dynamic_cast<AliAODMCParticle *>(arr->At(ind));
3068  if(!particle) continue ;
3069 
3070  if(label == particle->Label())
3071  {
3072  label = ind;
3073  //printf("AliAnalysisTaskEMCALClusterize::RemapMCLabelForAODs() - New Label Index %d \n",label);
3074  return;
3075  }
3076  }
3077 
3078  label = -1;
3079 
3080  //printf("AliCaloTrackReader::RemapMCLabelForAODs() - Label not found set to -1 \n");
3081 }
3082 
3083 //___________________________________
3085 //___________________________________
3087 {
3088  if(fCTSTracks) fCTSTracks -> Clear();
3089  if(fEMCALClusters) fEMCALClusters -> Clear("C");
3090  if(fPHOSClusters) fPHOSClusters -> Clear("C");
3091 
3092  fV0ADC[0] = 0; fV0ADC[1] = 0;
3093  fV0Mul[0] = 0; fV0Mul[1] = 0;
3094 
3095  if(fNonStandardJets) fNonStandardJets -> Clear("C");
3096  fBackgroundJets->Reset();
3097 }
3098 
3099 //___________________________________________
3103 //___________________________________________
3105 {
3106  fEventTrigMinBias = kFALSE;
3107  fEventTrigCentral = kFALSE;
3108  fEventTrigSemiCentral = kFALSE;
3109  fEventTrigEMCALL0 = kFALSE;
3110  fEventTrigEMCALL1Gamma1 = kFALSE;
3111  fEventTrigEMCALL1Gamma2 = kFALSE;
3112  fEventTrigEMCALL1Jet1 = kFALSE;
3113  fEventTrigEMCALL1Jet2 = kFALSE;
3114 
3115  AliDebug(1,Form("Select trigger mask bit %d - Trigger Event %s - Select <%s>",
3117 
3118  if(fEventTriggerMask <=0 )// in case no mask set
3119  {
3120  // EMC triggered event? Which type?
3121  if( GetFiredTriggerClasses().Contains("-B-") || GetFiredTriggerClasses().Contains("-S-") || GetFiredTriggerClasses().Contains("-I-") )
3122  {
3123  if ( GetFiredTriggerClasses().Contains("EGA" ) ||
3124  GetFiredTriggerClasses().Contains("EG1" ) )
3125  {
3126  fEventTrigEMCALL1Gamma1 = kTRUE;
3127  if( GetFiredTriggerClasses().Contains("EG1" ) && !fFiredTriggerClassName.Contains("EG1") ) fEventTrigEMCALL1Gamma1 = kFALSE;
3128  }
3129  else if( GetFiredTriggerClasses().Contains("EG2" ) )
3130  {
3131  fEventTrigEMCALL1Gamma2 = kTRUE;
3132  if( !fFiredTriggerClassName.Contains("EG2") ) fEventTrigEMCALL1Gamma2 = kFALSE;
3133  }
3134  else if( GetFiredTriggerClasses().Contains("EJE" ) ||
3135  GetFiredTriggerClasses().Contains("EJ1" ) )
3136  {
3137  fEventTrigEMCALL1Jet1 = kTRUE;
3138  if( GetFiredTriggerClasses().Contains("EJ1" ) && !fFiredTriggerClassName.Contains("EJ1") )
3139  fEventTrigEMCALL1Jet1 = kFALSE;
3140  }
3141  else if( GetFiredTriggerClasses().Contains("EJ2" ) )
3142  {
3143  fEventTrigEMCALL1Jet2 = kTRUE;
3144  if( !fFiredTriggerClassName.Contains("EJ2") ) fEventTrigEMCALL1Jet2 = kFALSE;
3145  }
3146  else if( GetFiredTriggerClasses().Contains("CEMC") &&
3147  !GetFiredTriggerClasses().Contains("EGA" ) &&
3148  !GetFiredTriggerClasses().Contains("EJE" ) &&
3149  !GetFiredTriggerClasses().Contains("EG1" ) &&
3150  !GetFiredTriggerClasses().Contains("EJ1" ) &&
3151  !GetFiredTriggerClasses().Contains("EG2" ) &&
3152  !GetFiredTriggerClasses().Contains("EJ2" ) )
3153  {
3154  fEventTrigEMCALL0 = kTRUE;
3155  }
3156 
3157  //Min bias event trigger?
3158  if (GetFiredTriggerClasses().Contains("CCENT_R2-B-NOPF-ALLNOTRD"))
3159  {
3160  fEventTrigCentral = kTRUE;
3161  }
3162  else if(GetFiredTriggerClasses().Contains("CSEMI_R1-B-NOPF-ALLNOTRD"))
3163  {
3164  fEventTrigSemiCentral = kTRUE;
3165  }
3166  else if((GetFiredTriggerClasses().Contains("CINT") || GetFiredTriggerClasses().Contains("CPBI2_B1") ) &&
3167  GetFiredTriggerClasses().Contains("-NOPF-ALLNOTRD") )
3168  {
3169  fEventTrigMinBias = kTRUE;
3170  }
3171  }
3172  }
3173  else
3174  {
3175  // EMC L1 Gamma
3176  if ( fEventTriggerMask & AliVEvent::kEMCEGA )
3177  {
3178  //printf("EGA trigger bit\n");
3179  if (GetFiredTriggerClasses().Contains("EG"))
3180  {
3181  if (GetFiredTriggerClasses().Contains("EGA")) fEventTrigEMCALL1Gamma1 = kTRUE;
3182  else
3183  {
3184  if(GetFiredTriggerClasses().Contains("EG1")) fEventTrigEMCALL1Gamma1 = kTRUE;
3185  if(GetFiredTriggerClasses().Contains("EG2")) fEventTrigEMCALL1Gamma2 = kTRUE;
3186  }
3187  }
3188  }
3189  // EMC L1 Jet
3190  else if( fEventTriggerMask & AliVEvent::kEMCEJE )
3191  {
3192  //printf("EGA trigger bit\n");
3193  if (GetFiredTriggerClasses().Contains("EJ"))
3194  {
3195  if (GetFiredTriggerClasses().Contains("EJE")) fEventTrigEMCALL1Jet1 = kTRUE;
3196  else
3197  {
3198  if(GetFiredTriggerClasses().Contains("EJ1")) fEventTrigEMCALL1Jet1 = kTRUE;
3199  if(GetFiredTriggerClasses().Contains("EJ2")) fEventTrigEMCALL1Jet2 = kTRUE;
3200  }
3201  }
3202  }
3203  // EMC L0
3204  else if((fEventTriggerMask & AliVEvent::kEMC7) ||
3205  (fEventTriggerMask & AliVEvent::kEMC1) )
3206  {
3207  //printf("L0 trigger bit\n");
3208  fEventTrigEMCALL0 = kTRUE;
3209  }
3210  // Min Bias Pb-Pb
3211  else if( fEventTriggerMask & AliVEvent::kCentral )
3212  {
3213  //printf("MB semi central trigger bit\n");
3214  fEventTrigSemiCentral = kTRUE;
3215  }
3216  // Min Bias Pb-Pb
3217  else if( fEventTriggerMask & AliVEvent::kSemiCentral )
3218  {
3219  //printf("MB central trigger bit\n");
3220  fEventTrigCentral = kTRUE;
3221  }
3222  // Min Bias pp, PbPb, pPb
3223  else if((fEventTriggerMask & AliVEvent::kMB ) ||
3224  (fEventTriggerMask & AliVEvent::kINT7) ||
3225  (fEventTriggerMask & AliVEvent::kINT8) ||
3226  (fEventTriggerMask & AliVEvent::kAnyINT) )
3227  {
3228  //printf("MB trigger bit\n");
3229  fEventTrigMinBias = kTRUE;
3230  }
3231  }
3232 
3233  AliDebug(1,Form("Event bits: \n \t MB %d, Cen %d, Sem %d, L0 %d, L1G1 %d, L1G2 %d, L1J1 %d, L1J2 %d",
3237 
3238  // L1 trigger bit
3239  if( fBitEGA == 0 && fBitEJE == 0 )
3240  {
3241  // Init the trigger bit once, correct depending on AliESD(AOD)CaloTrigger header version
3242 
3243  // Simpler way to do it ...
3244 // if( fInputEvent->GetRunNumber() < 172000 )
3245 // reader->SetEventTriggerL1Bit(4,5); // current LHC11 data
3246 // else
3247 // reader->SetEventTriggerL1Bit(6,8); // LHC12-13 data
3248 
3249  // Old values
3250  fBitEGA = 4;
3251  fBitEJE = 5;
3252 
3253  TFile* file = AliAnalysisManager::GetAnalysisManager()->GetTree()->GetCurrentFile();
3254 
3255  const TList *clist = file->GetStreamerInfoCache();
3256 
3257  if(clist)
3258  {
3259  TStreamerInfo *cinfo = (TStreamerInfo*)clist->FindObject("AliESDCaloTrigger");
3260  Int_t verid = 5; // newer ESD header version
3261  if(!cinfo)
3262  {
3263  cinfo = (TStreamerInfo*)clist->FindObject("AliAODCaloTrigger");
3264  verid = 3; // newer AOD header version
3265  }
3266 
3267  if(cinfo)
3268  {
3269  Int_t classversionid = cinfo->GetClassVersion();
3270  //printf("********* Header class version %d *********** \n",classversionid);
3271 
3272  if (classversionid >= verid)
3273  {
3274  fBitEGA = 6;
3275  fBitEJE = 8;
3276  }
3277  } else AliInfo("AliCaloTrackReader()::SetEventTriggerBit() - Streamer info for trigger class not available, bit not changed");
3278  } else AliInfo("AliCaloTrackReader::SetEventTriggerBit() - Streamer list not available!, bit not changed");
3279 
3280  } // set once the EJE, EGA trigger bit
3281 }
3282 
3283 //____________________________________________________________
3286 //____________________________________________________________
3287 void AliCaloTrackReader::SetInputEvent(AliVEvent* const input)
3288 {
3289  fInputEvent = input;
3290  fMixedEvent = dynamic_cast<AliMixedEvent*>(GetInputEvent()) ;
3291  if (fMixedEvent)
3292  fNMixedEvent = fMixedEvent->GetNumberOfEvents() ;
3293 
3294  //Delete previous vertex
3295  if(fVertex)
3296  {
3297  for (Int_t i = 0; i < fNMixedEvent; i++)
3298  {
3299  delete [] fVertex[i] ;
3300  }
3301  delete [] fVertex ;
3302  }
3303 
3304  fVertex = new Double_t*[fNMixedEvent] ;
3305  for (Int_t i = 0; i < fNMixedEvent; i++)
3306  {
3307  fVertex[i] = new Double_t[3] ;
3308  fVertex[i][0] = 0.0 ;
3309  fVertex[i][1] = 0.0 ;
3310  fVertex[i][2] = 0.0 ;
3311  }
3312 }
3313 
3314 
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
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