]>
Commit | Line | Data |
---|---|---|
50f3bbf4 | 1 | |
2 | #include "AliTRDtrackingAnalysis.h" | |
3 | ||
4 | ||
5 | #include "TFile.h" | |
6 | #include "TTree.h" | |
7 | #include "TH1D.h" | |
8 | #include "TH2D.h" | |
9 | #include "TObjArray.h" | |
10 | #include "TCanvas.h" | |
11 | #include "TGeoMatrix.h" | |
12 | #include "TStyle.h" | |
13 | #include "TGraphErrors.h" | |
14 | #include "TF1.h" | |
74cc0efb | 15 | #include "TMath.h" |
50f3bbf4 | 16 | |
17 | #include "AliRunLoader.h" | |
18 | #include "AliTRDgeometry.h" | |
19 | #include "AliRun.h" | |
20 | #include "AliESD.h" | |
21 | #include "AliESDtrack.h" | |
22 | #include "AliTrackReference.h" | |
23 | ||
24 | #include "AliTRDcluster.h" | |
25 | #include "AliTRDCommonParam.h" | |
26 | #include "AliTRDpadPlane.h" | |
27 | #include "AliTRDcalibDB.h" | |
28 | #include "AliTracker.h" | |
29 | #include "AliTRDtracker.h" | |
30 | //#include "AliTRDtracklet.h" | |
31 | ||
32 | #include "TGeoManager.h" | |
33 | ||
34 | ////////////////////////////////////////////////////////////////////////////////////////// | |
35 | ||
36 | AliTRDtrackingAnalysis::AliTRDtrackingAnalysis(): | |
37 | TObject(), | |
ae6817b3 | 38 | fPath(0), |
39 | fRefTPC(0), | |
40 | fRefTRD(0), | |
50f3bbf4 | 41 | fLoader(0), |
42 | fEsdTree(0), | |
43 | fESD(0), | |
ae6817b3 | 44 | fTracker(0), |
45 | fDeltaPt(0), | |
46 | fDeltaZ(0), | |
47 | fDeltaX(0), | |
48 | fDeltaYPos(0), | |
49 | fDeltaYNeg(0), | |
50 | fNPoints(0), | |
51 | fNGood(0), | |
52 | fRefSpace(0), | |
53 | fGeo(0), | |
54 | fClY2(0), | |
55 | fClY3(0), | |
56 | fTgPhi(0), | |
57 | fGrResTgPhi(0), | |
58 | fGrMeanTgPhi(0), | |
59 | fTrklY(0), | |
60 | fTrklZ(0), | |
61 | fClZ(0), | |
62 | fClZZ(0), | |
63 | fClYY(0), | |
64 | fClYX(0), | |
65 | fNLabels(0), | |
66 | fBits(0), | |
67 | fRefDx(0), | |
68 | fClZXref(0), | |
69 | fClZXcl(0), | |
70 | fClPos(0) | |
50f3bbf4 | 71 | { |
72 | ||
73 | fDeltaX = new TH1D("deltaX", ";delta X (cm)", 100, -1, 1); | |
74 | fDeltaZ = new TH1D("deltaZ", ";delta Z (cm)", 100, -2, 2); | |
75 | ||
76 | fDeltaYPos = new TH1D("deltaYpos", ";delta Y (mm)", 100, -1, 1); | |
77 | fDeltaYNeg = new TH1D("deltaYneg", ";delta Y (mm)", 100, -1, 1); | |
78 | ||
79 | fNPoints = new TH1D("nPoints", ";np", 40, -0.5, 39.5); | |
80 | fNGood = new TH1D("nGood", ";np", 40, -0.5, 39.5); | |
81 | ||
82 | fDeltaPt = new TH1D("deltaPt", ";delta Pt/Pt (%)", 100, -10, 10); | |
83 | fRefSpace = new TH2D("refSpace", ";y;x", 120, -60, 60, 200, -4, 1); | |
84 | ||
85 | fTrklY = new TH1D("trklY", ";delta Y (mm)", 100, -1, 1); | |
86 | fTrklZ = new TH1D("trklZ", ";delta Z (cm)", 100, -10, 10); | |
87 | ||
88 | ||
89 | // cluster studies | |
90 | fClY2 = new TH1D("clY2", ";delta Y (mm)", 100, -10, 10); | |
91 | fClY3 = new TH1D("clY3", ";delta Y (mm)", 100, -10, 10); | |
92 | ||
93 | for(int i=0; i<12; i++) // bewere hidden constants in the code | |
94 | fClYTgPhi[i] = new TH1D(Form("clYtgPhi%d", i), ";delta Y (mm)", 100, -3, 3); | |
95 | ||
96 | fTgPhi = new TH1D("tgPhi", ";Tg(#phi)", 100, -0.3, 0.3); | |
97 | fGrResTgPhi = new TGraphErrors(); | |
98 | fGrMeanTgPhi = new TGraphErrors(); | |
99 | ||
100 | //fPullY2 = new TH1D("pullY2", ";pulls Y", 100, -5, 5); | |
101 | //fPullY3 = new TH1D("pullY3", ";pulls Y", 100, -5, 5); | |
102 | ||
103 | ||
104 | fClZ = new TH1D("clZ", ";delta Z (cm)", 200, -20, 20); | |
105 | fClZZ = new TH2D("clZZ", ";z ref;z cl", 600, -300, 300, 600, -300, 300); | |
106 | ||
107 | fClYY = new TH2D("clYY", ";dY;dY", 100, -3, 3, 100, -3, 3); | |
108 | fClYX = new TH2D("clYX", ";Y;X", 250, -60, 60, 100, -4, 1); | |
109 | ||
110 | fNLabels = new TH1D("clLabels", ";n labels", 10, -0.5, 9.5); | |
111 | fBits = new TH1D("bits", ";bits", 10, -0.5, 9.5); | |
112 | ||
113 | fRefDx = new TH1D("refDX", ";delta X", 100, 0, 20); | |
114 | fClPos = new TH2D("clPos", ";z;y", 400, -400, 400, 120, -60, 60); | |
115 | ||
116 | fClZXref = new TH2D("clZXref", ";z;x", 36, -54, 54, 300, 280, 380); | |
117 | fClZXcl = new TH2D("clZXcl", ";z;x", 36, -54, 54, 300, 280, 380); | |
118 | ||
119 | //fGeo = new AliTRDgeometry(); | |
120 | } | |
121 | ||
122 | ////////////////////////////////////////////////////////////////////////////////////////// | |
123 | ||
124 | void AliTRDtrackingAnalysis::DrawResolutionPt(int startEvent, int stopEvent) { | |
125 | ||
126 | CheckFiles(); | |
127 | ||
128 | // loop over ESD events | |
129 | int nevents = fEsdTree->GetEntries(); | |
130 | ||
131 | for(int iEvent=startEvent; iEvent<nevents && iEvent < stopEvent; iEvent++) { | |
132 | ||
133 | Info("Draw", "Event = %d", iEvent); | |
134 | ||
135 | fEsdTree->GetEvent(iEvent); | |
136 | fLoader->GetEvent(iEvent); | |
137 | LoadRefs(); | |
138 | ||
139 | int nTracks = fESD->GetNumberOfTracks(); | |
140 | for(int iTrack=0; iTrack<nTracks; iTrack++) { | |
141 | ||
142 | //Info("Track", "Track = %d", iTrack); | |
143 | AliESDtrack *esdTrack = fESD->GetTrack(iTrack); | |
144 | if (!esdTrack->GetInnerParam()) continue; | |
145 | const AliExternalTrackParam *param = esdTrack->GetOuterParam(); | |
146 | int status = esdTrack->GetStatus(); | |
147 | ||
148 | if (!(status & AliESDtrack::kTRDout)) continue; | |
149 | if (!(status & AliESDtrack::kTRDrefit)) continue; | |
74cc0efb | 150 | if (TMath::Abs(esdTrack->GetOuterParam()->GetPt()) < 1.0) continue; |
50f3bbf4 | 151 | |
152 | int ch=0; | |
153 | while(param->GetX() > fGeo->GetTime0(ch)+2) ch++; | |
154 | fRefSpace->Fill(2.*ch+0.5, param->GetX() - fGeo->GetTime0(ch)); | |
155 | //if (ch < 5) continue; | |
156 | ||
157 | double lastX = 0; | |
158 | int label = abs(esdTrack->GetTRDLabel()); | |
159 | int ntr = 0; | |
160 | int ngood = 0; | |
161 | ||
162 | for(int iPoint=GetReference(label); iPoint<fRefTRD->GetEntries(); iPoint++) { | |
163 | ||
164 | AliTrackReference *aRef = (AliTrackReference*)(*fRefTRD)[iPoint]; | |
165 | if (aRef->GetTrack() != label) break; | |
166 | ntr++; | |
167 | ||
168 | lastX = (aRef->LocalX() < lastX)? lastX : aRef->LocalX(); | |
169 | double dx = aRef->LocalX() - param->GetX(); | |
74cc0efb | 170 | if (TMath::Abs(dx) > 1.) continue; |
50f3bbf4 | 171 | ngood++; |
172 | ||
173 | double bz=fESD->GetMagneticField(); | |
174 | AliExternalTrackParam out(*param); | |
175 | out.PropagateTo(aRef->LocalX(),bz); | |
176 | ||
177 | double dp = aRef->Pt() + out.GetPt(); | |
178 | double dy = 10. * (aRef->LocalY() - out.GetY()); // in mm | |
179 | ||
180 | fDeltaPt->Fill(100. * dp / aRef->Pt()); | |
181 | //fDeltaPt->Fill(out.GetPt() / aRef->Pt()); | |
182 | fDeltaX->Fill(dx); | |
183 | ||
184 | if (esdTrack->GetSign() > 0) fDeltaYPos->Fill(dy); | |
185 | else fDeltaYNeg->Fill(dy); | |
186 | ||
187 | fDeltaZ->Fill(aRef->Z() - out.GetZ()); | |
188 | } | |
189 | ||
190 | //if (ngood == 0) Info("X", "N = %d, X = %f, DX = %f", ntr, param->GetX(), param->GetX()-lastX); | |
191 | ||
192 | fNPoints->Fill(ntr); | |
193 | fNGood->Fill(ngood); | |
194 | } | |
195 | } | |
196 | ||
197 | new TCanvas(); | |
198 | fDeltaPt->Draw(); | |
199 | ||
200 | new TCanvas(); | |
201 | fDeltaX->Draw(); | |
202 | ||
203 | /* | |
204 | new TCanvas(); | |
205 | fNPoints->Draw(); | |
206 | ||
207 | new TCanvas(); | |
208 | fNGood->Draw(); | |
209 | */ | |
210 | ||
211 | new TCanvas(); | |
212 | fDeltaYPos->Draw(); | |
213 | ||
214 | new TCanvas(); | |
215 | fDeltaYNeg->Draw(); | |
216 | ||
217 | new TCanvas(); | |
218 | fDeltaZ->Draw(); | |
219 | ||
220 | new TCanvas(); | |
221 | fRefSpace->Draw(); | |
222 | ||
223 | } | |
224 | ||
225 | ////////////////////////////////////////////////////////////////////////////////////////// | |
226 | ||
227 | // void AliTRDtrackingAnalysis::DrawTrackletResolution(int startEvent, int stopEvent) { | |
228 | ||
229 | // LoadRecPointsFile(); | |
230 | ||
231 | // TFile *file = new TFile(Form("%s/TRD.Tracklets.root", fPath), "read"); | |
232 | // TTree *tree = (TTree*)file->Get("TRDtracklets"); | |
233 | ||
234 | // AliTRDtracklet *tracklet = new AliTRDtracklet(); | |
235 | // tree->SetBranchAddress("tracklets", &tracklet); | |
236 | ||
237 | // for(int ev=startEvent; ev<stopEvent; ev++) { | |
238 | ||
239 | // gAlice->GetEvent(ev); | |
240 | // LoadRefs(); | |
241 | ||
242 | // int N = tree->GetEntries(); | |
243 | // for(int i=0; i<N; i++) { | |
244 | ||
245 | // tree->GetEntry(i); | |
246 | ||
247 | // Double_t yref, zref, tgphi; | |
248 | // int stat = GetMCPosition(tracklet->GetLabel(), tracklet->GetX(), yref, zref, tgphi); | |
249 | // if (stat < 0) continue; | |
250 | ||
251 | // int plane = tracklet->GetPlane(); | |
252 | // Double_t h01 = tracklet->GetTilt(); | |
253 | ||
254 | // //printf("Tile = %f\tcorrection = %f um \n", h01, 1e4 * h01 * (tracklet->GetZ()-zref)); | |
255 | // //double dz = zref - tracklet->GetZ() - cls->GetZ(); | |
256 | ||
257 | // fTrklY->Fill(10 * (tracklet->GetY() - yref)); | |
258 | // fTrklZ->Fill(tracklet->GetZ() - zref); | |
259 | ||
260 | // int ch=0; | |
261 | // while(tracklet->GetX() > fGeo->GetTime0(ch)+2) ch++; | |
262 | // fRefSpace->Fill(tracklet->GetY(), tracklet->GetX() - fGeo->GetTime0(ch)); | |
263 | // } | |
264 | // } | |
265 | ||
266 | // new TCanvas(); | |
267 | // fTrklZ->Draw(); | |
268 | ||
269 | // new TCanvas(); | |
270 | // fRefSpace->Draw(); | |
271 | ||
272 | // gStyle->SetOptFit(1); | |
273 | // new TCanvas(); | |
274 | // fTrklY->Draw(); | |
275 | // fTrklY->Fit("gaus"); | |
276 | ||
277 | // //new TCanvas(); | |
278 | // //fClZ->Draw(); | |
279 | // } | |
280 | ||
281 | ////////////////////////////////////////////////////////////////////////////////////////// | |
282 | ||
283 | void AliTRDtrackingAnalysis::DrawRecPointResolution(int startEvent, int stopEvent) { | |
284 | ||
285 | LoadRecPointsFile(); | |
286 | TObjArray *module = new TObjArray(); | |
287 | ||
288 | int nEvents = gAlice->GetEventsPerRun(); | |
289 | ||
290 | for(int ev=startEvent; ev<nEvents && ev < stopEvent; ev++) { | |
291 | ||
292 | gAlice->GetEvent(ev); | |
293 | LoadRefs(); | |
294 | ||
295 | TTree *tree = fLoader->GetTreeR("TRD", 0); | |
296 | tree->SetBranchAddress("TRDcluster", &module); | |
297 | ||
298 | Info("Res", "Refs Loaded"); | |
299 | ||
300 | int N = tree->GetEntries(); | |
301 | for(int i=0; i<N; i++) { | |
302 | ||
303 | tree->GetEntry(i); | |
304 | int m = module->GetEntries(); | |
305 | ||
306 | for(int j=0; j<m; j++) { | |
307 | ||
308 | AliTRDcluster *cls = (AliTRDcluster*)module->At(j); | |
309 | if (cls->GetQ() < 10) continue; | |
310 | fTracker->Transform(cls); | |
311 | fClPos->Fill(cls->GetZ(), cls->GetY()); | |
312 | ||
313 | int plane = fGeo->GetPlane(cls->GetDetector()); | |
314 | ||
315 | int nl = 0; | |
316 | for(int k=0; k<3; k++) if (cls->GetLabel(k) > -1) nl++; | |
317 | fNLabels->Fill(nl); | |
318 | ||
319 | Double_t yref, zref, tgphi; | |
320 | int stat = GetMCPosition(cls->GetLabel(0), cls->GetX(), yref, zref, tgphi); | |
321 | if (stat < 0) continue; | |
322 | ||
323 | fClZXcl->Fill(cls->GetZ(), cls->GetX()); | |
324 | fClZXref->Fill(zref, cls->GetX()); | |
325 | ||
326 | AliTRDpadPlane *padPlane = AliTRDCommonParam::Instance()->GetPadPlane(plane,0); | |
327 | Double_t h01 = TMath::Tan(-TMath::Pi() / 180.0 * padPlane->GetTiltingAngle()); | |
328 | ||
329 | //double dz = zref - padPlane->GetRow0(); | |
330 | double dz = zref - cls->GetZ(); | |
331 | double dy = dz * h01; | |
332 | double yy = cls->GetY() - dy; | |
333 | ||
334 | if (cls->From2pad()) fClY2->Fill(10 * (yy - yref)); | |
335 | if (cls->From3pad()) fClY3->Fill(10 * (yy - yref)); | |
336 | ||
337 | int idx = GetPhiBin(tgphi); | |
338 | if (idx >= 0 && idx < 12) fClYTgPhi[idx]->Fill(10 * (yy - yref)); | |
339 | ||
340 | fClZZ->Fill(zref, cls->GetZ()); | |
341 | fClZ->Fill(dz); | |
342 | fTgPhi->Fill(tgphi); | |
343 | fClYX->Fill(cls->GetY(), cls->GetX() - fGeo->GetTime0(plane)); | |
344 | } | |
345 | } | |
346 | } | |
347 | ||
348 | new TCanvas(); | |
349 | fClYX->Draw(); | |
350 | ||
351 | //new TCanvas(); | |
352 | //fNLabels->Draw(); | |
353 | new TCanvas(); | |
354 | fClZ->Draw(); | |
355 | ||
356 | new TCanvas(); | |
357 | fTgPhi->Draw(); | |
358 | ||
359 | gStyle->SetOptFit(1); | |
360 | ||
361 | new TCanvas(); | |
362 | gPad->SetLogy(); | |
363 | fClY2->Draw(); | |
364 | fClY2->Fit("gaus", "", "", -2*fClY2->GetRMS(), 2*fClY2->GetRMS()); | |
365 | ||
366 | new TCanvas(); | |
367 | gPad->SetLogy(); | |
368 | fClY3->Draw(); | |
369 | fClY3->Fit("gaus", "", "", -2*fClY3->GetRMS(), 2*fClY3->GetRMS()); | |
370 | ||
371 | //new TCanvas(); | |
372 | //fRefSpace->Draw(); | |
373 | ||
374 | //new TCanvas(); | |
375 | //fClZXcl->Draw(); | |
376 | ||
377 | //new TCanvas(); | |
378 | //fClZXref->Draw(); | |
379 | ||
380 | /**/ | |
381 | TCanvas *c = new TCanvas(); | |
382 | c->Divide(4,3); | |
383 | ||
384 | for(int i=0; i<12; i++) { | |
385 | ||
386 | c->cd(i+1); | |
387 | fClYTgPhi[i]->Draw(); | |
388 | if (fClYTgPhi[i]->GetSum() < 100) continue; | |
389 | ||
390 | double mean = fClYTgPhi[i]->GetMean(); | |
391 | double rms = fClYTgPhi[i]->GetRMS(); | |
392 | ||
393 | fClYTgPhi[i]->Fit("gaus", "", "", mean-2*rms, mean+2*rms); | |
394 | TF1 *f = fClYTgPhi[i]->GetFunction("gaus"); | |
395 | ||
396 | int n = fGrResTgPhi->GetN(); | |
397 | fGrResTgPhi->SetPoint(n, GetPhi(i), f->GetParameter(2)); | |
398 | fGrResTgPhi->SetPointError(n, 0, f->GetParError(2)); | |
399 | ||
400 | fGrMeanTgPhi->SetPoint(n, GetPhi(i), f->GetParameter(1)); | |
401 | fGrMeanTgPhi->SetPointError(n, 0, f->GetParError(1)); | |
402 | } | |
403 | ||
404 | //gSystem->Sleep(1000); | |
405 | gStyle->SetOptStat(0); | |
406 | ||
407 | c = new TCanvas(); | |
408 | TH1D *dummy = new TH1D("dummy", "", 100, -0.3, 0.3); | |
409 | ||
410 | dummy->SetTitle(";tg(#phi);resolution (mm)"); | |
411 | dummy->SetMinimum(0); | |
412 | dummy->SetMaximum(1); | |
413 | ||
414 | //c->cd(); | |
415 | (dummy->Clone("dummy1"))->Draw(); | |
416 | ||
417 | fGrResTgPhi->Draw("PL"); | |
418 | fGrResTgPhi->GetHistogram()->SetTitle(";tg(#phi);resolution (mm)"); | |
419 | fGrResTgPhi->SetMarkerStyle(20); | |
420 | ||
421 | c = new TCanvas(); | |
422 | dummy->SetTitle(";tg(#phi);mean value (mm)"); | |
423 | dummy->SetMinimum(-0.3); | |
424 | dummy->SetMaximum(0.3); | |
425 | ||
426 | //c->cd(); | |
427 | (dummy->Clone("dummy2"))->Draw(); | |
428 | ||
429 | fGrMeanTgPhi->Draw("PL"); | |
430 | fGrMeanTgPhi->GetHistogram()->SetTitle(";tg(#phi);mean value (mm)"); | |
431 | fGrMeanTgPhi->SetMarkerStyle(20); | |
432 | /**/ | |
433 | ||
434 | ||
435 | //new TCanvas(); | |
436 | //fClZZ->Draw("colz"); | |
437 | ||
438 | //new TCanvas(); | |
439 | //fClPos->Draw("colz"); | |
440 | ||
441 | //new TCanvas(); | |
442 | //fBits->Draw(); | |
443 | ||
444 | //new TCanvas(); | |
445 | //fRefDx->Draw(); | |
446 | ||
447 | } | |
448 | ||
449 | ////////////////////////////////////////////////////////////////////////////////////////// | |
450 | ||
451 | void AliTRDtrackingAnalysis::LoadRecPointsFile() { | |
452 | ||
453 | char filename[256]; | |
454 | sprintf(filename, "%s/galice.root", fPath); | |
455 | ||
456 | fLoader = AliRunLoader::Open(filename); | |
457 | if (!fLoader) { | |
458 | Error("CheckFiles", "getting run loader from file %s/galice.root failed", filename); | |
459 | return; | |
460 | } | |
461 | ||
462 | fLoader->LoadgAlice(); | |
463 | gAlice = fLoader->GetAliRun(); | |
464 | ||
465 | if (!gAlice) { | |
466 | Error("CheckFiles", "no galice object found"); | |
467 | return; | |
468 | } | |
469 | ||
470 | fLoader->LoadKinematics(); | |
471 | fLoader->LoadHeader(); | |
472 | fLoader->LoadTrackRefs(); | |
473 | ||
474 | //TGeoManager::Import("/data/alice_u/radomski/condor/run_0/geometry.root"); | |
475 | TGeoManager::Import(Form("%s/geometry.root", fPath)); | |
476 | ||
477 | ||
478 | fLoader->CdGAFile(); | |
479 | fGeo = (AliTRDgeometry*)gDirectory->Get("TRDgeometry"); | |
480 | fTracker = new AliTRDtracker(gFile); | |
481 | ||
482 | fLoader->LoadRecPoints("TRD"); | |
483 | ||
484 | AliTracker::SetFieldMap(gAlice->Field(), 1); | |
485 | } | |
486 | ||
487 | ////////////////////////////////////////////////////////////////////////////////////////// | |
488 | ||
489 | void AliTRDtrackingAnalysis::CheckFiles() { | |
490 | ||
491 | // MC info | |
492 | ||
493 | char filename[256]; | |
494 | sprintf(filename, "%s/galice.root", fPath); | |
495 | ||
496 | fLoader = AliRunLoader::Open(filename); | |
497 | if (!fLoader) { | |
498 | Error("CheckFiles", "getting run loader from file %s/galice.root failed", filename); | |
499 | return; | |
500 | } | |
501 | ||
502 | fLoader->LoadgAlice(); | |
503 | gAlice = fLoader->GetAliRun(); | |
504 | ||
505 | if (!gAlice) { | |
506 | Error("CheckFiles", "no galice object found"); | |
507 | return; | |
508 | } | |
509 | ||
510 | fLoader->LoadKinematics(); | |
511 | fLoader->LoadHeader(); | |
512 | fLoader->LoadTrackRefs(); | |
513 | ||
514 | fLoader->CdGAFile(); | |
515 | fGeo = (AliTRDgeometry*)gDirectory->Get("TRDgeometry"); | |
516 | //fGeo->ReadGeoMatrices(); | |
517 | ||
518 | // ESD | |
519 | ||
520 | sprintf(filename,"%s/AliESDs.root", fPath); | |
521 | TFile *esdFile = new TFile(filename, "READ"); | |
522 | ||
523 | if (esdFile->IsZombie()) { | |
524 | Error("CheckFiles", "file not present: AliESDs.root"); | |
525 | return; | |
526 | } | |
527 | ||
528 | fEsdTree = (TTree*)esdFile->Get("esdTree"); | |
529 | fESD = new AliESD(); | |
530 | ||
531 | fEsdTree->SetBranchAddress("ESD", &fESD); | |
532 | } | |
533 | ||
534 | ////////////////////////////////////////////////////////////////////////////////////////// | |
535 | ||
536 | void AliTRDtrackingAnalysis::LoadRefs() { | |
537 | ||
538 | if (fRefTPC) delete fRefTPC; | |
539 | if (fRefTRD) delete fRefTRD; | |
540 | ||
541 | fRefTPC = new TObjArray(); | |
542 | fRefTRD = new TObjArray(); | |
543 | ||
544 | //fLoader->GetEvent(event); | |
545 | //AliStack* stack = gAlice->Stack(); | |
546 | TTree *refTree = fLoader->TreeTR(); | |
547 | ||
548 | const int nBranch = 2; | |
549 | const char *brName[] = {"TPC", "TRD"}; | |
550 | TClonesArray *clRefs = new TClonesArray("AliTrackReference"); | |
551 | ||
552 | for(int b=0; b<nBranch; b++) { | |
553 | ||
554 | TBranch *branch = refTree->GetBranch(brName[b]); | |
555 | refTree->SetBranchAddress(brName[b],&clRefs); | |
556 | ||
557 | int nEntries = branch->GetEntries(); | |
558 | for(int iTrack = 0; iTrack < nEntries; iTrack++) { | |
559 | ||
560 | refTree->GetEvent(iTrack); | |
561 | int nPoints = clRefs->GetEntries(); | |
562 | for(int iPoint=0; iPoint<nPoints; iPoint++) { | |
563 | AliTrackReference *ref = (AliTrackReference*)clRefs->At(iPoint); | |
564 | if (b == 0) fRefTPC->Add(new AliTrackReference(*ref)); | |
565 | if (b == 1) fRefTRD->Add(new AliTrackReference(*ref)); | |
566 | } | |
567 | } | |
568 | } | |
569 | ||
570 | fRefTPC->Sort(); | |
571 | fRefTRD->Sort(); | |
572 | ||
573 | for(int i=0; i<fRefTRD->GetEntries(); i++) { | |
574 | AliTrackReference *ref = (AliTrackReference*)(*fRefTRD)[i]; | |
575 | fLabels[i] = ref->GetTrack(); | |
576 | ||
577 | int p=0; | |
578 | while(ref->LocalX() > fGeo->GetTime0(p)+2) p++; | |
579 | fRefSpace->Fill(ref->LocalY(), ref->LocalX()-fGeo->GetTime0(p)); | |
580 | ||
581 | //for(int bit=0; bit<9; bit++) if (ref->TestBit(bit)) fBits->Fill(bit); | |
582 | } | |
583 | ||
584 | delete clRefs; | |
585 | Info("LoadRefs", "TPC = %d\t TRD = %d", fRefTPC->GetEntries(), fRefTRD->GetEntries()); | |
586 | } | |
587 | ||
588 | ////////////////////////////////////////////////////////////////////////////////////////// | |
589 | ||
590 | Int_t AliTRDtrackingAnalysis::GetReference(Int_t label) { | |
591 | ||
592 | int start = TMath::BinarySearch(fRefTRD->GetEntries(), fLabels, label); | |
593 | ||
594 | while (start >= 0) { | |
595 | AliTrackReference *ref = (AliTrackReference*)(*fRefTRD)[start]; | |
596 | if (ref->GetTrack() != label) return start+1; | |
597 | start--; | |
598 | } | |
599 | ||
600 | return 0; | |
601 | } | |
602 | ||
603 | ////////////////////////////////////////////////////////////////////////////////////////// | |
604 | ||
605 | int AliTRDtrackingAnalysis::GetMCPosition(Int_t label, Double_t x, Double_t &Y, Double_t &Z, Double_t &tgphi) { | |
606 | ||
607 | double lowX = 100.; | |
608 | double highX = 100.; | |
609 | int idLow = -1; | |
610 | int idHigh = -1; | |
611 | ||
612 | int nref= 0; | |
613 | int idx = GetReference(label); | |
614 | for(int i=idx; i<fRefTRD->GetEntries(); i++) { | |
615 | ||
616 | AliTrackReference *ref = (AliTrackReference*)(*fRefTRD)[i]; | |
617 | if (ref->GetTrack() != label) break; | |
618 | nref++; | |
619 | ||
620 | //int p=0; | |
621 | //while(ref->LocalX() > fGeo->GetTime0(p)+2) p++; | |
622 | //if (p != layer) continue; | |
623 | ||
624 | double dX = ref->LocalX()-x; | |
625 | if ( dX > 0 ) { | |
626 | if (dX < highX) { | |
627 | idHigh = i; | |
628 | highX = dX; | |
629 | } | |
630 | } else { | |
74cc0efb | 631 | dX = TMath::Abs(dX); |
50f3bbf4 | 632 | if (dX < lowX) { |
633 | idLow = i; | |
634 | lowX = dX; | |
635 | } | |
636 | } | |
637 | } | |
638 | ||
639 | if (idLow == -1 || idHigh == -1) return -1; | |
640 | ||
641 | AliTrackReference *refI = (AliTrackReference*)(*fRefTRD)[idLow]; | |
642 | AliTrackReference *refO = (AliTrackReference*)(*fRefTRD)[idHigh]; | |
643 | ||
644 | ||
645 | double dx = refO->LocalX() - refI->LocalX(); | |
646 | double dy = refO->LocalY() - refI->LocalY(); | |
647 | double dz = refO->Z() - refI->Z(); | |
648 | double ddx = (x - refI->LocalX())/dx; | |
649 | ||
650 | fRefDx->Fill(dx); | |
651 | ||
652 | Y = refI->LocalY() + ddx * dy; | |
653 | Z = refI->Z() + ddx * dz; | |
654 | ||
655 | tgphi = dy/dx; | |
656 | ||
657 | return 0; | |
658 | } | |
659 | ||
660 | ////////////////////////////////////////////////////////////////////////////////////////// | |
661 | Int_t AliTRDtrackingAnalysis::GetPhiBin(Double_t phi) { | |
662 | return (int)((phi+0.3)/0.05); | |
663 | } | |
664 | ||
665 | ////////////////////////////////////////////////////////////////////////////////////////// | |
666 | ||
667 | Double_t AliTRDtrackingAnalysis::GetPhi(Int_t bin) { | |
668 | return bin * 0.05 - 0.3 + 0.025; | |
669 | } | |
670 | ////////////////////////////////////////////////////////////////////////////////////////// |