AliPhysics  0e0bd91 (0e0bd91)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros
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 fUseEMCALTimeCut(1), fUseParamTimeCut(0),
71 fUseTrackTimeCut(0), fAccessTrackTOF(0),
72 fEMCALTimeCutMin(-10000), fEMCALTimeCutMax(10000),
73 fEMCALParamTimeCutMin(), fEMCALParamTimeCutMax(),
74 fTrackTimeCutMin(-10000), fTrackTimeCutMax(10000),
75 fUseTrackDCACut(0),
76 fAODBranchList(0x0),
77 fCTSTracks(0x0), fEMCALClusters(0x0),
78 fDCALClusters(0x0), fPHOSClusters(0x0),
79 fEMCALCells(0x0), fPHOSCells(0x0),
80 fInputEvent(0x0), fOutputEvent(0x0),fMC(0x0),
81 fFillCTS(0), fFillEMCAL(0),
82 fFillDCAL(0), fFillPHOS(0),
83 fFillEMCALCells(0), fFillPHOSCells(0),
84 fRecalculateClusters(kFALSE),fCorrectELinearity(kTRUE),
85 fSelectEmbeddedClusters(kFALSE),
86 fSmearShowerShape(0), fSmearShowerShapeWidth(0), fRandom(),
87 fTrackStatus(0), fSelectSPDHitTracks(0),
88 fTrackMult(0), fTrackMultEtaCut(0.9),
89 fReadStack(kFALSE), fReadAODMCParticles(kFALSE),
90 fDeltaAODFileName(""), fFiredTriggerClassName(""),
91 
92 fEventTriggerMask(0), fMixEventTriggerMask(0), fEventTriggerAtSE(0),
93 fEventTrigMinBias(0), fEventTrigCentral(0),
94 fEventTrigSemiCentral(0), fEventTrigEMCALL0(0),
95 fEventTrigEMCALL1Gamma1(0), fEventTrigEMCALL1Gamma2(0),
96 fEventTrigEMCALL1Jet1(0), fEventTrigEMCALL1Jet2(0),
97 fBitEGA(0), fBitEJE(0),
98 
99 fEventType(-1),
100 fTaskName(""), fCaloUtils(0x0),
101 fWeightUtils(0x0), fEventWeight(1),
102 fMixedEvent(NULL), fNMixedEvent(0), fVertex(NULL),
103 fListMixedTracksEvents(), fListMixedCaloEvents(),
104 fLastMixedTracksEvent(-1), fLastMixedCaloEvent(-1),
105 fWriteOutputDeltaAOD(kFALSE),
106 fEMCALClustersListName(""), fZvtxCut(0.),
107 fAcceptFastCluster(kFALSE), fRemoveLEDEvents(kFALSE),
108 //Trigger rejection
109 fRemoveBadTriggerEvents(0), fTriggerPatchClusterMatch(1),
110 fTriggerPatchTimeWindow(), fTriggerL0EventThreshold(0),
111 fTriggerL1EventThreshold(0), fTriggerL1EventThresholdFix(0),
112 fTriggerClusterBC(0), fTriggerClusterIndex(0), fTriggerClusterId(0),
113 fIsExoticEvent(0), fIsBadCellEvent(0), fIsBadMaxCellEvent(0),
114 fIsTriggerMatch(0), fIsTriggerMatchOpenCut(),
115 fTriggerClusterTimeRecal(kTRUE), fRemoveUnMatchedTriggers(kTRUE),
116 fDoPileUpEventRejection(kFALSE), fDoV0ANDEventSelection(kFALSE),
117 fDoVertexBCEventSelection(kFALSE),
118 fDoRejectNoTrackEvents(kFALSE),
119 fUseEventsWithPrimaryVertex(kFALSE),
120 //fTriggerAnalysis (0x0),
121 fTimeStampEventSelect(0),
122 fTimeStampEventFracMin(0), fTimeStampEventFracMax(0),
123 fTimeStampRunMin(0), fTimeStampRunMax(0),
124 fNPileUpClusters(-1), fNNonPileUpClusters(-1), fNPileUpClustersCut(3),
125 fVertexBC(-200), fRecalculateVertexBC(0),
126 fUseAliCentrality(0), fCentralityClass(""), fCentralityOpt(0),
127 fEventPlaneMethod(""),
128 fAcceptOnlyHIJINGLabels(0), fNMCProducedMin(0), fNMCProducedMax(0),
129 fFillInputNonStandardJetBranch(kFALSE),
130 fNonStandardJets(new TClonesArray("AliAODJet",100)),fInputNonStandardJetBranchName("jets"),
131 fFillInputBackgroundJetBranch(kFALSE),
132 fBackgroundJets(0x0),fInputBackgroundJetBranchName("jets"),
133 fAcceptEventsWithBit(0), fRejectEventsWithBit(0), fRejectEMCalTriggerEventsWith2Tresholds(0),
134 fMomentum()
135 {
136  InitParameters();
137 }
138 
139 //_______________________________________
141 //_______________________________________
143 {
144  DeletePointers();
145 }
146 
147 //_______________________________________
150 //_______________________________________
152 {
153  delete fFiducialCut ;
154 
155  if(fAODBranchList)
156  {
157  fAODBranchList->Delete();
158  delete fAODBranchList ;
159  }
160 
161  if(fCTSTracks)
162  {
163  if(fDataType!=kMC)fCTSTracks->Clear() ;
164  else fCTSTracks->Delete() ;
165  delete fCTSTracks ;
166  }
167 
168  if(fEMCALClusters)
169  {
170  if(fDataType!=kMC)fEMCALClusters->Clear("C") ;
171  else fEMCALClusters->Delete() ;
172  delete fEMCALClusters ;
173  }
174 
175  if(fDCALClusters)
176  {
177  if(fDataType!=kMC)fDCALClusters->Clear("C") ;
178  else fDCALClusters->Delete() ;
179  delete fDCALClusters ;
180  }
181 
182  if(fPHOSClusters)
183  {
184  if(fDataType!=kMC)fPHOSClusters->Clear("C") ;
185  else fPHOSClusters->Delete() ;
186  delete fPHOSClusters ;
187  }
188 
189  if(fVertex)
190  {
191  for (Int_t i = 0; i < fNMixedEvent; i++)
192  {
193  delete [] fVertex[i] ;
194 
195  }
196  delete [] fVertex ;
197  }
198 
199  //delete fTriggerAnalysis;
200 
201  if(fNonStandardJets)
202  {
203  if(fDataType!=kMC) fNonStandardJets->Clear("C") ;
204  else fNonStandardJets->Delete() ;
205  delete fNonStandardJets ;
206  }
207  delete fBackgroundJets ;
208 
209  fRejectEventsWithBit.Reset();
210  fAcceptEventsWithBit.Reset();
211 
212  if ( fWeightUtils ) delete fWeightUtils ;
213 
214  // Pointers not owned, done by the analysis frame
215  // if(fInputEvent) delete fInputEvent ;
216  // if(fOutputEvent) delete fOutputEvent ;
217  // if(fMC) delete fMC ;
218  // Pointer not owned, deleted by maker
219  // if (fCaloUtils) delete fCaloUtils ;
220 }
221 
222 //____________________________________________________________
226 //____________________________________________________________
227 Bool_t AliCaloTrackReader::AcceptDCA(Float_t pt, Float_t dca)
228 {
229  Float_t cut = fTrackDCACut[0]+fTrackDCACut[1]/TMath::Power(pt,fTrackDCACut[2]);
230 
231  if(TMath::Abs(dca) < cut)
232  return kTRUE;
233  else
234  return kFALSE;
235 }
236 
237 //_____________________________________________________
240 //_____________________________________________________
242 {
243  Int_t nAccept = fAcceptEventsWithBit.GetSize();
244 
245  //printf("N accept %d\n", nAccept);
246 
247  if( nAccept <= 0 )
248  return kTRUE ; // accept the event
249 
250  UInt_t trigFired = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
251 
252  for(Int_t ibit = 0; ibit < nAccept; ibit++)
253  {
254  Bool_t accept = (trigFired & fAcceptEventsWithBit.At(ibit));
255 
256  //printf("accept %d, ibit %d, bit %d \n",accept, ibit,fAcceptEventsWithBit.At(ibit));
257  if(accept) return kTRUE ; // accept the event
258  }
259 
260  return kFALSE ; // reject the event
261 }
262 
263 //_____________________________________________________
266 //_____________________________________________________
268 {
269  Int_t nReject = fRejectEventsWithBit.GetSize();
270 
271  //printf("N reject %d\n", nReject);
272 
273  if( nReject <= 0 )
274  return kTRUE ; // accept the event
275 
276  UInt_t trigFired = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
277 
278  for(Int_t ibit = 0; ibit < nReject; ibit++)
279  {
280  Bool_t reject = (trigFired & fRejectEventsWithBit.At(ibit));
281 
282  //printf("reject %d, ibit %d, bit %d \n",reject, ibit,fRejectEventsWithBit.At(ibit));
283  if(reject) return kFALSE ; // reject the event
284  }
285 
286  return kTRUE ; // accept the event
287 }
288 
289 //_____________________________________________
293 //_____________________________________________
295 {
296  //-----------------------------------------------------------
297  // Reject events depending on the trigger name
298  //-----------------------------------------------------------
299 
300  AliDebug(1,Form("FiredTriggerClass <%s>, selected class <%s>, compare name %d",
303 
304  if ( fFiredTriggerClassName != "" )
305  {
307  return kFALSE;
308  else
309  AliDebug(1,"Accepted triggered event");
310  }
311 
312  //-----------------------------------------------------------
313  // Reject events depending on the event species type
314  //-----------------------------------------------------------
315 
316  // Event types:
317  // kStartOfRun = 1, // START_OF_RUN
318  // kEndOfRun = 2, // END_OF_RUN
319  // kStartOfRunFiles = 3, // START_OF_RUN_FILES
320  // kEndOfRunFiles = 4, // END_OF_RUN_FILES
321  // kStartOfBurst = 5, // START_OF_BURST
322  // kEndOfBurst = 6, // END_OF_BURST
323  // kPhysicsEvent = 7, // PHYSICS_EVENT
324  // kCalibrationEvent = 8, // CALIBRATION_EVENT
325  // kFormatError = 9, // EVENT_FORMAT_ERROR
326  // kStartOfData = 10, // START_OF_DATA
327  // kEndOfData = 11, // END_OF_DATA
328  // kSystemSoftwareTriggerEvent = 12, // SYSTEM_SOFTWARE_TRIGGER_EVENT
329  // kDetectorSoftwareTriggerEvent = 13 // DETECTOR_SOFTWARE_TRIGGER_EVENT
330 
331  Int_t eventType = 0;
332  if(fInputEvent->GetHeader())
333  eventType = ((AliVHeader*)fInputEvent->GetHeader())->GetEventType();
334 
335  AliDebug(3,Form("Event type %d",eventType));
336 
337  // Select only Physics events in data, eventType = 7,
338  // usually done by physics selection
339  // MC not set, eventType =0, I think, so do not apply a value to fEventType by default
340  // LED events have eventType = 8, implemented a selection cut in the past
341  // but not really useful for EMCal data since LED are not reconstructed (unless wrongly assumed as physics)
342 
343  if ( fEventType >= 0 && eventType != fEventType ) return kFALSE ;
344 
345  AliDebug(1,"Pass event species selection");
346 
347  //-----------------------------------------------------------------
348  // In case of mixing analysis, select here the trigger of the event
349  //-----------------------------------------------------------------
350 
351  UInt_t isTrigger = kFALSE;
352  UInt_t isMB = kFALSE;
353 
354  if(!fEventTriggerAtSE)
355  {
356  // In case of mixing analysis, accept MB events, not only Trigger
357  // Track and cluster arrays filled for MB in order to create the pool in the corresponding analysis
358  // via de method in the base class FillMixedEventPool()
359 
360  AliAnalysisManager *manager = AliAnalysisManager::GetAnalysisManager();
361  AliInputEventHandler *inputHandler = dynamic_cast<AliInputEventHandler*>(manager->GetInputEventHandler());
362 
363  if(!inputHandler) return kFALSE ; // to content coverity
364 
365  isTrigger = inputHandler->IsEventSelected() & fEventTriggerMask;
366  isMB = inputHandler->IsEventSelected() & fMixEventTriggerMask;
367 
368  if(!isTrigger && !isMB) return kFALSE;
369 
370  //printf("Selected triggered event : %s\n",GetFiredTriggerClasses().Data());
371  AliDebug(0,"Pass uninteresting triggered events rejection in case of mixing analysis");
372  }
373 
374  //-------------------------------------------------------------------------------------
375  // Reject or accept events depending on the trigger bit
376  //-------------------------------------------------------------------------------------
377 
378  Bool_t okA = AcceptEventWithTriggerBit();
379  Bool_t okR = RejectEventWithTriggerBit();
380 
381  //printf("AliCaloTrackReader::FillInputEvent() - Accept event? %d, Reject event %d? \n",okA,okR);
382 
383  if(!okA || !okR) return kFALSE;
384 
385  AliDebug(1,"Pass event bit rejection");
386 
387  //----------------------------------------------------------------------
388  // Do not count events that were likely triggered by an exotic cluster
389  // or out BC cluster
390  //----------------------------------------------------------------------
391 
392  // Set a bit with the event kind, MB, L0, L1 ...
394 
395  // In case of Mixing, avoid checking the triggers in the min bias events
396  if(!fEventTriggerAtSE && (isMB && !isTrigger)) return kTRUE;
397 
399  {
401  {
402  // Reject triggered events when there is coincidence on both EMCal trigger thresholds,
403  // but the requested trigger is the low trigger threshold
404  if(IsEventEMCALL1Jet1 () && IsEventEMCALL1Jet2 () && fFiredTriggerClassName.Contains("EJ2")) return kFALSE;
405  if(IsEventEMCALL1Gamma1() && IsEventEMCALL1Gamma2() && fFiredTriggerClassName.Contains("EG2")) return kFALSE;
406  }
407 
408  //Get Patches that triggered
410 
411  MatchTriggerCluster(patches);
412 
413  patches.Reset();
414 
415  // If requested, remove badly triggeed events, but only when the EMCal trigger bit is set
417  {
418  AliDebug(1,Form("ACCEPT triggered event? \n exotic? %d - bad cell %d - bad Max cell %d - BC %d - Matched %d\n",
420  if (fIsExoticEvent) return kFALSE;
421  else if(fIsBadCellEvent) return kFALSE;
422  else if(fRemoveUnMatchedTriggers && !fIsTriggerMatch) return kFALSE ;
423  else if(fTriggerClusterBC != 0) return kFALSE;
424  AliDebug(1,Form("\t *** YES for %s",GetFiredTriggerClasses().Data()));
425  }
426 
427  AliDebug(1,"Pass EMCal triggered event rejection \n");
428  }
429 
430  //-------------------------------------------------------------------------------------
431  // Select events only fired by a certain trigger configuration if it is provided
432  //-------------------------------------------------------------------------------------
433 
434  if (GetFiredTriggerClasses().Contains("FAST") && !GetFiredTriggerClasses().Contains("ALL") && !fAcceptFastCluster)
435  {
436  AliDebug(1,Form("Do not count events from fast cluster, trigger name %s\n",fFiredTriggerClassName.Data()));
437  return kFALSE;
438  }
439 
440  //-------------------------------------------------------------------------------------
441  // Reject event if large clusters with large energy
442  // Use only for LHC11a data for the moment, and if input is clusterizer V1 or V1+unfolding
443  // If clusterzer NxN or V2 it does not help
444  //-------------------------------------------------------------------------------------
445 
446  //Int_t run = fInputEvent->GetRunNumber();
447 
448  if ( fRemoveLEDEvents )
449  {
450  Bool_t reject = RejectLEDEvents();
451 
452  if(reject) return kFALSE;
453 
454  AliDebug(1,"Pass LED event rejection");
455  } // Remove LED events
456 
457  // All selection criteria passed, accept the event
458  return kTRUE ;
459 }
460 
461 //________________________________________________
466 //________________________________________________
468 {
469  //printf("AliCaloTrackReader::ComparePtHardAndJetPt() - GenHeaderName : %s\n",GetGenEventHeader()->ClassName());
470 
471  if(!strcmp(GetGenEventHeader()->ClassName(), "AliGenPythiaEventHeader"))
472  {
473  TParticle * jet = 0;
474  AliGenPythiaEventHeader* pygeh= (AliGenPythiaEventHeader*) GetGenEventHeader();
475  Int_t nTriggerJets = pygeh->NTriggerJets();
476  Float_t ptHard = pygeh->GetPtHard();
477 
478  AliDebug(1,Form("Njets: %d, pT Hard %f",nTriggerJets, ptHard));
479 
480  Float_t tmpjet[]={0,0,0,0};
481  for(Int_t ijet = 0; ijet< nTriggerJets; ijet++)
482  {
483  pygeh->TriggerJet(ijet, tmpjet);
484  jet = new TParticle(94, 21, -1, -1, -1, -1, tmpjet[0],tmpjet[1],tmpjet[2],tmpjet[3], 0,0,0,0);
485 
486  AliDebug(1,Form("jet %d; pycell jet pT %f",ijet, jet->Pt()));
487 
488  //Compare jet pT and pt Hard
489  if(jet->Pt() > fPtHardAndJetPtFactor * ptHard)
490  {
491  AliInfo(Form("Reject jet event with : pT Hard %2.2f, pycell jet pT %2.2f, rejection factor %1.1f\n",
492  ptHard, jet->Pt(), fPtHardAndJetPtFactor));
493  return kFALSE;
494  }
495  }
496 
497  if(jet) delete jet;
498  }
499 
500  return kTRUE ;
501 }
502 
503 //____________________________________________________
508 //____________________________________________________
510 {
511  if(!strcmp(GetGenEventHeader()->ClassName(), "AliGenPythiaEventHeader"))
512  {
513  AliGenPythiaEventHeader* pygeh= (AliGenPythiaEventHeader*) GetGenEventHeader();
514  Float_t ptHard = pygeh->GetPtHard();
515 
516  Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
517  for (Int_t iclus = 0; iclus < nclusters; iclus++)
518  {
519  AliVCluster * clus = fInputEvent->GetCaloCluster(iclus) ;
520  Float_t ecluster = clus->E();
521 
522  if(ecluster > fPtHardAndClusterPtFactor*ptHard)
523  {
524  AliInfo(Form("Reject : ecluster %2.2f, calo %d, factor %2.2f, ptHard %f",ecluster,clus->GetType(),fPtHardAndClusterPtFactor,ptHard));
525 
526  return kFALSE;
527  }
528  }
529 
530  }
531 
532  return kTRUE ;
533 }
534 
535 //____________________________________________
537 //____________________________________________
539 {
540  if(fMC)
541  return fMC->Stack();
542  else
543  {
544  AliDebug(1,"Stack is not available");
545  return 0x0 ;
546  }
547 }
548 
549 //______________________________________________
551 //______________________________________________
553 {
554  if(fMC)
555  {
556  return fMC->Header();
557  }
558  else
559  {
560  AliInfo("Header is not available");
561  return 0x0 ;
562  }
563 }
564 
565 //____________________________________________________
569 //____________________________________________________
571 {
572  fNMCProducedMin = 0;
573  fNMCProducedMax = 0;
574 
575  if (ReadStack() && fMC)
576  {
577  AliGenEventHeader * eventHeader = fMC->GenEventHeader();
578 
579  if(!fAcceptOnlyHIJINGLabels) return ;
580 
581  // TODO Check if it works from here ...
582 
583  AliGenCocktailEventHeader *cocktail = dynamic_cast<AliGenCocktailEventHeader *>(eventHeader);
584 
585  if(!cocktail) return ;
586 
587  TList *genHeaders = cocktail->GetHeaders();
588 
589  Int_t nGenerators = genHeaders->GetEntries();
590  //printf("N generators %d \n", nGenerators);
591 
592  for(Int_t igen = 0; igen < nGenerators; igen++)
593  {
594  AliGenEventHeader * eventHeader2 = (AliGenEventHeader*)genHeaders->At(igen) ;
595  TString name = eventHeader2->GetName();
596 
597  //printf("Generator %d: Class Name %s, Name %s, title %s \n",igen, eventHeader2->ClassName(), name.Data(), eventHeader2->GetTitle());
598 
600  fNMCProducedMax+= eventHeader2->NProduced();
601 
602  if(name.Contains("Hijing",TString::kIgnoreCase)) return ;
603  }
604 
605  }
606  else if(ReadAODMCParticles() && GetAODMCHeader())
607  {
608  Int_t nGenerators = GetAODMCHeader()->GetNCocktailHeaders();
609  //printf("AliCaloTrackReader::GetGenEventHeader() - N headers %d\n",nGenerators);
610 
611  if( nGenerators <= 0) return ;
612 
613  if(!fAcceptOnlyHIJINGLabels) return ;
614 
615  for(Int_t igen = 0; igen < nGenerators; igen++)
616  {
617  AliGenEventHeader * eventHeader = GetAODMCHeader()->GetCocktailHeader(igen) ;
618  TString name = eventHeader->GetName();
619 
620  //printf("Generator %d: Class Name %s, Name %s, title %s \n",igen, eventHeader->ClassName(), name.Data(), eventHeader->GetTitle());
621 
623  fNMCProducedMax+= eventHeader->NProduced();
624 
625  if(name.Contains("Hijing",TString::kIgnoreCase)) return ;
626  }
627 
628  }
629 }
630 
631 //______________________________________________________________
634 //______________________________________________________________
635 AliGenEventHeader* AliCaloTrackReader::GetGenEventHeader() const
636 {
637  if (ReadStack() && fMC)
638  {
639  AliGenEventHeader * eventHeader = fMC->GenEventHeader();
640 
641  if(!fAcceptOnlyHIJINGLabels) return eventHeader ;
642 
643  // TODO Check if it works from here ...
644 
645  AliGenCocktailEventHeader *cocktail = dynamic_cast<AliGenCocktailEventHeader *>(eventHeader);
646 
647  if(!cocktail) return 0x0 ;
648 
649  TList *genHeaders = cocktail->GetHeaders();
650 
651  Int_t nGenerators = genHeaders->GetEntries();
652  //printf("N generators %d \n", nGenerators);
653 
654  for(Int_t igen = 0; igen < nGenerators; igen++)
655  {
656  AliGenEventHeader * eventHeader2 = (AliGenEventHeader*)genHeaders->At(igen) ;
657  TString name = eventHeader2->GetName();
658 
659  //printf("Generator %d: Class Name %s, Name %s, title %s \n",igen, eventHeader2->ClassName(), name.Data(), eventHeader2->GetTitle());
660 
661  if(name.Contains("Hijing",TString::kIgnoreCase)) return eventHeader2 ;
662  }
663 
664  return 0x0;
665 
666  }
667  else if(ReadAODMCParticles() && GetAODMCHeader())
668  {
669  Int_t nGenerators = GetAODMCHeader()->GetNCocktailHeaders();
670  //printf("AliCaloTrackReader::GetGenEventHeader() - N headers %d\n",nGenerators);
671 
672  if( nGenerators <= 0) return 0x0;
673 
674  if(!fAcceptOnlyHIJINGLabels) return GetAODMCHeader()->GetCocktailHeader(0);
675 
676  for(Int_t igen = 0; igen < nGenerators; igen++)
677  {
678  AliGenEventHeader * eventHeader = GetAODMCHeader()->GetCocktailHeader(igen) ;
679  TString name = eventHeader->GetName();
680 
681  //printf("Generator %d: Class Name %s, Name %s, title %s \n",igen, eventHeader->ClassName(), name.Data(), eventHeader->GetTitle());
682 
683  if(name.Contains("Hijing",TString::kIgnoreCase)) return eventHeader ;
684  }
685 
686  return 0x0;
687 
688  }
689  else
690  {
691  //printf("AliCaloTrackReader::GetGenEventHeader() - MC header not available! \n");
692  return 0x0;
693  }
694 }
695 
696 //_________________________________________________________
699 //_________________________________________________________
701 {
702  AliInfo("Input are not AODs");
703 
704  return NULL ;
705 }
706 
707 //________________________________________________________
710 //________________________________________________________
711 AliAODMCHeader* AliCaloTrackReader::GetAODMCHeader() const
712 {
713  AliInfo("Input are not AODs");
714 
715  return NULL ;
716 }
717 
718 //___________________________________________________________
726 //___________________________________________________________
727 Int_t AliCaloTrackReader::GetVertexBC(const AliVVertex * vtx)
728 {
729  Int_t vertexBC=vtx->GetBC();
730 
731  if(!fRecalculateVertexBC) return vertexBC;
732 
733  // Value not available, recalculate it.
734 
735  Double_t bz = fInputEvent->GetMagneticField();
736  Bool_t bc0 = kFALSE;
737  Int_t ntr = GetCTSTracks()->GetEntriesFast();
738  //printf("N Tracks %d\n",ntr);
739 
740  for(Int_t i = 0 ; i < ntr ; i++)
741  {
742  AliVTrack * track = (AliVTrack*) (GetCTSTracks()->At(i));
743 
744  //Check if has TOF info, if not skip
745  ULong_t status = track->GetStatus();
746  Bool_t okTOF = (status & AliVTrack::kTOFout) == AliVTrack::kTOFout ;
747  vertexBC = track->GetTOFBunchCrossing(bz);
748  Float_t pt = track->Pt();
749 
750  if(!okTOF) continue;
751 
752  // Get DCA x, y
753  Double_t dca[2] = {1e6,1e6};
754  Double_t covar[3] = {1e6,1e6,1e6};
755  track->PropagateToDCA(vtx,bz,100.,dca,covar);
756 
757  if(AcceptDCA(pt,dca[0]))
758  {
759  if (vertexBC !=0 && fVertexBC != AliVTrack::kTOFBCNA) return vertexBC;
760  else if(vertexBC == 0) bc0 = kTRUE;
761  }
762  }
763 
764  if( bc0 ) vertexBC = 0 ;
765  else vertexBC = AliVTrack::kTOFBCNA ;
766 
767  return vertexBC;
768 }
769 
770 //_____________________________
773 //_____________________________
775 {
777  {
778  AliInfo("Cannot access stack and mcparticles at the same time, change them");
779  fReadStack = kFALSE;
780  fReadAODMCParticles = kFALSE;
781  }
782 
783  // Activate debug level in AliAnaWeights
784  if( fWeightUtils->GetDebug() >= 0 )
785  (AliAnalysisManager::GetAnalysisManager())->AddClassDebug(fWeightUtils->ClassName(), fWeightUtils->GetDebug());
786 }
787 
788 //_______________________________________
790 //_______________________________________
792 {
793  fDataType = kESD ;
794  fCTSPtMin = 0.1 ;
795  fEMCALPtMin = 0.1 ;
796  fPHOSPtMin = 0.1 ;
797  fCTSPtMax = 1000. ;
798  fEMCALPtMax = 1000. ;
799  fPHOSPtMax = 1000. ;
800 
801  //Track DCA cuts
802  // dca_xy cut = 0.0105+0.0350/TMath::Power(pt,1.1);
803  fTrackDCACut[0] = 0.0105;
804  fTrackDCACut[1] = 0.0350;
805  fTrackDCACut[2] = 1.1;
806 
807  //Do not filter the detectors input by default.
808  fFillEMCAL = kFALSE;
809  fFillDCAL = kFALSE;
810  fFillPHOS = kFALSE;
811  fFillCTS = kFALSE;
812  fFillEMCALCells = kFALSE;
813  fFillPHOSCells = kFALSE;
814 
815  fReadStack = kFALSE; // Check in the constructor of the other readers if it was set or in the configuration file
816  fReadAODMCParticles = kFALSE; // Check in the constructor of the other readers if it was set or in the configuration file
817  fDeltaAODFileName = "deltaAODPartCorr.root";
819  fEventTriggerMask = AliVEvent::kAny;
820  fMixEventTriggerMask = AliVEvent::kAnyINT;
821  fEventTriggerAtSE = kTRUE; // Use only events that pass event selection at SE base class
822 
823  fAcceptFastCluster = kTRUE;
824  fEventType = -1;
825 
826  //We want tracks fitted in the detectors:
827  //fTrackStatus=AliESDtrack::kTPCrefit;
828  //fTrackStatus|=AliESDtrack::kITSrefit;
829  fTrackStatus = 0;
830 
831  fV0ADC[0] = 0; fV0ADC[1] = 0;
832  fV0Mul[0] = 0; fV0Mul[1] = 0;
833 
834  fZvtxCut = 10.;
835 
836  fNMixedEvent = 1;
837 
840 
841  //Centrality
842  fUseAliCentrality = kFALSE;
843  fCentralityClass = "V0M";
844  fCentralityOpt = 100;
845  fCentralityBin[0] = fCentralityBin[1]=-1;
846 
847  fEventPlaneMethod = "V0";
848 
849  // Allocate memory (not sure this is the right place)
850  fCTSTracks = new TObjArray();
851  fEMCALClusters = new TObjArray();
852  fDCALClusters = new TObjArray();
853  fPHOSClusters = new TObjArray();
854  //fTriggerAnalysis = new AliTriggerAnalysis;
855  fAODBranchList = new TList ;
856 
857  fPileUpParamSPD[0] = 3 ; fPileUpParamSPD[1] = 0.8 ;
858  fPileUpParamSPD[2] = 3.0 ; fPileUpParamSPD[3] = 2.0 ; fPileUpParamSPD[4] = 5.0;
859 
860  // Parametrized time cut (LHC11d)
863 
864  // Parametrized time cut (LHC11c)
865  //fEMCALParamTimeCutMin[0] =-5; fEMCALParamTimeCutMin[1] =-1 ; fEMCALParamTimeCutMin[2] = 1.87; fEMCALParamTimeCutMin[3] = 0.4;
866  //fEMCALParamTimeCutMax[0] = 3.5; fEMCALParamTimeCutMax[1] = 50; fEMCALParamTimeCutMax[2] = 0.15; fEMCALParamTimeCutMax[3] = 1.6;
867 
868  fTimeStampRunMin = -1;
869  fTimeStampRunMax = 1e12;
872 
873  for(Int_t i = 0; i < 19; i++)
874  {
875  fEMCalBCEvent [i] = 0;
876  fEMCalBCEventCut[i] = 0;
877  fTrackBCEvent [i] = 0;
878  fTrackBCEventCut[i] = 0;
879  }
880 
881  // Trigger match-rejection
884 
885  fTriggerClusterBC = -10000 ;
889  fTriggerClusterId = -1;
890 
891  //Jets
894  if(!fNonStandardJets) fNonStandardJets = new TClonesArray("AliAODJet",100);
897  if(!fBackgroundJets) fBackgroundJets = new AliAODJetEventBackground();
898 
899  fSmearShowerShapeWidth = 0.005;
900 
901  fWeightUtils = new AliAnaWeights() ;
902  fEventWeight = 1 ;
903 
904 }
905 
906 //___________________________________________________
911 //___________________________________________________
913 {
914  AliGenHijingEventHeader* hijingHeader = dynamic_cast<AliGenHijingEventHeader *> (GetGenEventHeader());
915 
916  //printf("header %p, label %d\n",hijingHeader,label);
917 
918  if(!hijingHeader || label < 0 ) return kFALSE;
919 
920 
921  //printf("pass a), N produced %d\n",nproduced);
922 
923  if(label >= fNMCProducedMin && label < fNMCProducedMax)
924  {
925  //printf(" accept!, label is smaller than produced, N %d\n",nproduced);
926 
927  return kTRUE;
928  }
929 
930  if(ReadStack())
931  {
932  if(!GetStack()) return kFALSE;
933 
934  Int_t nprimaries = GetStack()->GetNtrack();
935 
936  if(label > nprimaries) return kFALSE;
937 
938  TParticle * mom = GetStack()->Particle(label);
939 
940  Int_t iMom = label;
941  Int_t iParent = mom->GetFirstMother();
942  while(iParent!=-1)
943  {
944  if(iParent >= fNMCProducedMin && iParent < fNMCProducedMax)
945  {
946  //printf("\t accept, mother is %d \n",iParent)
947  return kTRUE;
948  }
949 
950  iMom = iParent;
951  mom = GetStack()->Particle(iMom);
952  iParent = mom->GetFirstMother();
953  }
954 
955  return kFALSE ;
956 
957  } // ESD
958  else
959  {
960  TClonesArray* mcparticles = GetAODMCParticles();
961 
962  if(!mcparticles) return kFALSE;
963 
964  Int_t nprimaries = mcparticles->GetEntriesFast();
965 
966  if(label > nprimaries) return kFALSE;
967 
968  //printf("pass b) N primaries %d \n",nprimaries);
969 
970  if(label >= fNMCProducedMin && label < fNMCProducedMax)
971  {
972  return kTRUE;
973  }
974 
975  // Find grand parent, check if produced in the good range
976  AliAODMCParticle * mom = (AliAODMCParticle *) mcparticles->At(label);
977 
978  Int_t iMom = label;
979  Int_t iParent = mom->GetMother();
980  while(iParent!=-1)
981  {
982  if(iParent >= fNMCProducedMin && iParent < fNMCProducedMax)
983  {
984  //printf("\t accept, mother is %d, with nProduced %d \n",iParent, nproduced);
985  return kTRUE;
986  }
987 
988  iMom = iParent;
989  mom = (AliAODMCParticle *) mcparticles->At(iMom);
990  iParent = mom->GetMother();
991 
992  }
993 
994  //printf("pass c), no match found \n");
995 
996  return kFALSE ;
997 
998  }//AOD
999 }
1000 
1001 //__________________________________________________________________________
1004 //__________________________________________________________________________
1005 Bool_t AliCaloTrackReader::IsInTimeWindow(Double_t tof, Float_t energy) const
1006 {
1007  // Parametrized cut depending on E
1008  if(fUseParamTimeCut)
1009  {
1010  Float_t minCut= fEMCALParamTimeCutMin[0]+fEMCALParamTimeCutMin[1]*TMath::Exp(-(energy-fEMCALParamTimeCutMin[2])/fEMCALParamTimeCutMin[3]);
1011  Float_t maxCut= fEMCALParamTimeCutMax[0]+fEMCALParamTimeCutMax[1]*TMath::Exp(-(energy-fEMCALParamTimeCutMax[2])/fEMCALParamTimeCutMax[3]);
1012  //printf("tof %f, minCut %f, maxCut %f\n",tof,minCut,maxCut);
1013  if( tof < minCut || tof > maxCut ) return kFALSE ;
1014  }
1015 
1016  //In any case, the time should to be larger than the fixed window ...
1017  if( tof < fEMCALTimeCutMin || tof > fEMCALTimeCutMax ) return kFALSE ;
1018 
1019  return kTRUE ;
1020 }
1021 
1022 //________________________________________________
1025 //________________________________________________
1027 {
1028  return fInputEvent->IsPileupFromSPD((Int_t) fPileUpParamSPD[0] , fPileUpParamSPD[1] ,
1029  fPileUpParamSPD[2] , fPileUpParamSPD[3] , fPileUpParamSPD[4] );
1030  //printf("Param : %d, %2.2f, %2.2f, %2.2f, %2.2f\n",(Int_t) fPileUpParamSPD[0], fPileUpParamSPD[1], fPileUpParamSPD[2], fPileUpParamSPD[3], fPileUpParamSPD[4]);
1031 }
1032 
1033 //__________________________________________________
1035 //__________________________________________________
1037 {
1038  if(fNPileUpClusters > fNPileUpClustersCut) return kTRUE ;
1039  else return kFALSE;
1040 }
1041 
1042 //________________________________________________________
1044 //________________________________________________________
1046 {
1047  if( IsPileUpFromSPD() && IsPileUpFromEMCal()) return kTRUE ;
1048  else return kFALSE;
1049 }
1050 
1051 //_______________________________________________________
1053 //_______________________________________________________
1055 {
1056  if( IsPileUpFromSPD() || IsPileUpFromEMCal()) return kTRUE ;
1057  else return kFALSE;
1058 }
1059 
1060 //___________________________________________________________
1062 //___________________________________________________________
1064 {
1065  if( IsPileUpFromSPD() && !IsPileUpFromEMCal()) return kTRUE ;
1066  else return kFALSE;
1067 }
1068 
1069 //___________________________________________________________
1071 //___________________________________________________________
1073 {
1074  if( !IsPileUpFromSPD() && IsPileUpFromEMCal()) return kTRUE ;
1075  else return kFALSE;
1076 }
1077 
1078 //______________________________________________________________
1080 //______________________________________________________________
1082 {
1083  if( !IsPileUpFromSPD() && !IsPileUpFromEMCal()) return kTRUE ;
1084  else return kFALSE;
1085 }
1086 
1087 //___________________________________________________________________________________
1095 //___________________________________________________________________________________
1096 Bool_t AliCaloTrackReader::FillInputEvent(Int_t iEntry, const char * /*curFileName*/)
1097 {
1098  fEventNumber = iEntry;
1099  fTriggerClusterIndex = -1;
1100  fTriggerClusterId = -1;
1101  fIsTriggerMatch = kFALSE;
1102  fTriggerClusterBC = -10000;
1103  fIsExoticEvent = kFALSE;
1104  fIsBadCellEvent = kFALSE;
1105  fIsBadMaxCellEvent = kFALSE;
1106 
1107  fIsTriggerMatchOpenCut[0] = kFALSE ;
1108  fIsTriggerMatchOpenCut[1] = kFALSE ;
1109  fIsTriggerMatchOpenCut[2] = kFALSE ;
1110 
1111  //fCurrentFileName = TString(currentFileName);
1112  if(!fInputEvent)
1113  {
1114  AliInfo("Input event not available, skip event analysis");
1115  return kFALSE;
1116  }
1117 
1118  //-----------------------------------------------
1119  // Select the event depending on the trigger type
1120  // and other event characteristics
1121  // like the goodness of the EMCal trigger
1122  //-----------------------------------------------
1123 
1124  Bool_t accept = CheckEventTriggers();
1125  if(!accept) return kFALSE;
1126 
1127  AliDebug(1,"Pass Event trigger selection");
1128 
1129  //---------------------------------------------------------------------------
1130  // In case of analysis of events with jets, skip those with jet pt > 5 pt hard
1131  // To be used on for MC data in pT hard bins
1132  //---------------------------------------------------------------------------
1133 
1135  {
1136  if(!ComparePtHardAndJetPt()) return kFALSE ;
1137  AliDebug(1,"Pass Pt Hard - Jet rejection");
1138  }
1139 
1141  {
1142  if(!ComparePtHardAndClusterPt()) return kFALSE ;
1143  AliDebug(1,"Pass Pt Hard - Cluster rejection");
1144  }
1145 
1146  //------------------------------------------------------
1147  // Event rejection depending on time stamp
1148  //------------------------------------------------------
1149 
1151  {
1152  AliESDEvent* esd = dynamic_cast<AliESDEvent*> (fInputEvent);
1153  if(esd)
1154  {
1155  Int_t timeStamp = esd->GetTimeStamp();
1156  Float_t timeStampFrac = 1.*(timeStamp-fTimeStampRunMin) / (fTimeStampRunMax-fTimeStampRunMin);
1157 
1158  //printf("stamp0 %d, max0 %d, frac %f\n", timeStamp-fTimeStampRunMin,fTimeStampRunMax-fTimeStampRunMin, timeStampFrac);
1159 
1160  if(timeStampFrac < fTimeStampEventFracMin || timeStampFrac > fTimeStampEventFracMax) return kFALSE;
1161  }
1162 
1163  AliDebug(1,"Pass Time Stamp rejection");
1164  }
1165 
1166  //------------------------------------------------------
1167  // Event rejection depending on vertex, pileup, v0and
1168  //------------------------------------------------------
1169 
1170  FillVertexArray();
1171 
1172  //Reject events with Z vertex too large, only for SE analysis, if not, cut on the analysis code
1173  if(!GetMixedEvent() && TMath::Abs(fVertex[0][2]) > fZvtxCut) return kFALSE;
1174 
1176  {
1177  if( !CheckForPrimaryVertex() ) return kFALSE; // algorithm in ESD/AOD Readers
1178  if( TMath::Abs(fVertex[0][0] ) < 1.e-6 &&
1179  TMath::Abs(fVertex[0][1] ) < 1.e-6 &&
1180  TMath::Abs(fVertex[0][2] ) < 1.e-6 ) return kFALSE;
1181  }
1182 
1183  AliDebug(1,"Pass primary vertex rejection");
1184 
1185  //printf("Reader : IsPileUp %d, Multi %d\n",IsPileUpFromSPD(),fInputEvent->IsPileupFromSPDInMultBins());
1186 
1188  {
1189  // Do not analyze events with pileup
1190  Bool_t bPileup = IsPileUpFromSPD();
1191  //IsPileupFromSPDInMultBins() // method to try
1192  //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]);
1193  if(bPileup) return kFALSE;
1194 
1195  AliDebug(1,"Pass Pile-Up event rejection");
1196  }
1197 
1199  {
1200  AliVVZERO* v0 = fInputEvent->GetVZEROData();
1201 
1202  Bool_t bV0AND = ((v0->GetV0ADecision()==1) && (v0->GetV0CDecision()==1));
1203  //bV0AND = fTriggerAnalysis->IsOfflineTriggerFired((AliESDEvent*)fInputEvent, AliTriggerAnalysis::kV0AND);
1204  //printf("V0AND event? %d\n",bV0AND);
1205 
1206  if(!bV0AND)
1207  {
1208  AliDebug(1,"Reject event by V0AND");
1209  return kFALSE;
1210  }
1211 
1212  AliDebug(1,"Pass V0AND event rejection");
1213  }
1214 
1215  //------------------------------------------------------
1216  // Check if there is a centrality value, PbPb analysis,
1217  // and if a centrality bin selection is requested
1218  //------------------------------------------------------
1219 
1220  //If we need a centrality bin, we select only those events in the corresponding bin.
1221  Int_t cen = -1;
1222  if ( fCentralityBin[0] >= 0 && fCentralityBin[1] >= 0 )
1223  {
1224  cen = GetEventCentrality();
1225 
1226  if(cen > fCentralityBin[1] || cen < fCentralityBin[0]) return kFALSE; //reject events out of bin.
1227 
1228  AliDebug(1,"Pass centrality rejection");
1229  }
1230 
1231  // Recover the weight assigned to the event, if provided
1232  // right now only for pT-hard bins and centrality depedent weights
1234  {
1236 
1238  }
1239  //---------------------------------------------------------------------------
1240  // In case of MC analysis, set the range of interest in the MC particles list
1241  // The particle selection is done inside the FillInputDetector() methods
1242  //---------------------------------------------------------------------------
1244 
1245  //printf("N min %d, N max %d\n",fNMCProducedMin,fNMCProducedMax);
1246 
1247  //-------------------------------------------------------
1248  // Get the main vertex BC, in case not available
1249  // it is calculated in FillCTS checking the BC of tracks
1250  //------------------------------------------------------
1251  fVertexBC = fInputEvent->GetPrimaryVertex()->GetBC();
1252 
1253  //-----------------------------------------------
1254  // Fill the arrays with cluster/tracks/cells data
1255  //-----------------------------------------------
1256 
1257  if(fFillCTS)
1258  {
1259  FillInputCTS();
1260  //Accept events with at least one track
1261  if(fTrackMult == 0 && fDoRejectNoTrackEvents) return kFALSE ;
1262 
1263  AliDebug(1,"Pass rejection of null track events");
1264  }
1265 
1267  {
1268  if(fVertexBC != 0 && fVertexBC != AliVTrack::kTOFBCNA) return kFALSE ;
1269 
1270  AliDebug(1,"Pass rejection of events with vertex at BC!=0");
1271  }
1272 
1273  if(fFillEMCALCells)
1275 
1276  if(fFillPHOSCells)
1278 
1279  if(fFillEMCAL || fFillDCAL)
1280  FillInputEMCAL();
1281 
1282  if(fFillPHOS)
1283  FillInputPHOS();
1284 
1285  FillInputVZERO();
1286 
1287  //one specified jet branch
1292 
1293  AliDebug(1,"Event accepted for analysis");
1294 
1295  return kTRUE ;
1296 }
1297 
1298 //__________________________________________________
1301 //__________________________________________________
1303 {
1304  if(fUseAliCentrality)
1305  {
1306  if ( !GetCentrality() ) return -1;
1307 
1308  if (fCentralityOpt == 100) return (Int_t) GetCentrality()->GetCentralityPercentile(fCentralityClass); // 100 bins max
1309  else if(fCentralityOpt == 10) return GetCentrality()->GetCentralityClass10(fCentralityClass);// 10 bins max
1310  else if(fCentralityOpt == 20) return GetCentrality()->GetCentralityClass5(fCentralityClass); // 20 bins max
1311  else
1312  {
1313  AliInfo(Form("Unknown centrality option %d, use 10, 20 or 100\n",fCentralityOpt));
1314  return -1;
1315  }
1316  }
1317  else
1318  {
1319  if ( !GetMultSelCen() ) return -1;
1320 
1321  return GetMultSelCen()->GetMultiplicityPercentile(fCentralityClass, kTRUE); // returns centrality only for events used in calibration
1322 
1323  // equivalent to
1324  //GetMultSelCen()->GetMultiplicityPercentile("V0M", kFALSE); // returns centrality for any event
1325  //Int_t qual = GetMultSelCen()->GetEvSelCode(); if (qual ! = 0) cent = qual;
1326  }
1327 }
1328 
1329 //_____________________________________________________
1332 //_____________________________________________________
1334 {
1335  if( !GetEventPlane() ) return -1000;
1336 
1337  Float_t ep = GetEventPlane()->GetEventplane(GetEventPlaneMethod(), GetInputEvent());
1338 
1339  if(GetEventPlaneMethod()=="Q" && (ep < 0 || ep > TMath::Pi()))
1340  {
1341  AliDebug(1,Form("Bad EP for <Q> method : %f\n",ep));
1342  return -1000;
1343  }
1344  else if(GetEventPlaneMethod().Contains("V0") )
1345  {
1346  if((ep > TMath::Pi()/2 || ep < -TMath::Pi()/2))
1347  {
1348  AliDebug(1,Form("Bad EP for <%s> method : %f\n",GetEventPlaneMethod().Data(), ep));
1349  return -1000;
1350  }
1351 
1352  ep+=TMath::Pi()/2; // put same range as for <Q> method
1353  }
1354 
1355  AliDebug(3,Form("Event plane angle %f",ep));
1356 
1357 // if(fDebug > 0 )
1358 // {
1359 // if (ep > TMath::Pi()) printf("AliCaloTrackReader::GetEventPlaneAngle() - Too large angle = %f\n",ep);
1360 // else if(ep < 0 ) printf("AliCaloTrackReader::GetEventPlaneAngle() - Negative angle = %f\n" ,ep);
1361 // }
1362 
1363  return ep;
1364 }
1365 
1366 //__________________________________________________________
1368 //__________________________________________________________
1369 void AliCaloTrackReader::GetVertex(Double_t vertex[3]) const
1370 {
1371  vertex[0] = fVertex[0][0];
1372  vertex[1] = fVertex[0][1];
1373  vertex[2] = fVertex[0][2];
1374 }
1375 
1376 //__________________________________________________________________________
1378 //__________________________________________________________________________
1379 void AliCaloTrackReader::GetVertex(Double_t vertex[3], Int_t evtIndex) const
1380 {
1381  vertex[0] = fVertex[evtIndex][0];
1382  vertex[1] = fVertex[evtIndex][1];
1383  vertex[2] = fVertex[evtIndex][2];
1384 }
1385 
1386 //________________________________________
1389 //________________________________________
1391 {
1392  // Delete previous vertex
1393  if(fVertex)
1394  {
1395  for (Int_t i = 0; i < fNMixedEvent; i++)
1396  {
1397  delete [] fVertex[i] ;
1398  }
1399  delete [] fVertex ;
1400  }
1401 
1402  fVertex = new Double_t*[fNMixedEvent] ;
1403  for (Int_t i = 0; i < fNMixedEvent; i++)
1404  {
1405  fVertex[i] = new Double_t[3] ;
1406  fVertex[i][0] = 0.0 ;
1407  fVertex[i][1] = 0.0 ;
1408  fVertex[i][2] = 0.0 ;
1409  }
1410 
1411  if (!fMixedEvent)
1412  { // Single event analysis
1413  if(fDataType!=kMC)
1414  {
1415 
1416  if(fInputEvent->GetPrimaryVertex())
1417  {
1418  fInputEvent->GetPrimaryVertex()->GetXYZ(fVertex[0]);
1419  }
1420  else
1421  {
1422  AliWarning("NULL primary vertex");
1423  fVertex[0][0]=0.; fVertex[0][1]=0.; fVertex[0][2]=0.;
1424  }//Primary vertex pointer do not exist
1425 
1426  } else
1427  {// MC read event
1428  fVertex[0][0]=0.; fVertex[0][1]=0.; fVertex[0][2]=0.;
1429  }
1430 
1431  AliDebug(1,Form("Single Event Vertex : %f,%f,%f\n",fVertex[0][0],fVertex[0][1],fVertex[0][2]));
1432 
1433  } else
1434  { // MultiEvent analysis
1435  for (Int_t iev = 0; iev < fNMixedEvent; iev++)
1436  {
1437  if (fMixedEvent->GetVertexOfEvent(iev))
1438  fMixedEvent->GetVertexOfEvent(iev)->GetXYZ(fVertex[iev]);
1439  else
1440  AliWarning("No vertex found");
1441 
1442  AliDebug(1,Form("Multi Event %d Vertex : %f,%f,%f",iev,fVertex[iev][0],fVertex[iev][1],fVertex[iev][2]));
1443  }
1444  }
1445 }
1446 
1447 //_____________________________________
1453 //_____________________________________
1455 {
1456  AliDebug(1,"Begin");
1457 
1458  Double_t pTrack[3] = {0,0,0};
1459 
1460  Int_t nTracks = fInputEvent->GetNumberOfTracks() ;
1461  fTrackMult = 0;
1462  Int_t nstatus = 0;
1463  Double_t bz = GetInputEvent()->GetMagneticField();
1464 
1465  for(Int_t i = 0; i < 19; i++)
1466  {
1467  fTrackBCEvent [i] = 0;
1468  fTrackBCEventCut[i] = 0;
1469  }
1470 
1471  Bool_t bc0 = kFALSE;
1472  if(fRecalculateVertexBC) fVertexBC = AliVTrack::kTOFBCNA;
1473 
1474  for (Int_t itrack = 0; itrack < nTracks; itrack++)
1475  {
1476  AliVTrack * track = (AliVTrack*)fInputEvent->GetTrack(itrack) ; // retrieve track from esd
1477 
1478  if(fAcceptOnlyHIJINGLabels && !IsHIJINGLabel(track->GetLabel())) continue ;
1479 
1480  //Select tracks under certain conditions, TPCrefit, ITSrefit ... check the set bits
1481  ULong_t status = track->GetStatus();
1482 
1483  if (fTrackStatus && !((status & fTrackStatus) == fTrackStatus))
1484  continue ;
1485 
1486  nstatus++;
1487 
1488  // Select the tracks depending on cuts of AOD or ESD
1489  if(!SelectTrack(track, pTrack)) continue ;
1490 
1491  // TOF cuts
1492  Bool_t okTOF = ( (status & AliVTrack::kTOFout) == AliVTrack::kTOFout ) ;
1493  Double_t tof = -1000;
1494  Int_t trackBC = -1000 ;
1495 
1496  if(fAccessTrackTOF)
1497  {
1498  if(okTOF)
1499  {
1500  trackBC = track->GetTOFBunchCrossing(bz);
1501  SetTrackEventBC(trackBC+9);
1502 
1503  tof = track->GetTOFsignal()*1e-3;
1504 
1505  //After selecting tracks with small DCA, pointing to vertex, set vertex BC depeding on tracks BC
1507  {
1508  if (trackBC != 0 && trackBC != AliVTrack::kTOFBCNA) fVertexBC = trackBC;
1509  else if(trackBC == 0) bc0 = kTRUE;
1510  }
1511 
1512  //In any case, the time should to be larger than the fixed window ...
1513  if( fUseTrackTimeCut && (trackBC !=0 || tof < fTrackTimeCutMin || tof > fTrackTimeCutMax) )
1514  {
1515  //printf("Remove track time %f and bc = %d\n",tof,trackBC);
1516  continue ;
1517  }
1518  //else printf("Accept track time %f and bc = %d\n",tof,trackBC);
1519  }
1520  }
1521 
1522  fMomentum.SetPxPyPzE(pTrack[0],pTrack[1],pTrack[2],0);
1523 
1524  // DCA cuts
1525 
1526  if(fUseTrackDCACut)
1527  {
1528  Float_t dcaTPC =-999;
1529  //In case of AODs, TPC tracks cannot be propagated back to primary vertex,
1530  if( fDataType == kAOD ) dcaTPC = ((AliAODTrack*) track)->DCA();
1531 
1532  //normal way to get the dca, cut on dca_xy
1533  if(dcaTPC==-999)
1534  {
1535  Double_t dca[2] = {1e6,1e6};
1536  Double_t covar[3] = {1e6,1e6,1e6};
1537  Bool_t okDCA = track->PropagateToDCA(fInputEvent->GetPrimaryVertex(),bz,100.,dca,covar);
1538  if( okDCA) okDCA = AcceptDCA(fMomentum.Pt(),dca[0]);
1539  if(!okDCA)
1540  {
1541  //printf("AliCaloTrackReader::FillInputCTS() - Reject track pt %2.2f, dca_xy %2.4f\n",fMomentum.Pt(),dca[0]);
1542  continue ;
1543  }
1544  }
1545  }// DCA cuts
1546 
1547  //Count the tracks in eta < 0.9
1548  if(TMath::Abs(track->Eta())< fTrackMultEtaCut) fTrackMult++;
1549 
1550  if(fCTSPtMin > fMomentum.Pt() || fCTSPtMax < fMomentum.Pt()) continue ;
1551 
1552  // Check effect of cuts on track BC
1553  if(fAccessTrackTOF && okTOF) SetTrackEventBCcut(trackBC+9);
1554 
1555  if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(fMomentum.Eta(),fMomentum.Phi(),kCTS)) continue;
1556 
1557  AliDebug(2,Form("Selected tracks pt %3.2f, phi %3.2f deg, eta %3.2f",
1558  fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta()));
1559 
1560  if (fMixedEvent) track->SetID(itrack);
1561 
1562 
1563  fCTSTracks->Add(track);
1564 
1565  }// track loop
1566 
1567  if( fRecalculateVertexBC && (fVertexBC == 0 || fVertexBC == AliVTrack::kTOFBCNA))
1568  {
1569  if( bc0 ) fVertexBC = 0 ;
1570  else fVertexBC = AliVTrack::kTOFBCNA ;
1571  }
1572 
1573  AliDebug(1,Form("AOD entries %d, input tracks %d, pass status %d, multipliticy %d", fCTSTracks->GetEntriesFast(), nTracks, nstatus, fTrackMult));//fCTSTracksNormalInputEntries);
1574 }
1575 
1576 //_______________________________________________________________________________
1583 //_______________________________________________________________________________
1584 void AliCaloTrackReader::FillInputEMCALAlgorithm(AliVCluster * clus, Int_t iclus)
1585 {
1586  // Accept clusters with the proper label, TO DO, use the new more General methods!!!
1587  if(fAcceptOnlyHIJINGLabels && !IsHIJINGLabel( clus->GetLabel() )) return ;
1588 
1589  Int_t vindex = 0 ;
1590  if (fMixedEvent)
1591  vindex = fMixedEvent->EventIndexForCaloCluster(iclus);
1592 
1593  clus->GetMomentum(fMomentum, fVertex[vindex]);
1594 
1595  //if( (fDebug > 2 && fMomentum.E() > 0.1) || fDebug > 10 )
1596  AliDebug(2,Form("Input cluster E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f",
1597  fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta()));
1598 
1600  {
1601  //Recalibrate the cluster energy
1602  if(GetCaloUtils()->IsRecalibrationOn())
1603  {
1605 
1606  clus->SetE(energy);
1607  //printf("Recalibrated Energy %f\n",clus->E());
1608 
1611 
1612  } // recalculate E
1613 
1614  //Recalculate distance to bad channels, if new list of bad channels provided
1616 
1617  //Recalculate cluster position
1618  if(GetCaloUtils()->IsRecalculationOfClusterPositionOn())
1619  {
1621  //clus->GetPosition(pos);
1622  //printf("After Corrections: e %f, x %f, y %f, z %f\n",clus->E(),pos[0],pos[1],pos[2]);
1623  }
1624 
1625  // Recalculate TOF
1626  if(GetCaloUtils()->GetEMCALRecoUtils()->IsTimeRecalibrationOn())
1627  {
1628  Double_t tof = clus->GetTOF();
1629  Float_t frac =-1;
1630  Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(fEMCALCells, clus,frac);
1631 
1633  {
1634  tof = fEMCALCells->GetCellTime(absIdMax);
1635  }
1636 
1637  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absIdMax,fInputEvent->GetBunchCrossNumber(),tof);
1638 
1639  //additional L1 phase shift
1640  if(GetCaloUtils()->GetEMCALRecoUtils()->IsL1PhaseInTimeRecalibrationOn()){
1641  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTimeL1Phase(GetCaloUtils()->GetEMCALGeometry()->GetSuperModuleNumber(absIdMax), fInputEvent->GetBunchCrossNumber(), tof);
1642  }
1643 
1644  clus->SetTOF(tof);
1645 
1646  }// Time recalibration
1647  }
1648 
1649  //Reject clusters with bad channels, close to borders and exotic;
1650  Bool_t goodCluster = GetCaloUtils()->GetEMCALRecoUtils()->IsGoodCluster(clus,
1651  GetCaloUtils()->GetEMCALGeometry(),
1652  GetEMCALCells(),fInputEvent->GetBunchCrossNumber());
1653 
1654  if(!goodCluster)
1655  {
1656  //if( (fDebug > 2 && fMomentum.E() > 0.1) || fDebug > 10 )
1657  AliDebug(2,Form("Bad cluster E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f",
1658  fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta()));
1659 
1660  return;
1661  }
1662 
1663  //Mask all cells in collumns facing ALICE thick material if requested
1664  if(GetCaloUtils()->GetNMaskCellColumns())
1665  {
1666  Int_t absId = -1;
1667  Int_t iSupMod = -1;
1668  Int_t iphi = -1;
1669  Int_t ieta = -1;
1670  Bool_t shared = kFALSE;
1671  GetCaloUtils()->GetEMCALRecoUtils()->GetMaxEnergyCell(GetCaloUtils()->GetEMCALGeometry(), GetEMCALCells(),clus,absId,iSupMod,ieta,iphi,shared);
1672 
1673  AliDebug(2,Form("Masked collumn: cluster E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f",
1674  fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta()));
1675 
1676 
1677  if(GetCaloUtils()->MaskFrameCluster(iSupMod, ieta)) return;
1678  }
1679 
1681  {
1682  if(clus->GetNLabels()==0 || clus->GetLabel() < 0) return;
1683  //else printf("Embedded cluster, %d, n label %d label %d \n",iclus,clus->GetNLabels(),clus->GetLabel());
1684  }
1685 
1686  //Float_t pos[3];
1687  //clus->GetPosition(pos);
1688  //printf("Before Corrections: e %f, x %f, y %f, z %f\n",clus->E(),pos[0],pos[1],pos[2]);
1689 
1690  //Correct non linearity or smear energy
1692  {
1694 
1695  //if( (fDebug > 5 && fMomentum.E() > 0.1) || fDebug > 10 )
1696  AliDebug(5,Form("Correct Non Lin: Old E %3.2f, New E %3.2f",
1697  fMomentum.E(),clus->E()));
1698 
1699  // In case of MC analysis, to match resolution/calibration in real data
1700  // Not needed anymore, just leave for MC studies on systematics
1701  if( GetCaloUtils()->GetEMCALRecoUtils()->IsClusterEnergySmeared() )
1702  {
1703  Float_t rdmEnergy = GetCaloUtils()->GetEMCALRecoUtils()->SmearClusterEnergy(clus);
1704 
1705  //if( (fDebug > 5 && fMomentum.E() > 0.1) || fDebug > 10 )
1706  AliDebug(5,Form("Smear energy: Old E %3.2f, New E %3.2f",clus->E(),rdmEnergy));
1707 
1708  clus->SetE(rdmEnergy);
1709  }
1710  }
1711 
1712  Double_t tof = clus->GetTOF()*1e9;
1713 
1714  Int_t bc = TMath::Nint(tof/50) + 9;
1715  //printf("tof %2.2f, bc+5=%d\n",tof,bc);
1716 
1717  SetEMCalEventBC(bc);
1718 
1719  if(fEMCALPtMin > clus->E() || fEMCALPtMax < clus->E())
1720  {
1721  AliDebug(2,Form("Too low energy cluster E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f",
1722  fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta()));
1723 
1724  return ;
1725  }
1726 
1727  clus->GetMomentum(fMomentum, fVertex[vindex]);
1728 
1729  SetEMCalEventBCcut(bc);
1730 
1731  if(!IsInTimeWindow(tof,clus->E()))
1732  {
1733  fNPileUpClusters++ ;
1734  if(fUseEMCALTimeCut)
1735  {
1736  AliDebug(2,Form("Out of time window E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f, time %e",
1737  fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta(),tof));
1738 
1739  return ;
1740  }
1741  }
1742  else
1744 
1745  if (fMixedEvent)
1746  clus->SetID(iclus) ;
1747 
1748  // Select cluster fiducial region
1749  Bool_t bEMCAL = kFALSE;
1750  Bool_t bDCAL = kFALSE;
1751  if(fCheckFidCut)
1752  {
1753  if(fFillEMCAL && fFiducialCut->IsInFiducialCut(fMomentum.Eta(),fMomentum.Phi(),kEMCAL)) bEMCAL = kTRUE ;
1754  if(fFillDCAL && fFiducialCut->IsInFiducialCut(fMomentum.Eta(),fMomentum.Phi(),kDCAL )) bDCAL = kTRUE ;
1755  }
1756  else
1757  {
1758  bEMCAL = kTRUE;
1759  }
1760 
1761  //if((fDebug > 2 && fMomentum.E() > 0.1) || fDebug > 10)
1762  AliDebug(2,Form("Selected clusters (EMCAL%d, DCAL%d), E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f",
1763  bEMCAL,bDCAL,fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta()));
1764 
1765 
1766  // Smear the SS to try to match data and simulations,
1767  // do it only for simulations.
1768  if( fSmearShowerShape && clus->GetNCells() > 2)
1769  {
1770  AliDebug(2,Form("Smear shower shape - Original: %2.4f\n", clus->GetM02()));
1771  clus->SetM02( clus->GetM02() + fRandom.Landau(0, fSmearShowerShapeWidth) );
1772  //clus->SetM02( fRandom.Landau(clus->GetM02(), fSmearShowerShapeWidth) );
1773  AliDebug(2,Form("Width %2.4f Smeared : %2.4f\n", fSmearShowerShapeWidth,clus->GetM02()));
1774  }
1775 
1776  // Fill the corresponding array. Usually just filling EMCal array with upper or lower clusters is enough, but maybe we want to do EMCal-DCal correlations.
1777  if (bEMCAL) fEMCALClusters->Add(clus);
1778  else if(bDCAL ) fDCALClusters ->Add(clus);
1779 }
1780 
1781 //_______________________________________
1788 //_______________________________________
1790 {
1791  AliDebug(1,"Begin");
1792 
1793  // First recalibrate cells, time or energy
1794  // if(GetCaloUtils()->IsRecalibrationOn())
1795  // GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCells(GetCaloUtils()->GetEMCALGeometry(),
1796  // GetEMCALCells(),
1797  // fInputEvent->GetBunchCrossNumber());
1798 
1799  fNPileUpClusters = 0; // Init counter
1800  fNNonPileUpClusters = 0; // Init counter
1801  for(Int_t i = 0; i < 19; i++)
1802  {
1803  fEMCalBCEvent [i] = 0;
1804  fEMCalBCEventCut[i] = 0;
1805  }
1806 
1807  //Loop to select clusters in fiducial cut and fill container with aodClusters
1808  if(fEMCALClustersListName=="")
1809  {
1810  Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
1811  for (Int_t iclus = 0; iclus < nclusters; iclus++)
1812  {
1813  AliVCluster * clus = 0;
1814  if ( (clus = fInputEvent->GetCaloCluster(iclus)) )
1815  {
1816  if (clus->IsEMCAL())
1817  {
1818  FillInputEMCALAlgorithm(clus, iclus);
1819  }//EMCAL cluster
1820  }// cluster exists
1821  }// cluster loop
1822 
1823  //Recalculate track matching
1825 
1826  }//Get the clusters from the input event
1827  else
1828  {
1829  TClonesArray * clusterList = 0x0;
1830 
1831  if (fInputEvent->FindListObject(fEMCALClustersListName))
1832  {
1833  clusterList = dynamic_cast<TClonesArray*> (fInputEvent->FindListObject(fEMCALClustersListName));
1834  }
1835  else if(fOutputEvent)
1836  {
1837  clusterList = dynamic_cast<TClonesArray*> (fOutputEvent->FindListObject(fEMCALClustersListName));
1838  }
1839 
1840  if(!clusterList)
1841  {
1842  AliWarning(Form("Wrong name of list with clusters? <%s>",fEMCALClustersListName.Data()));
1843  return;
1844  }
1845 
1846  Int_t nclusters = clusterList->GetEntriesFast();
1847  for (Int_t iclus = 0; iclus < nclusters; iclus++)
1848  {
1849  AliVCluster * clus = dynamic_cast<AliVCluster*> (clusterList->At(iclus));
1850  //printf("E %f\n",clus->E());
1851  if (clus) FillInputEMCALAlgorithm(clus, iclus);
1852  else AliWarning("Null cluster in list!");
1853  }// cluster loop
1854 
1855  // Recalculate the pile-up time, in case long time clusters removed during clusterization
1856  //printf("Input event INIT : Pile-up clusters %d, NO pile-up %d\n",fNPileUpClusters,fNNonPileUpClusters);
1857 
1858  fNPileUpClusters = 0; // Init counter
1859  fNNonPileUpClusters = 0; // Init counter
1860  for(Int_t i = 0; i < 19; i++)
1861  {
1862  fEMCalBCEvent [i] = 0;
1863  fEMCalBCEventCut[i] = 0;
1864  }
1865 
1866  for (Int_t iclus = 0; iclus < fInputEvent->GetNumberOfCaloClusters(); iclus++)
1867  {
1868  AliVCluster * clus = 0;
1869 
1870  if ( (clus = fInputEvent->GetCaloCluster(iclus)) )
1871  {
1872  if (clus->IsEMCAL())
1873  {
1874 
1875  Float_t frac =-1;
1876  Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(fEMCALCells, clus,frac);
1877  Double_t tof = clus->GetTOF();
1878  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absIdMax,fInputEvent->GetBunchCrossNumber(),tof);
1879  //additional L1 phase shift
1880  if(GetCaloUtils()->GetEMCALRecoUtils()->IsL1PhaseInTimeRecalibrationOn()){
1881  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTimeL1Phase(GetCaloUtils()->GetEMCALGeometry()->GetSuperModuleNumber(absIdMax), fInputEvent->GetBunchCrossNumber(), tof);
1882  }
1883 
1884  tof*=1e9;
1885 
1886  //printf("Input event cluster : AbsIdMax %d, E %2.2f, time %2.2f \n", absIdMax,clus->E(),tof);
1887 
1888  //Reject clusters with bad channels, close to borders and exotic;
1889  if(!GetCaloUtils()->GetEMCALRecoUtils()->IsGoodCluster(clus,GetCaloUtils()->GetEMCALGeometry(),GetEMCALCells(),fInputEvent->GetBunchCrossNumber())) continue;
1890 
1891  Int_t bc = TMath::Nint(tof/50) + 9;
1892  SetEMCalEventBC(bc);
1893 
1894  if(fEMCALPtMin > clus->E() || fEMCALPtMax < clus->E()) continue ;
1895 
1896  clus->GetMomentum(fMomentum, fVertex[0]);
1897 
1898  if(fCheckFidCut && !fFiducialCut->IsInFiducialCut(fMomentum.Eta(),fMomentum.Phi(),kEMCAL)) return ;
1899 
1900  SetEMCalEventBCcut(bc);
1901 
1902  if(!IsInTimeWindow(tof,clus->E()))
1903  fNPileUpClusters++ ;
1904  else
1906 
1907  }
1908  }
1909  }
1910 
1911  //printf("Input event : Pile-up clusters %d, NO pile-up %d\n",fNPileUpClusters,fNNonPileUpClusters);
1912 
1913  // Recalculate track matching, not necessary if already done in the reclusterization task.
1914  // in case it was not done ...
1916 
1917  }
1918 
1919  AliDebug(1,Form("AOD entries %d, n pile-up clusters %d, n non pile-up %d", fEMCALClusters->GetEntriesFast(),fNPileUpClusters,fNNonPileUpClusters));
1920 }
1921 
1922 //_______________________________________
1924 //_______________________________________
1926 {
1927  AliDebug(1,"Begin");
1928 
1929  // Loop to select clusters in fiducial cut and fill container with aodClusters
1930  Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
1931  for (Int_t iclus = 0; iclus < nclusters; iclus++)
1932  {
1933  AliVCluster * clus = fInputEvent->GetCaloCluster(iclus) ;
1934  if ( !clus ) continue ;
1935 
1936  if ( !clus->IsPHOS() ) continue ;
1937 
1938  // Skip CPV input
1939  if( clus->GetType() == AliVCluster::kPHOSCharged ) continue ;
1940 
1941  if(fAcceptOnlyHIJINGLabels && !IsHIJINGLabel(clus->GetLabel())) continue ;
1942 
1943  // Check if the cluster contains any bad channel and if close to calorimeter borders
1944  if( GetCaloUtils()->ClusterContainsBadChannel(kPHOS,clus->GetCellsAbsId(), clus->GetNCells()))
1945  continue;
1946 
1947  if(!GetCaloUtils()->CheckCellFiducialRegion(clus, fInputEvent->GetPHOSCells()))
1948  continue;
1949 
1951  {
1952  // Recalibrate the cluster energy
1954  {
1955  Float_t energy = GetCaloUtils()->RecalibrateClusterEnergy(clus, (AliAODCaloCells*)GetPHOSCells());
1956  clus->SetE(energy);
1957  }
1958  }
1959 
1960  // Dead code? remove?
1961  Int_t vindex = 0 ;
1962  if (fMixedEvent)
1963  vindex = fMixedEvent->EventIndexForCaloCluster(iclus);
1964 
1965  clus->GetMomentum(fMomentum, fVertex[vindex]);
1966 
1968  continue ;
1969 
1970  if (fPHOSPtMin > fMomentum.E() || fPHOSPtMax < fMomentum.E() )
1971  continue ;
1972 
1973  //if(fDebug > 2 && fMomentum.E() > 0.1)
1974  AliDebug(2,Form("Selected clusters E %3.2f, pt %3.2f, phi %3.2f deg, eta %3.2f",
1975  fMomentum.E(),fMomentum.Pt(),RadToDeg(GetPhi(fMomentum.Phi())),fMomentum.Eta()));
1976 
1977  // Dead code? remove?
1978  if (fMixedEvent)
1979  clus->SetID(iclus) ;
1980 
1981  fPHOSClusters->Add(clus);
1982 
1983  } // esd cluster loop
1984 
1985  AliDebug(1,Form("AOD entries %d",fPHOSClusters->GetEntriesFast()));
1986 
1987 }
1988 
1989 //____________________________________________
1991 //____________________________________________
1993 {
1994  fEMCALCells = fInputEvent->GetEMCALCells();
1995 }
1996 
1997 //___________________________________________
1999 //___________________________________________
2001 {
2002  fPHOSCells = fInputEvent->GetPHOSCells();
2003 }
2004 
2005 //_______________________________________
2008 //_______________________________________
2010 {
2011  AliVVZERO* v0 = fInputEvent->GetVZEROData();
2012  //printf("Init V0: ADC (%d,%d), Multiplicity (%d,%d) \n",fV0ADC[0],fV0ADC[1],fV0Mul[0],fV0Mul[1]);
2013 
2014  if (v0)
2015  {
2016  AliESDVZERO* esdV0 = dynamic_cast<AliESDVZERO*> (v0);
2017  for (Int_t i = 0; i < 32; i++)
2018  {
2019  if(esdV0)
2020  {//Only available in ESDs
2021  fV0ADC[0] += (Int_t)esdV0->GetAdcV0C(i);
2022  fV0ADC[1] += (Int_t)esdV0->GetAdcV0A(i);
2023  }
2024 
2025  fV0Mul[0] += (Int_t)v0->GetMultiplicityV0C(i);
2026  fV0Mul[1] += (Int_t)v0->GetMultiplicityV0A(i);
2027  }
2028 
2029  AliDebug(1,Form("ADC (%d,%d), Multiplicity (%d,%d)",fV0ADC[0],fV0ADC[1],fV0Mul[0],fV0Mul[1]));
2030  }
2031  else
2032  {
2033  AliDebug(1,"Cannot retrieve V0 ESD! Run w/ null V0 charges");
2034  }
2035 }
2036 
2037 //_________________________________________________
2041 //_________________________________________________
2043 {
2044  AliDebug(2,"Begin");
2045 
2046  //
2047  //check if branch name is given
2048  if(!fInputNonStandardJetBranchName.Length())
2049  {
2050  fInputEvent->Print();
2051  AliFatal("No non-standard jet branch name specified. Specify among existing ones.");
2052  }
2053 
2054  fNonStandardJets = dynamic_cast<TClonesArray*>(fInputEvent->FindListObject(fInputNonStandardJetBranchName.Data()));
2055 
2056  if(!fNonStandardJets)
2057  {
2058  //check if jet branch exist; exit if not
2059  fInputEvent->Print();
2060 
2061  AliFatal(Form("%s:%d no reconstructed jet array with name %s in AOD", (char*)__FILE__,__LINE__,fInputNonStandardJetBranchName.Data()));
2062  }
2063  else
2064  {
2065  AliDebug(1,Form("AOD input jets %d", fNonStandardJets->GetEntriesFast()));
2066  }
2067 }
2068 
2069 //_________________________________________________
2073 //_________________________________________________
2075 {
2076  AliDebug(1,"Begin");
2077  //
2078  //check if branch name is given
2079  if(!fInputBackgroundJetBranchName.Length())
2080  {
2081  fInputEvent->Print();
2082 
2083  AliFatal("No background jet branch name specified. Specify among existing ones.");
2084  }
2085 
2086  fBackgroundJets = (AliAODJetEventBackground*)(fInputEvent->FindListObject(fInputBackgroundJetBranchName.Data()));
2087 
2088  if(!fBackgroundJets)
2089  {
2090  //check if jet branch exist; exit if not
2091  fInputEvent->Print();
2092 
2093  AliFatal(Form("%s:%d no reconstructed jet array with name %s in AOD", (char*)__FILE__,__LINE__,fInputBackgroundJetBranchName.Data()));
2094  }
2095  else
2096  {
2097  AliDebug(1,"FillInputBackgroundJets");
2098  fBackgroundJets->Print("");
2099  }
2100 }
2101 
2102 //____________________________________________________________________
2108 //____________________________________________________________________
2109 TArrayI AliCaloTrackReader::GetTriggerPatches(Int_t tmin, Int_t tmax )
2110 {
2111  // init some variables
2112  Int_t trigtimes[30], globCol, globRow,ntimes, i;
2113  Int_t absId = -1; //[100];
2114  Int_t nPatch = 0;
2115 
2116  TArrayI patches(0);
2117 
2118  // get object pointer
2119  AliVCaloTrigger *caloTrigger = GetInputEvent()->GetCaloTrigger( "EMCAL" );
2120 
2121  if(!caloTrigger)
2122  {
2123  AliError("Trigger patches input (AliVCaloTrigger) not available in data!");
2124  return patches;
2125  }
2126 
2127  //printf("CaloTrigger Entries %d\n",caloTrigger->GetEntries() );
2128 
2129  // class is not empty
2130  if( caloTrigger->GetEntries() > 0 )
2131  {
2132  // must reset before usage, or the class will fail
2133  caloTrigger->Reset();
2134 
2135  // go throuth the trigger channels
2136  while( caloTrigger->Next() )
2137  {
2138  // get position in global 2x2 tower coordinates
2139  caloTrigger->GetPosition( globCol, globRow );
2140 
2141  //L0
2142  if(IsEventEMCALL0())
2143  {
2144  // get dimension of time arrays
2145  caloTrigger->GetNL0Times( ntimes );
2146 
2147  // no L0s in this channel
2148  // presence of the channel in the iterator still does not guarantee that L0 was produced!!
2149  if( ntimes < 1 )
2150  continue;
2151 
2152  // get timing array
2153  caloTrigger->GetL0Times( trigtimes );
2154  //printf("Get L0 patch : n times %d - trigger time window %d - %d\n",ntimes, tmin,tmax);
2155 
2156  // go through the array
2157  for( i = 0; i < ntimes; i++ )
2158  {
2159  // check if in cut - 8,9 shall be accepted in 2011
2160  if( trigtimes[i] >= tmin && trigtimes[i] <= tmax )
2161  {
2162  //printf("Accepted trigger time %d \n",trigtimes[i]);
2163  //if(nTrig > 99) continue;
2164  GetCaloUtils()->GetEMCALGeometry()->GetAbsFastORIndexFromPositionInEMCAL(globCol,globRow, absId);
2165  //printf("pass the time cut globCol %d, globRow %d absId %d\n",globCol,globRow, absId);
2166  patches.Set(nPatch+1);
2167  patches.AddAt(absId,nPatch++);
2168  }
2169  } // trigger time array
2170  }//L0
2171  else if(IsEventEMCALL1()) // L1
2172  {
2173  Int_t bit = 0;
2174  caloTrigger->GetTriggerBits(bit);
2175 
2176  Int_t sum = 0;
2177  caloTrigger->GetL1TimeSum(sum);
2178  //fBitEGA-=2;
2179  Bool_t isEGA1 = ((bit >> fBitEGA ) & 0x1) && IsEventEMCALL1Gamma1() ;
2180  Bool_t isEGA2 = ((bit >> (fBitEGA+1)) & 0x1) && IsEventEMCALL1Gamma2() ;
2181  Bool_t isEJE1 = ((bit >> fBitEJE ) & 0x1) && IsEventEMCALL1Jet1 () ;
2182  Bool_t isEJE2 = ((bit >> (fBitEJE+1)) & 0x1) && IsEventEMCALL1Jet2 () ;
2183 
2184  //if((bit>> fBitEGA )&0x1) printf("Trig Bit %d - bit %d - EG1 %d - EG2 %d\n",fBitEGA ,bit,IsEventEMCALL1Gamma1(),IsEventEMCALL1Gamma2());
2185  //if((bit>>(fBitEGA+1))&0x1) printf("Trig Bit %d - bit %d - EG1 %d - EG2 %d\n",fBitEGA+1,bit,IsEventEMCALL1Gamma1(),IsEventEMCALL1Gamma2());
2186 
2187  if(!isEGA1 && !isEJE1 && !isEGA2 && !isEJE2) continue;
2188 
2189  Int_t patchsize = -1;
2190  if (isEGA1 || isEGA2) patchsize = 2;
2191  else if (isEJE1 || isEJE2) patchsize = 16;
2192 
2193  //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",
2194  // bit,sum,patchsize,isEGA1,isEGA2,isEJE1,isEJE2,fBitEGA,fBitEJE,IsEventEMCALL1Gamma(),IsEventEMCALL1Jet());
2195 
2196 
2197  // add 2x2 (EGA) or 16x16 (EJE) patches
2198  for(Int_t irow=0; irow < patchsize; irow++)
2199  {
2200  for(Int_t icol=0; icol < patchsize; icol++)
2201  {
2202  GetCaloUtils()->GetEMCALGeometry()->GetAbsFastORIndexFromPositionInEMCAL(globCol+icol,globRow+irow, absId);
2203  //printf("pass the time cut globCol %d, globRow %d absId %d\n",globCol,globRow, absId);
2204  patches.Set(nPatch+1);
2205  patches.AddAt(absId,nPatch++);
2206  }
2207  }
2208 
2209  } // L1
2210 
2211  } // trigger iterator
2212  } // go through triggers
2213 
2214  if(patches.GetSize()<=0) AliInfo(Form("No patch found! for triggers: %s and selected <%s>",
2216  //else printf(">>>>> N patches %d, test %d,first %d, last %d\n",patches.GetSize(), nPatch, patches.At(0), patches.At(patches.GetSize()-1));
2217 
2218  return patches;
2219 }
2220 
2221 //____________________________________________________________
2225 //____________________________________________________________
2227 {
2228  // Init info from previous event
2229  fTriggerClusterIndex = -1;
2230  fTriggerClusterId = -1;
2231  fTriggerClusterBC = -10000;
2232  fIsExoticEvent = kFALSE;
2233  fIsBadCellEvent = kFALSE;
2234  fIsBadMaxCellEvent = kFALSE;
2235 
2236  fIsTriggerMatch = kFALSE;
2237  fIsTriggerMatchOpenCut[0] = kFALSE;
2238  fIsTriggerMatchOpenCut[1] = kFALSE;
2239  fIsTriggerMatchOpenCut[2] = kFALSE;
2240 
2241  // Do only analysis for triggered events
2242  if(!IsEventEMCALL1() && !IsEventEMCALL0())
2243  {
2244  fTriggerClusterBC = 0;
2245  return;
2246  }
2247 
2248  //printf("***** Try to match trigger to cluster %d **** L0 %d, L1 %d\n",fTriggerPatchClusterMatch,IsEventEMCALL0(),IsEventEMCALL1());
2249 
2250  //Recover the list of clusters
2251  TClonesArray * clusterList = 0;
2252  if (fInputEvent->FindListObject(fEMCALClustersListName))
2253  {
2254  clusterList = dynamic_cast<TClonesArray*> (fInputEvent->FindListObject(fEMCALClustersListName));
2255  }
2256  else if(fOutputEvent)
2257  {
2258  clusterList = dynamic_cast<TClonesArray*> (fOutputEvent->FindListObject(fEMCALClustersListName));
2259  }
2260 
2261  // Get number of clusters and of trigger patches
2262  Int_t nclusters = fInputEvent->GetNumberOfCaloClusters();
2263  if(clusterList)
2264  nclusters = clusterList->GetEntriesFast();
2265 
2266  Int_t nPatch = patches.GetSize();
2267  Float_t exoDiffTime = GetCaloUtils()->GetEMCALRecoUtils()->GetExoticCellDiffTimeCut();
2268 
2269  //Init some variables used in the cluster loop
2270  Float_t tofPatchMax = 100000;
2271  Float_t ePatchMax =-1;
2272 
2273  Float_t tofMax = 100000;
2274  Float_t eMax =-1;
2275 
2276  Int_t clusMax =-1;
2277  Int_t idclusMax =-1;
2278  Bool_t badClMax = kFALSE;
2279  Bool_t badCeMax = kFALSE;
2280  Bool_t exoMax = kFALSE;
2281 // Int_t absIdMaxTrig= -1;
2282  Int_t absIdMaxMax = -1;
2283 
2284  Int_t nOfHighECl = 0 ;
2285 
2286  //
2287  // Check what is the trigger threshold
2288  // set minimu energym of candidate for trigger cluster
2289  //
2291 
2292  Float_t triggerThreshold = fTriggerL1EventThreshold;
2293  if(IsEventEMCALL0()) triggerThreshold = fTriggerL0EventThreshold;
2294  Float_t minE = triggerThreshold / 2.;
2295 
2296  // This method is not really suitable for JET trigger
2297  // but in case, reduce the energy cut since we do not trigger on high energy particle
2298  if(IsEventEMCALL1Jet() || minE < 1) minE = 1;
2299 
2300  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));
2301 
2302  //
2303  // Loop on the clusters, check if there is any that falls into one of the patches
2304  //
2305  for (Int_t iclus = 0; iclus < nclusters; iclus++)
2306  {
2307  AliVCluster * clus = 0;
2308  if(clusterList) clus = (AliVCluster*) clusterList->At(iclus);
2309  else clus = fInputEvent->GetCaloCluster(iclus);
2310 
2311  if ( !clus ) continue ;
2312 
2313  if ( !clus->IsEMCAL() ) continue ;
2314 
2315  //Skip clusters with too low energy to be triggering
2316  if ( clus->E() < minE ) continue ;
2317 
2318  Float_t frac = -1;
2319  Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(fInputEvent->GetEMCALCells(), clus,frac);
2320 
2321  Bool_t badCluster = GetCaloUtils()->GetEMCALRecoUtils()->ClusterContainsBadChannel(GetCaloUtils()->GetEMCALGeometry(),
2322  clus->GetCellsAbsId(),clus->GetNCells());
2323  UShort_t cellMax[] = {(UShort_t) absIdMax};
2324  Bool_t badCell = GetCaloUtils()->GetEMCALRecoUtils()->ClusterContainsBadChannel(GetCaloUtils()->GetEMCALGeometry(),cellMax,1);
2325 
2326  // if cell is bad, it can happen that time calibration is not available,
2327  // when calculating if it is exotic, this can make it to be exotic by default
2328  // open it temporarily for this cluster
2329  if(badCell)
2330  GetCaloUtils()->GetEMCALRecoUtils()->SetExoticCellDiffTimeCut(10000000);
2331 
2332  Bool_t exotic = GetCaloUtils()->GetEMCALRecoUtils()->IsExoticCluster(clus, fInputEvent->GetEMCALCells());
2333 
2334  // Set back the time cut on exotics
2335  if(badCell)
2336  GetCaloUtils()->GetEMCALRecoUtils()->SetExoticCellDiffTimeCut(exoDiffTime);
2337 
2338  // Energy threshold for exotic Ecross typically at 4 GeV,
2339  // for lower energy, check that there are more than 1 cell in the cluster
2340  if(!exotic && clus->GetNCells() < 2) exotic = kTRUE;
2341 
2342  Float_t energy = clus->E();
2343  Int_t idclus = clus->GetID();
2344 
2345  Double_t tof = clus->GetTOF();
2346  if(GetCaloUtils()->GetEMCALRecoUtils()->IsTimeRecalibrationOn() && fTriggerClusterTimeRecal){
2347  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTime(absIdMax,fInputEvent->GetBunchCrossNumber(),tof);
2348  //additional L1 phase shift
2349  if(GetCaloUtils()->GetEMCALRecoUtils()->IsL1PhaseInTimeRecalibrationOn()){
2350  GetCaloUtils()->GetEMCALRecoUtils()->RecalibrateCellTimeL1Phase(GetCaloUtils()->GetEMCALGeometry()->GetSuperModuleNumber(absIdMax), fInputEvent->GetBunchCrossNumber(), tof);
2351  }
2352  }
2353  tof *=1.e9;
2354 
2355  //printf("cluster %d, ID %d, E %2.2f, tof %2.2f, AbsId max %d, exotic %d, bad Cluster %d, bad Cell %d\n",
2356  // iclus,idclus, energy,tof,absIdMax, exotic, badCluster,badCell);
2357 
2358  // Find the highest energy cluster, avobe trigger threshold
2359  // in the event in case no match to trigger is found
2360  if( energy > eMax )
2361  {
2362  tofMax = tof;
2363  eMax = energy;
2364  badClMax = badCluster;
2365  badCeMax = badCell;
2366  exoMax = exotic;
2367  clusMax = iclus;
2368  idclusMax = idclus;
2369  absIdMaxMax = absIdMax;
2370  }
2371 
2372  // count the good clusters in the event avobe the trigger threshold
2373  // to check the exotic events
2374  if(!badCluster && !exotic)
2375  nOfHighECl++;
2376 
2377  // Find match to trigger
2378  if(fTriggerPatchClusterMatch && nPatch>0)
2379  {
2380  for(Int_t iabsId =0; iabsId < nPatch; iabsId++)
2381  {
2382  Int_t absIDCell[4];
2383  GetCaloUtils()->GetEMCALGeometry()->GetCellIndexFromFastORIndex(patches.At(iabsId), absIDCell);
2384  //if(tof > 75 ) printf("E %2.2f TOF %2.2f Trigger patch %d, cells : %d, %d, %d, %d\n",
2385  // clus->E(),tof,patches.At(iabsId), absIDCell[0],absIDCell[1],absIDCell[2],absIDCell[3]);
2386 
2387  for(Int_t ipatch = 0; ipatch < 4; ipatch++)
2388  {
2389  if(absIdMax == absIDCell[ipatch])
2390  {
2391  //printf("*** Patches : absId %d, E %2.1f, tof %f \n",absIdMax,clus->E(), tof);
2392  if(energy > ePatchMax)
2393  {
2394  tofPatchMax = tof;
2395  ePatchMax = energy;
2396  fIsBadCellEvent = badCluster;
2397  fIsBadMaxCellEvent = badCell;
2398  fIsExoticEvent = exotic;
2399  fTriggerClusterIndex = iclus;
2400  fTriggerClusterId = idclus;
2401  fIsTriggerMatch = kTRUE;
2402 // absIdMaxTrig = absIdMax;
2403  }
2404  }
2405  }// cell patch loop
2406  }// trigger patch loop
2407  } // Do trigger patch matching
2408 
2409  }// Cluster loop
2410 
2411  // If there was no match, assign as trigger
2412  // the highest energy cluster in the event
2413  if(!fIsTriggerMatch)
2414  {
2415  tofPatchMax = tofMax;
2416  ePatchMax = eMax;
2417  fIsBadCellEvent = badClMax;
2418  fIsBadMaxCellEvent = badCeMax;
2419  fIsExoticEvent = exoMax;
2420  fTriggerClusterIndex = clusMax;
2421  fTriggerClusterId = idclusMax;
2422  }
2423 
2424  Double_t tofPatchMaxUS = TMath::Abs(tofPatchMax);
2425 
2426  if (tofPatchMaxUS < 28 ) fTriggerClusterBC = 0 ;
2427  else if(tofPatchMaxUS < 75 ) fTriggerClusterBC = 1 ;
2428  else if(tofPatchMaxUS < 125) fTriggerClusterBC = 2 ;
2429  else if(tofPatchMaxUS < 175) fTriggerClusterBC = 3 ;
2430  else if(tofPatchMaxUS < 225) fTriggerClusterBC = 4 ;
2431  else if(tofPatchMaxUS < 275) fTriggerClusterBC = 5 ;
2432  else
2433  {
2434  //printf("AliCaloTrackReader::MatchTriggerCluster() - Large BC - tof %2.3f - Index %d\n",tofPatchMaxUS,fTriggerClusterIndex);
2436  else
2437  {
2438  fTriggerClusterIndex = -2;
2439  fTriggerClusterId = -2;
2440  }
2441  }
2442 
2443  if(tofPatchMax < 0) fTriggerClusterBC*=-1;
2444 
2445 
2446  // 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",
2447  // fTriggerClusterIndex, fTriggerClusterId,ePatchMax, tofPatchMax,
2448  // fTriggerClusterBC, fIsBadCellEvent,fIsBadMaxCellEvent,fIsExoticEvent, fIsTriggerMatch, nOfHighECl,absIdMaxMax);
2449  //
2450  // 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",
2451  // clusMax, idclusMax, eMax,tofMax, badClMax, badCeMax,exoMax);
2452 
2453  //Redo matching but open cuts
2454  if(!fIsTriggerMatch && fTriggerClusterId >= 0)
2455  {
2456  // Open time patch time
2457  TArrayI patchOpen = GetTriggerPatches(7,10);
2458 
2459  Int_t patchAbsIdOpenTime = -1;
2460  for(Int_t iabsId =0; iabsId < patchOpen.GetSize(); iabsId++)
2461  {
2462  Int_t absIDCell[4];
2463  patchAbsIdOpenTime = patchOpen.At(iabsId);
2464  GetCaloUtils()->GetEMCALGeometry()->GetCellIndexFromFastORIndex(patchAbsIdOpenTime, absIDCell);
2465  //if(tof > 75 ) printf("E %2.2f TOF %2.2f Trigger patch %d, cells : %d, %d, %d, %d\n",
2466  // clus->E(),tof,patches.At(iabsId), absIDCell[0],absIDCell[1],absIDCell[2],absIDCell[3]);
2467 
2468  for(Int_t ipatch = 0; ipatch < 4; ipatch++)
2469  {
2470  if(absIdMaxMax == absIDCell[ipatch])
2471  {
2472  fIsTriggerMatchOpenCut[0] = kTRUE;
2473  break;
2474  }
2475  }// cell patch loop
2476  }// trigger patch loop
2477 
2478  // Check neighbour patches
2479  Int_t patchAbsId = -1;
2480  Int_t globalCol = -1;
2481  Int_t globalRow = -1;
2482  GetCaloUtils()->GetEMCALGeometry()->GetFastORIndexFromCellIndex(absIdMaxMax, patchAbsId);
2483  GetCaloUtils()->GetEMCALGeometry()->GetPositionInEMCALFromAbsFastORIndex(patchAbsId,globalCol,globalRow);
2484 
2485  // Check patches with strict time cut
2486  Int_t patchAbsIdNeigh = -1;
2487  for(Int_t icol = globalCol-1; icol <= globalCol+1; icol++)
2488  {
2489  if(icol < 0 || icol > 47) continue;
2490 
2491  for(Int_t irow = globalRow; irow <= globalRow+1; irow++)
2492  {
2493  if(irow < 0 || irow > 63) continue;
2494 
2495  GetCaloUtils()->GetEMCALGeometry()->GetAbsFastORIndexFromPositionInEMCAL(icol, irow, patchAbsIdNeigh);
2496 
2497  if ( patchAbsIdNeigh < 0 ) continue;
2498 
2499  for(Int_t iabsId =0; iabsId < patches.GetSize(); iabsId++)
2500  {
2501  if(patchAbsIdNeigh == patches.At(iabsId))
2502  {
2503  fIsTriggerMatchOpenCut[1] = kTRUE;
2504  break;
2505  }
2506  }// trigger patch loop
2507 
2508  }// row
2509  }// col
2510 
2511  // Check patches with open time cut
2512  Int_t patchAbsIdNeighOpenTime = -1;
2513  for(Int_t icol = globalCol-1; icol <= globalCol+1; icol++)
2514  {
2515  if(icol < 0 || icol > 47) continue;
2516 
2517  for(Int_t irow = globalRow; irow <= globalRow+1; irow++)
2518  {
2519  if(irow < 0 || irow > 63) continue;
2520 
2521  GetCaloUtils()->GetEMCALGeometry()->GetAbsFastORIndexFromPositionInEMCAL(icol, irow, patchAbsIdNeighOpenTime);
2522 
2523  if ( patchAbsIdNeighOpenTime < 0 ) continue;
2524 
2525  for(Int_t iabsId =0; iabsId < patchOpen.GetSize(); iabsId++)
2526  {
2527  if(patchAbsIdNeighOpenTime == patchOpen.At(iabsId))
2528  {
2529  fIsTriggerMatchOpenCut[2] = kTRUE;
2530  break;
2531  }
2532  }// trigger patch loop
2533 
2534  }// row
2535  }// col
2536 
2537  // printf("No match, new match: Open time %d-%d, open Neigh %d-%d, both open %d-%d\n",fIsTriggerMatchOpenCut[0],patchAbsIdOpenTime,
2538  // fIsTriggerMatchOpenCut[1],patchAbsIdNeigh,
2539  // fIsTriggerMatchOpenCut[2],patchAbsIdNeighOpenTime);
2540 
2541  patchOpen.Reset();
2542 
2543  }// No trigger match found
2544  //printf("Trigger BC %d, Id %d, Index %d\n",fTriggerClusterBC,fTriggerClusterId,fTriggerClusterIndex);
2545 }
2546 
2547 //_________________________________________________________
2551 //_________________________________________________________
2553 {
2555  {
2556  // get object pointer
2557  AliVCaloTrigger *caloTrigger = GetInputEvent()->GetCaloTrigger( "EMCAL" );
2558 
2559  if ( fBitEGA == 6 )
2560  {
2561  if (IsEventEMCALL1Gamma1()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(1);
2562  else if(IsEventEMCALL1Gamma2()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(3);
2563  else if(IsEventEMCALL1Jet1 ()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(0);
2564  else if(IsEventEMCALL1Jet2 ()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(2);
2565 
2566  // printf("L1 trigger Threshold Jet1 %f, Gamma1 %f, Jet2 %f, Gamma2 %f\n",
2567  // 0.07874*caloTrigger->GetL1Threshold(0),
2568  // 0.07874*caloTrigger->GetL1Threshold(1),
2569  // 0.07874*caloTrigger->GetL1Threshold(2),
2570  // 0.07874*caloTrigger->GetL1Threshold(3));
2571  }
2572  else
2573  {
2574  // Old AOD data format, in such case, order of thresholds different!!!
2575  if (IsEventEMCALL1Gamma1()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(0);
2576  else if(IsEventEMCALL1Gamma2()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(2);
2577  else if(IsEventEMCALL1Jet1 ()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(1);
2578  else if(IsEventEMCALL1Jet2 ()) fTriggerL1EventThreshold = 0.07874*caloTrigger->GetL1Threshold(3);
2579 
2580  // printf("L1 trigger Threshold Jet1 %f, Gamma1 %f, Jet2 %f, Gamma2 %f\n",
2581  // 0.07874*caloTrigger->GetL1Threshold(1),
2582  // 0.07874*caloTrigger->GetL1Threshold(0),
2583  // 0.07874*caloTrigger->GetL1Threshold(3),
2584  // 0.07874*caloTrigger->GetL1Threshold(2));
2585  }
2586  }
2587 
2588  // Set L0 threshold, if not set by user
2590  {
2591  // Revise for periods > LHC11d
2592  Int_t runNumber = fInputEvent->GetRunNumber();
2593  if (runNumber < 146861) fTriggerL0EventThreshold = 3. ;
2594  else if(runNumber < 154000) fTriggerL0EventThreshold = 4. ;
2595  else if(runNumber < 165000) fTriggerL0EventThreshold = 5.5;
2596  else fTriggerL0EventThreshold = 2 ;
2597  }
2598 }
2599 
2600 //________________________________________________________
2602 //________________________________________________________
2603 void AliCaloTrackReader::Print(const Option_t * opt) const
2604 {
2605  if(! opt)
2606  return;
2607 
2608  printf("***** Print: %s %s ******\n", GetName(), GetTitle() ) ;
2609  printf("Task name : %s\n", fTaskName.Data()) ;
2610  printf("Data type : %d\n", fDataType) ;
2611  printf("CTS Min pT : %2.1f GeV/c\n", fCTSPtMin) ;
2612  printf("EMCAL Min pT : %2.1f GeV/c\n", fEMCALPtMin) ;
2613  printf("PHOS Min pT : %2.1f GeV/c\n", fPHOSPtMin) ;
2614  printf("CTS Max pT : %2.1f GeV/c\n", fCTSPtMax) ;
2615  printf("EMCAL Max pT : %2.1f GeV/c\n", fEMCALPtMax) ;
2616  printf("PHOS Max pT : %2.1f GeV/c\n", fPHOSPtMax) ;
2617  printf("EMCAL Time Cut: %3.1f < TOF < %3.1f\n", fEMCALTimeCutMin, fEMCALTimeCutMax);
2618  printf("Use CTS = %d\n", fFillCTS) ;
2619  printf("Use EMCAL = %d\n", fFillEMCAL) ;
2620  printf("Use DCAL = %d\n", fFillDCAL) ;
2621  printf("Use PHOS = %d\n", fFillPHOS) ;
2622  printf("Use EMCAL Cells = %d\n", fFillEMCALCells) ;
2623  printf("Use PHOS Cells = %d\n", fFillPHOSCells) ;
2624  printf("Track status = %d\n", (Int_t) fTrackStatus) ;
2625  //printf("AODs Track filter mask = %d or hybrid %d (if filter bit comp %d), select : SPD hit %d, primary %d\n",
2626  // (Int_t) fTrackFilterMask, fSelectHybridTracks, (Int_t) fTrackFilterMaskComplementary, fSelectSPDHitTracks,fSelectPrimaryTracks) ;
2627  printf("Track Mult Eta Cut = %d\n", (Int_t) fTrackMultEtaCut) ;
2628  printf("Write delta AOD = %d\n", fWriteOutputDeltaAOD) ;
2629  printf("Recalculate Clusters = %d, E linearity = %d\n", fRecalculateClusters, fCorrectELinearity) ;
2630 
2631  printf("Use Triggers selected in SE base class %d; If not what Trigger Mask? %d; MB Trigger Mask for mixed %d \n",
2633 
2635  printf("Compare jet pt and pt hard to accept event, factor = %2.2f",fPtHardAndJetPtFactor);
2636 
2638  printf("Compare cluster pt and pt hard to accept event, factor = %2.2f",fPtHardAndClusterPtFactor);
2639 
2640  printf("Read Kine from, stack? %d, AOD ? %d \n", fReadStack, fReadAODMCParticles) ;
2641  printf("Delta AOD File Name = %s\n", fDeltaAODFileName.Data()) ;
2642  printf("Centrality: Class %s, Option %d, Bin [%d,%d] \n", fCentralityClass.Data(),fCentralityOpt,fCentralityBin[0], fCentralityBin[1]) ;
2643 
2644  printf(" \n") ;
2645 }
2646 
2647 //__________________________________________
2651 //__________________________________________
2653 {
2654  // Count number of cells with energy larger than 0.1 in SM3, cut on this number
2655  Int_t ncellsSM3 = 0;
2656  for(Int_t icell = 0; icell < fInputEvent->GetEMCALCells()->GetNumberOfCells(); icell++)
2657  {
2658  Int_t absID = fInputEvent->GetEMCALCells()->GetCellNumber(icell);
2659  Int_t sm = GetCaloUtils()->GetEMCALGeometry()->GetSuperModuleNumber(absID);
2660  if(fInputEvent->GetEMCALCells()->GetAmplitude(icell) > 0.1 && sm==3) ncellsSM3++;
2661  }
2662 
2663  Int_t ncellcut = 21;
2664  if(GetFiredTriggerClasses().Contains("EMC")) ncellcut = 35;
2665 
2666  if(ncellsSM3 >= ncellcut)
2667  {
2668  AliDebug(1,Form("Reject event with ncells in SM3 %d, cut %d, trig %s",
2669  ncellsSM3,ncellcut,GetFiredTriggerClasses().Data()));
2670  return kTRUE;
2671  }
2672 
2673  return kFALSE;
2674 }
2675 
2676 //_________________________________________________________
2680 //_________________________________________________________
2682 {
2683  if(label < 0) return ;
2684 
2685  AliAODEvent * evt = dynamic_cast<AliAODEvent*> (fInputEvent) ;
2686  if(!evt) return ;
2687 
2688  TClonesArray * arr = dynamic_cast<TClonesArray*>(evt->FindListObject("mcparticles")) ;
2689  if(!arr) return ;
2690 
2691  if(label < arr->GetEntriesFast())
2692  {
2693  AliAODMCParticle * particle = dynamic_cast<AliAODMCParticle *>(arr->At(label));
2694  if(!particle) return ;
2695 
2696  if(label == particle->Label()) return ; // label already OK
2697  //else printf("AliCaloTrackReader::RemapMCLabelForAODs() - Label %d - AOD stack %d \n",label, particle->Label());
2698  }
2699  //else printf("AliCaloTrackReader::RemapMCLabelForAODs() - Label %d > AOD labels %d \n",label, arr->GetEntriesFast());
2700 
2701  // loop on the particles list and check if there is one with the same label
2702  for(Int_t ind = 0; ind < arr->GetEntriesFast(); ind++ )
2703  {
2704  AliAODMCParticle * particle = dynamic_cast<AliAODMCParticle *>(arr->At(ind));
2705  if(!particle) continue ;
2706 
2707  if(label == particle->Label())
2708  {
2709  label = ind;
2710  //printf("AliAnalysisTaskEMCALClusterize::RemapMCLabelForAODs() - New Label Index %d \n",label);
2711  return;
2712  }
2713  }
2714 
2715  label = -1;
2716 
2717  //printf("AliCaloTrackReader::RemapMCLabelForAODs() - Label not found set to -1 \n");
2718 }
2719 
2720 //___________________________________
2722 //___________________________________
2724 {
2725  if(fCTSTracks) fCTSTracks -> Clear();
2726  if(fEMCALClusters) fEMCALClusters -> Clear("C");
2727  if(fPHOSClusters) fPHOSClusters -> Clear("C");
2728 
2729  fV0ADC[0] = 0; fV0ADC[1] = 0;
2730  fV0Mul[0] = 0; fV0Mul[1] = 0;
2731 
2732  if(fNonStandardJets) fNonStandardJets -> Clear("C");
2733  fBackgroundJets->Reset();
2734 }
2735 
2736 //___________________________________________
2740 //___________________________________________
2742 {
2743  fEventTrigMinBias = kFALSE;
2744  fEventTrigCentral = kFALSE;
2745  fEventTrigSemiCentral = kFALSE;
2746  fEventTrigEMCALL0 = kFALSE;
2747  fEventTrigEMCALL1Gamma1 = kFALSE;
2748  fEventTrigEMCALL1Gamma2 = kFALSE;
2749  fEventTrigEMCALL1Jet1 = kFALSE;
2750  fEventTrigEMCALL1Jet2 = kFALSE;
2751 
2752  AliDebug(1,Form("Select trigger mask bit %d - Trigger Event %s - Select <%s>",
2754 
2755  if(fEventTriggerMask <=0 )// in case no mask set
2756  {
2757  // EMC triggered event? Which type?
2758  if( GetFiredTriggerClasses().Contains("-B-") || GetFiredTriggerClasses().Contains("-S-") || GetFiredTriggerClasses().Contains("-I-") )
2759  {
2760  if ( GetFiredTriggerClasses().Contains("EGA" ) ||
2761  GetFiredTriggerClasses().Contains("EG1" ) )
2762  {
2763  fEventTrigEMCALL1Gamma1 = kTRUE;
2764  if( GetFiredTriggerClasses().Contains("EG1" ) && !fFiredTriggerClassName.Contains("EG1") ) fEventTrigEMCALL1Gamma1 = kFALSE;
2765  }
2766  else if( GetFiredTriggerClasses().Contains("EG2" ) )
2767  {
2768  fEventTrigEMCALL1Gamma2 = kTRUE;
2769  if( !fFiredTriggerClassName.Contains("EG2") ) fEventTrigEMCALL1Gamma2 = kFALSE;
2770  }
2771  else if( GetFiredTriggerClasses().Contains("EJE" ) ||
2772  GetFiredTriggerClasses().Contains("EJ1" ) )
2773  {
2774  fEventTrigEMCALL1Jet1 = kTRUE;
2775  if( GetFiredTriggerClasses().Contains("EJ1" ) && !fFiredTriggerClassName.Contains("EJ1") )
2776  fEventTrigEMCALL1Jet1 = kFALSE;
2777  }
2778  else if( GetFiredTriggerClasses().Contains("EJ2" ) )
2779  {
2780  fEventTrigEMCALL1Jet2 = kTRUE;
2781  if( !fFiredTriggerClassName.Contains("EJ2") ) fEventTrigEMCALL1Jet2 = kFALSE;
2782  }
2783  else if( GetFiredTriggerClasses().Contains("CEMC") &&
2784  !GetFiredTriggerClasses().Contains("EGA" ) &&
2785  !GetFiredTriggerClasses().Contains("EJE" ) &&
2786  !GetFiredTriggerClasses().Contains("EG1" ) &&
2787  !GetFiredTriggerClasses().Contains("EJ1" ) &&
2788  !GetFiredTriggerClasses().Contains("EG2" ) &&
2789  !GetFiredTriggerClasses().Contains("EJ2" ) )
2790  {
2791  fEventTrigEMCALL0 = kTRUE;
2792  }
2793 
2794  //Min bias event trigger?
2795  if (GetFiredTriggerClasses().Contains("CCENT_R2-B-NOPF-ALLNOTRD"))
2796  {
2797  fEventTrigCentral = kTRUE;
2798  }
2799  else if(GetFiredTriggerClasses().Contains("CSEMI_R1-B-NOPF-ALLNOTRD"))
2800  {
2801  fEventTrigSemiCentral = kTRUE;
2802  }
2803  else if((GetFiredTriggerClasses().Contains("CINT") || GetFiredTriggerClasses().Contains("CPBI2_B1") ) &&
2804  GetFiredTriggerClasses().Contains("-NOPF-ALLNOTRD") )
2805  {
2806  fEventTrigMinBias = kTRUE;
2807  }
2808  }
2809  }
2810  else
2811  {
2812  // EMC L1 Gamma
2813  if ( fEventTriggerMask & AliVEvent::kEMCEGA )
2814  {
2815  //printf("EGA trigger bit\n");
2816  if (GetFiredTriggerClasses().Contains("EG"))
2817  {
2818  if (GetFiredTriggerClasses().Contains("EGA")) fEventTrigEMCALL1Gamma1 = kTRUE;
2819  else
2820  {
2821  if(GetFiredTriggerClasses().Contains("EG1")) fEventTrigEMCALL1Gamma1 = kTRUE;
2822  if(GetFiredTriggerClasses().Contains("EG2")) fEventTrigEMCALL1Gamma2 = kTRUE;
2823  }
2824  }
2825  }
2826  // EMC L1 Jet
2827  else if( fEventTriggerMask & AliVEvent::kEMCEJE )
2828  {
2829  //printf("EGA trigger bit\n");
2830  if (GetFiredTriggerClasses().Contains("EJ"))
2831  {
2832  if (GetFiredTriggerClasses().Contains("EJE")) fEventTrigEMCALL1Jet1 = kTRUE;
2833  else
2834  {
2835  if(GetFiredTriggerClasses().Contains("EJ1")) fEventTrigEMCALL1Jet1 = kTRUE;
2836  if(GetFiredTriggerClasses().Contains("EJ2")) fEventTrigEMCALL1Jet2 = kTRUE;
2837  }
2838  }
2839  }
2840  // EMC L0
2841  else if((fEventTriggerMask & AliVEvent::kEMC7) ||
2842  (fEventTriggerMask & AliVEvent::kEMC1) )
2843  {
2844  //printf("L0 trigger bit\n");
2845  fEventTrigEMCALL0 = kTRUE;
2846  }
2847  // Min Bias Pb-Pb
2848  else if( fEventTriggerMask & AliVEvent::kCentral )
2849  {
2850  //printf("MB semi central trigger bit\n");
2851  fEventTrigSemiCentral = kTRUE;
2852  }
2853  // Min Bias Pb-Pb
2854  else if( fEventTriggerMask & AliVEvent::kSemiCentral )
2855  {
2856  //printf("MB central trigger bit\n");
2857  fEventTrigCentral = kTRUE;
2858  }
2859  // Min Bias pp, PbPb, pPb
2860  else if((fEventTriggerMask & AliVEvent::kMB ) ||
2861  (fEventTriggerMask & AliVEvent::kINT7) ||
2862  (fEventTriggerMask & AliVEvent::kINT8) ||
2863  (fEventTriggerMask & AliVEvent::kAnyINT) )
2864  {
2865  //printf("MB trigger bit\n");
2866  fEventTrigMinBias = kTRUE;
2867  }
2868  }
2869 
2870  AliDebug(1,Form("Event bits: \n \t MB %d, Cen %d, Sem %d, L0 %d, L1G1 %d, L1G2 %d, L1J1 %d, L1J2 %d",
2874 
2875  // L1 trigger bit
2876  if( fBitEGA == 0 && fBitEJE == 0 )
2877  {
2878  // Init the trigger bit once, correct depending on AliESD(AOD)CaloTrigger header version
2879 
2880  // Simpler way to do it ...
2881 // if( fInputEvent->GetRunNumber() < 172000 )
2882 // reader->SetEventTriggerL1Bit(4,5); // current LHC11 data
2883 // else
2884 // reader->SetEventTriggerL1Bit(6,8); // LHC12-13 data
2885 
2886  // Old values
2887  fBitEGA = 4;
2888  fBitEJE = 5;
2889 
2890  TFile* file = AliAnalysisManager::GetAnalysisManager()->GetTree()->GetCurrentFile();
2891 
2892  const TList *clist = file->GetStreamerInfoCache();
2893 
2894  if(clist)
2895  {
2896  TStreamerInfo *cinfo = (TStreamerInfo*)clist->FindObject("AliESDCaloTrigger");
2897  Int_t verid = 5; // newer ESD header version
2898  if(!cinfo)
2899  {
2900  cinfo = (TStreamerInfo*)clist->FindObject("AliAODCaloTrigger");
2901  verid = 3; // newer AOD header version
2902  }
2903 
2904  if(cinfo)
2905  {
2906  Int_t classversionid = cinfo->GetClassVersion();
2907  //printf("********* Header class version %d *********** \n",classversionid);
2908 
2909  if (classversionid >= verid)
2910  {
2911  fBitEGA = 6;
2912  fBitEJE = 8;
2913  }
2914  } else AliInfo("AliCaloTrackReader()::SetEventTriggerBit() - Streamer info for trigger class not available, bit not changed");
2915  } else AliInfo("AliCaloTrackReader::SetEventTriggerBit() - Streamer list not available!, bit not changed");
2916 
2917  } // set once the EJE, EGA trigger bit
2918 }
2919 
2920 //____________________________________________________________
2923 //____________________________________________________________
2924 void AliCaloTrackReader::SetInputEvent(AliVEvent* const input)
2925 {
2926  fInputEvent = input;
2927  fMixedEvent = dynamic_cast<AliMixedEvent*>(GetInputEvent()) ;
2928  if (fMixedEvent)
2929  fNMixedEvent = fMixedEvent->GetNumberOfEvents() ;
2930 
2931  //Delete previous vertex
2932  if(fVertex)
2933  {
2934  for (Int_t i = 0; i < fNMixedEvent; i++)
2935  {
2936  delete [] fVertex[i] ;
2937  }
2938  delete [] fVertex ;
2939  }
2940 
2941  fVertex = new Double_t*[fNMixedEvent] ;
2942  for (Int_t i = 0; i < fNMixedEvent; i++)
2943  {
2944  fVertex[i] = new Double_t[3] ;
2945  fVertex[i][0] = 0.0 ;
2946  fVertex[i][1] = 0.0 ;
2947  fVertex[i][2] = 0.0 ;
2948  }
2949 }
2950 
2951 
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)
ClassImp(AliAnalysisTaskTriggerRates) AliAnalysisTaskTriggerRates
virtual void FillInputVZERO()
Int_t fV0ADC[2]
Integrated V0 signal.
Bool_t fComparePtHardAndClusterPt
In MonteCarlo, jet events, reject events with too large cluster energy.
AliAnaWeights * fWeightUtils
Pointer to AliAnaWeights.
AliCalorimeterUtils * GetCaloUtils() const
Float_t fTimeStampEventFracMin
Minimum value of time stamp fraction event.
Bool_t fReadAODMCParticles
Access kine information from filtered AOD MC particles.
virtual void FillInputNonStandardJets()
Double_t fEMCALTimeCutMax
Remove clusters/cells with time larger than this value, in ns.
Int_t fBitEJE
Trigger bit on VCaloTrigger for EJE.
TLorentzVector fMomentum
! Temporal TLorentzVector container, avoid declaration of TLorentzVectors per event.
Double_t fTimeStampRunMin
Minimum value of time stamp in run.
Bool_t fDoPileUpEventRejection
Select pile-up events by SPD.
void SetCentrality(Float_t cen)
Definition: AliAnaWeights.h:50
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.
Bool_t ClusterContainsBadChannel(Int_t calo, UShort_t *cellList, Int_t nCells)
virtual AliHeader * GetHeader() const
Bool_t IsEventEMCALL1Jet1() const
Bool_t IsEventEMCALL1Gamma1() const
AliEMCALRecoUtils * GetEMCALRecoUtils() const
Bool_t ReadAODMCParticles() const
TString fEventPlaneMethod
Name of event plane method, by default "Q".
Int_t fTrackBCEventCut[19]
Fill one entry per event if there is a track in a given BC, depend on track pT, acceptance cut...
Bool_t fDoVertexBCEventSelection
Select events with vertex on BC=0 or -100.
AliMixedEvent * fMixedEvent
! Mixed event object. This class is not the owner.
virtual AliVEvent * GetInputEvent() const
Bool_t IsPileUpFromSPDAndNotEMCal() const
Check if event is from pile-up determined by SPD and not by EMCal.
Bool_t IsEventEMCALL1() const
Float_t fPtHardAndClusterPtFactor
Factor between ptHard and cluster pT to reject/accept event.
Bool_t fAcceptFastCluster
Accept events from fast cluster, exclude these events for LHC11a.
Bool_t fSmearShowerShape
Smear shower shape (use in MC).
AliVEvent * fInputEvent
! pointer to esd or aod input.
Calculate the weight to the event to be applied when filling histograms.
Definition: AliAnaWeights.h:31
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)
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.
Int_t GetDebug() const
Definition: AliAnaWeights.h:82
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.
Bool_t IsEventEMCALL1Gamma2() const
virtual AliMixedEvent * GetMixedEvent() const
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.
Int_t fTriggerClusterId
Id of trigger cluster (cluster->GetID()).
Float_t fCTSPtMax
pT Threshold on charged particles.
Double_t fEMCALParamTimeCutMax[4]
Remove clusters/cells with time larger than parametrized value, in ns.
virtual void FillInputEMCALCells()
Connects the array with EMCAL cells and the pointer.
Bool_t fFillPHOSCells
Use data from PHOS.
virtual AliVCaloCells * GetEMCALCells() const
Bool_t fRemoveLEDEvents
Remove events where LED was wrongly firing - EMCAL LHC11a.
AliEMCALGeometry * GetEMCALGeometry() const
Bool_t fIsBadMaxCellEvent
Bad cell triggered event flag, only max energy cell is bad.
void RecalculateClusterShowerShapeParameters(AliVCaloCells *cells, AliVCluster *clu)
TString fFiredTriggerClassName
Name of trigger event type used to do the analysis.
virtual TClonesArray * GetAODMCParticles() const
TString GetFiredTriggerClasses() const
Bool_t ReadStack() const
TObjArray * fEMCALClusters
Temporal array with EMCAL CaloClusters.
Bool_t IsEventEMCALL1Jet() const
Float_t fPHOSPtMax
pT Threshold on phos clusters.
Bool_t fAccessTrackTOF
Access the track TOF, in case of problems when accessing GetTOFBunchCrossing.
virtual Bool_t CheckForPrimaryVertex() const
Int_t fEMCalBCEvent[19]
Fill one entry per event if there is a cluster in a given BC.
TString fInputBackgroundJetBranchName
Name of background jet branch.
virtual AliEventplane * GetEventPlane() const
Bool_t fEventTriggerAtSE
Select triggered event at SE base task or here.
Float_t fEMCALPtMax
pT Threshold on emcal clusters.
Int_t fEventType
Set the event species: 7 physics, 0 MC, 8 LED (not useful now)
virtual Bool_t ComparePtHardAndJetPt()
Bool_t IsInFiducialCut(Float_t eta, Float_t phi, Int_t det) const
Bool_t fEventTrigCentral
Event is AliVEvent::kCentral on its name, it should correspond to PbPb.
Float_t fTrackMultEtaCut
Track multiplicity eta cut.
Int_t fCentralityBin[2]
Minimum and maximum value of the centrality for the analysis.
virtual void FillInputPHOSCells()
Connects the array with PHOS cells and the pointer.
Base class for event, clusters and tracks filtering and preparation for the analysis.
TObjArray * fDCALClusters
Temporal array with DCAL CaloClusters, not needed in the normal case, use just EMCal array with DCal ...
Bool_t fUseTrackTimeCut
Do time cut selection.
virtual TString GetEventPlaneMethod() const
TArrayI fAcceptEventsWithBit
Accept events if trigger bit is on.
virtual void Print(const Option_t *opt) const
Print some relevant parameters set for the analysis.
void MatchTriggerCluster(TArrayI patches)
Bool_t fTimeStampEventSelect
Select events within a fraction of data taking time.
Double_t fTimeStampRunMax
Maximum value of time stamp in run.
Double_t fEMCALTimeCutMin
Remove clusters/cells with time smaller than this value, in ns.
Bool_t fEventTrigEMCALL1Jet2
Event is L1-Gamma, threshold 2 on its name, it should correspond kEMCEGA.
void SetEMCalEventBCcut(Int_t bc)
Int_t fNPileUpClustersCut
Cut to select event as pile-up.
Bool_t fCheckFidCut
Do analysis for clusters in defined region.
virtual ~AliCaloTrackReader()
Destructor.
Bool_t fEventTrigEMCALL0
Event is EMCal L0 on its name, it should correspond to AliVEvent::kEMC7, AliVEvent::kEMC1.
Bool_t IsWeightSettingOn() const
Definition: AliAnaWeights.h:80
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.
AliFiducialCut * fFiducialCut
Acceptance cuts.
Bool_t fRejectEMCalTriggerEventsWith2Tresholds
Reject events EG2 also triggered by EG1 or EJ2 also triggered by EJ1.
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
void CorrectClusterEnergy(AliVCluster *cl)
Correct cluster energy non linearity.
Float_t fCTSPtMin
pT Threshold on charged particles.
virtual AliStack * GetStack() const
Float_t fTriggerL0EventThreshold
L0 Threshold to look for triggered events, set outside.
Float_t fTimeStampEventFracMax
Maximum value of time stamp fraction event.
virtual void ResetLists()
Reset lists, called in AliAnaCaloTrackCorrMaker.
Int_t fNPileUpClusters
Number of clusters with time avobe 20 ns.
Int_t fTrackMult
Track multiplicity.
AliAODEvent * fOutputEvent
! pointer to aod output.
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
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
TObjArray * fCTSTracks
Temporal array with tracks.
Int_t fTriggerPatchTimeWindow[2]
Trigger patch selection window.
AliMCEvent * fMC
! Monte Carlo Event Handler.
Float_t fPtHardAndJetPtFactor
Factor between ptHard and jet pT to reject/accept event.
TString fDeltaAODFileName
Delta AOD file name.
Int_t GetVertexBC() const
Bool_t IsPileUpFromEMCal() const
Check if event is from pile-up determined by EMCal.
virtual Double_t GetWeight()
Bool_t fFillInputNonStandardJetBranch
Flag to use data from non standard jets.
TClonesArray * fNonStandardJets
! Temporal array with jets.
Int_t fEMCalBCEventCut[19]
Fill one entry per event if there is a cluster in a given BC, depend on cluster E, acceptance cut.
Int_t fTriggerClusterIndex
Index in clusters array of trigger cluster.
Bool_t fTriggerL1EventThresholdFix
L1 Threshold is fix and set outside.
Float_t fSmearShowerShapeWidth
Smear shower shape landau function "width" (use in MC).
Int_t fBitEGA
Trigger bit on VCaloTrigger for EGA.
virtual void FillInputBackgroundJets()
Bool_t fDoV0ANDEventSelection
Select events depending on V0AND.
AliAODJetEventBackground * fBackgroundJets
! Background jets.
Bool_t fCorrectELinearity
Correct cluster linearity, always on.
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.
Bool_t fUseAliCentrality
Select as centrality estimator AliCentrality (Run1) or AliMultSelection (Run1 and Run2) ...
virtual TObjArray * GetCTSTracks() const
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