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