]>
Commit | Line | Data |
---|---|---|
3f83f224 | 1 | /**************************************************************************** |
2 | * This macro estimates efficiency of matching with the TOF. * | |
3 | * TOF "Good" tracks are those originating from the primary vertex, * | |
4 | * being "good" in the ITS and having at least one digit in the TOF. * | |
5 | * (To get the list of "good" tracks one should first run * | |
6 | * AliTPCComparison.C and AliITSComparisonV2.C macros) * | |
7 | * Origin: I.Belikov, CERN, Jouri.Belikov@cern.ch * | |
8 | ****************************************************************************/ | |
9 | ||
10 | #if !defined(__CINT__) || defined(__MAKECINT__) | |
11 | #include <Riostream.h> | |
12 | #include <fstream.h> | |
13 | ||
14 | #include "AliRun.h" | |
15 | #include "AliHeader.h" | |
16 | #include "AliStack.h" | |
17 | #include "AliRunLoader.h" | |
18 | #include "AliLoader.h" | |
19 | ||
20 | #include "AliESD.h" | |
21 | #include "AliTOFdigit.h" | |
22 | ||
23 | #include "TKey.h" | |
24 | #include "TFile.h" | |
25 | #include "TTree.h" | |
26 | #include "TH1.h" | |
27 | #include "TClonesArray.h" | |
28 | #include "TStyle.h" | |
29 | #include "TCanvas.h" | |
30 | #include "TLine.h" | |
31 | #include "TText.h" | |
32 | #include "TParticle.h" | |
33 | #endif | |
34 | ||
35 | struct GoodTrackTOF { | |
36 | Int_t lab; | |
37 | Int_t code; | |
38 | Float_t px,py,pz; | |
39 | Float_t x,y,z; | |
40 | }; | |
41 | ||
42 | extern AliRun *gAlice; | |
43 | ||
44 | Int_t AliTOFComparison() { | |
45 | Int_t good_tracks_tof(GoodTrackTOF *gt, const Int_t max); | |
46 | ||
47 | cerr<<"Doing comparison...\n"; | |
48 | ||
49 | if (gAlice) { | |
50 | delete gAlice->GetRunLoader(); | |
51 | delete gAlice;//if everything was OK here it is already NULL | |
52 | gAlice = 0x0; | |
53 | } | |
54 | ||
55 | TClonesArray dummy("AliTOFdigit",10000), *digits=&dummy; | |
56 | ||
57 | { | |
58 | AliRunLoader *rl = AliRunLoader::Open("galice.root"); | |
59 | if (rl == 0x0) { | |
60 | cerr<<"Can not open session"<<endl; | |
61 | return 1; | |
62 | } | |
63 | AliLoader* tofl = rl->GetLoader("TOFLoader"); | |
64 | if (tofl == 0x0) { | |
65 | cerr<<"Can not get the TOF loader"<<endl; | |
66 | return 2; | |
67 | } | |
68 | tofl->LoadDigits("read"); | |
69 | ||
70 | rl->GetEvent(0); | |
71 | ||
72 | TTree *tofTree=tofl->TreeD(); | |
73 | if (!tofTree) { | |
74 | cerr<<"Can't get the TOF cluster tree !\n"; | |
75 | return 3; | |
76 | } | |
77 | TBranch *branch=tofTree->GetBranch("TOF"); | |
78 | if (!branch) { | |
79 | cerr<<"Can't get the branch with the TOF digits !\n"; | |
80 | return 4; | |
81 | } | |
82 | branch->SetAddress(&digits); | |
83 | ||
84 | tofTree->GetEvent(0); | |
85 | ||
86 | delete rl; | |
87 | } | |
88 | ||
89 | Int_t nd=digits->GetEntriesFast(); | |
90 | cerr<<"Number of digits: "<<nd<<endl; | |
91 | ||
92 | const Int_t MAX=15000; | |
93 | GoodTrackTOF gt[MAX]; | |
94 | Int_t ngood=0; | |
95 | ifstream in("good_tracks_tof"); | |
96 | if (in) { | |
97 | cerr<<"Reading good tracks...\n"; | |
98 | while (in>>gt[ngood].lab>>gt[ngood].code>> | |
99 | gt[ngood].px>>gt[ngood].py>>gt[ngood].pz>> | |
100 | gt[ngood].x >>gt[ngood].y >>gt[ngood].z) { | |
101 | ngood++; | |
102 | cerr<<ngood<<'\r'; | |
103 | if (ngood==MAX) { | |
104 | cerr<<"Too many good tracks !\n"; | |
105 | break; | |
106 | } | |
107 | } | |
108 | if (!in.eof()) cerr<<"Read error (good_tracks_its) !\n"; | |
109 | } else { | |
110 | cerr<<"Marking good tracks (this will take a while)...\n"; | |
111 | ngood=good_tracks_tof(gt,MAX); | |
112 | ofstream out("good_tracks_tof"); | |
113 | if (out) { | |
114 | for (Int_t ngd=0; ngd<ngood; ngd++) | |
115 | out<<gt[ngd].lab<<' '<<gt[ngd].code<<' ' | |
116 | <<gt[ngd].px<<' '<<gt[ngd].py<<' '<<gt[ngd].pz<<' ' | |
117 | <<gt[ngd].x <<' '<<gt[ngd].y <<' '<<gt[ngd].z <<endl; | |
118 | } else cerr<<"Can not open file (good_tracks_tof) !\n"; | |
119 | out.close(); | |
120 | } | |
121 | ||
122 | Double_t pmin=0.2; | |
123 | Double_t pmax=4.0; | |
124 | ||
125 | TH1F *hgood=new TH1F("hgood","Good tracks",30,pmin,pmax); | |
126 | TH1F *hfound=new TH1F("hfound","Matched tracks",30,pmin,pmax); | |
127 | TH1F *hfake=new TH1F("hfake","Mismatched tracks",30,pmin,pmax); | |
128 | TH1F *hgp=new TH1F("hgp","",30,pmin,pmax); //efficiency for good tracks | |
129 | hgp->SetLineColor(4); hgp->SetLineWidth(2); | |
130 | TH1F *hfp=new TH1F("hfp","Probability of mismatching",30,pmin,pmax); | |
131 | hfp->SetFillColor(1); hfp->SetFillStyle(3013); hfp->SetLineWidth(2); | |
132 | ||
133 | TH1F *hgoo=new TH1F("hgoo","Good tracks",30,-1,1); | |
134 | TH1F *hfoun=new TH1F("hfoun","Matched tracks",30,-1,1); | |
135 | TH1F *hfak=new TH1F("hfak","Mismatched tracks",30,-1,1); | |
136 | TH1F *hgl=new TH1F("hgl","",30,-1,1); //efficiency for good tracks | |
137 | hgl->SetLineColor(4); hgl->SetLineWidth(2); | |
138 | TH1F *hfl=new TH1F("hfl","Probability of mismatching",30,-1,1); | |
139 | hfl->SetFillColor(1); hfl->SetFillStyle(3013); hfl->SetLineWidth(2); | |
140 | ||
141 | TFile *ef=TFile::Open("AliESDs.root"); | |
142 | if (!ef->IsOpen()) {cerr<<"Can't AliESDs.root !\n"; return 1;} | |
143 | ||
144 | TIter next(ef->GetListOfKeys()); | |
145 | TKey *key=0; | |
146 | Int_t nev=0; | |
147 | while ((key=(TKey*)next())!=0) { | |
148 | cerr<<"Processing event number : "<<nev++<<endl; | |
149 | ||
150 | AliESD *event=(AliESD*)key->ReadObj(); | |
151 | Int_t ntrk=event->GetNumberOfTracks(); | |
152 | cerr<<"Number of ESD tracks : "<<ntrk<<endl; | |
153 | ||
154 | Int_t matched=0; | |
155 | Int_t mismatched=0; | |
156 | for (Int_t i=0; i<ngood; i++) { | |
157 | Int_t lab=gt[i].lab; | |
158 | Double_t pxg=gt[i].px, pyg=gt[i].py, pzg=gt[i].pz; | |
159 | Double_t ptg=TMath::Sqrt(pxg*pxg+pyg*pyg); | |
160 | ||
161 | if (ptg<0.1) continue; | |
162 | ||
163 | Double_t tgl=pzg/ptg; //tan(lambda) | |
164 | ||
165 | if (ptg>pmin) { hgood->Fill(ptg); hgoo->Fill(tgl); } | |
166 | ||
167 | Int_t j; | |
168 | AliESDtrack *t=0; | |
169 | for (j=0; j<ntrk; j++) { | |
170 | AliESDtrack *tt=event->GetTrack(j); | |
171 | if (lab!=TMath::Abs(tt->GetLabel())) continue; | |
172 | t=tt; | |
173 | //if ((tt->GetStatus()&AliESDtrack::kTOFpid) == 0) continue; | |
174 | if (tt->GetTOFsignal() < 0) continue; | |
175 | UInt_t idx=tt->GetTOFcluster(); | |
176 | if ((Int_t)idx>=nd) { | |
177 | cerr<<"Wrong digit index ! "<<idx<<endl; | |
178 | return 5; | |
179 | } | |
180 | AliTOFdigit *dig=(AliTOFdigit*)digits->UncheckedAt(idx); | |
181 | Int_t *label=dig->GetTracks(); | |
182 | if (label[0]!=lab) | |
183 | if (label[1]!=lab) | |
184 | if (label[2]!=lab) { | |
185 | mismatched++; | |
186 | if (ptg>pmin) { hfake->Fill(ptg); hfak->Fill(tgl); } | |
187 | break; | |
188 | } | |
189 | if (ptg>pmin) { hfound->Fill(ptg); hfoun->Fill(tgl); } | |
190 | matched++; | |
191 | break; | |
192 | } | |
193 | if (j==ntrk) { | |
194 | cerr<<"Not matched: "<<lab<<" "; | |
195 | if (t) { | |
196 | cerr<<(t->GetStatus()&AliESDtrack::kITSout)<<' ' | |
197 | <<(t->GetStatus()&AliESDtrack::kTPCout)<<' ' | |
198 | <<(t->GetStatus()&AliESDtrack::kTRDout)<<' ' | |
199 | <<(t->GetStatus()&AliESDtrack::kTIME); | |
200 | } else cerr<<"No ESD track !"; | |
201 | cerr<<endl; | |
202 | } | |
203 | } | |
204 | ||
205 | cerr<<"Number of good tracks: "<<ngood<<endl; | |
206 | cerr<<"Number of matched tracks: "<<matched<<endl; | |
207 | cerr<<"Number of mismatched tracks: "<<mismatched<<endl; | |
208 | if (ngood!=0) cerr<<"Efficiency: "<<Float_t(matched)/ngood<<endl; | |
209 | ||
210 | hfound->Sumw2(); hgood->Sumw2(); hfake->Sumw2(); | |
211 | hgp->Divide(hfound,hgood,1,1.,"b"); | |
212 | hfp->Divide(hfake,hgood,1,1.,"b"); | |
213 | hgp->SetMaximum(1.4); | |
214 | hgp->SetYTitle("Matching efficiency"); | |
215 | hgp->SetXTitle("Pt (GeV/c)"); | |
216 | ||
217 | hfoun->Sumw2(); hgoo->Sumw2(); hfak->Sumw2(); | |
218 | hgl->Divide(hfoun,hgoo,1,1.,"b"); | |
219 | hfl->Divide(hfak,hgoo,1,1.,"b"); | |
220 | hgl->SetMaximum(1.4); | |
221 | hgl->SetYTitle("Matching efficiency"); | |
222 | hgl->SetXTitle("Tan(lambda)"); | |
223 | ||
224 | TCanvas *c1=new TCanvas("c1","",0,0,600,900); | |
225 | c1->Divide(1,2); | |
226 | ||
227 | c1->cd(1); | |
228 | ||
229 | hgp->Draw(); | |
230 | hfp->Draw("histsame"); | |
231 | TLine *line1 = new TLine(pmin,1.0,pmax,1.0); line1->SetLineStyle(4); | |
232 | line1->Draw("same"); | |
233 | TLine *line2 = new TLine(pmin,0.9,pmax,0.9); line2->SetLineStyle(4); | |
234 | line2->Draw("same"); | |
235 | ||
236 | c1->cd(2); | |
237 | ||
238 | hgl->Draw(); | |
239 | hfl->Draw("histsame"); | |
240 | TLine *line3 = new TLine(-1,1.0,1,1.0); line3->SetLineStyle(4); | |
241 | line3->Draw("same"); | |
242 | TLine *line4 = new TLine(-1,0.9,1,0.9); line4->SetLineStyle(4); | |
243 | line4->Draw("same"); | |
244 | ||
245 | break; | |
246 | } | |
247 | ||
248 | return 0; | |
249 | } | |
250 | ||
251 | Int_t good_tracks_tof(GoodTrackTOF *gt, const Int_t max) { | |
252 | ifstream in("good_tracks_its"); | |
253 | if (!in) { | |
254 | cerr<<"Can't get good_tracks_its !\n"; exit(11); | |
255 | } | |
256 | ||
257 | AliRunLoader *rl = AliRunLoader::Open("galice.root","COMPARISON"); | |
258 | if (!rl) { | |
259 | cerr<<"Can't start session !\n"; | |
260 | exit(1); | |
261 | } | |
262 | ||
263 | rl->GetEvent(0); | |
264 | ||
265 | ||
266 | rl->LoadgAlice(); | |
267 | rl->LoadHeader(); | |
268 | Int_t np = rl->GetHeader()->GetNtrack(); | |
269 | ||
270 | Int_t *good=new Int_t[np]; | |
271 | Int_t k; | |
272 | for (k=0; k<np; k++) good[k]=0; | |
273 | ||
274 | AliLoader* tofl = rl->GetLoader("TOFLoader"); | |
275 | if (tofl == 0x0) { | |
276 | cerr<<"Can not get the TOF loader"<<endl; | |
277 | exit(2); | |
278 | } | |
279 | tofl->LoadDigits("read"); | |
280 | ||
281 | TTree *dTree=tofl->TreeD(); | |
282 | if (!dTree) { | |
283 | cerr<<"Can't get the TOF cluster tree !\n"; | |
284 | exit(3); | |
285 | } | |
286 | ||
287 | TBranch *branch=dTree->GetBranch("TOF"); | |
288 | if (!branch) { | |
289 | cerr<<"Can't get the branch with the TOF digits !\n"; | |
290 | exit(4); | |
291 | } | |
292 | TClonesArray dummy("AliTOFdigit",10000), *digits=&dummy; | |
293 | branch->SetAddress(&digits); | |
294 | ||
295 | dTree->GetEvent(0); | |
296 | Int_t nd=digits->GetEntriesFast(); | |
297 | cerr<<"Number of digits: "<<nd<<endl; | |
298 | ||
299 | for (Int_t i=0; i<nd; i++) { | |
300 | AliTOFdigit *d=(AliTOFdigit*)digits->UncheckedAt(i); | |
301 | Int_t l0=d->GetTrack(0); | |
302 | if (l0>=np) {cerr<<"Wrong label: "<<l0<<endl; continue;} | |
303 | Int_t l1=d->GetTrack(1); | |
304 | if (l1>=np) {cerr<<"Wrong label: "<<l1<<endl; continue;} | |
305 | Int_t l2=d->GetTrack(2); | |
306 | if (l2>=np) {cerr<<"Wrong label: "<<l2<<endl; continue;} | |
307 | if (l0>=0) good[l0]++; | |
308 | if (l1>=0) good[l1]++; | |
309 | if (l2>=0) good[l2]++; | |
310 | } | |
311 | ||
312 | ||
313 | rl->LoadKinematics(); | |
314 | AliStack* stack = rl->Stack(); | |
315 | Int_t nt=0; | |
316 | Double_t px,py,pz,x,y,z; | |
317 | Int_t code,lab; | |
318 | while (in>>lab>>code>>px>>py>>pz>>x>>y>>z) { | |
319 | if (good[lab] == 0) continue; | |
320 | TParticle *p = (TParticle*)stack->Particle(lab); | |
321 | if (p == 0x0) { | |
322 | cerr<<"Can not get particle "<<lab<<endl; | |
323 | exit(1); | |
324 | } | |
325 | if (TMath::Abs(p->Vx())>0.1) continue; | |
326 | if (TMath::Abs(p->Vy())>0.1) continue; | |
327 | if (TMath::Abs(p->Vz())>0.1) continue; | |
328 | ||
329 | gt[nt].lab=lab; | |
330 | gt[nt].code=p->GetPdgCode(); | |
331 | //**** px py pz - in global coordinate system | |
332 | gt[nt].px=p->Px(); gt[nt].py=p->Py(); gt[nt].pz=p->Pz(); | |
333 | gt[nt].x=p->Vx(); gt[nt].y=p->Vy(); gt[nt].z=p->Vz(); | |
334 | nt++; | |
335 | if (nt==max) {cerr<<"Too many good tracks !\n"; break;} | |
336 | } | |
337 | ||
338 | delete[] good; | |
339 | ||
340 | rl->UnloadKinematics(); | |
341 | rl->UnloadHeader(); | |
342 | rl->UnloadgAlice(); | |
343 | delete rl; | |
344 | ||
345 | return nt; | |
346 | } |