AliPhysics  59e0e03 (59e0e03)
AliCaloTrackAODReader.cxx
Go to the documentation of this file.
1 
2 /**************************************************************************
3  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4  * *
5  * Author: The ALICE Off-line Project. *
6  * Contributors are mentioned in the code where appropriate. *
7  * *
8  * Permission to use, copy, modify and distribute this software and its *
9  * documentation strictly for non-commercial purposes is hereby granted *
10  * without fee, provided that the above copyright notice appears in all *
11  * copies and that both the copyright notice and this permission notice *
12  * appear in the supporting documentation. The authors make no claims *
13  * about the suitability of this software for any purpose. It is *
14  * provided "as is" without express or implied warranty. *
15  **************************************************************************/
16 
17 //---- ANALYSIS system ----
18 #include "AliCaloTrackAODReader.h"
19 #include "AliAODInputHandler.h"
20 #include "AliMultiEventInputHandler.h"
21 #include "AliAnalysisManager.h"
22 #include "AliMixedEvent.h"
23 #include "AliAODEvent.h"
24 #include "AliGenEventHeader.h"
25 #include "AliLog.h"
26 
28 ClassImp(AliCaloTrackAODReader) ;
30 
31 //______________________________________________
33 //______________________________________________
35  AliCaloTrackReader(), fOrgInputEvent(0x0),
36  fSelectHybridTracks(0), fSelectPrimaryTracks(0),
37  fTrackFilterMask(0), fTrackFilterMaskComplementary(0),
38  fSelectFractionTPCSharedClusters(0), fCutTPCSharedClustersFraction(0)
39 {
40  fDataType = kAOD;
41 
42  //fTrackFilterMask = 128;
43  //fTrackFilterMaskComplementary = 0; // in case of hybrid tracks, without using the standard method
44 
45  //fSelectFractionTPCSharedClusters = kTRUE;
47 }
48 
49 //_________________________________________________________
52 //_________________________________________________________
54 {
55  AliAODEvent * aodevent = dynamic_cast<AliAODEvent*>(fInputEvent);
56  if(!aodevent) return kFALSE;
57 
58  if (aodevent->GetPrimaryVertex() != NULL)
59  {
60  if(aodevent->GetPrimaryVertex()->GetNContributors() > 0)
61  {
62  return kTRUE;
63  }
64  }
65 
66  if(aodevent->GetPrimaryVertexSPD() != NULL)
67  {
68  if(aodevent->GetPrimaryVertexSPD()->GetNContributors() > 0)
69  {
70  return kTRUE;
71  }
72  else
73  {
74  AliDebug(1,Form("Null number of contributors from bad vertex type:: %s",
75  aodevent->GetPrimaryVertex()->GetName()));
76  return kFALSE;
77  }
78  }
79 
80  return kFALSE;
81 }
82 
83 
84 //___________________________________________________
89 //___________________________________________________
91 {
93 
94  if(fFillCTS)
95  {
96  for(Int_t i = 0; i < 4; i++)
97  {
98  TString names[] = {"FilterBit_Hybrid", "SPDHit", "SharedCluster", "Primary"};
99 
100  fhCTSAODTrackCutsPt[i] = new TH1F(Form("hCTSReaderAODClusterCuts_%d_%s",i,names[i].Data()),
101  Form("AOD CTS Cut %d, %s",i,names[i].Data()),
103  fhCTSAODTrackCutsPt[i]->SetYTitle("# tracks");
104  fhCTSAODTrackCutsPt[i]->SetXTitle("#it{p}_{T} (GeV)");
106  }
107  }
108 
109  return fOutputContainer ;
110 }
111 
112 //________________________________________________________
114 //________________________________________________________
116 {
117  // Recover the string from the mother class
118  TString parList = (AliCaloTrackReader::GetListOfParameters())->GetString();
119 
120  const Int_t buffersize = 255;
121  char onePar[buffersize] ;
122 
123  snprintf(onePar,buffersize,"AOD Track: Hybrid %d, Filter bit %d, Complementary bit %d, Primary %d; ",
125  parList+=onePar ;
126 
128  {
129  snprintf(onePar,buffersize,"Fraction of TPC shared clusters ON: %2.2f ", fCutTPCSharedClustersFraction) ;
130  parList+=onePar ;
131  }
132 
133  return new TObjString(parList) ;
134 
135 }
136 
137 //____________________________________________________________
139 //____________________________________________________________
141 {
142  TClonesArray * particles = NULL ;
143 
144  AliAODEvent * aod = dynamic_cast<AliAODEvent*> (fInputEvent) ;
145  if(aod) particles = (TClonesArray*) aod->FindListObject("mcparticles");
146 
147  return particles ;
148 }
149 
150 //___________________________________________________________
152 //___________________________________________________________
154 {
155  AliAODMCHeader *mch = NULL;
156 
157  AliAODEvent * aod = dynamic_cast<AliAODEvent*> (fInputEvent);
158  if(aod) mch = dynamic_cast<AliAODMCHeader*>(aod->FindListObject("mcHeader"));
159 
160  return mch;
161 }
162 
163 //______________________________________________________________
165 //______________________________________________________________
167 {
168  if ( !GetAODMCHeader() ) return 0x0;
169 
170  if ( fGenEventHeader ) return fGenEventHeader;
171 
172  Int_t nGenerators = GetAODMCHeader()->GetNCocktailHeaders();
173 
174  if ( nGenerators <= 0 ) return 0x0;
175 
176  if ( fMCGenerEventHeaderToAccept=="" )
177  return GetAODMCHeader()->GetCocktailHeader(0);
178 
179  //if(nGenerators < 1) printf("No generators %d\n",nGenerators);
180 
181  for(Int_t igen = 0; igen < nGenerators; igen++)
182  {
183  AliGenEventHeader * eventHeader = GetAODMCHeader()->GetCocktailHeader(igen) ;
184 
185  TString name = eventHeader->GetName();
186 
187  AliDebug(2,Form("AOD Event header %d name %s event header class %s: Select if contains <%s>",
188  igen,name.Data(),eventHeader->ClassName(),fMCGenerEventHeaderToAccept.Data()));
189 
190 // if(nGenerators < 2)
191 // printf("AOD Event header %d name %s event header class %s: Select if contains <%s>\n",
192 // igen,name.Data(),eventHeader->ClassName(),fMCGenerEventHeaderToAccept.Data());
193 
194  if(name.Contains(fMCGenerEventHeaderToAccept,TString::kIgnoreCase))
195  return eventHeader ;
196  }
197 
198  return 0x0;
199 }
200 
201 //_____________________________
212 //_____________________________
214 {
215  Int_t id = track->GetID();
216 
217  if( id < 0 ) id = TMath::Abs(id) - 1;
218 
219  return id;
220 }
221 
222 
223 //_____________________________________________________________________________
225 //_____________________________________________________________________________
227 {
228  AliAODTrack *aodtrack = dynamic_cast <AliAODTrack*>(track);
229 
230  if(!aodtrack) return kFALSE;
231 
232  track->GetPxPyPz(pTrack) ;
233 
234  AliDebug(2,Form("AOD track type: %d (primary %d), hybrid? %d",
235  aodtrack->GetType(),AliAODTrack::kPrimary,
236  aodtrack->IsHybridGlobalConstrainedGlobal()));
237 
238  // Hybrid?
240  {
241  if ( !aodtrack->IsHybridGlobalConstrainedGlobal() ) return kFALSE ;
242  }
243  else if ( fTrackFilterMask > 0 || fTrackFilterMaskComplementary > 0 ) // Filter Bit?
244  {
245  Bool_t accept = aodtrack->TestFilterBit(fTrackFilterMask);
246 
247  if ( !fSelectHybridTracks && !accept ) return kFALSE ;
248 
249  if ( fSelectHybridTracks ) // Second filter bit for hybrids?
250  {
251  Bool_t acceptcomplement = aodtrack->TestFilterBit(fTrackFilterMaskComplementary);
252  if ( !accept && !acceptcomplement ) return kFALSE ;
253  }
254  }
255 
256  fhCTSAODTrackCutsPt[0]->Fill(aodtrack->Pt());
257 
258  //
259  if ( fSelectSPDHitTracks )
260  { // Not much sense to use with TPC only or Hybrid tracks
261  if(!aodtrack->HasPointOnITSLayer(0) && !aodtrack->HasPointOnITSLayer(1)) return kFALSE ;
262  }
263 
264  fhCTSAODTrackCutsPt[1]->Fill(aodtrack->Pt());
265 
266  //
268  {
269  Double_t frac = 0;
270  Float_t ncls = Float_t(aodtrack->GetTPCncls ());
271  Float_t nclsS = Float_t(aodtrack->GetTPCnclsS());
272  if ( ncls> 0 ) frac = nclsS / ncls ;
273 
274  if ( frac > fCutTPCSharedClustersFraction )
275  {
276  AliDebug(2,Form("\t Reject track, shared cluster fraction %f > %f",frac, fCutTPCSharedClustersFraction));
277  return kFALSE ;
278  }
279  }
280 
281  fhCTSAODTrackCutsPt[2]->Fill(aodtrack->Pt());
282 
283  //
284  if ( fSelectPrimaryTracks )
285  {
286  if ( aodtrack->GetType()!= AliAODTrack::kPrimary )
287  {
288  AliDebug(2,"\t Remove not primary track");
289  return kFALSE ;
290  }
291  }
292 
293  AliDebug(2,"\t accepted track!");
294 
295  fhCTSAODTrackCutsPt[3]->Fill(aodtrack->Pt());
296 
297 
298  return kTRUE;
299 }
300 
301 //_________________________________________________________________
304 //_________________________________________________________________
306  AliAODEvent* aod,
307  AliMCEvent* mc)
308 {
309  //printf("AODInputHandler %p, MergeEvents %d \n",aodIH, aodIH->GetMergeEvents());
310 
311  Bool_t tesd = kFALSE ;
312  Bool_t taod = kTRUE ;
313  if ( strcmp(input->GetName(), "AliMixedEvent") == 0 )
314  {
315  AliMultiEventInputHandler* multiEH = dynamic_cast<AliMultiEventInputHandler*>((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
316  if(multiEH){
317  if (multiEH->GetFormat() == 0 )
318  {
319  tesd = kTRUE ;
320  } else if (multiEH->GetFormat() == 1)
321  {
322  taod = kTRUE ;
323  }
324  }
325  else
326  {
327  AliFatal("MultiEventHandler is NULL");
328  return;
329  }
330  }
331  if (strcmp(input->GetName(),"AliESDEvent") == 0)
332  {
333  tesd = kTRUE ;
334  } else if (strcmp(input->GetName(),"AliAODEvent") == 0)
335  {
336  taod = kTRUE ;
337  }
338 
339  if(tesd)
340  {
341  SetInputEvent(aod);
342  SetOutputEvent(aod);
343  fOrgInputEvent = input;
344  }
345  else if(taod)
346  {
347  AliAODInputHandler* aodIH = dynamic_cast<AliAODInputHandler*>((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
348 
349  if (aodIH && aodIH->GetMergeEvents())
350  {
351  //Merged events, use output AOD.
352  SetInputEvent(aod);
353  SetOutputEvent(aod);
354  fOrgInputEvent = input;
355  }
356  else
357  {
358  SetInputEvent(input);
359  SetOutputEvent(aod);
360  }
361  }
362  else
363  {
364  AliFatal(Form("STOP : Wrong data format: %s",input->GetName()));
365  }
366 
367  SetMC(mc);
368 }
369 
double Double_t
Definition: External.C:58
virtual void SetMC(AliMCEvent *const mc)
Bool_t fSelectHybridTracks
Select CTS tracks of type hybrid.
AliVEvent * fOrgInputEvent
! Original input event, not from filtering
AliVEvent * fInputEvent
! pointer to esd or aod input.
ULong_t fTrackFilterMaskComplementary
Complementary Track selection bit, for AODs in case hybrid option selected.
virtual void SetInputEvent(AliVEvent *input)
Bool_t SelectTrack(AliVTrack *track, Double_t *pTrack)
Select AOD track using the AOD filter bits or predefined selection methods.
TList * fOutputContainer
! Output container with cut control histograms.
Int_t fEnergyHistogramNbins
Binning of the control histograms, min and max window.
Bool_t fFillCTS
Use data from CTS.
AliAODMCHeader * GetAODMCHeader() const
ULong_t fTrackFilterMask
Track selection bit, for AODs (any difference with track status?)
Class for event, clusters and tracks filtering and preparation for the AOD analysis.
Float_t fEnergyHistogramLimit[2]
Binning of the control histograms, number of bins.
int Int_t
Definition: External.C:63
Bool_t fSelectSPDHitTracks
Ensure that track hits SPD layers.
AliGenEventHeader * fGenEventHeader
! Event header
float Float_t
Definition: External.C:68
virtual void SetOutputEvent(AliAODEvent *aod)
Bool_t fSelectPrimaryTracks
Select CTS tracks of type primary.
virtual TObjString * GetListOfParameters()
Save parameters used for analysis in a string.
Base class for event, clusters and tracks filtering and preparation for the analysis.
TClonesArray * GetAODMCParticles() const
Float_t fCutTPCSharedClustersFraction
Fraction of TPC shared clusters to be accepted.
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)
TObjString * GetListOfParameters()
Save parameters used for analysis in a string.
void SetInputOutputMCEvent(AliVEvent *esd, AliAODEvent *aod, AliMCEvent *mc)
AliGenEventHeader * GetGenEventHeader() const
Int_t GetTrackID(AliVTrack *track)
TString fMCGenerEventHeaderToAccept
Accept events that contain at least this event header name.
Bool_t fSelectFractionTPCSharedClusters
Accept only TPC tracks with over a given fraction of shared clusters.
AliCaloTrackAODReader()
Default constructor. Initialize parameters.
bool Bool_t
Definition: External.C:53
Bool_t CheckForPrimaryVertex() const
virtual TList * GetCreateControlHistograms()
Int_t fDataType
Select MC: Kinematics, Data: ESD/AOD, MCData: Both.
TH1F * fhCTSAODTrackCutsPt[4]
! control histogram on the different CTS tracks selection cuts, pT