]>
Commit | Line | Data |
---|---|---|
5707a463 | 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 | /* $Id$ */ | |
17 | ||
18 | ||
07d83718 | 19 | // Macro to produce trigger single muon efficiency versus pt plots for the |
20 | // 2 pt cuts. | |
21 | // see full description in the README file | |
22 | // Author: Fabien Guerin (LPC) | |
5707a463 | 23 | |
24 | // ROOT includes | |
25 | #include "TBranch.h" | |
26 | #include "TClonesArray.h" | |
27 | #include "TFile.h" | |
28 | #include "TF1.h" | |
29 | #include "TH1.h" | |
30 | #include "TH2.h" | |
31 | #include "TParticle.h" | |
32 | #include "TTree.h" | |
33 | #include "TMath.h" | |
34 | #include "TStyle.h" | |
35 | #include "TCanvas.h" | |
36 | #include "TLegend.h" | |
b19dcb40 | 37 | #include "Riostream.h" |
5707a463 | 38 | |
39 | // STEER includes | |
40 | #include "AliRun.h" | |
5707a463 | 41 | #include "AliRunLoader.h" |
42 | #include "AliHeader.h" | |
43 | #include "AliLoader.h" | |
44 | #include "AliStack.h" | |
5707a463 | 45 | |
46 | // MUON includes | |
b19dcb40 | 47 | #include "AliMUONSimData.h" |
48 | #include "AliMUONRecData.h" | |
5707a463 | 49 | #include "AliMUONHit.h" |
5707a463 | 50 | #include "AliMUONDigit.h" |
5707a463 | 51 | #include "AliMUONGlobalTrigger.h" |
52 | #include "AliMUONLocalTrigger.h" | |
53 | #include "AliMUONTrack.h" | |
54 | ||
55 | Double_t fitArc(Double_t *x,Double_t *par) | |
56 | { | |
57 | Double_t h1=par[1]*(x[0]+par[2]); | |
58 | return par[0]*TMath::ATan(TMath::Power(h1,6))+par[3]; | |
59 | } | |
60 | ||
61 | Double_t fitArch(Double_t *x,Double_t *par) | |
62 | { | |
63 | Double_t h0=2*par[0]/(TMath::Pi()); | |
64 | Double_t h1=par[1]*x[0]+par[2]*x[0]*x[0]; | |
65 | return h0*TMath::TanH(h1)+par[3]; | |
66 | } | |
67 | ||
b19dcb40 | 68 | void MUONTriggerEfficiencyPt(char filenameSim[10]="galice.root", |
69 | char filenameRec[10]="galice.root", | |
70 | Bool_t readFromRP = 0) | |
5707a463 | 71 | { |
72 | ||
73 | // define style | |
74 | TStyle *st1 = new TStyle("st1","My STYLE"); | |
75 | ||
76 | st1->SetOptStat(0); | |
77 | st1->SetOptFit(111); | |
78 | st1->SetFrameFillColor(10); | |
79 | st1->SetCanvasColor(10); | |
80 | st1->SetFillColor(10); | |
81 | st1->SetTitleBorderSize(0); | |
82 | st1->SetTitleOffset(1.2,"XY"); | |
83 | st1->SetTitleSize(0.05,"XY"); | |
84 | st1->SetLabelSize(0.045,"XY"); | |
85 | st1->SetLabelFont(22,"XY"); | |
86 | st1->SetTitleFont(22,"XY"); | |
87 | st1->SetOptTitle(0); | |
88 | st1->SetStatColor(10); | |
89 | st1->SetStatFont(62); | |
90 | st1->SetFitFormat("4.4g"); | |
91 | st1->SetPadLeftMargin(0.15); | |
92 | st1->SetPadRightMargin(0.06); | |
93 | st1->SetPadTopMargin(0.06); | |
94 | st1->SetPadBottomMargin(0.15); | |
95 | st1->cd(); | |
96 | ||
07d83718 | 97 | // gROOT->ForceStyle(); |
5707a463 | 98 | //TGaxis::SetMaxDigits(3); |
99 | ||
100 | // beginning of macro | |
5707a463 | 101 | TH1F *ptlpt = new TH1F("ptlpt","",50,0.15,5.); |
102 | TH1F *pthpt = new TH1F("pthpt","",50,0.15,5.); | |
103 | TH1F *ptcoinc = new TH1F("ptcoinc","",50,0.15,5.); | |
104 | ||
105 | TParticle *particle; | |
106 | AliStack* stack; | |
107 | ||
b19dcb40 | 108 | Int_t nevents; |
5707a463 | 109 | Double_t coincmuon=0; |
5707a463 | 110 | Double_t lptmuon=0; |
111 | Double_t hptmuon=0; | |
522b4724 | 112 | Double_t ptmu=0; |
5707a463 | 113 | |
114 | // output file | |
115 | char digitdat[100]; | |
b19dcb40 | 116 | char currentfile[100]; |
5707a463 | 117 | sprintf(digitdat,"MUONTriggerEfficiencyPt.out"); |
118 | FILE *fdat=fopen(digitdat,"w"); | |
b19dcb40 | 119 | |
120 | // Initialise AliRoot | |
5707a463 | 121 | // Creating Run Loader and openning file containing Hits |
b19dcb40 | 122 | AliRunLoader * RunLoaderSim = AliRunLoader::Open(filenameSim,"MUONFolderSim","READ"); |
123 | ||
124 | if (RunLoaderSim ==0x0) { | |
125 | printf(">>> Error : Error Opening %s file \n",currentfile); | |
126 | return; | |
127 | } | |
5707a463 | 128 | |
b19dcb40 | 129 | AliRunLoader * RunLoaderRec = AliRunLoader::Open(filenameRec,"MUONFolder","READ"); |
130 | ||
131 | if (RunLoaderRec ==0x0) { | |
132 | printf(">>> Error : Error Opening %s file \n",currentfile); | |
5707a463 | 133 | return; |
134 | } | |
135 | ||
b19dcb40 | 136 | nevents = RunLoaderSim->GetNumberOfEvents(); |
5707a463 | 137 | |
b19dcb40 | 138 | AliLoader * MUONLoaderSim = RunLoaderSim->GetLoader("MUONLoader"); |
139 | AliLoader * MUONLoaderRec = RunLoaderRec->GetLoader("MUONLoader"); | |
5707a463 | 140 | |
b19dcb40 | 141 | if (!readFromRP) { |
142 | cout << " reading from digits \n"; | |
143 | MUONLoaderSim->LoadDigits("READ"); | |
144 | } else { | |
145 | cout << " reading from RecPoints \n"; | |
146 | MUONLoaderRec->LoadRecPoints("READ"); | |
147 | } | |
148 | MUONLoaderSim->LoadHits("READ"); | |
149 | RunLoaderSim->LoadKinematics("READ"); | |
150 | ||
151 | ||
152 | // Creating MUON data containers | |
153 | AliMUONSimData muondataSim(MUONLoaderSim,"MUON","MUON"); | |
154 | AliMUONRecData muondataRec(MUONLoaderRec,"MUON","MUON"); | |
155 | ||
5707a463 | 156 | TClonesArray * globalTrigger; |
157 | AliMUONGlobalTrigger * gloTrg; | |
158 | ||
5707a463 | 159 | for (Int_t ievent=0; ievent<nevents; ievent++) { // Event loop |
160 | ||
b19dcb40 | 161 | RunLoaderSim->GetEvent(ievent); |
162 | RunLoaderRec->GetEvent(ievent); | |
5707a463 | 163 | |
164 | if (ievent%500==0) printf("ievent = %d \n",ievent); | |
165 | ||
5707a463 | 166 | // kine |
167 | Int_t iparticle, nparticles; | |
b19dcb40 | 168 | stack = RunLoaderSim->Stack(); |
5707a463 | 169 | nparticles = (Int_t) stack->GetNtrack(); |
170 | for (iparticle=0; iparticle<nparticles; iparticle++) { | |
171 | particle = stack->Particle(iparticle); | |
172 | Float_t pt = particle->Pt(); | |
173 | Int_t pdgcode = particle->GetPdgCode(); | |
174 | if (pdgcode==-13 || pdgcode==13) ptmu = pt; | |
175 | } | |
176 | ||
177 | // trigger | |
07d83718 | 178 | if (!readFromRP) { |
b19dcb40 | 179 | muondataSim.SetTreeAddress("D,GLT"); |
180 | muondataSim.GetTriggerD(); | |
181 | globalTrigger = muondataSim.GlobalTrigger(); | |
07d83718 | 182 | } else { |
b19dcb40 | 183 | muondataRec.SetTreeAddress("RC,TC"); |
184 | muondataRec.GetTrigger(); | |
185 | globalTrigger = muondataRec.GlobalTrigger(); | |
186 | ||
07d83718 | 187 | } |
5707a463 | 188 | |
189 | Int_t nglobals = (Int_t) globalTrigger->GetEntriesFast(); // should be 1 | |
190 | for (Int_t iglobal=0; iglobal<nglobals; iglobal++) { // Global Trigger | |
191 | gloTrg = static_cast<AliMUONGlobalTrigger*>(globalTrigger->At(iglobal)); | |
192 | ||
522b4724 | 193 | if (gloTrg->SingleLpt()>=1 ) { |
5707a463 | 194 | lptmuon++; |
195 | ptlpt->Fill(ptmu); | |
196 | } | |
522b4724 | 197 | if (gloTrg->SingleHpt()>=1 ) { |
5707a463 | 198 | hptmuon++; |
199 | pthpt->Fill(ptmu); | |
200 | } | |
201 | } // end of loop on Global Trigger | |
b19dcb40 | 202 | muondataSim.ResetTrigger(); |
203 | muondataRec.ResetTrigger(); | |
5707a463 | 204 | |
205 | // Hits | |
b19dcb40 | 206 | muondataSim.SetTreeAddress("H"); |
207 | ||
5707a463 | 208 | Int_t itrack, ntracks, NbHits[4]; |
209 | Int_t SumNbHits; | |
522b4724 | 210 | |
5707a463 | 211 | for (Int_t j=0; j<4; j++) NbHits[j]=0; |
212 | ||
b19dcb40 | 213 | ntracks = (Int_t) muondataSim.GetNtracks(); |
214 | ||
5707a463 | 215 | for (itrack=0; itrack<ntracks; itrack++) { // track loop |
b19dcb40 | 216 | muondataSim.GetTrack(itrack); |
5707a463 | 217 | |
218 | Int_t ihit, nhits; | |
b19dcb40 | 219 | nhits = (Int_t) muondataSim.Hits()->GetEntriesFast(); |
5707a463 | 220 | AliMUONHit* mHit; |
221 | ||
222 | for (ihit=0; ihit<nhits; ihit++) { | |
b19dcb40 | 223 | mHit = static_cast<AliMUONHit*>(muondataSim.Hits()->At(ihit)); |
5707a463 | 224 | Int_t Nch = mHit->Chamber(); |
225 | Int_t hittrack = mHit->Track(); | |
09f9a8f1 | 226 | Float_t IdPart = mHit->Particle(); |
522b4724 | 227 | |
5707a463 | 228 | for (Int_t j=0;j<4;j++) { |
229 | Int_t kch=11+j; | |
522b4724 | 230 | if (hittrack==0) { |
5707a463 | 231 | if (Nch==kch && (IdPart==-13 || IdPart==13) && NbHits[j]==0) NbHits[j]++; |
522b4724 | 232 | } |
5707a463 | 233 | } |
234 | } | |
b19dcb40 | 235 | muondataSim.ResetHits(); |
522b4724 | 236 | } // end track loop |
07d83718 | 237 | |
5707a463 | 238 | // 3/4 coincidence |
239 | SumNbHits=NbHits[0]+NbHits[1]+NbHits[2]+NbHits[3]; | |
522b4724 | 240 | |
5707a463 | 241 | if (SumNbHits==3 || SumNbHits==4) { |
242 | coincmuon++; | |
243 | ptcoinc->Fill(ptmu); | |
244 | } | |
245 | ||
b19dcb40 | 246 | } // end loop on event |
247 | ||
248 | if (coincmuon==0) { | |
249 | cout << " >>> <E> coincmuon = 0 after event loop " << "\n"; | |
250 | cout << " >>> this probably means that input does not contain one (and only one) " << "\n"; | |
251 | cout << " >>> muon track per event as it should " << "\n"; | |
252 | cout << " >>> see README for further information " << "\n"; | |
253 | cout << " >>> exiting now ! " << "\n"; | |
254 | return; | |
255 | } | |
5707a463 | 256 | |
b19dcb40 | 257 | |
258 | MUONLoaderSim->UnloadHits(); | |
259 | if (!readFromRP) { | |
260 | MUONLoaderSim->UnloadDigits(); | |
261 | } else { | |
262 | MUONLoaderRec->UnloadRecPoints(); | |
263 | } | |
264 | RunLoaderSim->UnloadKinematics(); | |
5707a463 | 265 | |
266 | fprintf(fdat,"\n"); | |
267 | fprintf(fdat,"\n"); | |
268 | fprintf(fdat," Number of events = %d \n",nevents); | |
269 | fprintf(fdat," Number of events with 3/4 coinc = %d \n",(Int_t)coincmuon); | |
5707a463 | 270 | fprintf(fdat," Nomber of dimuons with 3/4 coinc Lpt cut = %d \n",(Int_t)lptmuon); |
271 | fprintf(fdat," Number of dimuons with 3/4 coinc Hpt cut = %d \n",(Int_t)hptmuon); | |
272 | fprintf(fdat,"\n"); | |
273 | ||
522b4724 | 274 | Double_t efficiency,error; |
275 | ||
5707a463 | 276 | efficiency=lptmuon/coincmuon; |
277 | error=efficiency*TMath::Sqrt((lptmuon+coincmuon)/(lptmuon*coincmuon)); | |
278 | fprintf(fdat," Efficiency Lpt cut = %4.4f +/- %4.4f\n",efficiency,error); | |
279 | ||
280 | efficiency=hptmuon/coincmuon; | |
281 | error=efficiency*TMath::Sqrt((hptmuon+coincmuon)/(hptmuon*coincmuon)); | |
282 | fprintf(fdat," Efficiency Hpt cut = %4.4f +/- %4.4f\n",efficiency,error); | |
283 | ||
284 | fclose(fdat); | |
285 | ||
286 | Double_t x1,x2,xval,xerr; | |
287 | ||
288 | for (Int_t i=0;i<50;i++) { | |
289 | x1=ptcoinc->GetBinContent(i+1); | |
290 | ||
5707a463 | 291 | x2=ptlpt->GetBinContent(i+1); |
292 | if (x1!=0 && x2!=0) { | |
293 | xval=x2/x1; | |
294 | xerr=xval*TMath::Sqrt((x1+x2)/x1*x2); | |
295 | ptlpt->SetBinContent(i+1,xval); | |
296 | ptlpt->SetBinError(i+1,0); | |
297 | } else { | |
298 | ptlpt->SetBinContent(i+1,0.); | |
299 | ptlpt->SetBinError(i+1,0.); | |
300 | } | |
301 | ||
302 | x2=pthpt->GetBinContent(i+1); | |
303 | if (x1!=0 && x2!=0) { | |
304 | xval=x2/x1; | |
305 | xerr=xval*TMath::Sqrt((x1+x2)/x1*x2); | |
306 | pthpt->SetBinContent(i+1,xval); | |
307 | pthpt->SetBinError(i+1,0); | |
308 | } else { | |
309 | pthpt->SetBinContent(i+1,0.); | |
310 | pthpt->SetBinError(i+1,0.); | |
311 | } | |
b19dcb40 | 312 | } |
5707a463 | 313 | |
5707a463 | 314 | TF1 *fitlpt = new TF1("fitlpt",fitArc,0.,5.,4); |
315 | TF1 *fithpt = new TF1("fithpt",fitArc,0.,5.,4); | |
316 | ||
317 | TCanvas *c1 = new TCanvas("c1","Trigger efficiency vs pt and y for single muon",200,0,900,400); | |
522b4724 | 318 | c1->Divide(2,1); |
5707a463 | 319 | |
320 | c1->cd(1); | |
5707a463 | 321 | ptlpt->SetMinimum(0.); |
322 | ptlpt->SetMaximum(1.05); | |
323 | ptlpt->SetTitle(""); | |
324 | ptlpt->GetXaxis()->SetTitle("P_{T} (GeV/c)"); | |
325 | ptlpt->GetYaxis()->SetTitle("Efficiency"); | |
326 | //ptlpt->GetXaxis()->SetRange(3); | |
327 | ptlpt->Draw(""); | |
328 | fitlpt->SetParameters(0.602,0.774,0.273,0.048); | |
329 | fitlpt->SetLineColor(2); | |
330 | fitlpt->SetLineWidth(2); | |
331 | fitlpt->Draw("SAME"); | |
522b4724 | 332 | TLegend * leg = new TLegend(0.5,0.38,0.85,0.53); |
5707a463 | 333 | leg = new TLegend(0.5,0.38,0.85,0.53); |
334 | leg->SetBorderSize(0); | |
335 | leg->SetFillColor(0); | |
336 | leg->SetTextSize(0.05); | |
337 | leg->SetTextFont(22); | |
338 | leg->SetHeader("Lpt trigger pt cut"); | |
339 | leg->AddEntry(fitlpt,"Ref.","l"); | |
340 | leg->AddEntry(ptlpt,"New","l"); | |
341 | leg->Draw("SAME"); | |
342 | ||
522b4724 | 343 | c1->cd(2); |
5707a463 | 344 | pthpt->SetMinimum(0.); |
345 | pthpt->SetMaximum(1.05); | |
346 | pthpt->SetTitle(""); | |
347 | pthpt->GetXaxis()->SetTitle("P_{T} (GeV/c)"); | |
348 | pthpt->GetYaxis()->SetTitle("Efficiency"); | |
349 | //pthpt->GetXaxis()->SetRange(3); | |
350 | pthpt->Draw(""); | |
351 | fithpt->SetParameters(0.627,0.393,0.855,0.0081); | |
352 | fithpt->SetLineColor(2); | |
353 | fithpt->SetLineWidth(2); | |
354 | fithpt->Draw("SAME"); | |
355 | leg = new TLegend(0.5,0.38,0.85,0.53); | |
356 | leg->SetBorderSize(0); | |
357 | leg->SetFillColor(0); | |
358 | leg->SetTextSize(0.05); | |
359 | leg->SetTextFont(22); | |
360 | leg->SetHeader("Hpt trigger pt cut"); | |
361 | leg->AddEntry(fithpt,"Ref.","l"); | |
362 | leg->AddEntry(pthpt,"New","l"); | |
363 | leg->Draw("SAME"); | |
364 | ||
365 | c1->SaveAs("MUONTriggerEfficiencyPt.gif"); | |
07d83718 | 366 | c1->SaveAs("MUONTriggerEfficiencyPt.eps"); |
5707a463 | 367 | } |