In AliMUONTriggerQADataMakerRec:
[u/mrichter/AliRoot.git] / MUON / MUONTriggerEfficiency.C
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 /// \ingroup macros
19 /// \file MUONTriggerEfficiency.C
20 /// \brief Macro for checking trigger integrated efficiency for dimuons.
21 ///
22 /// The efficiency is calculated with respect to the 3/4 coincidence.
23 ///
24 /// \author Fabien Guerin, LPC Clermont-Ferrand, Jan. 2006
25
26 #if !defined(__CINT__) || defined(__MAKECINT__)
27 // ROOT includes
28 #include "TBranch.h"
29 #include "TClonesArray.h"
30 #include "TObject.h"
31 #include "TFile.h"
32 #include "TH1.h"
33 #include "TH2.h"
34 #include "TTree.h"
35 #include "TParticle.h"
36
37 // STEER includes
38 #include "AliRun.h"
39 #include "AliRunLoader.h"
40 #include "AliHeader.h"
41 #include "AliLoader.h"
42 #include "AliStack.h"
43 #include "AliCDBManager.h"
44
45 // MUON includes
46 #include "AliMUON.h"
47 #include "AliMUONHit.h"
48 #include "AliMUONConstants.h"
49 #include "AliMUONDigit.h"
50 #include "AliMUONGlobalTrigger.h"
51 #include "AliMUONTrack.h"
52
53 #include "AliMUONDataInterface.h"
54 #include "AliMUONMCDataInterface.h"
55 #include "AliMUONVHitStore.h"
56 #include "AliMUONVTriggerStore.h"
57
58 #endif
59
60 // Upsilon(1S)
61
62
63 void MUONTriggerEfficiency(const char* filenameSim="galice_sim.root", 
64                            const char* filenameRec="galice.root",  
65                            Bool_t readFromRP = 0)
66 {
67   
68   // Set default CDB storage
69   AliCDBManager* man = AliCDBManager::Instance();
70   man->SetDefaultStorage("local://$ALICE_ROOT/OCDB");
71
72   // output file
73   
74   AliMUONMCDataInterface diSim(filenameSim);
75   AliMUONDataInterface diRec(filenameRec);
76   
77   if (!diSim.IsValid())
78   {
79     cerr << "Cannot access " << filenameSim << endl;
80     return;
81   }
82   
83   if (!diRec.IsValid())
84   {
85     cerr << "Cannot access " << filenameRec << endl;
86     return;
87   }
88   
89   FILE* fdat = fopen("MUONTriggerEfficiency.out","w");          
90   if (!fdat)
91   {
92     cerr << "Cannot create output file" << endl;
93     return;
94   }
95   
96   Int_t coincmuon,muonlpt,muonhpt;
97   Int_t CoincMuPlus,CoincMuMinus;
98   coincmuon=0;
99   muonlpt=0;  
100   muonhpt=0;
101   
102   Int_t nevents = diSim.NumberOfEvents();          
103   
104   for (Int_t ievent=0; ievent<nevents; ++ievent) 
105   {   
106     CoincMuPlus=0;
107     CoincMuMinus=0;
108     
109     if (ievent%1000==0) printf("\t Event = %d\n",ievent);    
110     
111     // Hits
112     
113     Int_t NbHits[2][4];
114     for (Int_t j=0; j<2; j++) 
115     {
116       for (Int_t jj=0; jj<4; jj++) 
117       {
118         NbHits[j][jj]=0;
119       }
120     } 
121     
122     Int_t ntracks = (Int_t) diSim.NumberOfTracks(ievent);
123     
124     for ( Int_t itrack=0; itrack<ntracks; ++itrack ) 
125     {
126       AliMUONVHitStore* hitStore = diSim.HitStore(ievent,itrack);      
127       AliMUONHit* mHit;
128       TIter next(hitStore->CreateIterator());
129       
130       while ( ( mHit = static_cast<AliMUONHit*>(next()) ) )
131       {
132         Int_t Nch = mHit->Chamber(); 
133         Int_t hittrack = mHit->Track();
134         Float_t IdPart = mHit->Particle();
135         
136         for (Int_t j=0;j<4;j++) 
137         {
138           Int_t kch=11+j;
139           if (hittrack==1 || hittrack==2) 
140           {
141             if (Nch==kch && IdPart==-13 && NbHits[0][j]==0) NbHits[0][j]++; 
142             if (Nch==kch && IdPart==13 && NbHits[1][j]==0) NbHits[1][j]++; 
143           }       
144         }
145       }
146     } // end track loop     
147     
148     // 3/4 coincidence
149     Int_t SumNbHits=NbHits[0][0]+NbHits[0][1]+NbHits[0][2]+NbHits[0][3];
150     if (SumNbHits==3 || SumNbHits==4) CoincMuPlus=1;
151     
152     SumNbHits=NbHits[1][0]+NbHits[1][1]+NbHits[1][2]+NbHits[1][3];
153     if (SumNbHits==3 || SumNbHits==4) CoincMuMinus=1;
154     
155     if (CoincMuPlus==1 && CoincMuMinus==1) coincmuon++;         
156           
157     TString tree("D");
158     if ( readFromRP ) tree = "R";
159     
160     AliMUONVTriggerStore* triggerStore = diRec.TriggerStore(ievent,tree.Data());
161
162     AliMUONGlobalTrigger* gloTrg = triggerStore->Global();
163     
164     if (gloTrg->PairUnlikeLpt()>=1) muonlpt++;
165     if (gloTrg->PairUnlikeHpt()>=1) muonhpt++;
166
167   } // end loop on event  
168  
169   // calculate efficiency with as a ref. at least 3/4 planes fired
170   Float_t efficiencylpt,efficiencyhpt;
171   Double_t coincmu,lptmu,hptmu;  
172   Float_t error;  
173   coincmu=Double_t(coincmuon);
174   lptmu=Double_t(muonlpt); 
175   hptmu=Double_t(muonhpt); 
176
177   // output
178   fprintf(fdat,"\n"); 
179   fprintf(fdat," Number of events = %d \n",nevents);  
180   fprintf(fdat," Number of events with 3/4 coinc = %d \n",coincmuon);  
181   fprintf(fdat," Nomber of dimuons with 3/4 coinc Lpt cut = %d \n",muonlpt);  
182   fprintf(fdat," Number of dimuons with 3/4 coinc Hpt cut = %d \n",muonhpt);  
183   fprintf(fdat,"\n");
184   
185   efficiencylpt=lptmu/coincmu;  
186   error=efficiencylpt*TMath::Sqrt((lptmu+coincmu)/(lptmu*coincmu));  
187   fprintf(fdat," Efficiency Lpt cut = %4.4f +/- %4.4f\n",efficiencylpt,error); 
188  
189   efficiencyhpt=hptmu/coincmu; 
190   error=efficiencyhpt*TMath::Sqrt((hptmu+coincmu)/(hptmu*coincmu));          
191   fprintf(fdat," Efficiency Hpt cut = %4.4f +/- %4.4f\n",efficiencyhpt,error);
192       
193
194   fclose(fdat);  
195 }
196
197
198
199
200