]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TRD/AliTRDtrackerDebug.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Remove obsolete classes
[u/mrichter/AliRoot.git] / TRD / AliTRDtrackerDebug.cxx
CommitLineData
eb38ed55 1/**************************************************************************
e3cf3d02 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**************************************************************************/
eb38ed55 15
16/* $Id: AliTRDtrackerDebug.cxx 23810 2008-02-08 09:00:27Z hristov $ */
17
18///////////////////////////////////////////////////////////////////////////////
19// //
20// Tracker debug streamer //
21// //
22// Authors: //
23// Alex Bercuci <A.Bercuci@gsi.de> //
24// Markus Fasel <M.Fasel@gsi.de> //
25// //
26///////////////////////////////////////////////////////////////////////////////
27
28#include "AliTRDtrackerDebug.h"
29
30#include "TFile.h"
31#include "TTree.h"
32#include "TTreeStream.h"
33#include "TLinearFitter.h"
bb56afff 34#include "TGraph.h"
35#include "TCanvas.h"
36#include "TMath.h"
eb38ed55 37
38#include "AliLog.h"
39#include "AliTRDgeometry.h"
40#include "AliTRDtrackV1.h"
41#include "AliTRDseedV1.h"
42#include "AliTRDseed.h"
43#include "AliTRDcluster.h"
bb56afff 44#include "AliTRDgeometry.h"
eb38ed55 45
46ClassImp(AliTRDtrackerDebug)
47
bb56afff 48Int_t AliTRDtrackerDebug::fgEventNumber = 0;
49Int_t AliTRDtrackerDebug::fgTrackNumber = 0;
50Int_t AliTRDtrackerDebug::fgCandidateNumber = 0;
51
eb38ed55 52//____________________________________________________
53AliTRDtrackerDebug::AliTRDtrackerDebug() : AliTRDtrackerV1()
e3cf3d02 54 ,fOutputStreamer(0x0)
55 ,fTree(0x0)
56 ,fTracklet(0x0)
57 ,fTrack(0x0)
58 ,fNClusters(0)
59 ,fAlpha(0.)
eb38ed55 60{
61 //
e3cf3d02 62 // Default constructor
63 //
64 fOutputStreamer = new TTreeSRedirector("TRD.Debug.root");
eb38ed55 65}
66
67//____________________________________________________
68AliTRDtrackerDebug::~AliTRDtrackerDebug()
69{
e3cf3d02 70 // destructor
71
72 delete fOutputStreamer;
eb38ed55 73}
74
75
76//____________________________________________________
c202e52e 77void AliTRDtrackerDebug::Draw(Option_t *)
eb38ed55 78{
79// steer draw function
80}
81
82
83//____________________________________________________
84Bool_t AliTRDtrackerDebug::Init()
85{
86// steer linking data for various debug streams
e3cf3d02 87 fTrack = new AliTRDtrackV1();
88 fTree->SetBranchAddress("ncl", &fNClusters);
89 fTree->SetBranchAddress("track.", &fTrack);
90 return kTRUE;
eb38ed55 91}
92
93//____________________________________________________
94Bool_t AliTRDtrackerDebug::Open(const char *method)
95{
e3cf3d02 96 // Connect to the tracker debug file
97
98 TDirectory *savedir = gDirectory;
99 TFile::Open("TRD.TrackerDebugger.root");
100 fTree = (TTree*)gFile->Get(method);
101 if(!fTree){
102 AliInfo(Form("Can not find debug stream for the %s method.\n", method));
103 savedir->cd();
104 return kFALSE;
105 }
106 savedir->cd();
107 return kTRUE;
eb38ed55 108}
109
110//____________________________________________________
111Int_t AliTRDtrackerDebug::Process()
112{
113// steer debug process threads
e3cf3d02 114
115 for(int it = 0; it<fTree->GetEntries(); it++){
116 if(!fTree->GetEntry(it)) continue;
117 if(!fNClusters) continue;
118 fAlpha = fTrack->GetAlpha();
119 //printf("Processing track %d [%d] ...\n", it, fNClusters);
120 ResidualsTrackletsTrack();
121
122 const AliTRDseedV1 *tracklet = 0x0;
123 for(int ip = 5; ip>=0; ip--){
124 if(!(tracklet = fTrack->GetTracklet(ip))) continue;
125 if(!tracklet->GetN()) continue;
126
127 ResidualsClustersTrack(tracklet);
128 ResidualsClustersTracklet(tracklet);
129 ResidualsClustersParametrisation(tracklet);
130 }
131 }
132 return kTRUE;
eb38ed55 133}
134
135
136//____________________________________________________
137void AliTRDtrackerDebug::ResidualsClustersTrack(const AliTRDseedV1 *tracklet)
138{
139// Calculate averange distances from clusters to the TRD track
e3cf3d02 140
141 Double_t x[3];
142 AliTRDcluster *c = 0x0;
143 for(int ic=0; ic<35/*AliTRDseed:knTimebins*/; ic++){
144 if(!(c = tracklet->GetClusters(ic))) continue;
145 Double_t xc = c->GetX(), yc = c->GetY(), zc = c->GetZ();
146
147 // propagate track to cluster
148 PropagateToX(*fTrack, xc, 2.);
149 fTrack->GetXYZ(x);
150
151 // transform to local tracking coordinates
152 //Double_t xg = x[0] * TMath::Cos(fAlpha) + x[1] * TMath::Sin(fAlpha);
153 Double_t yg = -x[0] * TMath::Sin(fAlpha) + x[1] * TMath::Cos(fAlpha);
154
155 // apply tilt pad correction
156 yc+= (zc - x[2]) * tracklet->GetTilt();
157
158 Double_t dy = yc-yg;
159
160 TTreeSRedirector &cstreamer = *fOutputStreamer;
161 cstreamer << "ResidualsClustersTrack"
162 << "c.=" << c
163 << "dy=" << dy
164 << "\n";
165 }
eb38ed55 166}
167
168//____________________________________________________
169void AliTRDtrackerDebug::ResidualsClustersTracklet(const AliTRDseedV1 *tracklet) const
170{
171// Calculates distances from clusters to tracklets
e3cf3d02 172
173 Double_t x0 = tracklet->GetX0(),
174 y0 = tracklet->GetYfit(0),
175 ys = tracklet->GetYfit(1);
176 //z0 = tracklet->GetZfit(0),
177 //zs = tracklet->GetZfit(1);
178
179 AliTRDcluster *c = 0x0;
180 for(int ic=0; ic<35/*AliTRDseed:knTimebins*/; ic++){
181 if(!(c = tracklet->GetClusters(ic))) continue;
182 Double_t xc = c->GetX(), yc = c->GetY()/*, zc = c->GetZ()*/;
183 Double_t dy = yc- (y0-(x0-xc)*ys);
184
185 //To draw use :
eb38ed55 186 //ResidualsClustersTracklet->Draw("TMath::Abs(10.*dy):TMath::ATan(ys)*TMath::RadToDeg()>>h(20, -40, 40)", "", "prof");
e3cf3d02 187 TTreeSRedirector &cstreamer = *fOutputStreamer;
188 cstreamer << "ResidualsClustersTracklet"
189 << "c.=" << c
190 << "ys=" << ys
191 << "dy=" << dy
192 << "\n";
193 }
eb38ed55 194}
195
196//____________________________________________________
197void AliTRDtrackerDebug::ResidualsClustersParametrisation(const AliTRDseedV1 *tracklet) const
198{
199// Calculates distances from clusters to Rieman fit.
e3cf3d02 200
201 // store cluster positions
202 Double_t x0 = tracklet->GetX0();
203 AliTRDcluster *c = 0x0;
204
205 Double_t x[2]; Int_t ncl, mcl, jc;
206 TLinearFitter fitter(3, "hyp2");
207 for(int ic=0; ic<35/*AliTRDseed:knTimebins*/; ic++){
208 if(!(c = tracklet->GetClusters(ic))) continue;
209 Double_t xc = c->GetX(), yc = c->GetY()/*, zc = c->GetZ()*/;
210
211 jc = ic; ncl = 0; mcl=0; fitter.ClearPoints();
212 while(ncl<6){
213 // update index
214 mcl++;
215 jc = ic + ((mcl&1)?-1:1)*(mcl>>1);
216
217 if(jc<0 || jc>=35) continue;
218 if(!(c = tracklet->GetClusters(jc))) continue;
219
220 x[0] = c->GetX()-x0;
221 x[1] = x[0]*x[0];
222 fitter.AddPoint(x, c->GetY(), c->GetSigmaY2());
223 ncl++;
224 }
225 fitter.Eval();
226 Double_t dy = yc - fitter.GetParameter(0) -fitter.GetParameter(1) * (xc-x0) - fitter.GetParameter(2)* (xc-x0)*(xc-x0);
227
228 TTreeSRedirector &cstreamer = *fOutputStreamer;
229 cstreamer << "ResidualsClustersParametrisation"
230 << "dy=" << dy
231 << "\n";
232 }
eb38ed55 233}
234
235
236//____________________________________________________
237void AliTRDtrackerDebug::ResidualsTrackletsTrack() const
238{
239// Calculates distances from tracklets to the TRD track.
e3cf3d02 240
241 if(fTrack->GetNumberOfTracklets() < 6) return;
242
243 // build a working copy of the tracklets attached to the track
244 // and initialize working variables fX, fY and fZ
245 AliTRDseedV1 tracklet[6] = {0x0, 0x0, 0x0, 0x0, 0x0, 0x0};
246 const AliTRDseedV1 *ctracklet = 0x0;
247 for(int ip = 0; ip<6; ip++){
248 if(!(ctracklet = fTrack->GetTracklet(ip))) continue;
249 tracklet[ip] = (*ctracklet);
250// Double_t x0 = tracklet[ip].GetX0();
251// for(int ic=0; ic<AliTRDseedV1:knTimebins; ic++){
252// if(!(c = tracklet[ip].GetClusters(ic))) continue;
253// Double_t xc = c->GetX(), yc = c->GetY(), zc = c->GetZ();
254// tracklet[ip].SetX(ic, xc-x0);
255// tracklet[ip].SetY(ic, yc);
256// tracklet[ip].SetZ(ic, zc);
257// }
258 }
259
260 // Do a Rieman fit (with tilt correction) for all tracklets
261 // except the one which is tested.
262 // (Based on AliTRDseed::IsOK() return false)
263 for(int ip=0; ip<6; ip++){
264 // reset tracklet to be tested
265 Double_t x0 = tracklet[ip].GetX0();
266 new(&tracklet[ip]) AliTRDseedV1();
267 tracklet[ip].SetX0(x0);
268
269 // fit Rieman with tilt correction
270 AliTRDtrackerV1::FitRiemanTilt(0x0, &tracklet[0], kTRUE);
271
272 // make a copy of the fit result
273 Double_t
274 y0 = tracklet[ip].GetYref(0),
275 dydx = tracklet[ip].GetYref(1),
276 z0 = tracklet[ip].GetZref(0),
277 dzdx = tracklet[ip].GetZref(1);
278
279 // restore tracklet
280 tracklet[ip] = (*fTrack->GetTracklet(ip));
281// for(int ic=0; ic<AliTRDseedV1:knTimebins; ic++){
282// if(!(c = tracklet[ip].GetClusters(ic))) continue;
283// Double_t xc = c->GetX(), yc = c->GetY(), zc = c->GetZ();
284// tracklet[ip].SetX(ic, xc-x0);
285// tracklet[ip].SetY(ic, yc);
286// tracklet[ip].SetZ(ic, zc);
287// }
288
289 // fit clusters
290 AliTRDseedV1 ts(tracklet[ip]);
291 ts.SetYref(0, y0); ts.SetYref(1, dydx);
292 ts.SetZref(0, z0); ts.SetZref(1, dzdx);
f29f13a6 293 ts.Fit(kTRUE);
e3cf3d02 294
295 // save results for plotting
296 Int_t plane = tracklet[ip].GetPlane();
297 Double_t dy = tracklet[ip].GetYfit(0) - ts.GetYfit(0);
298 Double_t tgy = tracklet[ip].GetYfit(1);
299 Double_t dtgy = (tracklet[ip].GetYfit(1) - ts.GetYfit(1))/(1. + tracklet[ip].GetYfit(1) * ts.GetYfit(1));
300 Double_t dz = tracklet[ip].GetZfit(0) - ts.GetZfit(0);
301 Double_t tgz = tracklet[ip].GetZfit(1);
302 Double_t dtgz = (tracklet[ip].GetZfit(1) - ts.GetZfit(1))/(1. + tracklet[ip].GetZfit(1) * ts.GetZfit(1));
303 TTreeSRedirector &cstreamer = *fOutputStreamer;
304 cstreamer << "ResidualsTrackletsTrack"
305 << "ref.=" << &tracklet[ip]
306 << "fit.=" << &ts
307 << "plane=" << plane
308 << "dy=" << dy
309 << "tgy=" << tgy
310 << "dtgy=" << dtgy
311 << "dz=" << dz
312 << "tgz=" << tgz
313 << "dtgz=" << dtgz
314 << "\n";
315 }
eb38ed55 316}
317
bb56afff 318//____________________________________________________
319void AliTRDtrackerDebug::AnalyseFindable(Char_t *treename){
320//
321// Calculates the number of findable tracklets defined as the number of tracklets
322// per track candidate where the tan phi_tracklet is below 0.15 (maximum inclination
323// in y-direction.
324//
325// Parameters: -the treename (this method can be used for all trees which store the
326// tracklets
327// Output: -void
328//
329// A new TTree containing the number of findable tracklets and the number of clusters
330// attached to the full track is stored to disk
331//
e3cf3d02 332 // Link the File
333 TFile *debfile = TFile::Open("TRD.TrackerDebug.root");
334 fTree = (TTree *)(debfile->Get(treename));
335 if(!fTree){
336 AliError(Form("Tree %s not found in file TRDdebug.root. Abborting!", treename));
337 debfile->Close();
338 return;
339 }
340
341 AliTRDseedV1 *tracklets[kNPlanes];
342 for(Int_t iPlane = 0; iPlane < AliTRDtrackerV1::kNPlanes; iPlane++)
343 tracklets[iPlane] = 0x0;
344 for(Int_t iPlane = 0; iPlane < kNPlanes; iPlane++)
345 fTree->SetBranchAddress(Form("S%d.", iPlane), &tracklets[iPlane]);
346 fTree->SetBranchAddress("EventNumber", &fgEventNumber);
347 fTree->SetBranchAddress("CandidateNumber", &fgCandidateNumber);
348
349 Int_t findable = 0, nClusters = 0;
350 Int_t nEntries = fTree->GetEntriesFast();
351 for(Int_t iEntry = 0; iEntry < nEntries; iEntry++){
352 printf("Entry %d\n", iEntry);
353 fTree->GetEntry(iEntry);
354 findable = 0;
355 nClusters = 0;
356 // Calculate Findable
357 for(Int_t iPlane = 0; iPlane < kNPlanes; iPlane++){
358 if (TMath::Abs(tracklets[iPlane]->GetYref(0) / tracklets[iPlane]->GetX0()) < 0.15) findable++;
359 if (!tracklets[iPlane]->IsOK()) continue;
360 nClusters += tracklets[iPlane]->GetN2();
361 }
362
363 // Fill Histogramms
364 TTreeSRedirector &cstreamer = *fOutputStreamer;
365 cstreamer << "AnalyseFindable"
366 << "EventNumber=" << fgEventNumber
367 << "CandidateNumber=" << fgCandidateNumber
368 << "Findable=" << findable
369 << "NClusters=" << nClusters
370 << "\n";
371 }
bb56afff 372}
373//____________________________________________________
374void AliTRDtrackerDebug::AnalyseTiltedRiemanFit(){
375//
376// Creating a Data Set for the method FitTiltedRieman containing usefull variables
377// Each variable can be addressed to tracks later. Data can be processed later.
378//
379// Parameters: -
380// Output: -
381//
382// TODO: Match everything with Init and Process
383//
e3cf3d02 384 TFile *debfile = TFile::Open("TRD.TrackerDebug.root");
385 fTree = (TTree *)(debfile->Get("MakeSeeds2"));
386 if(!fTree) return;
387 Int_t nEntries = fTree->GetEntries();
388 TLinearFitter *tiltedRiemanFitter = 0x0;
389 fTree->SetBranchAddress("FitterT.", &tiltedRiemanFitter);
390 fTree->SetBranchAddress("EventNumber", &fgEventNumber);
391 fTree->SetBranchAddress("CandidateNumber", &fgCandidateNumber);
392 for(Int_t entry = 0; entry < nEntries; entry++){
393 fTree->GetEntry(entry);
394 Double_t a = tiltedRiemanFitter->GetParameter(0);
395 Double_t b = tiltedRiemanFitter->GetParameter(1);
396 Double_t c = tiltedRiemanFitter->GetParameter(2);
397 Double_t offset = tiltedRiemanFitter->GetParameter(3);
398 Double_t slope = tiltedRiemanFitter->GetParameter(4);
399 Float_t radius = GetTrackRadius(a, b, c);
400 Float_t curvature = GetTrackCurvature(a, b, c);
401 Float_t dca = GetDCA(a, b, c);
402 TTreeSRedirector &cstreamer = *fOutputStreamer;
403 cstreamer << "AnalyseTiltedRiemanFit"
404 << "EventNumber=" << fgEventNumber
405 << "CandidateNumber=" << fgCandidateNumber
406 << "Radius=" << radius
407 << "Curvature=" << curvature
408 << "DCA=" << dca
409 << "Offset=" << offset
410 << "Slope=" << slope
411 << "\n";
412 }
bb56afff 413}
414
415//____________________________________________________
416Float_t AliTRDtrackerDebug::GetTrackRadius(Float_t a, Float_t b, Float_t c) const {
417//
418// Calculates the track radius using the parameters given by the tilted Rieman fit
419//
420// Parameters: The three parameters from the Rieman fit
421// Output: The track radius
422//
e3cf3d02 423 Float_t radius = 0;
424 if(1.0 + b*b - c*a > 0.0)
425 radius = TMath::Sqrt(1.0 + b*b - c*a )/a;
426 return radius;
bb56afff 427}
428
429//____________________________________________________
430Float_t AliTRDtrackerDebug::GetTrackCurvature(Float_t a, Float_t b, Float_t c) const {
431//
432// Calculates the track curvature using the parameters given by the linear fitter
433//
434// Parameters: the three parameters from the tilted Rieman fitter
435// Output: the full track curvature
436//
e3cf3d02 437 Float_t curvature = 1.0 + b*b - c*a;
438 if (curvature > 0.0)
439 curvature = a / TMath::Sqrt(curvature);
440 return curvature;
bb56afff 441}
442
443//____________________________________________________
444Float_t AliTRDtrackerDebug::GetDCA(Float_t a, Float_t b, Float_t c) const {
445//
446// Calculates the Distance to Clostest Approach for the Vertex using the paramters
447// given by the tilted Rieman fit
448//
449// Parameters: the three parameters from the tilted Rieman fitter
450// Output: the Distance to Closest Approach
451//
e3cf3d02 452 Float_t dca = 0.0;
453 if (1.0 + b*b - c*a > 0.0) {
454 dca = -c / (TMath::Sqrt(1.0 + b*b - c*a) + TMath::Sqrt(1.0 + b*b));
455 }
456 return dca;
bb56afff 457}
458
459//____________________________________________________
460void AliTRDtrackerDebug::AnalyseMinMax()
461{
462//
e3cf3d02 463 TFile *debfile = TFile::Open("TRD.TrackerDebug.root");
464 if(!debfile){
465 AliError("File TRD.TrackerDebug.root not found!");
466 return;
467 }
468 fTree = (TTree *)(debfile->Get("MakeSeeds0"));
469 if(!fTree){
470 AliError("Tree MakeSeeds0 not found in File TRD.TrackerDebug.root.");
471 return;
472 }
473 AliTRDseedV1 *cseed[4] = {0x0, 0x0, 0x0, 0x0};
474 AliTRDcluster *c[4] = {0x0, 0x0, 0x0, 0x0};
475 for(Int_t il = 0; il < 4; il++){
476 fTree->SetBranchAddress(Form("Seed%d.", il), &cseed[il]);
477 fTree->SetBranchAddress(Form("c%d.",il), &c[il]);
478 }
479 fTree->SetBranchAddress("CandidateNumber", &fgCandidateNumber);
480 fTree->SetBranchAddress("EventNumber", &fgEventNumber);
481 Int_t entries = fTree->GetEntries();
482 for(Int_t ientry = 0; ientry < entries; ientry++){
483 fTree->GetEntry(ientry);
484 Float_t minmax[2] = { -100.0, 100.0 };
485 for (Int_t iLayer = 0; iLayer < 4; iLayer++) {
486 Float_t max = c[iLayer]->GetZ() + cseed[iLayer]->GetPadLength() * 0.5 + 1.0 - cseed[iLayer]->GetZref(0);
487 if (max < minmax[1]) minmax[1] = max;
488 Float_t min = c[iLayer]->GetZ()-cseed[iLayer]->GetPadLength() * 0.5 - 1.0 - cseed[iLayer]->GetZref(0);
489 if (min > minmax[0]) minmax[0] = min;
490 }
491 TTreeSRedirector &cstreamer = *fOutputStreamer;
492 cstreamer << "AnalyseMinMaxLayer"
493 << "EventNumber=" << fgEventNumber
494 << "CandidateNumber=" << fgCandidateNumber
495 << "Min=" << minmax[0]
496 << "Max=" << minmax[1]
497 << "\n";
498 }
bb56afff 499}
500
501//____________________________________________________
c202e52e 502TCanvas* AliTRDtrackerDebug::PlotSeedingConfiguration(const Char_t *direction, Int_t event, Int_t candidate){
bb56afff 503//
504// Plots the four seeding clusters, the helix fit and the reference Points for
505// a given combination consisting of eventnr. and candidatenr.
506//
507// Parameters: -direction (y or z)
508// -Event Nr
509// -Candidate that has to be plotted
510//
e3cf3d02 511 const Float_t kxmin = 280;
512 const Float_t kxmax = 380;
513 const Float_t kxdelta = (kxmax - kxmin)/1000;
514
515 if((strcmp(direction, "y") != 0) && (strcmp(direction, "z") != 0)){
516 AliError(Form("Direction %s does not exist. Abborting!", direction));
517 return 0x0;
518 }
519
520 TFile *debfile = TFile::Open("TRD.TrackerDebug.root");
521 if(!debfile){
522 AliError("File TRD.TrackerDebug.root not found!");
523 return 0x0;
524 }
525 fTree = (TTree *)(debfile->Get("MakeSeeds0"));
526 if(!fTree){
527 AliError("Tree MakeSeeds0 not found in File TRD.TrackerDebug.root.");
528 return 0x0;
529 }
530
531 TGraph *seedcl = new TGraph(4);
532 TGraph *seedRef = new TGraph(4);
533 TGraph *riemanFit = new TGraph(1000);
534 seedcl->SetMarkerStyle(20);
535 seedcl->SetMarkerColor(kRed);
536 seedRef->SetMarkerStyle(2);
537
538 AliTRDcluster *c[4] = {0x0, 0x0, 0x0, 0x0};
539 AliRieman *rim = 0x0;
540 Bool_t found = kFALSE;
541 for(Int_t il = 0; il < 4; il++) fTree->SetBranchAddress(Form("c%d.",il), &c[il]);
542 fTree->SetBranchAddress("EventNumber", &fgEventNumber);
543 fTree->SetBranchAddress("CandidateNumber", &fgCandidateNumber);
544 fTree->SetBranchAddress("RiemanFitter.", &rim);
545 Int_t entries = fTree->GetEntries();
546 for(Int_t entry = 0; entry < entries; entry++){
547 fTree->GetEntry(entry);
548 if(fgEventNumber < event) continue;
549 if(fgEventNumber > event) break;
550 // EventNumber matching: Do the same for the candidate number
551 if(fgCandidateNumber < candidate) continue;
552 if(fgCandidateNumber > candidate) break;
553 found = kTRUE;
554 Int_t nPoints = 0;
555 for(Int_t il = 0; il < 4; il++){
556 Float_t cluster = 0.0, reference = 0.0;
557 if(!strcmp(direction, "y")){
558 cluster = c[il]->GetY();
559 reference = rim->GetYat(c[il]->GetX());
560 }
561 else{
562 cluster = c[il]->GetZ();
563 reference = rim->GetZat(c[il]->GetX());
564 }
565 seedcl->SetPoint(nPoints, cluster, c[il]->GetX());
566 seedRef->SetPoint(nPoints, reference , c[il]->GetX());
567 nPoints++;
568 }
569 // evaluate the fitter Function numerically
570 nPoints = 0;
571 for(Int_t ipt = 0; ipt < 1000; ipt++){
572 Float_t x = kxmin + ipt * kxdelta;
573 Float_t point = 0.0;
574 if(!strcmp(direction, "y"))
575 point = rim->GetYat(x);
576 else
577 point = rim->GetZat(x);
578 riemanFit->SetPoint(nPoints++, point, x);
579 }
580 // We reached the End: break
581 break;
582 }
583 if(found){
584 seedcl->SetTitle(Form("Event %d, Candidate %d\n", fgEventNumber, fgCandidateNumber));
585 seedRef->SetTitle(Form("Event %d, Candidate %d\n", fgEventNumber, fgCandidateNumber));
586 riemanFit->SetTitle(Form("Event %d, Candidate %d\n", fgEventNumber, fgCandidateNumber));
587 TCanvas *c1 = new TCanvas();
588 seedcl->Draw("ap");
589 seedRef->Draw("psame");
590 riemanFit->Draw("lpsame");
591 return c1;
592 }
593 else{
594 AliError(Form("Combination consisting of event %d and candidate %d not found", event, candidate));
595 delete seedcl;
596 delete seedRef;
597 delete riemanFit;
598 return 0x0;
599 }
bb56afff 600}
601
602//____________________________________________________
c202e52e 603TCanvas* AliTRDtrackerDebug::PlotFullTrackFit(Int_t event, Int_t candidate, Int_t iteration, const Char_t *direction){
bb56afff 604//
605// Plots the tracklets (clusters and reference in y direction) and the fitted function for several iterations
606// in the function ImproveSeedQuality (default is before ImproveSeedQuality)
607//
608// Parameters: -Event Number
609// -Candidate Number
610// -Iteration Number in ImproveSeedQuality (default: -1 = before ImproveSeedQuality)
611// -direction (default: y)
612// Output: -TCanvas (containing the Picture);
613//
e3cf3d02 614 const Float_t kxmin = 280;
615 const Float_t kxmax = 380;
616 const Float_t kxdelta = (kxmax - kxmin)/1000;
617
618 if(strcmp(direction, "y") && strcmp(direction, "z")){
619 AliError(Form("Direction %s does not exist. Abborting!", direction));
620 return 0x0;
621 }
622
623 TFile *debfile = TFile::Open("TRD.TrackerDebug.root");
624 if(!debfile){
625 AliError("File TRD.TrackerDebug.root not found.");
626 return 0x0;
627 }
628 TString *treename = 0x0;
629 if(iteration > -1)
630 treename = new TString("ImproveSeedQuality");
631 else
632 treename = new TString("MakeSeeds1");
633 fTree = (TTree *)(debfile->Get(treename->Data()));
634 if(!fTree){
635 AliError(Form("Tree %s not found in File TRD.TrackerDebug.root.", treename->Data()));
636 delete treename;
637 return 0x0;
638 }
639 delete treename;
640
641 TGraph *fitfun = new TGraph(1000);
642 // Prepare containers
643 Float_t x0[AliTRDtrackerV1::kNPlanes],
644 refP[AliTRDtrackerV1::kNPlanes],
645 clx[AliTRDtrackerV1::kNPlanes * AliTRDtrackerV1::kNTimeBins],
646 clp[AliTRDtrackerV1::kNPlanes * AliTRDtrackerV1::kNTimeBins];
647 Int_t nLayers = 0, ncls = 0;
648
649 TLinearFitter *fitter = 0x0;
650 AliTRDseedV1 *tracklet[6] = {0x0, 0x0, 0x0, 0x0, 0x0, 0x0};
651 for(Int_t iLayer = 0; iLayer < 6; iLayer++)
652 fTree->SetBranchAddress(Form("S%d.", iLayer), &tracklet[iLayer]);
653 fTree->SetBranchAddress("FitterT.", &fitter);
654 fTree->SetBranchAddress("EventNumber", &fgEventNumber);
655 fTree->SetBranchAddress("CandidateNumber", &fgCandidateNumber);
656
657 Int_t nEntries = fTree->GetEntriesFast();
658 Bool_t found = kFALSE;
659 for(Int_t entry = 0; entry < nEntries; entry++){
660 fTree->GetEntry(entry);
661 if(fgEventNumber < event) continue;
662 if(fgEventNumber > event) break;
663 // EventNumber matching: Do the same for the candidate number
664 if(fgCandidateNumber < candidate) continue;
665 if(fgCandidateNumber > candidate) break;
666 found = kTRUE;
667
668 for(Int_t iLayer = 0; iLayer < 6; iLayer++){
669 if(!tracklet[iLayer]->IsOK()) continue;
670 x0[nLayers] = tracklet[iLayer]->GetX0();
671 if(!strcmp(direction, "y"))
672 refP[nLayers] = tracklet[iLayer]->GetYref(0);
673 else
674 refP[nLayers] = tracklet[iLayer]->GetZref(0);
675 nLayers++;
676 for(Int_t itb = 0; itb < 30; itb++){
677 if(!tracklet[iLayer]->IsUsable(itb)) continue;
678 AliTRDcluster *cl = tracklet[iLayer]->GetClusters(itb);
679 if(!strcmp(direction, "y"))
680 clp[ncls] = cl->GetY();
681 else
682 clp[ncls] = cl->GetZ();
683 clx[ncls] = cl->GetX();
684 ncls++;
685 }
686 }
687 // Add function derived by the tilted Rieman fit (Defined by the curvature)
688 Int_t nPoints = 0;
689 if(!strcmp(direction, "y")){
690 Double_t a = fitter->GetParameter(0);
691 Double_t b = fitter->GetParameter(1);
692 Double_t c = fitter->GetParameter(2);
693 Double_t curvature = 1.0 + b*b - c*a;
694 if (curvature > 0.0) {
695 curvature = a / TMath::Sqrt(curvature);
696 }
697 // Numerical evaluation of the function:
698 for(Int_t ipt = 0; ipt < 1000; ipt++){
699 Float_t x = kxmin + ipt * kxdelta;
700 Double_t res = (x * a + b); // = (x - x0)/y0
701 res *= res;
702 res = 1.0 - c * a + b * b - res; // = (R^2 - (x - x0)^2)/y0^2
703 Double_t y = 0.;
704 if (res >= 0) {
705 res = TMath::Sqrt(res);
706 y = (1.0 - res) / a;
707 }
708 fitfun->SetPoint(nPoints++, y, x);
709 }
710 }
711 else{
712 Double_t offset = fitter->GetParameter(3);
713 Double_t slope = fitter->GetParameter(4);
714 // calculate the reference x (defined as medium between layer 2 and layer 3)
715 // same procedure as in the tracker code
716 Float_t medx = 0, xref = 0;
717 Int_t startIndex = 5, nDistances = 0;
718 for(Int_t iLayer = 5; iLayer > 0; iLayer--){
719 if(tracklet[iLayer]->IsOK() && tracklet[iLayer - 1]->IsOK()){
720 medx += tracklet[iLayer]->GetX0() - tracklet[iLayer - 1]->GetX0();
721 startIndex = iLayer - 1;
722 nDistances++;
723 }
724 }
725 if(nDistances){
726 medx /= nDistances;
727 }
728 else{
729 Float_t xpos[2]; memset(xpos, 0, sizeof(Float_t) * 2);
730 Int_t ien = 0, idiff = 0;
731 for(Int_t iLayer = 5; iLayer > 0; iLayer--){
732 if(tracklet[iLayer]->IsOK()){
733 xpos[ien++] = tracklet[iLayer]->GetX0();
734 startIndex = iLayer;
735 }
736 if(ien)
737 idiff++;
738 if(ien >=2)
739 break;
740 }
741 medx = (xpos[0] - xpos[1])/idiff;
742 }
743 xref = tracklet[startIndex]->GetX0() + medx * (2.5 - startIndex) - 0.5 * (AliTRDgeometry::AmThick() + AliTRDgeometry::DrThick());
744
745 for(Int_t ipt = 0; ipt < 1000; ipt++){
746 Float_t x = kxmin + ipt * kxdelta;
747 Float_t z = offset + slope * (x - xref);
748 fitfun->SetPoint(nPoints++, z, x);
749 }
750 }
751 break;
752 }
753 if(found){
754 TGraph *trGraph = new TGraph(ncls);
755 TGraph *refPoints = new TGraph(nLayers);
756 trGraph->SetMarkerStyle(20);
757 trGraph->SetMarkerColor(kRed);
758 refPoints->SetMarkerStyle(21);
759 refPoints->SetMarkerColor(kBlue);
760 // fill the graphs
761 for(Int_t iLayer = 0; iLayer < nLayers; iLayer++)
762 refPoints->SetPoint(iLayer, refP[iLayer], x0[iLayer]);
763 for(Int_t icls = 0; icls < ncls; icls++)
764 trGraph->SetPoint(icls, clp[icls], clx[icls]);
765 TCanvas *c1 = new TCanvas();
766 trGraph->SetTitle(Form("Event %d, Candidate %d\n", fgEventNumber, fgCandidateNumber));
767 refPoints->SetTitle(Form("Event %d, Candidate %d\n", fgEventNumber, fgCandidateNumber));
768 fitfun->SetTitle(Form("Event %d, Candidate %d\n", fgEventNumber, fgCandidateNumber));
769 trGraph->Draw("ap");
770 refPoints->Draw("psame");
771 fitfun->Draw("lpsame");
772 return c1;
773 }
774 else{
775 AliError(Form("Combination consisting of event %d and candidate %d not found", event, candidate));
776 delete fitfun;
777 return 0x0;
778 }
bb56afff 779}
780