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