AliPhysics  32b88a8 (32b88a8)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliCaloTrackReader.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Author: The ALICE Off-line Project. *
5  * Contributors are mentioned in the code where appropriate. *
6  * *
7  * Permission to use, copy, modify and distribute this software and its *
8  * documentation strictly for non-commercial purposes is hereby granted *
9  * without fee, provided that the above copyright notice appears in all *
10  * copies and that both the copyright notice and this permission notice *
11  * appear in the supporting documentation. The authors make no claims *
12  * about the suitability of this software for any purpose. It is *
13  * provided "as is" without express or implied warranty. *
14  **************************************************************************/
15 
16 // --- ROOT system ---
17 #include <TFile.h>
18 #include <TGeoManager.h>
19 #include <TStreamerInfo.h>
20 
21 // ---- ANALYSIS system ----
22 #include "AliMCEvent.h"
23 #include "AliAODMCHeader.h"
24 #include "AliGenPythiaEventHeader.h"
25 #include "AliGenCocktailEventHeader.h"
26 #include "AliGenHijingEventHeader.h"
27 #include "AliESDEvent.h"
28 #include "AliAODEvent.h"
29 #include "AliVTrack.h"
30 #include "AliVParticle.h"
31 #include "AliMixedEvent.h"
32 //#include "AliTriggerAnalysis.h"
33 #include "AliESDVZERO.h"
34 #include "AliVCaloCells.h"
35 #include "AliAnalysisManager.h"
36 #include "AliInputEventHandler.h"
37 #include "AliAODMCParticle.h"
38 #include "AliStack.h"
39 #include "AliLog.h"
40 #include "AliMultSelection.h"
41 
42 // ---- Detectors ----
43 #include "AliPHOSGeoUtils.h"
44 #include "AliEMCALGeometry.h"
45 #include "AliEMCALRecoUtils.h"
46 
47 // ---- CaloTrackCorr ---
48 #include "AliCalorimeterUtils.h"
49 #include "AliCaloTrackReader.h"
50 
51 // ---- Jets ----
52 #include "AliAODJet.h"
53 #include "AliAODJetEventBackground.h"
54 
58 
59 //________________________________________
61 //________________________________________
63 TObject(), fEventNumber(-1), //fCurrentFileName(""),
64 fDataType(0), fDebug(0),
65 fFiducialCut(0x0), fCheckFidCut(kFALSE),
66 fComparePtHardAndJetPt(0), fPtHardAndJetPtFactor(0),
67 fComparePtHardAndClusterPt(0),fPtHardAndClusterPtFactor(0),
68 fCTSPtMin(0), fEMCALPtMin(0), fPHOSPtMin(0),
69 fCTSPtMax(0), fEMCALPtMax(0), fPHOSPtMax(0),
70 fEMCALBadChMinDist(0), fPHOSBadChMinDist (0),
71 fEMCALNCellsCut(0), fPHOSNCellsCut(0),
72 fUseEMCALTimeCut(1), fUseParamTimeCut(0),
73 fUseTrackTimeCut(0), fAccessTrackTOF(0),
74 fEMCALTimeCutMin(-10000), fEMCALTimeCutMax(10000),
75 fEMCALParamTimeCutMin(), fEMCALParamTimeCutMax(),
76 fTrackTimeCutMin(-10000), fTrackTimeCutMax(10000),
77 fUseTrackDCACut(0),
78 fAODBranchList(0x0),
79 fCTSTracks(0x0), fEMCALClusters(0x0),
80 fDCALClusters(0x0), fPHOSClusters(0x0),
81 fEMCALCells(0x0), fPHOSCells(0x0),
82 fInputEvent(0x0), fOutputEvent(0x0),fMC(0x0),
83 fFillCTS(0), fFillEMCAL(0),
84 fFillDCAL(0), fFillPHOS(0),
85 fFillEMCALCells(0), fFillPHOSCells(0),
86 fRecalculateClusters(kFALSE),fCorrectELinearity(kTRUE),
87 fSelectEmbeddedClusters(kFALSE),
88 fSmearShowerShape(0), fSmearShowerShapeWidth(0), fRandom(),
89 fSmearingFunction(0), fSmearNLMMin(0), fSmearNLMMax(0),
90 fTrackStatus(0), fSelectSPDHitTracks(0),
91 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 //____________________________________________________________
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 
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  {
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, %d<=NLM<=%d; ",
695  parList+=onePar ;
696  }
697 
699  {
700  snprintf(onePar,buffersize,"jet pt / pt hard < %2.1f; ",fPtHardAndJetPtFactor);
701  parList+=onePar ;
702  }
703 
705  {
706  snprintf(onePar,buffersize,"cluster pt / pt hard < %2.2f",fPtHardAndClusterPtFactor);
707  parList+=onePar ;
708  }
709 
710  snprintf(onePar,buffersize,"Centrality: Class %s, Option %d, Bin [%d,%d]; New centrality %d; Event plane method %s; ",
712  parList+=onePar ;
713 
714  return new TObjString(parList) ;
715 }
716 
717 
718 //____________________________________________
720 //____________________________________________
722 {
723  if(fMC)
724  return fMC->Stack();
725  else
726  {
727  AliDebug(1,"Stack is not available");
728  return 0x0 ;
729  }
730 }
731 
732 //______________________________________________
734 //______________________________________________
736 {
737  if(fMC)
738  {
739  return fMC->Header();
740  }
741  else
742  {
743  AliInfo("Header is not available");
744  return 0x0 ;
745  }
746 }
747 
748 //____________________________________________________
752 //____________________________________________________
754 {
755  fNMCProducedMin = 0;
756  fNMCProducedMax = 0;
757 
758  if (ReadStack() && fMC)
759  {
760  AliGenEventHeader * eventHeader = fMC->GenEventHeader();
761 
762  if(!fAcceptOnlyHIJINGLabels) return ;
763 
764  // TODO Check if it works from here ...
765 
766  AliGenCocktailEventHeader *cocktail = dynamic_cast<AliGenCocktailEventHeader *>(eventHeader);
767 
768  if(!cocktail) return ;
769 
770  TList *genHeaders = cocktail->GetHeaders();
771 
772  Int_t nGenerators = genHeaders->GetEntries();
773  //printf("N generators %d \n", nGenerators);
774 
775  for(Int_t igen = 0; igen < nGenerators; igen++)
776  {
777  AliGenEventHeader * eventHeader2 = (AliGenEventHeader*)genHeaders->At(igen) ;
778  TString name = eventHeader2->GetName();
779 
780  //printf("Generator %d: Class Name %s, Name %s, title %s \n",igen, eventHeader2->ClassName(), name.Data(), eventHeader2->GetTitle());
781 
783  fNMCProducedMax+= eventHeader2->NProduced();
784 
785  if(name.Contains("Hijing",TString::kIgnoreCase)) return ;
786  }
787 
788  }
789  else if(ReadAODMCParticles() && GetAODMCHeader())
790  {
791  Int_t nGenerators = GetAODMCHeader()->GetNCocktailHeaders();
792  //printf("AliCaloTrackReader::GetGenEventHeader() - N headers %d\n",nGenerators);
793 
794  if( nGenerators <= 0) return ;
795 
796  if(!fAcceptOnlyHIJINGLabels) return ;
797 
798  for(Int_t igen = 0; igen < nGenerators; igen++)
799  {
800  AliGenEventHeader * eventHeader = GetAODMCHeader()->GetCocktailHeader(igen) ;
801  TString name = eventHeader->GetName();
802 
803  //printf("Generator %d: Class Name %s, Name %s, title %s \n",igen, eventHeader->ClassName(), name.Data(), eventHeader->GetTitle());
804 
806  fNMCProducedMax+= eventHeader->NProduced();
807 
808  if(name.Contains("Hijing",TString::kIgnoreCase)) return ;
809  }
810 
811  }
812 }
813 
814 //______________________________________________________________
817 //______________________________________________________________
818 AliGenEventHeader* AliCaloTrackReader::GetGenEventHeader() const
819 {
820  if (ReadStack() && fMC)
821  {
822  AliGenEventHeader * eventHeader = fMC->GenEventHeader();
823 
824  if(!fAcceptOnlyHIJINGLabels) return eventHeader ;
825 
826  // TODO Check if it works from here ...
827 
828  AliGenCocktailEventHeader *cocktail = dynamic_cast<AliGenCocktailEventHeader *>(eventHeader);
829 
830  if(!cocktail) return 0x0 ;
831 
832  TList *genHeaders = cocktail->GetHeaders();
833 
834  Int_t nGenerators = genHeaders->GetEntries();
835  //printf("N generators %d \n", nGenerators);
836 
837  for(Int_t igen = 0; igen < nGenerators; igen++)
838  {
839  AliGenEventHeader * eventHeader2 = (AliGenEventHeader*)genHeaders->At(igen) ;
840  TString name = eventHeader2->GetName();
841 
842  //printf("Generator %d: Class Name %s, Name %s, title %s \n",igen, eventHeader2->ClassName(), name.Data(), eventHeader2->GetTitle());
843 
844  if(name.Contains("Hijing",TString::kIgnoreCase)) return eventHeader2 ;
845  }
846 
847  return 0x0;
848 
849  }
850  else if(ReadAODMCParticles() && GetAODMCHeader())
851  {
852  Int_t nGenerators = GetAODMCHeader()->GetNCocktailHeaders();
853  //printf("AliCaloTrackReader::GetGenEventHeader() - N headers %d\n",nGenerators);
854 
855  if( nGenerators <= 0) return 0x0;
856 
857  if(!fAcceptOnlyHIJINGLabels) return GetAODMCHeader()->GetCocktailHeader(0);
858 
859  for(Int_t igen = 0; igen < nGenerators; igen++)
860  {
861  AliGenEventHeader * eventHeader = GetAODMCHeader()->GetCocktailHeader(igen) ;
862  TString name = eventHeader->GetName();
863 
864  //printf("Generator %d: Class Name %s, Name %s, title %s \n",igen, eventHeader->ClassName(), name.Data(), eventHeader->GetTitle());
865 
866  if(name.Contains("Hijing",TString::kIgnoreCase)) return eventHeader ;
867  }
868 
869  return 0x0;
870 
871  }
872  else
873  {
874  //printf("AliCaloTrackReader::GetGenEventHeader() - MC header not available! \n");
875  return 0x0;
876  }
877 }
878 
879 //_________________________________________________________
882 //_________________________________________________________
884 {
885  AliInfo("Input are not AODs");
886 
887  return NULL ;
888 }
889 
890 //________________________________________________________
893 //________________________________________________________
894 AliAODMCHeader* AliCaloTrackReader::GetAODMCHeader() const
895 {
896  AliInfo("Input are not AODs");
897 
898  return NULL ;
899 }
900 
901 //___________________________________________________________
909 //___________________________________________________________
910 Int_t AliCaloTrackReader::GetVertexBC(const AliVVertex * vtx)
911 {
912  Int_t vertexBC=vtx->GetBC();
913 
914  if(!fRecalculateVertexBC) return vertexBC;
915 
916  // Value not available, recalculate it.
917 
918  Double_t bz = fInputEvent->GetMagneticField();
919  Bool_t bc0 = kFALSE;
920  Int_t ntr = GetCTSTracks()->GetEntriesFast();
921  //printf("N Tracks %d\n",ntr);
922 
923  for(Int_t i = 0 ; i < ntr ; i++)
924  {
925  AliVTrack * track = (AliVTrack*) (GetCTSTracks()->At(i));
926 
927  //Check if has TOF info, if not skip
928  ULong_t status = track->GetStatus();
929  Bool_t okTOF = (status & AliVTrack::kTOFout) == AliVTrack::kTOFout ;
930  vertexBC = track->GetTOFBunchCrossing(bz);
931  Float_t pt = track->Pt();
932 
933  if(!okTOF) continue;
934 
935  // Get DCA x, y
936  Double_t dca[2] = {1e6,1e6};
937  Double_t covar[3] = {1e6,1e6,1e6};
938  track->PropagateToDCA(vtx,bz,100.,dca,covar);
939 
940  if(AcceptDCA(pt,dca[0]))
941  {
942  if (vertexBC !=0 && fVertexBC != AliVTrack::kTOFBCNA) return vertexBC;
943  else if(vertexBC == 0) bc0 = kTRUE;
944  }
945  }
946 
947  if( bc0 ) vertexBC = 0 ;
948  else vertexBC = AliVTrack::kTOFBCNA ;
949 
950  return vertexBC;
951 }
952 
953 //_____________________________
960 //_____________________________
962 {
963  return track->GetID();
964 }
965 
966 //_____________________________
969 //_____________________________
971 {
973  {
974  AliInfo("Cannot access stack and mcparticles at the same time, change them");
975  fReadStack = kFALSE;
976  fReadAODMCParticles = kFALSE;
977  }
978 
979  // Activate debug level in AliAnaWeights
980  if( fWeightUtils->GetDebug() >= 0 )
981  (AliAnalysisManager::GetAnalysisManager())->AddClassDebug(fWeightUtils->ClassName(), fWeightUtils->GetDebug());
982 }
983 
984 //_______________________________________
986 //_______________________________________
988 {
989  fDataType = kESD ;
990  fCTSPtMin = 0.1 ;
991  fEMCALPtMin = 0.1 ;
992  fPHOSPtMin = 0.1 ;
993  fCTSPtMax = 1000. ;
994  fEMCALPtMax = 1000. ;
995  fPHOSPtMax = 1000. ;
996 
997  fEMCALBadChMinDist = 0; // open, 2; // standard
998  fPHOSBadChMinDist = 0; // open, 2; // standard
999 
1000  fEMCALNCellsCut = 0; // open, 1; // standard
1001  fPHOSNCellsCut = 0; // open, 2; // standard
1002 
1003  //Track DCA cuts
1004  // dca_xy cut = 0.0105+0.0350/TMath::Power(pt,1.1);
1005  fTrackDCACut[0] = 0.0105;
1006  fTrackDCACut[1] = 0.0350;
1007  fTrackDCACut[2] = 1.1;
1008 
1009  //Do not filter the detectors input by default.
1010  fFillEMCAL = kFALSE;
1011  fFillDCAL = kFALSE;
1012  fFillPHOS = kFALSE;
1013  fFillCTS = kFALSE;
1014  fFillEMCALCells = kFALSE;
1015  fFillPHOSCells = kFALSE;
1016 
1017  fReadStack = kFALSE; // Check in the constructor of the other readers if it was set or in the configuration file
1018  fReadAODMCParticles = kFALSE; // Check in the constructor of the other readers if it was set or in the configuration file
1019  fDeltaAODFileName = "deltaAODPartCorr.root";
1021  fEventTriggerMask = AliVEvent::kAny;
1022  fMixEventTriggerMask = AliVEvent::kAnyINT;
1023  fEventTriggerAtSE = kTRUE; // Use only events that pass event selection at SE base class
1024 
1025  fAcceptFastCluster = kTRUE;
1026  fEventType = -1;
1027 
1028  //We want tracks fitted in the detectors:
1029  //fTrackStatus=AliESDtrack::kTPCrefit;
1030  //fTrackStatus|=AliESDtrack::kITSrefit;
1031  fTrackStatus = 0;
1032 
1033  fV0ADC[0] = 0; fV0ADC[1] = 0;
1034  fV0Mul[0] = 0; fV0Mul[1] = 0;
1035 
1036  fZvtxCut = 10.;
1037 
1038  fNMixedEvent = 1;
1039 
1040  fPtHardAndJetPtFactor = 7.;
1042 
1043  //Centrality
1044  fUseAliCentrality = kFALSE;
1045  fCentralityClass = "V0M";
1046  fCentralityOpt = 100;
1047  fCentralityBin[0] = fCentralityBin[1]=-1;
1048 
1049  fEventPlaneMethod = "V0";
1050 
1051  // Allocate memory (not sure this is the right place)
1052  fCTSTracks = new TObjArray();
1053  fEMCALClusters = new TObjArray();
1054  fDCALClusters = new TObjArray();
1055  fPHOSClusters = new TObjArray();
1056  //fTriggerAnalysis = new AliTriggerAnalysis;
1057  fAODBranchList = new TList ;
1058  fOutputContainer = new TList ;
1059 
1060  fEnergyHistogramNbins = 200;
1061  fEnergyHistogramLimit[0] = 0 ;
1062  fEnergyHistogramLimit[1] = 100;
1063 
1064  fPileUpParamSPD[0] = 3 ; fPileUpParamSPD[1] = 0.8 ;
1065  fPileUpParamSPD[2] = 3.0 ; fPileUpParamSPD[3] = 2.0 ; fPileUpParamSPD[4] = 5.0;
1066 
1067  // Parametrized time cut (LHC11d)
1070 
1071  // Parametrized time cut (LHC11c)
1072  //fEMCALParamTimeCutMin[0] =-5; fEMCALParamTimeCutMin[1] =-1 ; fEMCALParamTimeCutMin[2] = 1.87; fEMCALParamTimeCutMin[3] = 0.4;
1073  //fEMCALParamTimeCutMax[0] = 3.5; fEMCALParamTimeCutMax[1] = 50; fEMCALParamTimeCutMax[2] = 0.15; fEMCALParamTimeCutMax[3] = 1.6;
1074 
1075  fTimeStampRunMin = -1;
1076  fTimeStampRunMax = 1e12;
1079 
1080  for(Int_t i = 0; i < 19; i++)
1081  {
1082  fEMCalBCEvent [i] = 0;
1083  fEMCalBCEventCut[i] = 0;
1084  fTrackBCEvent [i] = 0;
1085  fTrackBCEventCut[i] = 0;
1086  }
1087 
1088  // Trigger match-rejection
1089  fTriggerPatchTimeWindow[0] = 8;
1090  fTriggerPatchTimeWindow[1] = 9;
1091 
1092  fTriggerClusterBC = -10000 ;
1095  fTriggerClusterIndex = -1;
1096  fTriggerClusterId = -1;
1097 
1098  //Jets
1101  if(!fNonStandardJets) fNonStandardJets = new TClonesArray("AliAODJet",100);
1104  if(!fBackgroundJets) fBackgroundJets = new AliAODJetEventBackground();
1105 
1106  fSmearShowerShapeWidth = 0.005;
1107  fSmearNLMMin = 1;
1108  fSmearNLMMax = 1;
1109 
1110  fWeightUtils = new AliAnaWeights() ;
1111  fEventWeight = 1 ;
1112 
1113 }
1114 
1115 //___________________________________________________
1120 //___________________________________________________
1122 {
1123  AliGenHijingEventHeader* hijingHeader = dynamic_cast<AliGenHijingEventHeader *> (GetGenEventHeader());
1124 
1125  //printf("header %p, label %d\n",hijingHeader,label);
1126 
1127  if(!hijingHeader || label < 0 ) return kFALSE;
1128 
1129 
1130  //printf("pass a), N produced %d\n",nproduced);
1131 
1132  if(label >= fNMCProducedMin && label < fNMCProducedMax)
1133  {
1134  //printf(" accept!, label is smaller than produced, N %d\n",nproduced);
1135 
1136  return kTRUE;
1137  }
1138 
1139  if(ReadStack())
1140  {
1141  if(!GetStack()) return kFALSE;
1142 
1143  Int_t nprimaries = GetStack()->GetNtrack();
1144 
1145  if(label > nprimaries) return kFALSE;
1146 
1147  TParticle * mom = GetStack()->Particle(label);
1148 
1149  Int_t iMom = label;
1150  Int_t iParent = mom->GetFirstMother();
1151  while(iParent!=-1)
1152  {
1153  if(iParent >= fNMCProducedMin && iParent < fNMCProducedMax)
1154  {
1155  //printf("\t accept, mother is %d \n",iParent)
1156  return kTRUE;
1157  }
1158 
1159  iMom = iParent;
1160  mom = GetStack()->Particle(iMom);
1161  iParent = mom->GetFirstMother();
1162  }
1163 
1164  return kFALSE ;
1165 
1166  } // ESD
1167  else
1168  {
1169  TClonesArray* mcparticles = GetAODMCParticles();
1170 
1171  if(!mcparticles) return kFALSE;
1172 
1173  Int_t nprimaries = mcparticles->GetEntriesFast();
1174 
1175  if(label > nprimaries) return kFALSE;
1176 
1177  //printf("pass b) N primaries %d \n",nprimaries);
1178 
1179  if(label >= fNMCProducedMin && label < fNMCProducedMax)
1180  {
1181  return kTRUE;
1182  }
1183 
1184  // Find grand parent, check if produced in the good range
1185  AliAODMCParticle * mom = (AliAODMCParticle *) mcparticles->At(label);
1186 
1187  Int_t iMom = label;
1188  Int_t iParent = mom->GetMother();
1189  while(iParent!=-1)
1190  {
1191  if(iParent >= fNMCProducedMin && iParent < fNMCProducedMax)
1192  {
1193  //printf("\t accept, mother is %d, with nProduced %d \n",iParent, nproduced);
1194  return kTRUE;
1195  }
1196 
1197  iMom = iParent;
1198  mom = (AliAODMCParticle *) mcparticles->At(iMom);
1199  iParent = mom->GetMother();
1200 
1201  }
1202 
1203  //printf("pass c), no match found \n");
1204 
1205  return kFALSE ;
1206 
1207  }//AOD
1208 }
1209 
1210 //__________________________________________________________________________
1213 //__________________________________________________________________________
1215 {
1216  // Parametrized cut depending on E
1217  if(fUseParamTimeCut)
1218  {
1221  //printf("tof %f, minCut %f, maxCut %f\n",tof,minCut,maxCut);
1222  if( tof < minCut || tof > maxCut ) return kFALSE ;
1223  }
1224 
1225  //In any case, the time should to be larger than the fixed window ...
1226  if( tof < fEMCALTimeCutMin || tof > fEMCALTimeCutMax ) return kFALSE ;
1227 
1228  return kTRUE ;
1229 }
1230 
1231 //________________________________________________
1234 //________________________________________________
1236 {
1237  return fInputEvent->IsPileupFromSPD((Int_t) fPileUpParamSPD[0] , fPileUpParamSPD[1] ,
1238  fPileUpParamSPD[2] , fPileUpParamSPD[3] , fPileUpParamSPD[4] );
1239  //printf("Param : %d, %2.2f, %2.2f, %2.2f, %2.2f\n",(Int_t) fPileUpParamSPD[0], fPileUpParamSPD[1], fPileUpParamSPD[2], fPileUpParamSPD[3], fPileUpParamSPD[4]);
1240 }
1241 
1242 //__________________________________________________
1244 //__________________________________________________
1246 {
1247  if(fNPileUpClusters > fNPileUpClustersCut) return kTRUE ;
1248  else return kFALSE;
1249 }
1250 
1251 //________________________________________________________
1253 //________________________________________________________
1255 {
1256  if( IsPileUpFromSPD() && IsPileUpFromEMCal()) return kTRUE ;
1257  else return kFALSE;
1258 }
1259 
1260 //_______________________________________________________
1262 //_______________________________________________________
1264 {
1265  if( IsPileUpFromSPD() || IsPileUpFromEMCal()) return kTRUE ;
1266  else return kFALSE;
1267 }
1268 
1269 //___________________________________________________________
1271 //___________________________________________________________
1273 {
1274  if( IsPileUpFromSPD() && !IsPileUpFromEMCal()) return kTRUE ;
1275  else return kFALSE;
1276 }
1277 
1278 //___________________________________________________________
1280 //___________________________________________________________
1282 {
1283  if( !IsPileUpFromSPD() && IsPileUpFromEMCal()) return kTRUE ;
1284  else return kFALSE;
1285 }
1286 
1287 //______________________________________________________________
1289 //______________________________________________________________
1291 {
1292  if( !IsPileUpFromSPD() && !IsPileUpFromEMCal()) return kTRUE ;
1293  else return kFALSE;
1294 }
1295 
1296 //___________________________________________________________________________________
1304 //___________________________________________________________________________________
1305 Bool_t AliCaloTrackReader::FillInputEvent(Int_t iEntry, const char * /*curFileName*/)
1306 {
1307  fEventNumber = iEntry;
1308  fTriggerClusterIndex = -1;
1309  fTriggerClusterId = -1;
1310  fIsTriggerMatch = kFALSE;
1311  fTriggerClusterBC = -10000;
1312  fIsExoticEvent = kFALSE;
1313  fIsBadCellEvent = kFALSE;
1314  fIsBadMaxCellEvent = kFALSE;
1315 
1316  fIsTriggerMatchOpenCut[0] = kFALSE ;
1317  fIsTriggerMatchOpenCut[1] = kFALSE ;
1318  fIsTriggerMatchOpenCut[2] = kFALSE ;
1319 
1320  //fCurrentFileName = TString(currentFileName);
1321  if(!fInputEvent)
1322  {
1323  AliInfo("Input event not available, skip event analysis");
1324  return kFALSE;
1325  }
1326 
1327  fhNEventsAfterCut->Fill(0.5);
1328 
1329  //-----------------------------------------------
1330  // Select the event depending on the trigger type
1331  // and other event characteristics
1332  // like the goodness of the EMCal trigger
1333  //-----------------------------------------------
1334 
1335  Bool_t accept = CheckEventTriggers();
1336  if(!accept) return kFALSE;
1337 
1338  AliDebug(1,"Pass Event trigger selection");
1339 
1340  //---------------------------------------------------------------------------
1341  // In case of analysis of events with jets, skip those with jet pt > 5 pt hard
1342  // To be used on for MC data in pT hard bins
1343  //---------------------------------------------------------------------------
1344 
1346  {
1347  if(!ComparePtHardAndJetPt()) return kFALSE ;
1348  AliDebug(1,"Pass Pt Hard - Jet rejection");
1349  fhNEventsAfterCut->Fill(8.5);
1350  }
1351 
1353  {
1354  if(!ComparePtHardAndClusterPt()) return kFALSE ;
1355  AliDebug(1,"Pass Pt Hard - Cluster rejection");
1356  fhNEventsAfterCut->Fill(9.5);
1357  }
1358 
1359  //------------------------------------------------------
1360  // Event rejection depending on time stamp
1361  //------------------------------------------------------
1362 
1364  {
1365  AliESDEvent* esd = dynamic_cast<AliESDEvent*> (fInputEvent);
1366  if(esd)
1367  {
1368  Int_t timeStamp = esd->GetTimeStamp();
1369  Float_t timeStampFrac = 1.*(timeStamp-fTimeStampRunMin) / (fTimeStampRunMax-fTimeStampRunMin);
1370 
1371  //printf("stamp0 %d, max0 %d, frac %f\n", timeStamp-fTimeStampRunMin,fTimeStampRunMax-fTimeStampRunMin, timeStampFrac);
1372 
1373  if(timeStampFrac < fTimeStampEventFracMin || timeStampFrac > fTimeStampEventFracMax) return kFALSE;
1374  }
1375 
1376  AliDebug(1,"Pass Time Stamp rejection");
1377  fhNEventsAfterCut->Fill(10.5);
1378  }
1379 
1380  //------------------------------------------------------
1381  // Event rejection depending on vertex, pileup, v0and
1382  //------------------------------------------------------
1383 
1384  FillVertexArray();
1385 
1386  //Reject events with Z vertex too large, only for SE analysis, if not, cut on the analysis code
1387  if(!GetMixedEvent() && TMath::Abs(fVertex[0][2]) > fZvtxCut) return kFALSE;
1388 
1389  fhNEventsAfterCut->Fill(11.5);
1390 
1392  {
1393  if( !CheckForPrimaryVertex() ) return kFALSE; // algorithm in ESD/AOD Readers
1394  if( TMath::Abs(fVertex[0][0] ) < 1.e-6 &&
1395  TMath::Abs(fVertex[0][1] ) < 1.e-6 &&
1396  TMath::Abs(fVertex[0][2] ) < 1.e-6 ) return kFALSE;
1397  }
1398 
1399  AliDebug(1,"Pass primary vertex rejection");
1400 
1401  fhNEventsAfterCut->Fill(12.5);
1402 
1403  //printf("Reader : IsPileUp %d, Multi %d\n",IsPileUpFromSPD(),fInputEvent->IsPileupFromSPDInMultBins());
1404 
1406  {
1407  // Do not analyze events with pileup
1408  Bool_t bPileup = IsPileUpFromSPD();
1409  //IsPileupFromSPDInMultBins() // method to try
1410  //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]);
1411  if(bPileup) return kFALSE;
1412 
1413  AliDebug(1,"Pass Pile-Up event rejection");
1414 
1415  fhNEventsAfterCut->Fill(13.5);
1416  }
1417 
1419  {
1420  AliVVZERO* v0 = fInputEvent->GetVZEROData();
1421 
1422  Bool_t bV0AND = ((v0->GetV0ADecision()==1) && (v0->GetV0CDecision()==1));
1423  //bV0AND = fTriggerAnalysis->IsOfflineTriggerFired((AliESDEvent*)fInputEvent, AliTriggerAnalysis::kV0AND);
1424  //printf("V0AND event? %d\n",bV0AND);
1425 
1426  if(!bV0AND)
1427  {
1428  AliDebug(1,"Reject event by V0AND");
1429  return kFALSE;
1430  }
1431 
1432  AliDebug(1,"Pass V0AND event rejection");
1433 
1434  fhNEventsAfterCut->Fill(14.5);
1435  }
1436 
1437  //------------------------------------------------------
1438  // Check if there is a centrality value, PbPb analysis,
1439  // and if a centrality bin selection is requested
1440  //------------------------------------------------------
1441 
1442  //If we need a centrality bin, we select only those events in the corresponding bin.
1443  Int_t cen = -1;
1444  if ( fCentralityBin[0] >= 0 && fCentralityBin[1] >= 0 )
1445  {
1446  cen = GetEventCentrality();
1447 
1448  if(cen > fCentralityBin[1] || cen < fCentralityBin[0]) return kFALSE; //reject events out of bin.
1449 
1450  AliDebug(1,"Pass centrality rejection");
1451 
1452  fhNEventsAfterCut->Fill(15.5);
1453  }
1454 
1455  // Recover the weight assigned to the event, if provided
1456  // right now only for pT-hard bins and centrality depedent weights
1458  {
1460 
1461  fWeightUtils->SetPythiaEventHeader(((AliGenPythiaEventHeader*)GetGenEventHeader()));
1462 
1464  }
1465  //---------------------------------------------------------------------------
1466  // In case of MC analysis, set the range of interest in the MC particles list
1467  // The particle selection is done inside the FillInputDetector() methods
1468  //---------------------------------------------------------------------------
1470 
1471  //printf("N min %d, N max %d\n",fNMCProducedMin,fNMCProducedMax);
1472 
1473  //-------------------------------------------------------
1474  // Get the main vertex BC, in case not available
1475  // it is calculated in FillCTS checking the BC of tracks
1476  //------------------------------------------------------
1477  fVertexBC = fInputEvent->GetPrimaryVertex()->GetBC();
1478 
1479  //-----------------------------------------------
1480  // Fill the arrays with cluster/tracks/cells data
1481  //-----------------------------------------------
1482 
1483  if(fFillCTS)
1484  {
1485  FillInputCTS();
1486  //Accept events with at least one track
1487  if(fTrackMult == 0 && fDoRejectNoTrackEvents) return kFALSE ;
1488 
1489  fhNEventsAfterCut->Fill(16.5);
1490 
1491  AliDebug(1,"Pass rejection of null track events");
1492  }
1493 
1495  {
1496  if(fVertexBC != 0 && fVertexBC != AliVTrack::kTOFBCNA) return kFALSE ;
1497 
1498  AliDebug(1,"Pass rejection of events with vertex at BC!=0");
1499 
1500  fhNEventsAfterCut->Fill(17.5);
1501  }
1502 
1503  if(fFillEMCALCells)
1505 
1506  if(fFillPHOSCells)
1508 
1509  if(fFillEMCAL || fFillDCAL)
1510  FillInputEMCAL();
1511 
1512  if(fFillPHOS)
1513  FillInputPHOS();
1514 
1515  FillInputVZERO();
1516 
1517  //one specified jet branch
1522 
1523  AliDebug(1,"Event accepted for analysis");
1524 
1525  return kTRUE ;
1526 }
1527 
1528 //__________________________________________________
1531 //__________________________________________________
1533 {
1534  if(fUseAliCentrality)
1535  {
1536  if ( !GetCentrality() ) return -1;
1537 
1538  AliDebug(1,Form("Cent. Percentile: V0M %2.2f, CL0 %2.2f, CL1 %2.2f; selected class %s",
1539  GetCentrality()->GetCentralityPercentile("V0M"),
1540  GetCentrality()->GetCentralityPercentile("CL0"),
1541  GetCentrality()->GetCentralityPercentile("CL1"),
1542  fCentralityClass.Data()));
1543 
1544  if (fCentralityOpt == 100) return (Int_t) GetCentrality()->GetCentralityPercentile(fCentralityClass); // 100 bins max
1545  else if(fCentralityOpt == 10) return GetCentrality()->GetCentralityClass10(fCentralityClass);// 10 bins max
1546  else if(fCentralityOpt == 20) return GetCentrality()->GetCentralityClass5(fCentralityClass); // 20 bins max
1547  else
1548  {
1549  AliInfo(Form("Unknown centrality option %d, use 10, 20 or 100\n",fCentralityOpt));
1550  return -1;
1551  }
1552  }
1553  else
1554  {
1555  if ( !GetMultSelCen() ) return -1;
1556 
1557  AliDebug(1,Form("Mult. Percentile: V0M %2.2f, CL0 %2.2f, CL1 %2.2f; selected class %s",
1558  GetMultSelCen()->GetMultiplicityPercentile("V0M",1),
1559  GetMultSelCen()->GetMultiplicityPercentile("CL0",1),
1560  GetMultSelCen()->GetMultiplicityPercentile("CL1",1),
1561  fCentralityClass.Data()));
1562 
1563  return GetMultSelCen()->GetMultiplicityPercentile(fCentralityClass, kTRUE); // returns centrality only for events used in calibration
1564 
1565  // equivalent to
1566  //GetMultSelCen()->GetMultiplicityPercentile("V0M", kFALSE); // returns centrality for any event
1567  //Int_t qual = GetMultSelCen()->GetEvSelCode(); if (qual ! = 0) cent = qual;
1568  }
1569 }
1570 
1571 //_____________________________________________________
1574 //_____________________________________________________
1576 {
1577  if( !GetEventPlane() ) return -1000;
1578 
1579  Float_t ep = GetEventPlane()->GetEventplane(GetEventPlaneMethod(), GetInputEvent());
1580 
1581  if(GetEventPlaneMethod()=="Q" && (ep < 0 || ep > TMath::Pi()))
1582  {
1583  AliDebug(1,Form("Bad EP for <Q> method : %f\n",ep));
1584  return -1000;
1585  }
1586  else if(GetEventPlaneMethod().Contains("V0") )
1587  {
1588  if((ep > TMath::Pi()/2 || ep < -TMath::Pi()/2))
1589  {
1590  AliDebug(1,Form("Bad EP for <%s> method : %f\n",GetEventPlaneMethod().Data(), ep));
1591  return -1000;
1592  }
1593 
1594  ep+=TMath::Pi()/2; // put same range as for <Q> method
1595  }
1596 
1597  AliDebug(3,Form("Event plane angle %f",ep));
1598 
1599 // if(fDebug > 0 )
1600 // {
1601 // if (ep > TMath::Pi()) printf("AliCaloTrackReader::GetEventPlaneAngle() - Too large angle = %f\n",ep);
1602 // else if(ep < 0 ) printf("AliCaloTrackReader::GetEventPlaneAngle() - Negative angle = %f\n" ,ep);
1603 // }
1604 
1605  return ep;
1606 }
1607 
1608 //__________________________________________________________
1610 //__________________________________________________________
1612 {
1613  vertex[0] = fVertex[0][0];
1614  vertex[1] = fVertex[0][1];
1615  vertex[2] = fVertex[0][2];
1616 }
1617 
1618 //__________________________________________________________________________
1620 //__________________________________________________________________________
1621 void AliCaloTrackReader::GetVertex(Double_t vertex[3], Int_t evtIndex) const
1622 {
1623  vertex[0] = fVertex[evtIndex][0];
1624  vertex[1] = fVertex[evtIndex][1];
1625  vertex[2] = fVertex[evtIndex][2];
1626 }
1627 
1628 //________________________________________
1631 //________________________________________
1633 {
1634  // Delete previous vertex
1635  if(fVertex)
1636  {
1637  for (Int_t i = 0; i < fNMixedEvent; i++)
1638  {
1639  delete [] fVertex[i] ;
1640  }
1641  delete [] fVertex ;
1642  }
1643 
1644  fVertex = new Double_t*[fNMixedEvent] ;
1645  for (Int_t i = 0; i < fNMixedEvent; i++)
1646  {
1647  fVertex[i] = new Double_t[3] ;
1648  fVertex[i][0] = 0.0 ;
1649  fVertex[i][1] = 0.0 ;
1650  fVertex[i][2] = 0.0 ;
1651  }
1652 
1653  if (!fMixedEvent)
1654  { // Single event analysis
1655  if(fDataType!=kMC)
1656  {
1657 
1658  if(fInputEvent->GetPrimaryVertex())
1659  {
1660  fInputEvent->GetPrimaryVertex()->GetXYZ(fVertex[0]);
1661  }
1662  else
1663  {
1664  AliWarning("NULL primary vertex");
1665  fVertex[0][0]=0.; fVertex[0][1]=0.; fVertex[0][2]=0.;
1666  }//Primary vertex pointer do not exist
1667 
1668  } else
1669  {// MC read event
1670  fVertex[0][0]=0.; fVertex[0][1]=0.; fVertex[0][2]=0.;
1671  }
1672 
1673  AliDebug(1,Form("Single Event Vertex : %f,%f,%f\n",fVertex[0][0],fVertex[0][1],fVertex[0][2]));
1674 
1675  } else
1676  { // MultiEvent analysis
1677  for (Int_t iev = 0; iev < fNMixedEvent; iev++)
1678  {
1679  if (fMixedEvent->GetVertexOfEvent(iev))
1680  fMixedEvent->GetVertexOfEvent(iev)->GetXYZ(fVertex[iev]);
1681  else
1682  AliWarning("No vertex found");
1683 
1684  AliDebug(1,Form("Multi Event %d Vertex : %f,%f,%f",iev,fVertex[iev][0],fVertex[iev][1],fVertex[iev][2]));
1685  }
1686  }
1687 }
1688 
1689 //_____________________________________
1695 //_____________________________________
1697 {
1698  AliDebug(1,"Begin");
1699 
1700  Double_t pTrack[3] = {0,0,0};
1701 
1702  Int_t nTracks = fInputEvent->GetNumberOfTracks() ;
1703  fTrackMult = 0;
1704  Int_t nstatus = 0;
1705  Double_t bz = GetInputEvent()->GetMagneticField();
1706 
1707  for(Int_t i = 0; i < 19; i++)
1708  {
1709  fTrackBCEvent [i] = 0;
1710  fTrackBCEventCut[i] = 0;
1711  }
1712 
1713  Bool_t bc0 = kFALSE;
1714  if(fRecalculateVertexBC) fVertexBC = AliVTrack::kTOFBCNA;
1715 
1716  for (Int_t itrack = 0; itrack < nTracks; itrack++)
1717  {
1718  AliVTrack * track = (AliVTrack*)fInputEvent->GetTrack(itrack) ; // retrieve track from esd
1719 
1720  if(fAcceptOnlyHIJINGLabels && !IsHIJINGLabel(track->GetLabel())) continue ;
1721 
1722  fhCTSTrackCutsPt[0]->Fill(track->Pt());
1723 
1724  //Select tracks under certain conditions, TPCrefit, ITSrefit ... check the set bits
1725  ULong_t status = track->GetStatus();
1726 
1727  if (fTrackStatus && !((status & fTrackStatus) == fTrackStatus))
1728  continue ;
1729 
1730  fhCTSTrackCutsPt[1]->Fill(track->Pt());
1731 
1732  nstatus++;
1733 
1734  //-------------------------
1735  // Select the tracks depending on cuts of AOD or ESD
1736  if(!SelectTrack(track, pTrack)) continue ;
1737 
1738  fhCTSTrackCutsPt[2]->Fill(track->Pt());
1739 
1740  //-------------------------
1741  // TOF cuts
1742  Bool_t okTOF = ( (status & AliVTrack::kTOFout) == AliVTrack::kTOFout ) ;
1743  Double_t tof = -1000;
1744  Int_t trackBC = -1000 ;
1745 
1746  if(fAccessTrackTOF)
1747  {
1748  if(okTOF)
1749  {
1750  trackBC = track->GetTOFBunchCrossing(bz);
1751  SetTrackEventBC(trackBC+9);
1752 
1753  tof = track->GetTOFsignal()*1e-3;
1754 
1755  //After selecting tracks with small DCA, pointing to vertex, set vertex BC depeding on tracks BC
1757  {
1758  if (trackBC != 0 && trackBC != AliVTrack::kTOFBCNA) fVertexBC = trackBC;
1759  else if(trackBC == 0) bc0 = kTRUE;
1760  }
1761 
1762  //In any case, the time should to be larger than the fixed window ...
1763  if( fUseTrackTimeCut && (trackBC !=0 || tof < fTrackTimeCutMin || tof > fTrackTimeCutMax) )
1764  {
1765  //printf("Remove track time %f and bc = %d\n",tof,trackBC);
1766  continue ;
1767  }
1768  //else printf("Accept track time %f and bc = %d\n",tof,trackBC);
1769  }
1770  }
1771 
1772  fhCTSTrackCutsPt[3]->Fill(track->Pt());
1773 
1774  //---------------------
1775  // DCA cuts
1776  //
1777  fMomentum.SetPxPyPzE(pTrack[0],pTrack[1],pTrack[2],0);
1778 
1779  if(fUseTrackDCACut)
1780  {
1781  Float_t dcaTPC =-999;
1782  //In case of AODs, TPC tracks cannot be propagated back to primary vertex,
1783  if( fDataType == kAOD ) dcaTPC = ((AliAODTrack*) track)->DCA();
1784 
1785  //normal way to get the dca, cut on dca_xy
1786  if(dcaTPC==-999)
1787  {
1788  Double_t dca[2] = {1e6,1e6};
1789  Double_t covar[3] = {1e6,1e6,1e6};
1790  Bool_t okDCA = track->PropagateToDCA(fInputEvent->GetPrimaryVertex(),bz,100.,dca,covar);
1791  if( okDCA) okDCA = AcceptDCA(fMomentum.Pt(),dca[0]);
1792  if(!okDCA)
1793  {
1794  //printf("AliCaloTrackReader::FillInputCTS() - Reject track pt %2.2f, dca_xy %2.4f\n",fMomentum.Pt(),dca[0]);
1795  continue ;
1796  }
1797  }
1798  }// DCA cuts
1799 
1800  fhCTSTrackCutsPt[4]->Fill(track->Pt());
1801 
1802  //-------------------------
1803  // Kinematic/acceptance cuts
1804  //
1805  // Count the tracks in eta < 0.9
1806  if(TMath::Abs(track->Eta())< fTrackMultEtaCut) fTrackMult++;
1807 
1808  if(fCTSPtMin > fMomentum.Pt() || fCTSPtMax < fMomentum.Pt()) continue ;
1809 
1810  // Check effect of cuts on track BC
1811  if(fAccessTrackTOF && okTOF) SetTrackEventBCcut(trackBC+9);
1812 
1813  if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(fMomentum.Eta(),fMomentum.Phi(),kCTS)) continue;
1814 
1815  fhCTSTrackCutsPt[5]->Fill(track->Pt());
1816 
1817  // ------------------------------
1818  // Add selected tracks to array
1819  AliDebug(2,Form("Selected tracks pt %3.2f, phi %3.2f deg, eta %3.2f",
1820  fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta()));
1821 
1822  fCTSTracks->Add(track);
1823 
1824  // TODO, check if remove
1825  if (fMixedEvent) track->SetID(itrack);
1826 
1827  }// track loop
1828 
1829  if( fRecalculateVertexBC && (fVertexBC == 0 || fVertexBC == AliVTrack::kTOFBCNA))
1830  {
1831  if( bc0 ) fVertexBC = 0 ;
1832  else fVertexBC = AliVTrack::kTOFBCNA ;
1833  }
1834 
1835  AliDebug(1,Form("AOD entries %d, input tracks %d, pass status %d, multipliticy %d", fCTSTracks->GetEntriesFast(), nTracks, nstatus, fTrackMult));//fCTSTracksNormalInputEntries);
1836 }
1837 
1838 //_______________________________________________________________________________
1856 //_______________________________________________________________________________
1858 {
1859  // Accept clusters with the proper label, TO DO, use the new more General methods!!!
1860  if(fAcceptOnlyHIJINGLabels && !IsHIJINGLabel( clus->GetLabel() )) return ;
1861 
1862  // TODO, not sure if needed anymore
1863  Int_t vindex = 0 ;
1864  if (fMixedEvent)
1865  vindex = fMixedEvent->EventIndexForCaloCluster(iclus);
1866 
1867  clus->GetMomentum(fMomentum, fVertex[vindex]);
1868 
1869  // No correction/cut applied yet
1870  fhEMCALClusterCutsE[0]->Fill(clus->E());
1871 
1872  //if( (fDebug > 2 && fMomentum.E() > 0.1) || fDebug > 10 )
1873  AliDebug(2,Form("Input cluster E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f",
1874  fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta()));
1875 
1876  //---------------------------
1877  // Embedding case
1878  // TO BE REVISED
1880  {
1881  if(clus->GetNLabels()==0 || clus->GetLabel() < 0) return;
1882  //else printf("Embedded cluster, %d, n label %d label %d \n",iclus,clus->GetNLabels(),clus->GetLabel());
1883  }
1884 
1885  //--------------------------------------
1886  // Apply some corrections in the cluster
1887  //
1889  {
1890  //Recalibrate the cluster energy
1892  {
1894 
1895  clus->SetE(energy);
1896  //printf("Recalibrated Energy %f\n",clus->E());
1897 
1900 
1901  } // recalculate E
1902 
1903  //Recalculate distance to bad channels, if new list of bad channels provided
1905 
1906  //Recalculate cluster position
1907  if(GetCaloUtils()->IsRecalculationOfClusterPositionOn())
1908  {
1910  //clus->GetPosition(pos);
1911  //printf("After Corrections: e %f, x %f, y %f, z %f\n",clus->E(),pos[0],pos[1],pos[2]);
1912  }
1913 
1914  // Recalculate TOF
1915  if(GetCaloUtils()->GetEMCALRecoUtils()->IsTimeRecalibrationOn())
1916  {
1917  Double_t tof = clus->GetTOF();
1918  Float_t frac =-1;
1919  Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(fEMCALCells, clus,frac);
1920 
1921  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absIdMax,fInputEvent->GetBunchCrossNumber(),tof);
1922 
1923  //additional L1 phase shift
1924  if(GetCaloUtils()->GetEMCALRecoUtils()->IsL1PhaseInTimeRecalibrationOn())
1925  {
1926  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTimeL1Phase(GetCaloUtils()->GetEMCALGeometry()->GetSuperModuleNumber(absIdMax),
1927  fInputEvent->GetBunchCrossNumber(), tof);
1928  }
1929 
1930  clus->SetTOF(tof);
1931 
1932  }// Time recalibration
1933  }
1934 
1935  // Check effect of corrections
1936  fhEMCALClusterCutsE[1]->Fill(clus->E());
1937 
1938  //-----------------------------------------------------------------
1939  // Reject clusters with bad channels, close to borders and exotic
1940  //
1941  Bool_t goodCluster = GetCaloUtils()->GetEMCALRecoUtils()->IsGoodCluster(clus,
1942  GetCaloUtils()->GetEMCALGeometry(),
1943  GetEMCALCells(),fInputEvent->GetBunchCrossNumber());
1944 
1945  if(!goodCluster)
1946  {
1947  //if( (fDebug > 2 && fMomentum.E() > 0.1) || fDebug > 10 )
1948  AliDebug(2,Form("Bad cluster E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f",
1949  fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta()));
1950 
1951  return;
1952  }
1953 
1954  // Check effect of bad cluster removal
1955  fhEMCALClusterCutsE[2]->Fill(clus->E());
1956 
1957  //Float_t pos[3];
1958  //clus->GetPosition(pos);
1959  //printf("Before Corrections: e %f, x %f, y %f, z %f\n",clus->E(),pos[0],pos[1],pos[2]);
1960 
1961  //--------------------------------------
1962  // Correct non linearity or smear energy
1963  //
1965  {
1967 
1968  //if( (fDebug > 5 && fMomentum.E() > 0.1) || fDebug > 10 )
1969  AliDebug(5,Form("Correct Non Lin: Old E %3.2f, New E %3.2f",
1970  fMomentum.E(),clus->E()));
1971 
1972  // In case of MC analysis, to match resolution/calibration in real data
1973  // Not needed anymore, just leave for MC studies on systematics
1974  if( GetCaloUtils()->GetEMCALRecoUtils()->IsClusterEnergySmeared() )
1975  {
1976  Float_t rdmEnergy = GetCaloUtils()->GetEMCALRecoUtils()->SmearClusterEnergy(clus);
1977 
1978  //if( (fDebug > 5 && fMomentum.E() > 0.1) || fDebug > 10 )
1979  AliDebug(5,Form("Smear energy: Old E %3.2f, New E %3.2f",clus->E(),rdmEnergy));
1980 
1981  clus->SetE(rdmEnergy);
1982  }
1983  }
1984 
1985  clus->GetMomentum(fMomentum, fVertex[vindex]);
1986 
1987  // Check effect linearity correction, energy smearing
1988  fhEMCALClusterCutsE[3]->Fill(clus->E());
1989 
1990  // Check the event BC depending on EMCal clustr before final cuts
1991  Double_t tof = clus->GetTOF()*1e9;
1992 
1993  Int_t bc = TMath::Nint(tof/50) + 9;
1994  //printf("tof %2.2f, bc+5=%d\n",tof,bc);
1995 
1996  SetEMCalEventBC(bc);
1997 
1998  //--------------------------------------
1999  // Apply some kinematical/acceptance cuts
2000  //
2001  if(fEMCALPtMin > clus->E() || fEMCALPtMax < clus->E()) return ;
2002 
2003  // Select cluster fiducial region
2004  //
2005  Bool_t bEMCAL = kFALSE;
2006  Bool_t bDCAL = kFALSE;
2007  if(fCheckFidCut)
2008  {
2009  if(fFillEMCAL && fFiducialCut->IsInFiducialCut(fMomentum.Eta(),fMomentum.Phi(),kEMCAL)) bEMCAL = kTRUE ;
2010  if(fFillDCAL && fFiducialCut->IsInFiducialCut(fMomentum.Eta(),fMomentum.Phi(),kDCAL )) bDCAL = kTRUE ;
2011  }
2012  else
2013  {
2014  bEMCAL = kTRUE;
2015  }
2016 
2017  //---------------------------------------------------------------------
2018  // Mask all cells in collumns facing ALICE thick material if requested
2019  //
2021  {
2022  Int_t absId = -1;
2023  Int_t iSupMod = -1;
2024  Int_t iphi = -1;
2025  Int_t ieta = -1;
2026  Bool_t shared = kFALSE;
2027  GetCaloUtils()->GetEMCALRecoUtils()->GetMaxEnergyCell(GetCaloUtils()->GetEMCALGeometry(), GetEMCALCells(),clus,absId,iSupMod,ieta,iphi,shared);
2028 
2029  AliDebug(2,Form("Masked collumn: cluster E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f",
2030  fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta()));
2031 
2032 
2033  if(GetCaloUtils()->MaskFrameCluster(iSupMod, ieta)) return;
2034  }
2035 
2036  // Check effect of energy and fiducial cuts
2037  fhEMCALClusterCutsE[4]->Fill(clus->E());
2038 
2039 
2040  //------------------------------------------
2041  // Apply time cut, count EMCal BC before cut
2042  //
2043  SetEMCalEventBCcut(bc);
2044 
2045  if(!IsInTimeWindow(tof,clus->E()))
2046  {
2047  fNPileUpClusters++ ;
2048  if(fUseEMCALTimeCut)
2049  {
2050  AliDebug(2,Form("Out of time window E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f, time %e",
2051  fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta(),tof));
2052 
2053  return ;
2054  }
2055  }
2056  else
2058 
2059  // Check effect of time cut
2060  fhEMCALClusterCutsE[5]->Fill(clus->E());
2061 
2062 
2063  //----------------------------------------------------
2064  // Apply N cells cut
2065  //
2066  if(clus->GetNCells() <= fEMCALNCellsCut && fDataType != AliCaloTrackReader::kMC) return ;
2067 
2068  // Check effect of n cells cut
2069  fhEMCALClusterCutsE[6]->Fill(clus->E());
2070 
2071  //----------------------------------------------------
2072  // Apply distance to bad channel cut
2073  //
2074  Double_t distBad = clus->GetDistanceToBadChannel() ; //Distance to bad channel
2075 
2076  if(distBad < 0.) distBad=9999. ; //workout strange convension dist = -1. ;
2077 
2078  if(distBad < fEMCALBadChMinDist) return ;
2079 
2080  // Check effect distance to bad channel cut
2081  fhEMCALClusterCutsE[7]->Fill(clus->E());
2082 
2083  //----------------------------------------------------
2084  // Smear the SS to try to match data and simulations,
2085  // do it only for simulations.
2086  //
2087  Int_t nMaxima = GetCaloUtils()->GetNumberOfLocalMaxima(clus, GetEMCALCells());
2088  // Int_t nMaxima = clus->GetNExMax(); // For Run2
2089  if( fSmearShowerShape && clus->GetNCells() > 2 &&
2090  nMaxima >= fSmearNLMMin && nMaxima <= fSmearNLMMax )
2091  {
2092  AliDebug(2,Form("Smear shower shape - Original: %2.4f\n", clus->GetM02()));
2094  {
2095  clus->SetM02( clus->GetM02() + fRandom.Landau(0, fSmearShowerShapeWidth) );
2096  }
2098  {
2099  if(iclus%3 == 0 && clus->GetM02() > 0.1) clus->SetM02( clus->GetM02() + fRandom.Landau(0.05, fSmearShowerShapeWidth) ); //fSmearShowerShapeWidth = 0.035
2100  }
2101  //clus->SetM02( fRandom.Landau(clus->GetM02(), fSmearShowerShapeWidth) );
2102  AliDebug(2,Form("Width %2.4f Smeared : %2.4f\n", fSmearShowerShapeWidth,clus->GetM02()));
2103  }
2104 
2105  //--------------------------------------------------------
2106  // Fill the corresponding array with the selected clusters
2107  // Usually just filling EMCal array with upper or lower clusters is enough,
2108  // but maybe we want to do EMCal-DCal correlations.
2109 
2110  //if((fDebug > 2 && fMomentum.E() > 0.1) || fDebug > 10)
2111  AliDebug(2,Form("Selected clusters (EMCAL%d, DCAL%d), E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f",
2112  bEMCAL,bDCAL,fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta()));
2113 
2114 
2115  if (bEMCAL) fEMCALClusters->Add(clus);
2116  else if(bDCAL ) fDCALClusters ->Add(clus);
2117 
2118  // TODO, not sure if needed anymore
2119  if (fMixedEvent)
2120  clus->SetID(iclus) ;
2121 }
2122 
2123 //_______________________________________
2130 //_______________________________________
2132 {
2133  AliDebug(1,"Begin");
2134 
2135  // First recalibrate cells, time or energy
2136  // if(GetCaloUtils()->IsRecalibrationOn())
2137  // GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCells(GetCaloUtils()->GetEMCALGeometry(),
2138  // GetEMCALCells(),
2139  // fInputEvent->GetBunchCrossNumber());
2140 
2141  fNPileUpClusters = 0; // Init counter
2142  fNNonPileUpClusters = 0; // Init counter
2143  for(Int_t i = 0; i < 19; i++)
2144  {
2145  fEMCalBCEvent [i] = 0;
2146  fEMCalBCEventCut[i] = 0;
2147  }
2148 
2149  //Loop to select clusters in fiducial cut and fill container with aodClusters
2150  if(fEMCALClustersListName=="")
2151  {
2152  Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
2153  for (Int_t iclus = 0; iclus < nclusters; iclus++)
2154  {
2155  AliVCluster * clus = 0;
2156  if ( (clus = fInputEvent->GetCaloCluster(iclus)) )
2157  {
2158  if (clus->IsEMCAL())
2159  {
2160  FillInputEMCALAlgorithm(clus, iclus);
2161  }//EMCAL cluster
2162  }// cluster exists
2163  }// cluster loop
2164 
2165  //Recalculate track matching
2167 
2168  }//Get the clusters from the input event
2169  else
2170  {
2171  TClonesArray * clusterList = 0x0;
2172 
2173  if (fInputEvent->FindListObject(fEMCALClustersListName))
2174  {
2175  clusterList = dynamic_cast<TClonesArray*> (fInputEvent->FindListObject(fEMCALClustersListName));
2176  }
2177  else if(fOutputEvent)
2178  {
2179  clusterList = dynamic_cast<TClonesArray*> (fOutputEvent->FindListObject(fEMCALClustersListName));
2180  }
2181 
2182  if(!clusterList)
2183  {
2184  AliWarning(Form("Wrong name of list with clusters? <%s>",fEMCALClustersListName.Data()));
2185  return;
2186  }
2187 
2188  Int_t nclusters = clusterList->GetEntriesFast();
2189  for (Int_t iclus = 0; iclus < nclusters; iclus++)
2190  {
2191  AliVCluster * clus = dynamic_cast<AliVCluster*> (clusterList->At(iclus));
2192  //printf("E %f\n",clus->E());
2193  if (clus) FillInputEMCALAlgorithm(clus, iclus);
2194  else AliWarning("Null cluster in list!");
2195  }// cluster loop
2196 
2197  // Recalculate the pile-up time, in case long time clusters removed during clusterization
2198  //printf("Input event INIT : Pile-up clusters %d, NO pile-up %d\n",fNPileUpClusters,fNNonPileUpClusters);
2199 
2200  fNPileUpClusters = 0; // Init counter
2201  fNNonPileUpClusters = 0; // Init counter
2202  for(Int_t i = 0; i < 19; i++)
2203  {
2204  fEMCalBCEvent [i] = 0;
2205  fEMCalBCEventCut[i] = 0;
2206  }
2207 
2208  for (Int_t iclus = 0; iclus < fInputEvent->GetNumberOfCaloClusters(); iclus++)
2209  {
2210  AliVCluster * clus = 0;
2211 
2212  if ( (clus = fInputEvent->GetCaloCluster(iclus)) )
2213  {
2214  if (clus->IsEMCAL())
2215  {
2216 
2217  Float_t frac =-1;
2218  Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(fEMCALCells, clus,frac);
2219  Double_t tof = clus->GetTOF();
2220  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absIdMax,fInputEvent->GetBunchCrossNumber(),tof);
2221  //additional L1 phase shift
2222  if(GetCaloUtils()->GetEMCALRecoUtils()->IsL1PhaseInTimeRecalibrationOn())
2223  {
2224  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTimeL1Phase(GetCaloUtils()->GetEMCALGeometry()->GetSuperModuleNumber(absIdMax), fInputEvent->GetBunchCrossNumber(), tof);
2225  }
2226 
2227  tof*=1e9;
2228 
2229  //printf("Input event cluster : AbsIdMax %d, E %2.2f, time %2.2f \n", absIdMax,clus->E(),tof);
2230 
2231  //Reject clusters with bad channels, close to borders and exotic;
2232  if(!GetCaloUtils()->GetEMCALRecoUtils()->IsGoodCluster(clus,GetCaloUtils()->GetEMCALGeometry(),GetEMCALCells(),fInputEvent->GetBunchCrossNumber())) continue;
2233 
2234  Int_t bc = TMath::Nint(tof/50) + 9;
2235  SetEMCalEventBC(bc);
2236 
2237  if(fEMCALPtMin > clus->E() || fEMCALPtMax < clus->E()) continue ;
2238 
2239  clus->GetMomentum(fMomentum, fVertex[0]);
2240 
2241  if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(fMomentum.Eta(),fMomentum.Phi(),kEMCAL)) return ;
2242 
2243  SetEMCalEventBCcut(bc);
2244 
2245  if(!IsInTimeWindow(tof,clus->E()))
2246  fNPileUpClusters++ ;
2247  else
2249 
2250  }
2251  }
2252  }
2253 
2254  //printf("Input event : Pile-up clusters %d, NO pile-up %d\n",fNPileUpClusters,fNNonPileUpClusters);
2255 
2256  // Recalculate track matching, not necessary if already done in the reclusterization task.
2257  // in case it was not done ...
2259 
2260  }
2261 
2262  AliDebug(1,Form("AOD entries %d, n pile-up clusters %d, n non pile-up %d", fEMCALClusters->GetEntriesFast(),fNPileUpClusters,fNNonPileUpClusters));
2263 }
2264 
2265 //_______________________________________
2267 //_______________________________________
2269 {
2270  AliDebug(1,"Begin");
2271 
2272  // Loop to select clusters in fiducial cut and fill container with aodClusters
2273  Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
2274  for (Int_t iclus = 0; iclus < nclusters; iclus++)
2275  {
2276  AliVCluster * clus = fInputEvent->GetCaloCluster(iclus) ;
2277  if ( !clus ) continue ;
2278 
2279  if ( !clus->IsPHOS() ) continue ;
2280 
2281  if(fAcceptOnlyHIJINGLabels && !IsHIJINGLabel(clus->GetLabel())) continue ;
2282 
2283  fhPHOSClusterCutsE[0]->Fill(clus->E());
2284 
2285  // Skip CPV input
2286  if( clus->GetType() == AliVCluster::kPHOSCharged ) continue ;
2287 
2288  fhPHOSClusterCutsE[1]->Fill(clus->E());
2289 
2290  //---------------------------------------------
2291  // Old feature, try to rely on PHOS tender
2292  //
2294  {
2295  // Recalibrate the cluster energy
2297  {
2298  Float_t energy = GetCaloUtils()->RecalibrateClusterEnergy(clus, (AliAODCaloCells*)GetPHOSCells());
2299  clus->SetE(energy);
2300  }
2301  }
2302 
2303  //----------------------------------------------------------------------------------
2304  // Check if the cluster contains any bad channel and if close to calorimeter borders
2305  //
2306  // Old feature, try to rely on PHOS tender
2307  if( GetCaloUtils()->ClusterContainsBadChannel(kPHOS,clus->GetCellsAbsId(), clus->GetNCells()))
2308  continue;
2309 
2310  if(!GetCaloUtils()->CheckCellFiducialRegion(clus, fInputEvent->GetPHOSCells()))
2311  continue;
2312 
2313  // TODO, add exotic cut???
2314 
2315  fhPHOSClusterCutsE[2]->Fill(clus->E());
2316 
2317  // TODO Dead code? remove?
2318  Int_t vindex = 0 ;
2319  if (fMixedEvent)
2320  vindex = fMixedEvent->EventIndexForCaloCluster(iclus);
2321 
2322  clus->GetMomentum(fMomentum, fVertex[vindex]);
2323 
2324  //----------------------------------------------------------------------------------
2325  // Remove clusters close to borders
2326  //
2328  continue ;
2329 
2330  fhPHOSClusterCutsE[3]->Fill(clus->E());
2331 
2332  //----------------------------------------------------------------------------------
2333  // Remove clusters with too low energy
2334  //
2335  if (fPHOSPtMin > fMomentum.E() || fPHOSPtMax < fMomentum.E() )
2336  continue ;
2337 
2338  fhPHOSClusterCutsE[4]->Fill(clus->E());
2339 
2340  //----------------------------------------------------
2341  // Apply N cells cut
2342  //
2343  if(clus->GetNCells() <= fPHOSNCellsCut && fDataType != AliCaloTrackReader::kMC) return ;
2344 
2345  // Check effect of n cells cut
2346  fhPHOSClusterCutsE[5]->Fill(clus->E());
2347 
2348  //----------------------------------------------------
2349  // Apply distance to bad channel cut
2350  //
2351  Double_t distBad = clus->GetDistanceToBadChannel() ; //Distance to bad channel
2352 
2353  if(distBad < 0.) distBad=9999. ; //workout strange convension dist = -1. ;
2354 
2355  if(distBad < fPHOSBadChMinDist) return ;
2356 
2357  // Check effect distance to bad channel cut
2358  fhPHOSClusterCutsE[6]->Fill(clus->E());
2359 
2360  // TODO, add time cut
2361 
2362  //----------------------------------------------------------------------------------
2363  // Add selected clusters to array
2364  //
2365  //if(fDebug > 2 && fMomentum.E() > 0.1)
2366  AliDebug(2,Form("Selected clusters E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f",
2367  fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta()));
2368 
2369  fPHOSClusters->Add(clus);
2370 
2371  // TODO Dead code? remove?
2372  if (fMixedEvent)
2373  clus->SetID(iclus) ;
2374 
2375  } // esd/aod cluster loop
2376 
2377  AliDebug(1,Form("AOD entries %d",fPHOSClusters->GetEntriesFast())) ;
2378 }
2379 
2380 //____________________________________________
2382 //____________________________________________
2384 {
2385  fEMCALCells = fInputEvent->GetEMCALCells();
2386 }
2387 
2388 //___________________________________________
2390 //___________________________________________
2392 {
2393  fPHOSCells = fInputEvent->GetPHOSCells();
2394 }
2395 
2396 //_______________________________________
2399 //_______________________________________
2401 {
2402  AliVVZERO* v0 = fInputEvent->GetVZEROData();
2403  //printf("Init V0: ADC (%d,%d), Multiplicity (%d,%d) \n",fV0ADC[0],fV0ADC[1],fV0Mul[0],fV0Mul[1]);
2404 
2405  if (v0)
2406  {
2407  AliESDVZERO* esdV0 = dynamic_cast<AliESDVZERO*> (v0);
2408  for (Int_t i = 0; i < 32; i++)
2409  {
2410  if(esdV0)
2411  {//Only available in ESDs
2412  fV0ADC[0] += (Int_t)esdV0->GetAdcV0C(i);
2413  fV0ADC[1] += (Int_t)esdV0->GetAdcV0A(i);
2414  }
2415 
2416  fV0Mul[0] += (Int_t)v0->GetMultiplicityV0C(i);
2417  fV0Mul[1] += (Int_t)v0->GetMultiplicityV0A(i);
2418  }
2419 
2420  AliDebug(1,Form("ADC (%d,%d), Multiplicity (%d,%d)",fV0ADC[0],fV0ADC[1],fV0Mul[0],fV0Mul[1]));
2421  }
2422  else
2423  {
2424  AliDebug(1,"Cannot retrieve V0 ESD! Run w/ null V0 charges");
2425  }
2426 }
2427 
2428 //_________________________________________________
2432 //_________________________________________________
2434 {
2435  AliDebug(2,"Begin");
2436 
2437  //
2438  //check if branch name is given
2439  if(!fInputNonStandardJetBranchName.Length())
2440  {
2441  fInputEvent->Print();
2442  AliFatal("No non-standard jet branch name specified. Specify among existing ones.");
2443  }
2444 
2445  fNonStandardJets = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(fInputNonStandardJetBranchName.Data()));
2446 
2447  if(!fNonStandardJets)
2448  {
2449  //check if jet branch exist; exit if not
2450  fInputEvent->Print();
2451 
2452  AliFatal(Form("%s:%d no reconstructed jet array with name %s in AOD", (char*)__FILE__,__LINE__,fInputNonStandardJetBranchName.Data()));
2453  }
2454  else
2455  {
2456  AliDebug(1,Form("AOD input jets %d", fNonStandardJets->GetEntriesFast()));
2457  }
2458 }
2459 
2460 //_________________________________________________
2464 //_________________________________________________
2466 {
2467  AliDebug(1,"Begin");
2468  //
2469  //check if branch name is given
2470  if(!fInputBackgroundJetBranchName.Length())
2471  {
2472  fInputEvent->Print();
2473 
2474  AliFatal("No background jet branch name specified. Specify among existing ones.");
2475  }
2476 
2477  fBackgroundJets = (AliAODJetEventBackground*)(fInputEvent->FindListObject(fInputBackgroundJetBranchName.Data()));
2478 
2479  if(!fBackgroundJets)
2480  {
2481  //check if jet branch exist; exit if not
2482  fInputEvent->Print();
2483 
2484  AliFatal(Form("%s:%d no reconstructed jet array with name %s in AOD", (char*)__FILE__,__LINE__,fInputBackgroundJetBranchName.Data()));
2485  }
2486  else
2487  {
2488  AliDebug(1,"FillInputBackgroundJets");
2489  fBackgroundJets->Print("");
2490  }
2491 }
2492 
2493 //____________________________________________________________________
2499 //____________________________________________________________________
2501 {
2502  // init some variables
2503  Int_t trigtimes[30], globCol, globRow,ntimes, i;
2504  Int_t absId = -1; //[100];
2505  Int_t nPatch = 0;
2506 
2507  TArrayI patches(0);
2508 
2509  // get object pointer
2510  AliVCaloTrigger *caloTrigger = GetInputEvent()->GetCaloTrigger( "EMCAL" );
2511 
2512  if(!caloTrigger)
2513  {
2514  AliError("Trigger patches input (AliVCaloTrigger) not available in data!");
2515  return patches;
2516  }
2517 
2518  //printf("CaloTrigger Entries %d\n",caloTrigger->GetEntries() );
2519 
2520  // class is not empty
2521  if( caloTrigger->GetEntries() > 0 )
2522  {
2523  // must reset before usage, or the class will fail
2524  caloTrigger->Reset();
2525 
2526  // go throuth the trigger channels
2527  while( caloTrigger->Next() )
2528  {
2529  // get position in global 2x2 tower coordinates
2530  caloTrigger->GetPosition( globCol, globRow );
2531 
2532  //L0
2533  if(IsEventEMCALL0())
2534  {
2535  // get dimension of time arrays
2536  caloTrigger->GetNL0Times( ntimes );
2537 
2538  // no L0s in this channel
2539  // presence of the channel in the iterator still does not guarantee that L0 was produced!!
2540  if( ntimes < 1 )
2541  continue;
2542 
2543  // get timing array
2544  caloTrigger->GetL0Times( trigtimes );
2545  //printf("Get L0 patch : n times %d - trigger time window %d - %d\n",ntimes, tmin,tmax);
2546 
2547  // go through the array
2548  for( i = 0; i < ntimes; i++ )
2549  {
2550  // check if in cut - 8,9 shall be accepted in 2011
2551  if( trigtimes[i] >= tmin && trigtimes[i] <= tmax )
2552  {
2553  //printf("Accepted trigger time %d \n",trigtimes[i]);
2554  //if(nTrig > 99) continue;
2555  GetCaloUtils()->GetEMCALGeometry()->GetAbsFastORIndexFromPositionInEMCAL(globCol,globRow, absId);
2556  //printf("pass the time cut globCol %d, globRow %d absId %d\n",globCol,globRow, absId);
2557  patches.Set(nPatch+1);
2558  patches.AddAt(absId,nPatch++);
2559  }
2560  } // trigger time array
2561  }//L0
2562  else if(IsEventEMCALL1()) // L1
2563  {
2564  Int_t bit = 0;
2565  caloTrigger->GetTriggerBits(bit);
2566 
2567  Int_t sum = 0;
2568  caloTrigger->GetL1TimeSum(sum);
2569  //fBitEGA-=2;
2570  Bool_t isEGA1 = ((bit >> fBitEGA ) & 0x1) && IsEventEMCALL1Gamma1() ;
2571  Bool_t isEGA2 = ((bit >> (fBitEGA+1)) & 0x1) && IsEventEMCALL1Gamma2() ;
2572  Bool_t isEJE1 = ((bit >> fBitEJE ) & 0x1) && IsEventEMCALL1Jet1 () ;
2573  Bool_t isEJE2 = ((bit >> (fBitEJE+1)) & 0x1) && IsEventEMCALL1Jet2 () ;
2574 
2575  //if((bit>> fBitEGA )&0x1) printf("Trig Bit %d - bit %d - EG1 %d - EG2 %d\n",fBitEGA ,bit,IsEventEMCALL1Gamma1(),IsEventEMCALL1Gamma2());
2576  //if((bit>>(fBitEGA+1))&0x1) printf("Trig Bit %d - bit %d - EG1 %d - EG2 %d\n",fBitEGA+1,bit,IsEventEMCALL1Gamma1(),IsEventEMCALL1Gamma2());
2577 
2578  if(!isEGA1 && !isEJE1 && !isEGA2 && !isEJE2) continue;
2579 
2580  Int_t patchsize = -1;
2581  if (isEGA1 || isEGA2) patchsize = 2;
2582  else if (isEJE1 || isEJE2) patchsize = 16;
2583 
2584  //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",
2585  // bit,sum,patchsize,isEGA1,isEGA2,isEJE1,isEJE2,fBitEGA,fBitEJE,IsEventEMCALL1Gamma(),IsEventEMCALL1Jet());
2586 
2587 
2588  // add 2x2 (EGA) or 16x16 (EJE) patches
2589  for(Int_t irow=0; irow < patchsize; irow++)
2590  {
2591  for(Int_t icol=0; icol < patchsize; icol++)
2592  {
2593  GetCaloUtils()->GetEMCALGeometry()->GetAbsFastORIndexFromPositionInEMCAL(globCol+icol,globRow+irow, absId);
2594  //printf("pass the time cut globCol %d, globRow %d absId %d\n",globCol,globRow, absId);
2595  patches.Set(nPatch+1);
2596  patches.AddAt(absId,nPatch++);
2597  }
2598  }
2599 
2600  } // L1
2601 
2602  } // trigger iterator
2603  } // go through triggers
2604 
2605  if(patches.GetSize()<=0) AliInfo(Form("No patch found! for triggers: %s and selected <%s>",
2607  //else printf(">>>>> N patches %d, test %d,first %d, last %d\n",patches.GetSize(), nPatch, patches.At(0), patches.At(patches.GetSize()-1));
2608 
2609  return patches;
2610 }
2611 
2612 //____________________________________________________________
2616 //____________________________________________________________
2618 {
2619  // Init info from previous event
2620  fTriggerClusterIndex = -1;
2621  fTriggerClusterId = -1;
2622  fTriggerClusterBC = -10000;
2623  fIsExoticEvent = kFALSE;
2624  fIsBadCellEvent = kFALSE;
2625  fIsBadMaxCellEvent = kFALSE;
2626 
2627  fIsTriggerMatch = kFALSE;
2628  fIsTriggerMatchOpenCut[0] = kFALSE;
2629  fIsTriggerMatchOpenCut[1] = kFALSE;
2630  fIsTriggerMatchOpenCut[2] = kFALSE;
2631 
2632  // Do only analysis for triggered events
2633  if(!IsEventEMCALL1() && !IsEventEMCALL0())
2634  {
2635  fTriggerClusterBC = 0;
2636  return;
2637  }
2638 
2639  //printf("***** Try to match trigger to cluster %d **** L0 %d, L1 %d\n",fTriggerPatchClusterMatch,IsEventEMCALL0(),IsEventEMCALL1());
2640 
2641  //Recover the list of clusters
2642  TClonesArray * clusterList = 0;
2643  if (fInputEvent->FindListObject(fEMCALClustersListName))
2644  {
2645  clusterList = dynamic_cast<TClonesArray*> (fInputEvent->FindListObject(fEMCALClustersListName));
2646  }
2647  else if(fOutputEvent)
2648  {
2649  clusterList = dynamic_cast<TClonesArray*> (fOutputEvent->FindListObject(fEMCALClustersListName));
2650  }
2651 
2652  // Get number of clusters and of trigger patches
2653  Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
2654  if(clusterList)
2655  nclusters = clusterList->GetEntriesFast();
2656 
2657  Int_t nPatch = patches.GetSize();
2658  Float_t exoDiffTime = GetCaloUtils()->GetEMCALRecoUtils()->GetExoticCellDiffTimeCut();
2659 
2660  //Init some variables used in the cluster loop
2661  Float_t tofPatchMax = 100000;
2662  Float_t ePatchMax =-1;
2663 
2664  Float_t tofMax = 100000;
2665  Float_t eMax =-1;
2666 
2667  Int_t clusMax =-1;
2668  Int_t idclusMax =-1;
2669  Bool_t badClMax = kFALSE;
2670  Bool_t badCeMax = kFALSE;
2671  Bool_t exoMax = kFALSE;
2672 // Int_t absIdMaxTrig= -1;
2673  Int_t absIdMaxMax = -1;
2674 
2675  Int_t nOfHighECl = 0 ;
2676 
2677  //
2678  // Check what is the trigger threshold
2679  // set minimu energym of candidate for trigger cluster
2680  //
2682 
2683  Float_t triggerThreshold = fTriggerL1EventThreshold;
2684  if(IsEventEMCALL0()) triggerThreshold = fTriggerL0EventThreshold;
2685  Float_t minE = triggerThreshold / 2.;
2686 
2687  // This method is not really suitable for JET trigger
2688  // but in case, reduce the energy cut since we do not trigger on high energy particle
2689  if(IsEventEMCALL1Jet() || minE < 1) minE = 1;
2690 
2691  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));
2692 
2693  //
2694  // Loop on the clusters, check if there is any that falls into one of the patches
2695  //
2696  for (Int_t iclus = 0; iclus < nclusters; iclus++)
2697  {
2698  AliVCluster * clus = 0;
2699  if(clusterList) clus = (AliVCluster*) clusterList->At(iclus);
2700  else clus = fInputEvent->GetCaloCluster(iclus);
2701 
2702  if ( !clus ) continue ;
2703 
2704  if ( !clus->IsEMCAL() ) continue ;
2705 
2706  //Skip clusters with too low energy to be triggering
2707  if ( clus->E() < minE ) continue ;
2708 
2709  Float_t frac = -1;
2710  Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(fInputEvent->GetEMCALCells(), clus,frac);
2711 
2712  Bool_t badCluster = GetCaloUtils()->GetEMCALRecoUtils()->ClusterContainsBadChannel(GetCaloUtils()->GetEMCALGeometry(),
2713  clus->GetCellsAbsId(),clus->GetNCells());
2714  UShort_t cellMax[] = {(UShort_t) absIdMax};
2715  Bool_t badCell = GetCaloUtils()->GetEMCALRecoUtils()->ClusterContainsBadChannel(GetCaloUtils()->GetEMCALGeometry(),cellMax,1);
2716 
2717  // if cell is bad, it can happen that time calibration is not available,
2718  // when calculating if it is exotic, this can make it to be exotic by default
2719  // open it temporarily for this cluster
2720  if(badCell)
2721  GetCaloUtils()->GetEMCALRecoUtils()->SetExoticCellDiffTimeCut(10000000);
2722 
2723  Bool_t exotic = GetCaloUtils()->GetEMCALRecoUtils()->IsExoticCluster(clus, fInputEvent->GetEMCALCells());
2724 
2725  // Set back the time cut on exotics
2726  if(badCell)
2727  GetCaloUtils()->GetEMCALRecoUtils()->SetExoticCellDiffTimeCut(exoDiffTime);
2728 
2729  // Energy threshold for exotic Ecross typically at 4 GeV,
2730  // for lower energy, check that there are more than 1 cell in the cluster
2731  if(!exotic && clus->GetNCells() < 2) exotic = kTRUE;
2732 
2733  Float_t energy = clus->E();
2734  Int_t idclus = clus->GetID();
2735 
2736  Double_t tof = clus->GetTOF();
2737  if(GetCaloUtils()->GetEMCALRecoUtils()->IsTimeRecalibrationOn() && fTriggerClusterTimeRecal){
2738  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absIdMax,fInputEvent->GetBunchCrossNumber(),tof);
2739  //additional L1 phase shift
2740  if(GetCaloUtils()->GetEMCALRecoUtils()->IsL1PhaseInTimeRecalibrationOn()){
2741  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTimeL1Phase(GetCaloUtils()->GetEMCALGeometry()->GetSuperModuleNumber(absIdMax), fInputEvent->GetBunchCrossNumber(), tof);
2742  }
2743  }
2744  tof *=1.e9;
2745 
2746  //printf("cluster %d, ID %d, E %2.2f, tof %2.2f, AbsId max %d, exotic %d, bad Cluster %d, bad Cell %d\n",
2747  // iclus,idclus, energy,tof,absIdMax, exotic, badCluster,badCell);
2748 
2749  // Find the highest energy cluster, avobe trigger threshold
2750  // in the event in case no match to trigger is found
2751  if( energy > eMax )
2752  {
2753  tofMax = tof;
2754  eMax = energy;
2755  badClMax = badCluster;
2756  badCeMax = badCell;
2757  exoMax = exotic;
2758  clusMax = iclus;
2759  idclusMax = idclus;
2760  absIdMaxMax = absIdMax;
2761  }
2762 
2763  // count the good clusters in the event avobe the trigger threshold
2764  // to check the exotic events
2765  if(!badCluster && !exotic)
2766  nOfHighECl++;
2767 
2768  // Find match to trigger
2769  if(fTriggerPatchClusterMatch && nPatch>0)
2770  {
2771  for(Int_t iabsId =0; iabsId < nPatch; iabsId++)
2772  {
2773  Int_t absIDCell[4];
2774  GetCaloUtils()->GetEMCALGeometry()->GetCellIndexFromFastORIndex(patches.At(iabsId), absIDCell);
2775  //if(tof > 75 ) printf("E %2.2f TOF %2.2f Trigger patch %d, cells : %d, %d, %d, %d\n",
2776  // clus->E(),tof,patches.At(iabsId), absIDCell[0],absIDCell[1],absIDCell[2],absIDCell[3]);
2777 
2778  for(Int_t ipatch = 0; ipatch < 4; ipatch++)
2779  {
2780  if(absIdMax == absIDCell[ipatch])
2781  {
2782  //printf("*** Patches : absId %d, E %2.1f, tof %f \n",absIdMax,clus->E(), tof);
2783  if(energy > ePatchMax)
2784  {
2785  tofPatchMax = tof;
2786  ePatchMax = energy;
2787  fIsBadCellEvent = badCluster;
2788  fIsBadMaxCellEvent = badCell;
2789  fIsExoticEvent = exotic;
2790  fTriggerClusterIndex = iclus;
2791  fTriggerClusterId = idclus;
2792  fIsTriggerMatch = kTRUE;
2793 // absIdMaxTrig = absIdMax;
2794  }
2795  }
2796  }// cell patch loop
2797  }// trigger patch loop
2798  } // Do trigger patch matching
2799 
2800  }// Cluster loop
2801 
2802  // If there was no match, assign as trigger
2803  // the highest energy cluster in the event
2804  if(!fIsTriggerMatch)
2805  {
2806  tofPatchMax = tofMax;
2807  ePatchMax = eMax;
2808  fIsBadCellEvent = badClMax;
2809  fIsBadMaxCellEvent = badCeMax;
2810  fIsExoticEvent = exoMax;
2811  fTriggerClusterIndex = clusMax;
2812  fTriggerClusterId = idclusMax;
2813  }
2814 
2815  Double_t tofPatchMaxUS = TMath::Abs(tofPatchMax);
2816 
2817  if (tofPatchMaxUS < 28 ) fTriggerClusterBC = 0 ;
2818  else if(tofPatchMaxUS < 75 ) fTriggerClusterBC = 1 ;
2819  else if(tofPatchMaxUS < 125) fTriggerClusterBC = 2 ;
2820  else if(tofPatchMaxUS < 175) fTriggerClusterBC = 3 ;
2821  else if(tofPatchMaxUS < 225) fTriggerClusterBC = 4 ;
2822  else if(tofPatchMaxUS < 275) fTriggerClusterBC = 5 ;
2823  else
2824  {
2825  //printf("AliCaloTrackReader::MatchTriggerCluster() - Large BC - tof %2.3f - Index %d\n",tofPatchMaxUS,fTriggerClusterIndex);
2827  else
2828  {
2829  fTriggerClusterIndex = -2;
2830  fTriggerClusterId = -2;
2831  }
2832  }
2833 
2834  if(tofPatchMax < 0) fTriggerClusterBC*=-1;
2835 
2836 
2837  // 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",
2838  // fTriggerClusterIndex, fTriggerClusterId,ePatchMax, tofPatchMax,
2839  // fTriggerClusterBC, fIsBadCellEvent,fIsBadMaxCellEvent,fIsExoticEvent, fIsTriggerMatch, nOfHighECl,absIdMaxMax);
2840  //
2841  // 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",
2842  // clusMax, idclusMax, eMax,tofMax, badClMax, badCeMax,exoMax);
2843 
2844  //Redo matching but open cuts
2845  if(!fIsTriggerMatch && fTriggerClusterId >= 0)
2846  {
2847  // Open time patch time
2848  TArrayI patchOpen = GetTriggerPatches(7,10);
2849 
2850  Int_t patchAbsIdOpenTime = -1;
2851  for(Int_t iabsId =0; iabsId < patchOpen.GetSize(); iabsId++)
2852  {
2853  Int_t absIDCell[4];
2854  patchAbsIdOpenTime = patchOpen.At(iabsId);
2855  GetCaloUtils()->GetEMCALGeometry()->GetCellIndexFromFastORIndex(patchAbsIdOpenTime, absIDCell);
2856  //if(tof > 75 ) printf("E %2.2f TOF %2.2f Trigger patch %d, cells : %d, %d, %d, %d\n",
2857  // clus->E(),tof,patches.At(iabsId), absIDCell[0],absIDCell[1],absIDCell[2],absIDCell[3]);
2858 
2859  for(Int_t ipatch = 0; ipatch < 4; ipatch++)
2860  {
2861  if(absIdMaxMax == absIDCell[ipatch])
2862  {
2863  fIsTriggerMatchOpenCut[0] = kTRUE;
2864  break;
2865  }
2866  }// cell patch loop
2867  }// trigger patch loop
2868 
2869  // Check neighbour patches
2870  Int_t patchAbsId = -1;
2871  Int_t globalCol = -1;
2872  Int_t globalRow = -1;
2873  GetCaloUtils()->GetEMCALGeometry()->GetFastORIndexFromCellIndex(absIdMaxMax, patchAbsId);
2874  GetCaloUtils()->GetEMCALGeometry()->GetPositionInEMCALFromAbsFastORIndex(patchAbsId,globalCol,globalRow);
2875 
2876  // Check patches with strict time cut
2877  Int_t patchAbsIdNeigh = -1;
2878  for(Int_t icol = globalCol-1; icol <= globalCol+1; icol++)
2879  {
2880  if(icol < 0 || icol > 47) continue;
2881 
2882  for(Int_t irow = globalRow; irow <= globalRow+1; irow++)
2883  {
2884  if(irow < 0 || irow > 63) continue;
2885 
2886  GetCaloUtils()->GetEMCALGeometry()->GetAbsFastORIndexFromPositionInEMCAL(icol, irow, patchAbsIdNeigh);
2887 
2888  if ( patchAbsIdNeigh < 0 ) continue;
2889 
2890  for(Int_t iabsId =0; iabsId < patches.GetSize(); iabsId++)
2891  {
2892  if(patchAbsIdNeigh == patches.At(iabsId))
2893  {
2894  fIsTriggerMatchOpenCut[1] = kTRUE;
2895  break;
2896  }
2897  }// trigger patch loop
2898 
2899  }// row
2900  }// col
2901 
2902  // Check patches with open time cut
2903  Int_t patchAbsIdNeighOpenTime = -1;
2904  for(Int_t icol = globalCol-1; icol <= globalCol+1; icol++)
2905  {
2906  if(icol < 0 || icol > 47) continue;
2907 
2908  for(Int_t irow = globalRow; irow <= globalRow+1; irow++)
2909  {
2910  if(irow < 0 || irow > 63) continue;
2911 
2912  GetCaloUtils()->GetEMCALGeometry()->GetAbsFastORIndexFromPositionInEMCAL(icol, irow, patchAbsIdNeighOpenTime);
2913 
2914  if ( patchAbsIdNeighOpenTime < 0 ) continue;
2915 
2916  for(Int_t iabsId =0; iabsId < patchOpen.GetSize(); iabsId++)
2917  {
2918  if(patchAbsIdNeighOpenTime == patchOpen.At(iabsId))
2919  {
2920  fIsTriggerMatchOpenCut[2] = kTRUE;
2921  break;
2922  }
2923  }// trigger patch loop
2924 
2925  }// row
2926  }// col
2927 
2928  // printf("No match, new match: Open time %d-%d, open Neigh %d-%d, both open %d-%d\n",fIsTriggerMatchOpenCut[0],patchAbsIdOpenTime,
2929  // fIsTriggerMatchOpenCut[1],patchAbsIdNeigh,
2930  // fIsTriggerMatchOpenCut[2],patchAbsIdNeighOpenTime);
2931 
2932  patchOpen.Reset();
2933 
2934  }// No trigger match found
2935  //printf("Trigger BC %d, Id %d, Index %d\n",fTriggerClusterBC,fTriggerClusterId,fTriggerClusterIndex);
2936 }
2937 
2938 //_________________________________________________________
2942 //_________________________________________________________
2944 {
2946  {
2947  // get object pointer
2948  AliVCaloTrigger *caloTrigger = GetInputEvent()->GetCaloTrigger( "EMCAL" );
2949 
2950  if ( fBitEGA == 6 )
2951  {
2952  if (IsEventEMCALL1Gamma1()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(1);
2953  else if(IsEventEMCALL1Gamma2()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(3);
2954  else if(IsEventEMCALL1Jet1 ()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(0);
2955  else if(IsEventEMCALL1Jet2 ()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(2);
2956 
2957  // printf("L1 trigger Threshold Jet1 %f, Gamma1 %f, Jet2 %f, Gamma2 %f\n",
2958  // 0.07874*caloTrigger->GetL1Threshold(0),
2959  // 0.07874*caloTrigger->GetL1Threshold(1),
2960  // 0.07874*caloTrigger->GetL1Threshold(2),
2961  // 0.07874*caloTrigger->GetL1Threshold(3));
2962  }
2963  else
2964  {
2965  // Old AOD data format, in such case, order of thresholds different!!!
2966  if (IsEventEMCALL1Gamma1()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(0);
2967  else if(IsEventEMCALL1Gamma2()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(2);
2968  else if(IsEventEMCALL1Jet1 ()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(1);
2969  else if(IsEventEMCALL1Jet2 ()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(3);
2970 
2971  // printf("L1 trigger Threshold Jet1 %f, Gamma1 %f, Jet2 %f, Gamma2 %f\n",
2972  // 0.07874*caloTrigger->GetL1Threshold(1),
2973  // 0.07874*caloTrigger->GetL1Threshold(0),
2974  // 0.07874*caloTrigger->GetL1Threshold(3),
2975  // 0.07874*caloTrigger->GetL1Threshold(2));
2976  }
2977  }
2978 
2979  // Set L0 threshold, if not set by user
2981  {
2982  // Revise for periods > LHC11d
2983  Int_t runNumber = fInputEvent->GetRunNumber();
2984  if (runNumber < 146861) fTriggerL0EventThreshold = 3. ;
2985  else if(runNumber < 154000) fTriggerL0EventThreshold = 4. ;
2986  else if(runNumber < 165000) fTriggerL0EventThreshold = 5.5;
2987  else fTriggerL0EventThreshold = 2 ;
2988  }
2989 }
2990 
2991 //________________________________________________________
2993 //________________________________________________________
2994 void AliCaloTrackReader::Print(const Option_t * opt) const
2995 {
2996  if(! opt)
2997  return;
2998 
2999  printf("***** Print: %s %s ******\n", GetName(), GetTitle() ) ;
3000  printf("Task name : %s\n", fTaskName.Data()) ;
3001  printf("Data type : %d\n", fDataType) ;
3002  printf("CTS Min pT : %2.1f GeV/c\n", fCTSPtMin) ;
3003  printf("EMCAL Min pT : %2.1f GeV/c\n", fEMCALPtMin) ;
3004  printf("PHOS Min pT : %2.1f GeV/c\n", fPHOSPtMin) ;
3005  printf("CTS Max pT : %2.1f GeV/c\n", fCTSPtMax) ;
3006  printf("EMCAL Max pT : %2.1f GeV/c\n", fEMCALPtMax) ;
3007  printf("PHOS Max pT : %2.1f GeV/c\n", fPHOSPtMax) ;
3008  printf("EMCAL Bad Dist > %2.1f \n" , fEMCALBadChMinDist) ;
3009  printf("PHOS Bad Dist > %2.1f \n" , fPHOSBadChMinDist) ;
3010  printf("EMCAL N cells > %d \n" , fEMCALNCellsCut) ;
3011  printf("PHOS N cells > %d \n" , fPHOSNCellsCut) ;
3012  printf("EMCAL Time Cut: %3.1f < TOF < %3.1f\n", fEMCALTimeCutMin, fEMCALTimeCutMax);
3013  printf("Use CTS = %d\n", fFillCTS) ;
3014  printf("Use EMCAL = %d\n", fFillEMCAL) ;
3015  printf("Use DCAL = %d\n", fFillDCAL) ;
3016  printf("Use PHOS = %d\n", fFillPHOS) ;
3017  printf("Use EMCAL Cells = %d\n", fFillEMCALCells) ;
3018  printf("Use PHOS Cells = %d\n", fFillPHOSCells) ;
3019  printf("Track status = %d\n", (Int_t) fTrackStatus) ;
3020 
3021  printf("Track Mult Eta Cut = %2.2f\n", fTrackMultEtaCut) ;
3022  printf("Write delta AOD = %d\n", fWriteOutputDeltaAOD) ;
3023  printf("Recalculate Clusters = %d, E linearity = %d\n", fRecalculateClusters, fCorrectELinearity) ;
3024 
3025  printf("Use Triggers selected in SE base class %d; If not what Trigger Mask? %d; MB Trigger Mask for mixed %d \n",
3027 
3029  printf("Compare jet pt and pt hard to accept event, factor = %2.2f",fPtHardAndJetPtFactor);
3030 
3032  printf("Compare cluster pt and pt hard to accept event, factor = %2.2f",fPtHardAndClusterPtFactor);
3033 
3034  printf("Read Kine from, stack? %d, AOD ? %d \n", fReadStack, fReadAODMCParticles) ;
3035  printf("Delta AOD File Name = %s\n", fDeltaAODFileName.Data()) ;
3036  printf("Centrality: Class %s, Option %d, Bin [%d,%d] \n", fCentralityClass.Data(),fCentralityOpt,fCentralityBin[0], fCentralityBin[1]) ;
3037 
3038  printf(" \n") ;
3039 }
3040 
3041 //__________________________________________
3045 //__________________________________________
3047 {
3048  // Count number of cells with energy larger than 0.1 in SM3, cut on this number
3049  Int_t ncellsSM3 = 0;
3050  for(Int_t icell = 0; icell < fInputEvent->GetEMCALCells()->GetNumberOfCells(); icell++)
3051  {
3052  Int_t absID = fInputEvent->GetEMCALCells()->GetCellNumber(icell);
3053  Int_t sm = GetCaloUtils()->GetEMCALGeometry()->GetSuperModuleNumber(absID);
3054  if(fInputEvent->GetEMCALCells()->GetAmplitude(icell) > 0.1 && sm==3) ncellsSM3++;
3055  }
3056 
3057  Int_t ncellcut = 21;
3058  if(GetFiredTriggerClasses().Contains("EMC")) ncellcut = 35;
3059 
3060  if(ncellsSM3 >= ncellcut)
3061  {
3062  AliDebug(1,Form("Reject event with ncells in SM3 %d, cut %d, trig %s",
3063  ncellsSM3,ncellcut,GetFiredTriggerClasses().Data()));
3064  return kTRUE;
3065  }
3066 
3067  return kFALSE;
3068 }
3069 
3070 //_________________________________________________________
3074 //_________________________________________________________
3076 {
3077  if(label < 0) return ;
3078 
3079  AliAODEvent * evt = dynamic_cast<AliAODEvent*> (fInputEvent) ;
3080  if(!evt) return ;
3081 
3082  TClonesArray * arr = dynamic_cast<TClonesArray*>(evt->FindListObject("mcparticles")) ;
3083  if(!arr) return ;
3084 
3085  if(label < arr->GetEntriesFast())
3086  {
3087  AliAODMCParticle * particle = dynamic_cast<AliAODMCParticle *>(arr->At(label));
3088  if(!particle) return ;
3089 
3090  if(label == particle->Label()) return ; // label already OK
3091  //else printf("AliCaloTrackReader::RemapMCLabelForAODs() - Label %d - AOD stack %d \n",label, particle->Label());
3092  }
3093  //else printf("AliCaloTrackReader::RemapMCLabelForAODs() - Label %d > AOD labels %d \n",label, arr->GetEntriesFast());
3094 
3095  // loop on the particles list and check if there is one with the same label
3096  for(Int_t ind = 0; ind < arr->GetEntriesFast(); ind++ )
3097  {
3098  AliAODMCParticle * particle = dynamic_cast<AliAODMCParticle *>(arr->At(ind));
3099  if(!particle) continue ;
3100 
3101  if(label == particle->Label())
3102  {
3103  label = ind;
3104  //printf("AliAnalysisTaskEMCALClusterize::RemapMCLabelForAODs() - New Label Index %d \n",label);
3105  return;
3106  }
3107  }
3108 
3109  label = -1;
3110 
3111  //printf("AliCaloTrackReader::RemapMCLabelForAODs() - Label not found set to -1 \n");
3112 }
3113 
3114 //___________________________________
3116 //___________________________________
3118 {
3119  if(fCTSTracks) fCTSTracks -> Clear();
3120  if(fEMCALClusters) fEMCALClusters -> Clear("C");
3121  if(fPHOSClusters) fPHOSClusters -> Clear("C");
3122 
3123  fV0ADC[0] = 0; fV0ADC[1] = 0;
3124  fV0Mul[0] = 0; fV0Mul[1] = 0;
3125 
3126  if(fNonStandardJets) fNonStandardJets -> Clear("C");
3127  fBackgroundJets->Reset();
3128 }
3129 
3130 //___________________________________________
3134 //___________________________________________
3136 {
3137  fEventTrigMinBias = kFALSE;
3138  fEventTrigCentral = kFALSE;
3139  fEventTrigSemiCentral = kFALSE;
3140  fEventTrigEMCALL0 = kFALSE;
3141  fEventTrigEMCALL1Gamma1 = kFALSE;
3142  fEventTrigEMCALL1Gamma2 = kFALSE;
3143  fEventTrigEMCALL1Jet1 = kFALSE;
3144  fEventTrigEMCALL1Jet2 = kFALSE;
3145 
3146  AliDebug(1,Form("Select trigger mask bit %d - Trigger Event %s - Select <%s>",
3148 
3149  if(fEventTriggerMask <=0 )// in case no mask set
3150  {
3151  // EMC triggered event? Which type?
3152  if( GetFiredTriggerClasses().Contains("-B-") || GetFiredTriggerClasses().Contains("-S-") || GetFiredTriggerClasses().Contains("-I-") )
3153  {
3154  if ( GetFiredTriggerClasses().Contains("EGA" ) ||
3155  GetFiredTriggerClasses().Contains("EG1" ) )
3156  {
3157  fEventTrigEMCALL1Gamma1 = kTRUE;
3158  if( GetFiredTriggerClasses().Contains("EG1" ) && !fFiredTriggerClassName.Contains("EG1") ) fEventTrigEMCALL1Gamma1 = kFALSE;
3159  }
3160  else if( GetFiredTriggerClasses().Contains("EG2" ) )
3161  {
3162  fEventTrigEMCALL1Gamma2 = kTRUE;
3163  if( !fFiredTriggerClassName.Contains("EG2") ) fEventTrigEMCALL1Gamma2 = kFALSE;
3164  }
3165  else if( GetFiredTriggerClasses().Contains("EJE" ) ||
3166  GetFiredTriggerClasses().Contains("EJ1" ) )
3167  {
3168  fEventTrigEMCALL1Jet1 = kTRUE;
3169  if( GetFiredTriggerClasses().Contains("EJ1" ) && !fFiredTriggerClassName.Contains("EJ1") )
3170  fEventTrigEMCALL1Jet1 = kFALSE;
3171  }
3172  else if( GetFiredTriggerClasses().Contains("EJ2" ) )
3173  {
3174  fEventTrigEMCALL1Jet2 = kTRUE;
3175  if( !fFiredTriggerClassName.Contains("EJ2") ) fEventTrigEMCALL1Jet2 = kFALSE;
3176  }
3177  else if( GetFiredTriggerClasses().Contains("CEMC") &&
3178  !GetFiredTriggerClasses().Contains("EGA" ) &&
3179  !GetFiredTriggerClasses().Contains("EJE" ) &&
3180  !GetFiredTriggerClasses().Contains("EG1" ) &&
3181  !GetFiredTriggerClasses().Contains("EJ1" ) &&
3182  !GetFiredTriggerClasses().Contains("EG2" ) &&
3183  !GetFiredTriggerClasses().Contains("EJ2" ) )
3184  {
3185  fEventTrigEMCALL0 = kTRUE;
3186  }
3187 
3188  //Min bias event trigger?
3189  if (GetFiredTriggerClasses().Contains("CCENT_R2-B-NOPF-ALLNOTRD"))
3190  {
3191  fEventTrigCentral = kTRUE;
3192  }
3193  else if(GetFiredTriggerClasses().Contains("CSEMI_R1-B-NOPF-ALLNOTRD"))
3194  {
3195  fEventTrigSemiCentral = kTRUE;
3196  }
3197  else if((GetFiredTriggerClasses().Contains("CINT") || GetFiredTriggerClasses().Contains("CPBI2_B1") ) &&
3198  GetFiredTriggerClasses().Contains("-NOPF-ALLNOTRD") )
3199  {
3200  fEventTrigMinBias = kTRUE;
3201  }
3202  }
3203  }
3204  else
3205  {
3206  // EMC L1 Gamma
3207  if ( fEventTriggerMask & AliVEvent::kEMCEGA )
3208  {
3209  //printf("EGA trigger bit\n");
3210  if (GetFiredTriggerClasses().Contains("EG"))
3211  {
3212  if (GetFiredTriggerClasses().Contains("EGA")) fEventTrigEMCALL1Gamma1 = kTRUE;
3213  else
3214  {
3215  if(GetFiredTriggerClasses().Contains("EG1")) fEventTrigEMCALL1Gamma1 = kTRUE;
3216  if(GetFiredTriggerClasses().Contains("EG2")) fEventTrigEMCALL1Gamma2 = kTRUE;
3217  }
3218  }
3219  }
3220  // EMC L1 Jet
3221  else if( fEventTriggerMask & AliVEvent::kEMCEJE )
3222  {
3223  //printf("EGA trigger bit\n");
3224  if (GetFiredTriggerClasses().Contains("EJ"))
3225  {
3226  if (GetFiredTriggerClasses().Contains("EJE")) fEventTrigEMCALL1Jet1 = kTRUE;
3227  else
3228  {
3229  if(GetFiredTriggerClasses().Contains("EJ1")) fEventTrigEMCALL1Jet1 = kTRUE;
3230  if(GetFiredTriggerClasses().Contains("EJ2")) fEventTrigEMCALL1Jet2 = kTRUE;
3231  }
3232  }
3233  }
3234  // EMC L0
3235  else if((fEventTriggerMask & AliVEvent::kEMC7) ||
3236  (fEventTriggerMask & AliVEvent::kEMC1) )
3237  {
3238  //printf("L0 trigger bit\n");
3239  fEventTrigEMCALL0 = kTRUE;
3240  }
3241  // Min Bias Pb-Pb
3242  else if( fEventTriggerMask & AliVEvent::kCentral )
3243  {
3244  //printf("MB semi central trigger bit\n");
3245  fEventTrigSemiCentral = kTRUE;
3246  }
3247  // Min Bias Pb-Pb
3248  else if( fEventTriggerMask & AliVEvent::kSemiCentral )
3249  {
3250  //printf("MB central trigger bit\n");
3251  fEventTrigCentral = kTRUE;
3252  }
3253  // Min Bias pp, PbPb, pPb
3254  else if((fEventTriggerMask & AliVEvent::kMB ) ||
3255  (fEventTriggerMask & AliVEvent::kINT7) ||
3256  (fEventTriggerMask & AliVEvent::kINT8) ||
3257  (fEventTriggerMask & AliVEvent::kAnyINT) )
3258  {
3259  //printf("MB trigger bit\n");
3260  fEventTrigMinBias = kTRUE;
3261  }
3262  }
3263 
3264  AliDebug(1,Form("Event bits: \n \t MB %d, Cen %d, Sem %d, L0 %d, L1G1 %d, L1G2 %d, L1J1 %d, L1J2 %d",
3268 
3269  // L1 trigger bit
3270  if( fBitEGA == 0 && fBitEJE == 0 )
3271  {
3272  // Init the trigger bit once, correct depending on AliESD(AOD)CaloTrigger header version
3273 
3274  // Simpler way to do it ...
3275 // if( fInputEvent->GetRunNumber() < 172000 )
3276 // reader->SetEventTriggerL1Bit(4,5); // current LHC11 data
3277 // else
3278 // reader->SetEventTriggerL1Bit(6,8); // LHC12-13 data
3279 
3280  // Old values
3281  fBitEGA = 4;
3282  fBitEJE = 5;
3283 
3284  TFile* file = AliAnalysisManager::GetAnalysisManager()->GetTree()->GetCurrentFile();
3285 
3286  const TList *clist = file->GetStreamerInfoCache();
3287 
3288  if(clist)
3289  {
3290  TStreamerInfo *cinfo = (TStreamerInfo*)clist->FindObject("AliESDCaloTrigger");
3291  Int_t verid = 5; // newer ESD header version
3292  if(!cinfo)
3293  {
3294  cinfo = (TStreamerInfo*)clist->FindObject("AliAODCaloTrigger");
3295  verid = 3; // newer AOD header version
3296  }
3297 
3298  if(cinfo)
3299  {
3300  Int_t classversionid = cinfo->GetClassVersion();
3301  //printf("********* Header class version %d *********** \n",classversionid);
3302 
3303  if (classversionid >= verid)
3304  {
3305  fBitEGA = 6;
3306  fBitEJE = 8;
3307  }
3308  } else AliInfo("AliCaloTrackReader()::SetEventTriggerBit() - Streamer info for trigger class not available, bit not changed");
3309  } else AliInfo("AliCaloTrackReader::SetEventTriggerBit() - Streamer list not available!, bit not changed");
3310 
3311  } // set once the EJE, EGA trigger bit
3312 }
3313 
3314 //____________________________________________________________
3317 //____________________________________________________________
3318 void AliCaloTrackReader::SetInputEvent(AliVEvent* const input)
3319 {
3320  fInputEvent = input;
3321  fMixedEvent = dynamic_cast<AliMixedEvent*>(GetInputEvent()) ;
3322  if (fMixedEvent)
3323  fNMixedEvent = fMixedEvent->GetNumberOfEvents() ;
3324 
3325  //Delete previous vertex
3326  if(fVertex)
3327  {
3328  for (Int_t i = 0; i < fNMixedEvent; i++)
3329  {
3330  delete [] fVertex[i] ;
3331  }
3332  delete [] fVertex ;
3333  }
3334 
3335  fVertex = new Double_t*[fNMixedEvent] ;
3336  for (Int_t i = 0; i < fNMixedEvent; i++)
3337  {
3338  fVertex[i] = new Double_t[3] ;
3339  fVertex[i][0] = 0.0 ;
3340  fVertex[i][1] = 0.0 ;
3341  fVertex[i][2] = 0.0 ;
3342  }
3343 }
3344 
3345 
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.
double Double_t
Definition: External.C:58
virtual void FillInputNonStandardJets()
Double_t fEMCALTimeCutMax
Remove clusters/cells with time larger than this value, in ns.
Int_t fBitEJE
Trigger bit on VCaloTrigger for EJE.
TLorentzVector fMomentum
! Temporal TLorentzVector container, avoid declaration of TLorentzVectors per event.
Double_t fTimeStampRunMin
Minimum value of time stamp in run.
Bool_t fDoPileUpEventRejection
Select pile-up events by SPD.
void SetCentrality(Float_t cen)
Definition: AliAnaWeights.h: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
Int_t fSmearNLMMin
Do smearing for clusters with at least this value.
TString fEventPlaneMethod
Name of event plane method, by default "Q".
Int_t fTrackBCEventCut[19]
Fill one entry per event if there is a track in a given BC, depend on track pT, acceptance cut...
Bool_t fDoVertexBCEventSelection
Select events with vertex on BC=0 or -100.
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.
int Int_t
Definition: External.C:63
void RecalculateClusterShowerShapeParameters(AliVCaloCells *cells, AliVCluster *clu)
Bool_t fSelectSPDHitTracks
Ensure that track hits SPD layers.
TString fFiredTriggerClassName
Name of trigger event type used to do the analysis.
virtual TClonesArray * GetAODMCParticles() const
Definition: External.C:204
TString GetFiredTriggerClasses() const
Bool_t ReadStack() const
unsigned int UInt_t
Definition: External.C:33
TObjArray * fEMCALClusters
Temporal array with EMCAL CaloClusters.
Bool_t IsEventEMCALL1Jet() const
Float_t fPHOSPtMax
pT Threshold on phos clusters.
Bool_t fAccessTrackTOF
Access the track TOF, in case of problems when accessing GetTOFBunchCrossing.
float Float_t
Definition: External.C:68
virtual Bool_t CheckForPrimaryVertex() const
Int_t fEMCalBCEvent[19]
Fill one entry per event if there is a cluster in a given BC.
TString fInputBackgroundJetBranchName
Name of background jet branch.
virtual AliEventplane * GetEventPlane() const
Bool_t fEventTriggerAtSE
Select triggered event at SE base task or here.
Float_t fEMCALPtMax
pT Threshold on emcal clusters.
Int_t fEventType
Set the event species: 7 physics, 0 MC, 8 LED (not useful now)
virtual Bool_t ComparePtHardAndJetPt()
Bool_t IsInFiducialCut(Float_t eta, Float_t phi, Int_t det) const
Bool_t fEventTrigCentral
Event is AliVEvent::kCentral on its name, it should correspond to PbPb.
virtual TObjString * GetListOfParameters()
Save parameters used for analysis in a string.
Float_t fTrackMultEtaCut
Track multiplicity eta cut.
Int_t fCentralityBin[2]
Minimum and maximum value of the centrality for the analysis.
virtual void FillInputPHOSCells()
Connects the array with PHOS cells and the pointer.
Base class for event, clusters and tracks filtering and preparation for the analysis.
TObjArray * fDCALClusters
Temporal array with DCAL CaloClusters, not needed in the normal case, use just EMCal array with DCal ...
Bool_t fUseTrackTimeCut
Do time cut selection.
virtual TString GetEventPlaneMethod() const
virtual Int_t GetTrackID(AliVTrack *track)
TArrayI fAcceptEventsWithBit
Accept events if trigger bit is on.
unsigned long ULong_t
Definition: External.C:38
TH1F * fhEMCALClusterCutsE[8]
! Control histogram on the different EMCal cluster selection cuts, E
virtual void Print(const Option_t *opt) const
Print some relevant parameters set for the analysis.
void MatchTriggerCluster(TArrayI patches)
Bool_t fTimeStampEventSelect
Select events within a fraction of data taking time.
Double_t fTimeStampRunMax
Maximum value of time stamp in run.
Double_t fEMCALTimeCutMin
Remove clusters/cells with time smaller than this value, in ns.
Bool_t fEventTrigEMCALL1Jet2
Event is L1-Gamma, threshold 2 on its name, it should correspond kEMCEGA.
void SetEMCalEventBCcut(Int_t bc)
Int_t fNPileUpClustersCut
Cut to select event as pile-up.
Bool_t fCheckFidCut
Do analysis for clusters in defined region.
virtual ~AliCaloTrackReader()
Destructor.
Bool_t fEventTrigEMCALL0
Event is EMCal L0 on its name, it should correspond to AliVEvent::kEMC7, AliVEvent::kEMC1.
Int_t GetNumberOfLocalMaxima(AliVCluster *cluster, AliVCaloCells *cells)
Find the number of local maxima in cluster.
Bool_t IsWeightSettingOn() const
Definition: AliAnaWeights.h: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
Int_t fSmearNLMMax
Do smearing for clusters with at maximum this value.
Float_t fTriggerL0EventThreshold
L0 Threshold to look for triggered events, set outside.
Float_t fTimeStampEventFracMax
Maximum value of time stamp fraction event.
virtual void ResetLists()
Reset lists, called in AliAnaCaloTrackCorrMaker.
Int_t fNPileUpClusters
Number of clusters with time avobe 20 ns.
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.
unsigned short UShort_t
Definition: External.C:28
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.
const char Option_t
Definition: External.C:48
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.
bool Bool_t
Definition: External.C:53
Float_t fSmearShowerShapeWidth
Smear shower shape landau function "width" (use in MC).
Int_t fBitEGA
Trigger bit on VCaloTrigger for EGA.
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