AliPhysics  e59a9ba (e59a9ba)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliAnaRandomTrigger.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 <TH2F.h>
18 #include <TClonesArray.h>
19 
20 //---- AliRoot system ----
21 #include "AliAnaRandomTrigger.h"
22 #include "AliAODPWG4ParticleCorrelation.h"
23 #include "AliEMCALGeometry.h"
24 
28 
29 //__________________________________________
31 //__________________________________________
34  fTriggerDetector(kEMCAL),
35  fTriggerDetectorString("EMCAL"),
36  fRandom(0), fNRandom(0),
37  fMomentum(),
38  fhE(0), fhPt(0),
39  fhPhi(0), fhEta(0),
40  fhEtaPhi(0)
41 {
43 }
44 
45 //_________________________________________________________________________
48 //_________________________________________________________________________
49 Bool_t AliAnaRandomTrigger::ExcludeDeadBadRegions(Float_t eta, Float_t phi)
50 {
51  if(fTriggerDetector!=kEMCAL) return kFALSE;
52 
53  //-------------------------------------
54  // Get the corresponding cell in EMCAL, check if it exists in acceptance (phi gaps, borders)
55  //-------------------------------------
56 
57  Int_t absId = -1;
58  if(!GetEMCALGeometry()->GetAbsCellIdFromEtaPhi(eta,phi, absId)) return kTRUE; // remove if out of EMCAL acceptance, phi gaps
59 
60  Int_t icol = -1, irow = -1, iRCU = -1;
61  Int_t sm = GetCaloUtils()->GetModuleNumberCellIndexes(absId,kEMCAL, icol, irow, iRCU);
62 
63  //printf("eta %f, phi %f, ieta %d, iphi %d, sm %d\n",eta,phi,icol,irow,sm);
64 
65  //-------------------------------------
66  // Remove in case of close to border, by default always 1 but check what was set in reco utils
67  //-------------------------------------
68 
69  Bool_t okrow = kFALSE;
70  Bool_t okcol = kFALSE;
71  Int_t nborder = GetCaloUtils()->GetEMCALRecoUtils()->GetNumberOfCellsFromEMCALBorder();
72  if (nborder<1) nborder = 1;
73 
74  // Rows
75  if(sm < 10)
76  {
77  if(irow >= nborder && irow < 24-nborder) okrow =kTRUE;
78  }
79  else
80  {
81  if((GetCaloUtils()->EMCALGeometryName()).Contains("12SM")) // 1/3 SM
82  {
83  if(irow >= nborder && irow < 8-nborder) okrow =kTRUE;
84  }
85  else // 1/2 SM
86  {
87  if(irow >= nborder && irow <12-nborder) okrow =kTRUE;
88  }
89  }
90 
91  // Columns
92  if(sm%2==0)
93  {
94  if(icol >= nborder) okcol = kTRUE;
95  }
96  else
97  {
98  if(icol < 48-nborder) okcol = kTRUE;
99  }
100 
101  //printf("okcol %d, okrow %d\n",okcol,okrow);
102  if (!okcol || !okrow) return kTRUE;
103 
104  //-------------------------------------
105  // Check if the cell or those around are bad
106  //-------------------------------------
107 
108  if(GetCaloUtils()->GetEMCALChannelStatus(sm,icol, irow) > 0) return kTRUE ; // trigger falls into a bad channel
109 
110  // Check if close there was a bad channel
111 // for(Int_t i = -1; i <= 1; i++)
112 // {
113 // for(Int_t j = -1; j <= 1; j++)
114 // {
115 // //printf("\t check icol %d, irow %d \n",icol+i, irow+j);
116 // if(GetCaloUtils()->GetEMCALChannelStatus(sm,icol+i, irow+j) > 0) return kTRUE ; // trigger falls into a bad channel
117 // //printf("\t ok\n");
118 // }
119 // }
120 
121  //printf("\t OK\n");
122 
123  return kFALSE;
124 }
125 
126 //__________________________________________________
128 //__________________________________________________
130 {
131  TString parList ; //this will be list of parameters used for this analysis.
132  const Int_t buffersize = 255;
133  char onePar[buffersize] ;
134 
135  snprintf(onePar,buffersize,"--- AliAnaRandomTrigger ---:") ;
136  parList+=onePar ;
137  snprintf(onePar,buffersize,"Detector: %s;" , fTriggerDetectorString.Data()) ;
138  parList+=onePar ;
139  snprintf(onePar,buffersize,"N per event = %d;", fNRandom ) ;
140  parList+=onePar ;
141  snprintf(onePar,buffersize,"Min E = %3.2f - Max E = %3.2f;", GetMinPt(), GetMaxPt()) ;
142  parList+=onePar ;
143  snprintf(onePar,buffersize,"Min Eta = %3.2f - Max Eta = %3.2f;", fEtaCut[0], fEtaCut[1]) ;
144  parList+=onePar ;
145  snprintf(onePar,buffersize,"Min Phi = %3.2f - Max Phi = %3.2f;", fPhiCut[0], fPhiCut[1]) ;
146  parList+=onePar ;
147 
148  return new TObjString(parList) ;
149 }
150 
151 //____________________________________________________
154 //____________________________________________________
156 {
157  TList * outputContainer = new TList() ;
158  outputContainer->SetName("RandomTrigger") ;
159 
160  Int_t nptbins = GetHistogramRanges()->GetHistoPtBins(); Int_t nphibins = GetHistogramRanges()->GetHistoPhiBins(); Int_t netabins = GetHistogramRanges()->GetHistoEtaBins();
163 
164  fhE = new TH1F ("hE","Random E distribution", nptbins,ptmin,ptmax);
165  fhE->SetXTitle("E (GeV)");
166  outputContainer->Add(fhE);
167 
168  fhPt = new TH1F ("hPt","Random p_{T} distribution", nptbins,ptmin,ptmax);
169  fhPt->SetXTitle("p_{T} (GeV/c)");
170  outputContainer->Add(fhPt);
171 
172  fhPhi = new TH2F ("hPhi","Random #phi distribution",
173  nptbins,ptmin,ptmax, nphibins,phimin,phimax);
174  fhPhi->SetYTitle("#phi (rad)");
175  fhPhi->SetXTitle("p_{T} (GeV/c)");
176  outputContainer->Add(fhPhi);
177 
178  fhEta = new TH2F ("hEta","Random #eta distribution",
179  nptbins,ptmin,ptmax, netabins,etamin,etamax);
180  fhEta->SetYTitle("#eta ");
181  fhEta->SetXTitle("p_{T} (GeV/c)");
182  outputContainer->Add(fhEta);
183 
184  fhEtaPhi = new TH2F ("hEtaPhi","Random #eta vs #phi ",netabins,etamin,etamax, nphibins,phimin,phimax);
185  fhEtaPhi->SetXTitle("#eta ");
186  fhEtaPhi->SetYTitle("#phi (rad)");
187  outputContainer->Add(fhEtaPhi);
188 
189  return outputContainer;
190 }
191 
192 //________________________________________
194 //________________________________________
196 {
197  SetOutputAODClassName("AliAODPWG4ParticleCorrelation");
198  SetOutputAODName("RandomTrigger");
199 
200  AddToHistogramsName("AnaRandomTrigger_");
201 
202  fNRandom = 1 ;
203  fPhiCut[0] = 0. ;
204  fPhiCut[1] = TMath::TwoPi() ;
205  fEtaCut[0] =-1. ;
206  fEtaCut[1] = 1. ;
207 }
208 
209 //_________________________________________________________
211 //_________________________________________________________
212 void AliAnaRandomTrigger::Print(const Option_t * opt) const
213 {
214  if(! opt)
215  return;
216 
217  printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
219 
220  printf("Detector = %s\n", fTriggerDetectorString.Data());
221  printf("Min E = %3.2f - Max E = %3.2f\n", GetMinPt(), GetMaxPt());
222  printf("Min Eta = %3.2f - Max Eta = %3.2f\n", fEtaCut[0], fEtaCut[1]);
223  printf("Min Phi = %3.2f - Max Phi = %3.2f\n", fPhiCut[0], fPhiCut[1]);
224 }
225 
226 //______________________________________________
230 //______________________________________________
232 {
233  for(Int_t irandom = 0; irandom < fNRandom; irandom++)
234  {
235  // Get the random variables of the trigger
236  Float_t pt = fRandom.Uniform(GetMinPt(), GetMaxPt());
237  Float_t eta = fRandom.Uniform(fEtaCut[0], fEtaCut[1]);
238  Float_t phi = fRandom.Uniform(fPhiCut[0], fPhiCut[1]);
239 
240  // Check if particle falls into a dead region, if inside, get new
241  Bool_t excluded = ExcludeDeadBadRegions(eta,phi);
242 
243  // If excluded, generate a new trigger until accepted
244  while (excluded)
245  {
246  pt = fRandom.Uniform(GetMinPt(), GetMaxPt());
247  eta = fRandom.Uniform(fEtaCut[0], fEtaCut[1]);
248  phi = fRandom.Uniform(fPhiCut[0], fPhiCut[1]);
249 
250  excluded = ExcludeDeadBadRegions(eta,phi);
251  }
252 
253  // Create the AOD trigger object
254  fMomentum.SetPtEtaPhiM(pt,eta,phi,0);
255 
256  AliAODPWG4Particle trigger = AliAODPWG4Particle(fMomentum);
257  trigger.SetDetectorTag(fTriggerDetector);
258  trigger.SetSModNumber(GetModuleNumber(&trigger));
259 
260  AliDebug(1,Form("iRandom %d, Trigger e %2.2f pt %2.2f, phi %2.2f, eta %2.2f, SM %d",
261  irandom, trigger.E(), trigger.Pt(), trigger.Phi(), trigger.Eta(), trigger.GetSModNumber()));
262 
263  AddAODParticle(trigger);
264  }
265 
266  AliDebug(1,Form("Final aod branch entries %d", GetOutputAODBranch()->GetEntriesFast()));
267 }
268 
269 //_____________________________________________________
270 // Fill control histograms with generated trigger kinematics.
271 //_____________________________________________________
273 {
274  // Loop on stored AODParticles
275  Int_t naod = GetOutputAODBranch()->GetEntriesFast();
276 
277  AliDebug(1,Form("AOD branch entries %d, fNRandom %d", naod, fNRandom));
278 
279  for(Int_t iaod = 0; iaod < naod ; iaod++)
280  {
281  AliAODPWG4Particle* trigger = (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));
282 
283  fhPt ->Fill(trigger->Pt (), GetEventWeight());
284  fhE ->Fill(trigger->E (), GetEventWeight());
285  fhPhi ->Fill(trigger->Pt (), trigger->Phi(), GetEventWeight());
286  fhEta ->Fill(trigger->Pt (), trigger->Eta(), GetEventWeight());
287  fhEtaPhi->Fill(trigger->Eta(), trigger->Phi(), GetEventWeight());
288  }// aod branch loop
289 }
290 
291 //_________________________________________________________
293 //_________________________________________________________
295 {
297 
298  if (det=="EMCAL") fTriggerDetector = kEMCAL;
299  else if(det=="PHOS" ) fTriggerDetector = kPHOS;
300  else if(det=="CTS") fTriggerDetector = kCTS;
301  else if(det=="DCAL") fTriggerDetector = kDCAL;
302  else if(det.Contains("DCAL") && det.Contains("PHOS")) fTriggerDetector = kDCALPHOS;
303  else AliFatal(Form("Detector < %s > not known!", det.Data()));
304 }
305 
306 //______________________________________________________
307 // Set the calorimeter for the analysis.
308 //______________________________________________________
310 {
311  fTriggerDetector = det;
312 
313  if ( det == kEMCAL ) fTriggerDetectorString = "EMCAL";
314  else if( det == kPHOS ) fTriggerDetectorString = "PHOS";
315  else if( det == kCTS ) fTriggerDetectorString = "CTS";
316  else if( det == kDCAL ) fTriggerDetectorString = "DCAL";
317  else if( det == kDCALPHOS) fTriggerDetectorString = "DCAL_PHOS";
318  else AliFatal(Form("Detector < %d > not known!", det));
319 }
320 
321 
Float_t GetHistoPtMax() const
Float_t GetHistoPtMin() const
Int_t fNRandom
Number of random particles per event.
virtual void AddToHistogramsName(TString add)
Float_t fPhiCut[2]
Phi acceptance, radians.
Bool_t ExcludeDeadBadRegions(Float_t eta, Float_t phi)
TH2F * fhEta
! eta distribution vs pT, negative
AliEMCALRecoUtils * GetEMCALRecoUtils() const
void InitParameters()
Initialize the parameters of the analysis.
TObjString * GetAnalysisCuts()
Save parameters used for analysis.
virtual void SetOutputAODClassName(TString name)
Int_t GetHistoPhiBins() const
virtual TClonesArray * GetOutputAODBranch() const
Int_t fTriggerDetector
Detector : EMCAL, PHOS, CTS.
Float_t GetHistoPhiMin() const
TLorentzVector fMomentum
! Avoid generating TLorentzVectors per event.
void Print(const Option_t *opt) const
Print some relevant parameters set for the analysis.
virtual void SetOutputAODName(TString name)
void SetTriggerDetector(TString &det)
Set the detrimeter for the analysis.
const Double_t etamin
Base class for CaloTrackCorr analysis algorithms.
TString fTriggerDetectorString
Detector : EMCAL, PHOS, CTS.
virtual AliHistogramRanges * GetHistogramRanges()
const Double_t ptmax
virtual AliEMCALGeometry * GetEMCALGeometry() const
virtual AliCalorimeterUtils * GetCaloUtils() const
TRandom3 fRandom
Random generator.
const Double_t ptmin
virtual Double_t GetEventWeight() const
TH2F * fhEtaPhi
! eta vs phi distribution of positive charge
Float_t GetHistoEtaMin() const
virtual Int_t GetModuleNumber(AliAODPWG4Particle *part) const
TH1F * fhPt
! pT distribution
Float_t GetHistoEtaMax() const
virtual void AddAODParticle(AliAODPWG4Particle part)
Int_t GetHistoPtBins() const
ClassImp(AliAnalysisTaskCRC) AliAnalysisTaskCRC
TH1F * fhE
! E distribution
Gerenate a random trigger.
const Double_t etamax
virtual void Print(const Option_t *) const
Print some relevant parameters set for the analysis.
TH2F * fhPhi
! phi distribution vs pT, negative
AliAnaRandomTrigger()
Default Constructor. Initialize parameters.
Float_t GetHistoPhiMax() const
Int_t GetHistoEtaBins() const
DCal, not used so far, just in case.
Float_t fEtaCut[2]
Eta acceptance.
Int_t nptbins
Int_t GetModuleNumberCellIndexes(Int_t absId, Int_t calo, Int_t &icol, Int_t &irow, Int_t &iRCU) const
Get the EMCAL/PHOS module, columns, row and RCU/DDL number that corresponds to this absId...
const Double_t phimin