First commit.
[u/mrichter/AliRoot.git] / FASTSIM / AliFastMuonTriggerEff.cxx
CommitLineData
6255180c 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/*
17$Log$
18*/
19
20
21#include "AliFastMuonTriggerEff.h"
22#include <TFile.h>
23#include <TTree.h>
24#include <TROOT.h>
25#include <stdlib.h>
26
27// Debugging flag
28//#define MYTRIGDEBUG
29
30ClassImp(AliFastMuonTriggerEff);
31
32AliFastMuonTriggerEff::AliFastMuonTriggerEff():
33 AliFastResponse("Efficiency", "Muon Trigger Efficiency")
34{
35//
36// Default constructor
37//
38 SetCut(kLow);
39 fZones=0;
40}
41
42
43void AliFastMuonTriggerEff::Init()
44{
45//
46// Initialization
47//
48 printf("Initializing %s / %s", fName.Data(), fTitle.Data());
49//
50
51 fPhiMin = -90. ;
52 fPhiMax = 90. ;
53 fDphi = 9. ;
54 fThetaMin = 2. ;
55 fThetaMax = 9. ;
56 fDtheta = 0.7;
57 fDpt = 0.1;
58
59 InitTree();
60}
61
62void AliFastMuonTriggerEff::InitTree()
63{
64//
65// Initialize tables from tree
66//
67 TTree *chain; // pointer to the analyzed TTree or TChain
68 Int_t nSim, nZona;
69 Double_t pmin, pmax, tmin, tmax, phimin, phimax, bkg;
70 Double_t len50, hen50, leff, heff;
71 Double_t vLPt[50];
72 Double_t vHPt[50];
73 Char_t file[100]="$(ALICE_ROOT)/data/vettorpara.root";
74//
75// Avoid memory leak in case of reinitialization
76 if(fZones!=0) {
77 printf("\nWarning: reinitialization of an object of class: AliFastMuonTriggerEff\n");
78 for (Int_t i=0; i<fZones; i++) {
79 if(fEffLow [i])delete[] fEffLow[i];
80 if(fEffHigh[i])delete[] fEffHigh[i];
81 }
82 if(fEffLow) {
83 delete[] fEffLow;
84 fEffLow=0;
85 }
86 if(fEffHigh) {
87 delete[] fEffHigh;
88 fEffHigh=0;
89 }
90 }
91 printf("AliFastMuonTriggerEff: Initialization\n");
92 TFile *f = new TFile(file);
93 if(f->IsZombie()) {
94 printf("Cannot open file: %s\n",file);
95 return;
96 }
97 f->ls();
98
99 TTree* tree = (TTree*)gDirectory->Get("fitPar");
100
101//
102//
103 chain = tree;
104 chain->SetMakeClass(1);
105
106 chain->SetBranchAddress("nSim",&nSim);
107 chain->SetBranchAddress("nZona",&nZona);
108 chain->SetBranchAddress("ptmin",&pmin);
109 chain->SetBranchAddress("ptmax",&pmax);
110 chain->SetBranchAddress("Thetamin",&tmin);
111 chain->SetBranchAddress("Thetamax",&tmax);
112 chain->SetBranchAddress("Phimin",&phimin);
113 chain->SetBranchAddress("Phimax",&phimax);
114 chain->SetBranchAddress("Bkg",&bkg);
115 chain->SetBranchAddress("EffLPt",vLPt);
116 chain->SetBranchAddress("EffHPt",vHPt);
117 chain->SetBranchAddress("Pt0.5Low",&len50);
118 chain->SetBranchAddress("Pt0.5High",&hen50);
119 chain->SetBranchAddress("EffLowMax",&leff);
120 chain->SetBranchAddress("EffHighMax",&heff);
121//
122//
123 Int_t nentries = Int_t(chain->GetEntries());
124 Int_t nbytes = 0, nb = 0;
125
126// Count the number of zones of the parametrization
127 Int_t nzone0=0, nzone1=0;
128 for (Int_t jentry=0; jentry<nentries; jentry++) {
129 nb = chain->GetEntry(jentry);
130 if(bkg==0.) {
131 if(nSim==0)nzone0++;
132 if(nSim==1)nzone1++;
133 }
134 }
135
136 printf("Trigger parametrisation for %d zones for Pt: 0. 5. GeV/c\n",nzone0);
137 printf("and %d zones extended to 10. GeV/c\n",nzone1);
138 fZones=nzone0+nzone1;
139// printf("Ciao\n");
140 if(fZones<=0){
141 printf("Number of zones must be greater than 0\n");
142 exit(6);
143 }
144
145 fEffLow =new Float_t*[fZones];
146 fEffHigh=new Float_t*[fZones];
147 for (Int_t i=0; i<fZones; i++) {
148 fEffLow [i]=new Float_t[50];
149 fEffHigh[i]=new Float_t[50];
150 }
151
152// Initialize look-up table to standard values
153 Int_t isim, itheta, iphi;
154 for (isim=0; isim<2; isim++) {
155 for (itheta=0; itheta<10; itheta++) {
156 for (iphi=0; iphi<20; iphi++) {
157 fLook[isim][itheta][iphi]=0;
158 }
159 }
160 }
161
162// Loading Trigger efficiencies
163 Int_t myzone=0;
164#ifdef MYTRIGDEBUG
165 printf("Filling nSim nZona pmin pmax tmin tmax phimin phimax: ....\n");
166#endif
167 for (Int_t jentry=0; jentry<nentries; jentry++) {
168// Int_t ientry = LoadTree(jentry);
169 nb = chain->GetEntry(jentry);
170 nbytes += nb;
171#ifdef MYTRIGDEBUG
172 printf("Getting entry %d... ",jentry);
173#endif
174// For the time being it works with background 0
175 if ((nSim == 0 || nSim == 1)&&bkg==0.) {
176#ifdef MYTRIGDEBUG
177 printf("Filling %d %d %6.2f %6.2f %6.2f %6.2f %6.2f %6.2f: ",
178 nSim,nZona,pmin,pmax,tmin,tmax,phimin,phimax);
179#endif
180 for (Int_t k = 0; k < 50; k++) {
181 fEffLow [myzone][k] = vLPt[k];
182 fEffHigh[myzone][k] = vHPt[k];
183#ifdef MYTRIGDEBUG
184 if(k<15)printf(" %5.3f",vLPt[k]);
185#endif
186 }
187#ifdef MYTRIGDEBUG
188 printf("\n");
189#endif
190 myzone++;
191 iphi=Int_t(nZona/10)-10;
192 itheta=nZona%10;
193 if(iphi<0||iphi>19||itheta<0||itheta>9) {
194 printf("The format of data file is not consistent\nwith this version of the code\n");
195 printf("This should never happen: iphi %d, itheta: %d\n",iphi,itheta);
196 exit(7);
197 }
198 fLook[nSim][itheta][iphi]=myzone;
199 } else {
200 printf("Skipping entry with nSim=%d, bkg=%f\n",nSim,bkg);
201 }
202 }
203#ifdef MYTRIGDEBUG
204 printf("This is the content of the LUT after first step of initialization\n");
205 for(isim=0; isim<2; isim++) {
206 printf("isim=%d\n",isim);
207 for(iphi=0; iphi<20; iphi++) {
208 printf("iphi: %2d:",iphi);
209 for(itheta=0; itheta<10; itheta++) {
210 printf(" %4d",fLook[isim][itheta][iphi]);
211 }
212 printf("\n");
213 }
214 }
215#endif
216// Filling look=up table for the zones where the extended simulation does
217// not exists
218 for(iphi=0; iphi<20; iphi++) {
219 for(itheta=0; itheta<10; itheta++) {
220 if(fLook[0][itheta][iphi]==0) {
221 printf("Missing entry isim=%d itheta=%d iphi=%d\n",isim,itheta,iphi);
222 exit(8);
223 }
224 if(fLook[0][itheta][iphi]<0||fLook[0][itheta][iphi]>fZones) {
225 printf("Problem with entry isim=%d itheta=%d iphi=%d\n",isim,itheta,iphi);
226 exit(9);
227 }
228 if(fLook[1][itheta][iphi]==0) {
229 fLook[1][itheta][iphi]=-fLook[0][itheta][iphi];
230 }
231 }
232 }
233#ifdef MYTRIGDEBUG
234 for(isim=0; isim<2; isim++) {
235 printf("isim=%d\n",isim);
236 for(iphi=0; iphi<20; iphi++) {
237 printf("iphi: %2d:",iphi);
238 for(itheta=0; itheta<10; itheta++) {
239 printf(" %4d",fLook[isim][itheta][iphi]);
240 }
241 printf("\n");
242 }
243 }
244 for(iphi=0; iphi<20; iphi++) {
245 for(itheta=0; itheta<10; itheta++) {
246 for(isim=0; isim<2; isim++) {
247 printf("%d %2d %2d %4d:",isim,iphi,itheta,fLook[isim][itheta][iphi]);
248 if(fLook[isim][itheta][iphi]>0) {
249 myzone=fLook[isim][itheta][iphi]-1;
250 for(Int_t ipt=0; ipt<20; ipt++) {
251 //printf(" %5.3f",fEffLow[myzone][ipt]);
252 printf(" %5.3f",fEffHigh[myzone][ipt]);
253 }
254 printf(" ...");
255 for(Int_t ipt=40; ipt<50; ipt++) {
256 //printf(" %5.3f",fEffLow[myzone][ipt]);
257 printf(" %5.3f",fEffHigh[myzone][ipt]);
258 }
259 }
260 printf("\n");
261 }
262 }
263 }
264#endif
265 f->Close();
266}
267
268void AliFastMuonTriggerEff::Evaluate(Float_t charge, Float_t pt,
269 Float_t theta, Float_t phi,Float_t& effLow, Float_t& effHigh)
270{
271//
272// Trigger efficiency for pt, theta, phi (low and high cut)
273//
274 effLow=0.;
275 effHigh=0.;
276 if(fZones==0) {
277 printf("Call to uninitialized object of class: AliFastMuonTriggerEff\n");
278 return;
279 }
280 if(pt<0) {
281 printf("Warning: pt: %f < 0. GeV/c\n",pt);
282 return;
283 }
284 Int_t iPt = Int_t(pt/fDpt)%50;
285 Int_t iSim = Int_t(pt/fDpt)/50;
286 Int_t iPhi = Int_t((phi-fPhiMin)/fDphi);
287 if(phi<fPhiMin)iPhi=iPhi-1;
288 Int_t iTheta = Int_t((theta-fThetaMin)/fDtheta);
289#ifdef MYTRIGDEBUG
290 printf("iSim iPt iTheta iPhi: %d %d %d %d\n",iSim,iPt,iTheta,iPhi);
291#endif
292 iPhi=iPhi-40*(iPhi/40);
293#ifdef MYTRIGDEBUG
294 printf("1: iPhi converted to: %d for angle equivalence\n",iPhi);
295#endif
296 if(iPhi<0)iPhi=-iPhi-1;
297 if(iPhi>19)iPhi=39-iPhi;
298#ifdef MYTRIGDEBUG
299 printf("2: iPhi converted to: %d for the symmetry of the spectrometer\n",iPhi);
300#endif
301 if(charge==1.){
302 } else if(charge==-1.) {
303 iPhi=19-iPhi;
304#ifdef MYTRIGDEBUG
305 printf("3: iPhi converted to: %d for charge symmetry\n",iPhi);
306#endif
307 } else {
308 printf("Warning: not understand charge: %f\n",charge);
309 return;
310 }
311 if(iTheta<0||iTheta>9) {
312 printf("Warning: theta: %f outside acceptance\n",theta);
313 return;
314 }
315 if(iPt<0) {
316 printf("Warning: what do you mean with pt: %f <0?\n",pt);
317 return;
318 }
319 if(iSim>=fSim) {
320 iSim=fSim-1;
321 iPt=49;
322#ifdef MYTRIGDEBUG
323 printf("4: iSim iPt converted to: %d %d (last zone)\n",iSim,iPt);
324#endif
325 }
326 Int_t iLook=fLook[iSim][iTheta][iPhi];
327 if(iLook<0) {
328#ifdef MYTRIGDEBUG
329 printf("5: iLook iPt: %d %d converted to: ",iLook,iPt);
330#endif
331 iLook=-iLook-1;
332 iPt=49;
333#ifdef MYTRIGDEBUG
334 printf("%d %d from look up table contents\n",iLook,iPt);
335#endif
336 } else {
337 iLook=iLook-1;
338 }
339 effLow=fEffLow[iLook][iPt];
340 effHigh=fEffHigh[iLook][iPt];
341#ifdef MYTRIGDEBUG
342 printf("6: charge, iSim, iTheta, iPhi, iPt: %f %d %d %d %d effLow: %f, effHigh: %f\n",
343 charge,iSim,iTheta,iPhi,iPt,effLow,effHigh);
344#endif
345
346 //fEffLow [iPhi][iTheta][iPt];
347 //fEffHigh[iPhi][iTheta][iPt];
348}
349
350Float_t AliFastMuonTriggerEff::Evaluate(Float_t charge, Float_t pt,
351 Float_t theta, Float_t phi)
352{
353 //
354 // Trigger efficiency for pt, theta, phi depending of fCut
355 //
356 if(fZones==0) {
357 printf("Call to uninitialized object of class: AliFastMuonTriggerEff\n");
358 return 0.;
359 }
360 Float_t eff;
361 Float_t effLow, effHigh;
362
363 Evaluate(charge,pt,theta,phi,effLow,effHigh);
364 if (fCut == kLow)
365 eff = effLow;
366 else
367 eff = effHigh;
368
369 return eff;
370}