]>
Commit | Line | Data |
---|---|---|
1 | // $Id$ | |
2 | // Main authors: Matevz Tadel & Alja Mrak-Tadel: 2006, 2007 | |
3 | ||
4 | /************************************************************************** | |
5 | * Copyright(c) 1998-2008, ALICE Experiment at CERN, all rights reserved. * | |
6 | * See http://aliceinfo.cern.ch/Offline/AliRoot/License.html for * | |
7 | * full copyright notice. * | |
8 | **************************************************************************/ | |
9 | ||
10 | #include "TROOT.h" | |
11 | #include "TStyle.h" | |
12 | #include "TVector.h" | |
13 | #include "TLinearFitter.h" | |
14 | #include "TCanvas.h" | |
15 | #include "TGeoMatrix.h" | |
16 | #include "TGeoManager.h" | |
17 | ||
18 | #include "TEveTrans.h" | |
19 | #include "TEveManager.h" | |
20 | ||
21 | #include "EveBase/AliEveEventManager.h" | |
22 | ||
23 | #include "AliEveTRDData.h" | |
24 | #include "AliEveTRDModuleImp.h" | |
25 | #include "AliEveTRDLoader.h" | |
26 | #include "AliEveTRDLoaderImp.h" | |
27 | ||
28 | #include "AliGeomManager.h" | |
29 | #include "AliESDtrack.h" | |
30 | #include "AliLog.h" | |
31 | #include "AliPID.h" | |
32 | #include "AliTrackPointArray.h" | |
33 | #include "AliRieman.h" | |
34 | ||
35 | #include "AliTRDhit.h" | |
36 | #include "AliTRDcluster.h" | |
37 | #include "AliTRDseedV1.h" | |
38 | #include "AliTRDtrackletMCM.h" | |
39 | #include "AliTRDtrackletWord.h" | |
40 | #include "AliTRDmcmSim.h" | |
41 | #include "AliTRDtrackV1.h" | |
42 | #include "AliTRDtrackerV1.h" | |
43 | #include "AliTRDpadPlane.h" | |
44 | #include "AliTRDdigitsManager.h" | |
45 | #include "AliTRDmcmSim.h" | |
46 | #include "AliTRDarrayADC.h" | |
47 | #include "AliTRDSignalIndex.h" | |
48 | #include "AliTRDgeometry.h" | |
49 | #include "AliTRDtransform.h" | |
50 | #include "AliTRDReconstructor.h" | |
51 | #include "AliTRDrecoParam.h" | |
52 | ||
53 | ClassImp(AliEveTRDHits) | |
54 | ClassImp(AliEveTRDDigits) | |
55 | ClassImp(AliEveTRDClusters) | |
56 | ClassImp(AliEveTRDTracklet) | |
57 | ClassImp(AliEveTRDTrack) | |
58 | ClassImp(AliEveTRDTrackletOnline) | |
59 | ClassImp(AliEveTRDmcm) | |
60 | ||
61 | /////////////////////////////////////////////////////////// | |
62 | ///////////// AliEveTRDDigits ///////////////////// | |
63 | /////////////////////////////////////////////////////////// | |
64 | ||
65 | //______________________________________________________________________________ | |
66 | AliEveTRDDigits::AliEveTRDDigits(AliEveTRDChamber *p) | |
67 | :TEveQuadSet("digits", "") | |
68 | ,fParent(p) | |
69 | { | |
70 | // Constructor. | |
71 | SetOwnIds(kTRUE); | |
72 | gStyle->SetPalette(1, 0); | |
73 | SetPalette(new TEveRGBAPalette(0, 512)); | |
74 | Reset(TEveQuadSet::kQT_RectangleYZ, kFALSE, 32); | |
75 | } | |
76 | ||
77 | //______________________________________________________________________________ | |
78 | AliEveTRDDigits::~AliEveTRDDigits() | |
79 | { | |
80 | // AliInfo(GetTitle()); | |
81 | } | |
82 | ||
83 | //______________________________________________________________________________ | |
84 | void AliEveTRDDigits::SetData(AliTRDdigitsManager *digits) | |
85 | { | |
86 | // Set data source. | |
87 | ||
88 | Int_t det(fParent->GetID()); | |
89 | AliTRDarrayADC *data = digits->GetDigits(det); | |
90 | if(!data->GetDim()) return; | |
91 | data->Expand(); | |
92 | ||
93 | AliTRDSignalIndex *indexes = digits->GetIndexes(det); | |
94 | if(!indexes->IsAllocated()) digits->BuildIndexes(det); | |
95 | ||
96 | Double_t scale, dy, dz; | |
97 | Int_t ly = AliTRDgeometry::GetLayer(det), | |
98 | stk = AliTRDgeometry::GetStack(det), | |
99 | sec = AliTRDgeometry::GetSector(det), | |
100 | vid = AliGeomManager::LayerToVolUID(AliGeomManager::kTRD1 + ly, stk + AliTRDgeometry::Nstack() * sec); | |
101 | SetNameTitle(Form("digits%03d", det), Form("D-%03d [%02d_%d_%d]", det, sec, stk, ly)); | |
102 | Short_t sig[7]={0,0,0,10,0,0,0}; | |
103 | ||
104 | AliTRDtransform transform(det); | |
105 | AliTRDpadPlane *pp(fParent->fGeo->GetPadPlane(ly, stk)); | |
106 | ||
107 | Int_t row, col; | |
108 | AliTRDcluster c; | |
109 | indexes->ResetCounters(); | |
110 | while (indexes->NextRCIndex(row, col)){ | |
111 | dz = pp->GetRowSize(row); | |
112 | dy = pp->GetColSize(col); | |
113 | Short_t *const adc = data->GetDataAddress(row,col); | |
114 | for (Int_t time(0); time<30; time++){ | |
115 | if(data->IsPadCorrupted(row, col, time)) break; | |
116 | if(adc[time] <= 1) continue; | |
117 | new (&c) AliTRDcluster(det, col, row, time, sig, vid); | |
118 | transform.Transform(&c); | |
119 | ||
120 | scale = adc[time] < 512 ? adc[time]/512. : 1.; | |
121 | AddQuad(c.GetY()-0.5*dy, c.GetZ()-0.5*dz*scale, c.GetX(), dy*0.95, dz*scale); | |
122 | QuadValue(Float_t(adc[time])); | |
123 | QuadId(new TNamed(Form("ADC %d", adc[time]), Form("det[%3d(%02d_%d_%d)] col[%3d] row[%2d] tb[%2d]", det, sec, stk, ly, col, row, time))); | |
124 | } | |
125 | } | |
126 | ||
127 | // rotate to global coordinates | |
128 | RefitPlex(); | |
129 | TEveTrans& t = RefMainTrans(); | |
130 | t.SetRotByAngles((sec+.5)*AliTRDgeometry::GetAlpha(), 0.,0.); | |
131 | } | |
132 | ||
133 | ||
134 | // //______________________________________________________________________________ | |
135 | // void AliEveTRDDigits::Paint(Option_t *option) | |
136 | // { | |
137 | // // Paint the object. | |
138 | // | |
139 | // if(fParent->GetDigitsBox()) fBoxes.Paint(option); | |
140 | // else TEveQuadSet::Paint(option); | |
141 | // } | |
142 | ||
143 | // //______________________________________________________________________________ | |
144 | // void AliEveTRDDigits::Reset() | |
145 | // { | |
146 | // // Reset raw and visual data. | |
147 | // | |
148 | // TEveQuadSet::Reset(TEveQuadSet::kQT_RectangleYZ, kTRUE, 64); | |
149 | // // MT fBoxes.fBoxes.clear(); | |
150 | // fData.Reset(); | |
151 | // } | |
152 | ||
153 | /////////////////////////////////////////////////////////// | |
154 | ///////////// AliEveTRDHits ///////////////////// | |
155 | /////////////////////////////////////////////////////////// | |
156 | ||
157 | //______________________________________________________________________________ | |
158 | AliEveTRDHits::AliEveTRDHits() : TEvePointSet("hits", 20) | |
159 | { | |
160 | // Constructor. | |
161 | SetMarkerSize(.1); | |
162 | SetMarkerColor(kGreen); | |
163 | SetOwnIds(kTRUE); | |
164 | } | |
165 | ||
166 | //______________________________________________________________________________ | |
167 | AliEveTRDHits::~AliEveTRDHits() | |
168 | { | |
169 | //AliInfo(GetTitle()); | |
170 | } | |
171 | ||
172 | //______________________________________________________________________________ | |
173 | void AliEveTRDHits::PointSelected(Int_t n) | |
174 | { | |
175 | // Handle an individual point selection from GL. | |
176 | ||
177 | AliTRDhit *h = NULL; | |
178 | if(!(h = dynamic_cast<AliTRDhit*>(GetPointId(n)))) return; | |
179 | printf("Id[%3d] Det[%3d] Reg[%c] TR[%c] Q[%3d] MC[%d] t[%f]\n", | |
180 | n, h->GetDetector(), | |
181 | h->FromAmplification() ? 'A' : 'D', | |
182 | h->FromTRphoton() ? 'y' : 'n', | |
183 | h->GetCharge(), h->GetTrack(), h->GetTime()); | |
184 | } | |
185 | ||
186 | ||
187 | /////////////////////////////////////////////////////////// | |
188 | ///////////// AliEveTRDClusters ///////////////////// | |
189 | /////////////////////////////////////////////////////////// | |
190 | ||
191 | //______________________________________________________________________________ | |
192 | AliEveTRDClusters::AliEveTRDClusters():AliEveTRDHits() | |
193 | { | |
194 | // Constructor. | |
195 | SetName("clusters"); | |
196 | ||
197 | SetMarkerSize(.4); | |
198 | SetMarkerStyle(24); | |
199 | SetMarkerColor(kGray); | |
200 | SetOwnIds(kTRUE); | |
201 | } | |
202 | ||
203 | //______________________________________________________________________________ | |
204 | void AliEveTRDClusters::PointSelected(Int_t n) | |
205 | { | |
206 | // Handle an individual point selection from GL. | |
207 | ||
208 | AliTRDcluster *c = dynamic_cast<AliTRDcluster*>(GetPointId(n)); | |
209 | if(!c) return; | |
210 | c->Print(); | |
211 | Emit("PointSelected(Int_t)", n); | |
212 | // Bool_t AliCluster::GetGlobalCov(Float_t* cov) const | |
213 | // Bool_t AliCluster::GetGlobalXYZ(Float_t* xyz) const | |
214 | // Float_t AliCluster::GetSigmaY2() const | |
215 | // Float_t AliCluster::GetSigmaYZ() const | |
216 | // Float_t AliCluster::GetSigmaZ2() const | |
217 | } | |
218 | ||
219 | //______________________________________________________________________________ | |
220 | void AliEveTRDClusters::Print(Option_t *o) const | |
221 | { | |
222 | AliTRDcluster *c = NULL; | |
223 | ||
224 | for(Int_t n = GetN(); n--;){ | |
225 | if(!(c = dynamic_cast<AliTRDcluster*>(GetPointId(n)))) continue; | |
226 | c->Print(o); | |
227 | } | |
228 | } | |
229 | ||
230 | //______________________________________________________________________________ | |
231 | void AliEveTRDClusters::Load(Char_t *w) const | |
232 | { | |
233 | Int_t typ = -1; | |
234 | if(strcmp(w, "hit")==0) typ = 0; | |
235 | else if(strcmp(w, "dig")==0) typ = 1; | |
236 | else if(strcmp(w, "cls")==0) typ = 2; | |
237 | else if(strcmp(w, "all")==0) typ = 3; | |
238 | else{ | |
239 | AliInfo("The following arguments are accepted:"); | |
240 | AliInfo(" \"hit\" : loading of MC hits"); | |
241 | AliInfo(" \"dig\" : loading of digits"); | |
242 | AliInfo(" \"cls\" : loading of reconstructed clusters"); | |
243 | AliInfo(" \"all\" : loading of MC hits+digits+clusters"); | |
244 | return; | |
245 | } | |
246 | ||
247 | AliTRDcluster *c = NULL; | |
248 | Int_t n = 0; | |
249 | while((n = GetN() && !(c = dynamic_cast<AliTRDcluster*>(GetPointId(n))))) n++; | |
250 | if(!c) return; | |
251 | ||
252 | Int_t det = c->GetDetector(); | |
253 | AliEveTRDLoader *loader = NULL; | |
254 | switch(typ){ | |
255 | case 0: | |
256 | loader = new AliEveTRDLoader("Hits"); | |
257 | if(!loader->Open("TRD.Hits.root")){ | |
258 | delete loader; | |
259 | return; | |
260 | } | |
261 | loader->SetDataType(AliEveTRDLoader::kTRDHits); | |
262 | break; | |
263 | case 1: | |
264 | loader = new AliEveTRDLoader("Digits"); | |
265 | if(!loader->Open("TRD.Digits.root")){ | |
266 | delete loader; | |
267 | return; | |
268 | } | |
269 | loader->SetDataType(AliEveTRDLoader::kTRDDigits); | |
270 | break; | |
271 | case 2: | |
272 | loader = new AliEveTRDLoader("Clusters"); | |
273 | if(!loader->Open("TRD.RecPoints.root")){ | |
274 | delete loader; | |
275 | return; | |
276 | } | |
277 | loader->SetDataType(AliEveTRDLoader::kTRDClusters); | |
278 | break; | |
279 | case 3: | |
280 | loader = new AliEveTRDLoaderSim("MC"); | |
281 | if(!loader->Open("galice.root")){ | |
282 | delete loader; | |
283 | return; | |
284 | } | |
285 | loader->SetDataType(AliEveTRDLoader::kTRDHits | AliEveTRDLoader::kTRDDigits | AliEveTRDLoader::kTRDClusters); | |
286 | break; | |
287 | default: return; | |
288 | } | |
289 | ||
290 | loader->AddChambers(AliTRDgeometry::GetSector(det),AliTRDgeometry::GetStack(det), AliTRDgeometry::GetLayer(det)); | |
291 | // load first event | |
292 | loader->GoToEvent(AliEveEventManager::GetCurrent()->GetEventId()); | |
293 | ||
294 | // register loader with alieve | |
295 | gEve->AddElement(loader->GetChamber(det), *(BeginParents())); | |
296 | //loader->SpawnEditor(); | |
297 | gEve->Redraw3D(); | |
298 | } | |
299 | ||
300 | /////////////////////////////////////////////////////////// | |
301 | ///////////// AliEveTRDTracklet ///////////////////// | |
302 | /////////////////////////////////////////////////////////// | |
303 | ||
304 | //______________________________________________________________________________ | |
305 | AliEveTRDTracklet::AliEveTRDTracklet(AliTRDseedV1 *trklt):TEveLine() | |
306 | ,fClusters(NULL) | |
307 | { | |
308 | // Constructor. | |
309 | SetName("tracklet"); | |
310 | ||
311 | if(!gGeoManager){ | |
312 | AliEveEventManager::AssertGeometry(); | |
313 | AliInfo(Form("gGeo[%p] Closed[%c]", (void*)gGeoManager, gGeoManager->IsClosed()?'y':'n')); | |
314 | } | |
315 | SetUserData(trklt); | |
316 | Float_t dx; | |
317 | Float_t x0 = trklt->GetX0(); | |
318 | Float_t y0 = trklt->GetYref(0); | |
319 | Float_t z0 = trklt->GetZref(0); | |
320 | Float_t dydx = trklt->GetYref(1); | |
321 | Float_t dzdx = trklt->GetZref(1); | |
322 | Float_t tilt = trklt->GetTilt(); | |
323 | Float_t g[3]; | |
324 | AliTRDcluster *c = NULL; | |
325 | for(Int_t ic=0; ic<AliTRDseedV1::kNclusters; ic++){ | |
326 | if(!(c = trklt->GetClusters(ic))) continue; | |
327 | if(!fClusters) AddElement(fClusters = new AliEveTRDClusters()); | |
328 | dx = x0 - c->GetX(); | |
329 | //Float_t yt = y0 - dx*dydx; | |
330 | Float_t zt = z0 - dx*dzdx; | |
331 | // backup yc - for testing purposes | |
332 | Float_t yc = c->GetY(); | |
333 | c->SetY(yc-tilt*(c->GetZ()-zt)); | |
334 | c->GetGlobalXYZ(g); | |
335 | Int_t id = fClusters->SetNextPoint(g[0], g[1], g[2]); | |
336 | c->SetY(yc); | |
337 | fClusters->SetPointId(id, new AliTRDcluster(*c)); | |
338 | } | |
339 | if(fClusters){ | |
340 | fClusters->SetName("TC clusters"); | |
341 | fClusters->SetTitle(Form("N[%d]", trklt->GetN())); | |
342 | fClusters->SetMarkerColor(kMagenta); | |
343 | } | |
344 | ||
345 | SetTitle(Form("Det[%d] RC[%c] Layer[%d] P[%7.3f]", trklt->GetDetector(), trklt->IsRowCross()?'y':'n', trklt->GetPlane(), trklt->GetMomentum())); | |
346 | SetLineColor(kRed); | |
347 | //SetOwnIds(kTRUE); | |
348 | ||
349 | // init tracklet line | |
350 | Int_t sec = AliTRDgeometry::GetSector(trklt->GetDetector()); | |
351 | Double_t alpha = AliTRDgeometry::GetAlpha() * (sec<9 ? sec + .5 : sec - 17.5); | |
352 | ||
353 | //trklt->Fit(kTRUE); | |
354 | y0 = trklt->GetYfit(0); | |
355 | dydx = trklt->GetYfit(1); | |
356 | Double_t xg = x0 * TMath::Cos(alpha) - y0 * TMath::Sin(alpha); | |
357 | Double_t yg = x0 * TMath::Sin(alpha) + y0 * TMath::Cos(alpha); | |
358 | SetPoint(0, xg, yg, z0); | |
359 | //SetPoint(0, x0, y0, z0); | |
360 | ||
361 | ||
362 | dx = .5*AliTRDgeometry::CamHght()+AliTRDgeometry::CdrHght(); | |
363 | x0 -= dx; | |
364 | y0 -= dydx*dx, | |
365 | z0 -= dzdx*dx; | |
366 | xg = x0 * TMath::Cos(alpha) - y0 * TMath::Sin(alpha); | |
367 | yg = x0 * TMath::Sin(alpha) + y0 * TMath::Cos(alpha); | |
368 | SetPoint(1, xg, yg, z0); | |
369 | //SetPoint(1, x0, y0, z0); | |
370 | } | |
371 | ||
372 | ||
373 | AliEveTRDTracklet::~AliEveTRDTracklet() | |
374 | { | |
375 | ||
376 | } | |
377 | ||
378 | //______________________________________________________________________________ | |
379 | void AliEveTRDTracklet::Print(Option_t *o) const | |
380 | { | |
381 | AliTRDseedV1 *tracklet = (AliTRDseedV1*)GetUserData(); | |
382 | if(!tracklet) return; | |
383 | tracklet->Print(o); | |
384 | } | |
385 | ||
386 | /////////////////////////////////////////////////////////// | |
387 | ///////////// AliEveTRDTrack ///////////////////// | |
388 | /////////////////////////////////////////////////////////// | |
389 | ||
390 | //______________________________________________________________________________ | |
391 | AliEveTRDTrack::AliEveTRDTrack(AliTRDtrackV1 *trk) | |
392 | :TEveLine() | |
393 | ,fTrackState(0) | |
394 | ,fESDStatus(0) | |
395 | ,fAlpha(0.) | |
396 | ,fPoints(NULL) | |
397 | ,fRim(NULL) | |
398 | { | |
399 | // Constructor. | |
400 | SetUserData(trk); | |
401 | SetName(""); | |
402 | ||
403 | // AliTRDtrackerV1::SetNTimeBins(24); | |
404 | ||
405 | fRim = new AliRieman(trk->GetNumberOfClusters()); | |
406 | AliTRDseedV1 *tracklet = NULL; | |
407 | for(Int_t il=0; il<AliTRDgeometry::kNlayer; il++){ | |
408 | if(!(tracklet = trk->GetTracklet(il))) continue; | |
409 | if(!tracklet->IsOK()) continue; | |
410 | AddElement(new AliEveTRDTracklet(tracklet)); | |
411 | ||
412 | // tracklet->ResetClusterIter(kFALSE); | |
413 | // while((c = tracklet->PrevCluster())){ | |
414 | // AliTRDcluster *c(NULL); | |
415 | /* for(Int_t ic=AliTRDseedV1::kNtb; ic--;){ | |
416 | if(!(c=tracklet->GetClusters(ic))) continue; | |
417 | Float_t xc = c->GetX(); | |
418 | Float_t yc = c->GetY(); | |
419 | Float_t zc = c->GetZ(); | |
420 | Float_t zt = tracklet->GetZref(0) - (tracklet->GetX0()-xc)*tracklet->GetZref(1); | |
421 | yc -= tracklet->GetTilt()*(zc-zt); | |
422 | fRim->AddPoint(xc, yc, zc, .05, 2.3); | |
423 | }*/ | |
424 | } | |
425 | if(trk->GetNumberOfTracklets()>1) fRim->Update(); | |
426 | SetStatus(fTrackState); | |
427 | } | |
428 | ||
429 | //______________________________________________________________________________ | |
430 | AliEveTRDTrack::~AliEveTRDTrack() | |
431 | { | |
432 | if(fPoints) delete [] fPoints; fPoints = NULL; | |
433 | //delete dynamic_cast<AliTRDtrackV1*>(GetUserData()); | |
434 | } | |
435 | ||
436 | //______________________________________________________________________________ | |
437 | void AliEveTRDTrack::Print(Option_t *o) const | |
438 | { | |
439 | AliTRDtrackV1 *track = (AliTRDtrackV1*)GetUserData(); | |
440 | if(!track) return; | |
441 | track->Print(o); | |
442 | } | |
443 | ||
444 | //______________________________________________________________________________ | |
445 | void AliEveTRDTrack::SetStatus(UChar_t s) | |
446 | { | |
447 | // nothing to be done | |
448 | if(fPoints && fTrackState == s) return; | |
449 | //return; | |
450 | const Int_t nc = AliTRDtrackV1::kMAXCLUSTERSPERTRACK; | |
451 | AliTRDtrackV1 *trk(NULL); | |
452 | if(!(trk=static_cast<AliTRDtrackV1*>(GetUserData()))) { | |
453 | AliError("Failed casting data to TRD track."); | |
454 | return; | |
455 | } | |
456 | ||
457 | Bool_t BUILD = kFALSE; | |
458 | if(!fPoints){ | |
459 | fPoints = new AliTrackPoint[nc]; | |
460 | ||
461 | // define the radial span of the track in the TRD | |
462 | Double_t xmin = -1., xmax = -1.; | |
463 | Int_t det = 0; | |
464 | AliTRDseedV1 *trklt = NULL; | |
465 | for(Int_t ily=0; ily<AliTRDgeometry::kNlayer; ily++){ | |
466 | if(!(trklt = trk->GetTracklet(ily))) continue; | |
467 | if(xmin<0.) xmin = trklt->GetX0() - AliTRDgeometry::CamHght() - AliTRDgeometry::CdrHght(); | |
468 | if(trklt->GetX0()>xmax) xmax = trklt->GetX0(); | |
469 | det = trklt->GetDetector(); | |
470 | } | |
471 | Int_t sec = det/AliTRDgeometry::kNdets; | |
472 | fAlpha = AliTRDgeometry::GetAlpha() * (sec<9 ? sec + .5 : sec - 17.5); //trk->GetAlpha() | |
473 | ||
474 | Double_t dx =(xmax - xmin)/nc; | |
475 | for(Int_t ip=0; ip<nc; ip++){ | |
476 | fPoints[ip].SetXYZ(xmin, 0., 0.); | |
477 | xmin+=dx; | |
478 | } | |
479 | BUILD = kTRUE; | |
480 | } | |
481 | ||
482 | // select track model | |
483 | if(BUILD || ((s&12) != (fTrackState&12))){ | |
484 | if(TESTBIT(s, kTrackCosmics)){ | |
485 | //printf("Straight track\n"); | |
486 | AliTRDtrackerV1::FitLine(trk, NULL, kFALSE, nc, fPoints); | |
487 | } else { | |
488 | if(TESTBIT(s, kTrackModel)){ | |
489 | //printf("Kalman track\n"); | |
490 | if(trk->GetNumberOfTracklets() >=4) AliTRDtrackerV1::FitKalman(trk, NULL, kFALSE, nc, fPoints); | |
491 | } else { | |
492 | //printf("Rieman track rim[%p] nc[%d]\n", (void*)fRim, nc); | |
493 | //if(trk->GetNumberOfTracklets() >=4) AliTRDtrackerV1::FitRiemanTilt(trk, NULL, kTRUE, nc, fPoints); | |
494 | Float_t x = 0.; | |
495 | for(Int_t ip = nc; ip--;){ | |
496 | x = fPoints[ip].GetX(); | |
497 | //printf("%2d x[%f] y[%f] z[%f]\n", ip, x, fRim->GetYat(x), fRim->GetZat(x)); | |
498 | fPoints[ip].SetXYZ(x, fRim->GetYat(x), fRim->GetZat(x)); | |
499 | } | |
500 | } | |
501 | } | |
502 | ||
503 | Float_t global[3]; | |
504 | for(Int_t ip=0; ip<nc; ip++){ | |
505 | fPoints[ip].Rotate(-fAlpha).GetXYZ(global); | |
506 | SetPoint(ip, global[0], global[1], global[2]); | |
507 | //printf("*** %2d x[%f] y[%f] z[%f]\n", ip, global[0], global[1], global[2]); | |
508 | ||
509 | } | |
510 | SetSmooth(kTRUE); | |
511 | } | |
512 | ||
513 | // set color | |
514 | if(BUILD || ((s&3) != (fTrackState&3))){ | |
515 | if(TESTBIT(s, kSource)){ | |
516 | //printf("Source color\n"); | |
517 | if(fESDStatus&AliESDtrack::kTRDin){ | |
518 | SetMarkerColor(kGreen); | |
519 | SetLineColor(kGreen); | |
520 | } else { | |
521 | SetMarkerColor(kMagenta); | |
522 | SetLineColor(kMagenta); | |
523 | } | |
524 | } else { | |
525 | if(TESTBIT(s, kPID) == AliTRDpidUtil::kLQ){ | |
526 | //printf("PID color kLQPID\n"); | |
527 | //trk->GetReconstructor()->SetOption("!nn"); | |
528 | } else { | |
529 | //printf("PID color kNNPID\n"); | |
530 | //trk->GetReconstructor()->SetOption("nn"); | |
531 | } | |
532 | //trk->CookPID(); | |
533 | ||
534 | Int_t species = 0; Float_t pid = 0.; | |
535 | for(Int_t is=0; is<AliPID::kSPECIES; is++) | |
536 | if(trk->GetPID(is) > pid){ | |
537 | pid = trk->GetPID(is); | |
538 | species = is; | |
539 | } | |
540 | switch(species){ | |
541 | case AliPID::kElectron: | |
542 | SetMarkerColor(kRed); | |
543 | SetLineColor(kRed); | |
544 | break; | |
545 | default: | |
546 | SetMarkerColor(kBlue); | |
547 | SetLineColor(kBlue); | |
548 | break; | |
549 | } | |
550 | } | |
551 | SetLineWidth(2); | |
552 | } | |
553 | ||
554 | const Char_t *model = "line"; | |
555 | if(!TESTBIT(s, kTrackCosmics)){ | |
556 | if(TESTBIT(s, kTrackModel)) model = "kalman"; | |
557 | else model = "rieman"; | |
558 | } | |
559 | Int_t species = 0; Float_t pid = 0.; | |
560 | for(Int_t is=0; is<AliPID::kSPECIES; is++) | |
561 | if(trk->GetPID(is) > pid){ | |
562 | pid = trk->GetPID(is); | |
563 | species = is; | |
564 | } | |
565 | ||
566 | SetTitle(Form( | |
567 | "Tracklets[%d] Clusters[%d]\n" | |
568 | "Reconstruction Source[%s]\n" | |
569 | "PID[%4.1f %4.1f %4.1f %4.1f %4.1f]\n" | |
570 | "MC[%d]", trk->GetNumberOfTracklets(), trk->GetNumberOfClusters(), fESDStatus&AliESDtrack::kTRDin ? "barrel" : "sa", | |
571 | 1.E2*trk->GetPID(0), 1.E2*trk->GetPID(1), | |
572 | 1.E2*trk->GetPID(2), 1.E2*trk->GetPID(3), 1.E2*trk->GetPID(4), trk->GetLabel())); | |
573 | ||
574 | if(GetName()){ | |
575 | char id[6]; snprintf(id, 6, "%s", GetName()); | |
576 | SetName(Form("%s %s", id, AliPID::ParticleName(species))); | |
577 | } | |
578 | ||
579 | // save track status | |
580 | fTrackState = s; | |
581 | } | |
582 | ||
583 | //______________________________________________________________________________ | |
584 | void AliEveTRDTrack::Load(Char_t *what) const | |
585 | { | |
586 | // Spread downwards to tracklets the command "what" | |
587 | ||
588 | const AliEveTRDTracklet* trklt(NULL); | |
589 | TEveElement::List_ci itrklt=BeginChildren(); | |
590 | while(itrklt!=EndChildren()){ | |
591 | if((trklt = dynamic_cast<const AliEveTRDTracklet*>(*itrklt))) trklt->Load(what); | |
592 | itrklt++; | |
593 | } | |
594 | } | |
595 | ||
596 | //______________________________________________________________________________ | |
597 | AliEveTRDTrackletOnline::AliEveTRDTrackletOnline(AliTRDtrackletMCM *tracklet) : | |
598 | TEveLine(), | |
599 | fDetector(-1), | |
600 | fROB(-1), | |
601 | fMCM(-1) | |
602 | { | |
603 | AliTRDtrackletMCM *trkl = new AliTRDtrackletMCM(*tracklet); | |
604 | SetUserData(trkl); | |
605 | ||
606 | fDetector = trkl->GetDetector(); | |
607 | fROB = trkl->GetROB(); | |
608 | fMCM = trkl->GetMCM(); | |
609 | ||
610 | SetName("sim. tracklet"); | |
611 | SetTitle(Form("Det: %i, ROB: %i, MCM: %i, Label: %i\n0x%08x", | |
612 | trkl->GetDetector(), trkl->GetROB(), trkl->GetMCM(), trkl->GetLabel(), | |
613 | trkl->GetTrackletWord())); | |
614 | SetLineColor(kGreen); | |
615 | SetLineWidth(3); | |
616 | ||
617 | AliTRDgeometry geo; | |
618 | TGeoHMatrix *matrix = geo.GetClusterMatrix(trkl->GetDetector()); | |
619 | ||
620 | fDetector = trkl->GetDetector(); | |
621 | fROB = trkl->GetROB(); | |
622 | fMCM = trkl->GetMCM(); | |
623 | ||
624 | Float_t length = 3.; | |
625 | Double_t x[3]; | |
626 | Double_t p[3]; | |
627 | Double_t p2[3]; | |
628 | x[0] = AliTRDgeometry::AnodePos(); | |
629 | x[1] = trkl->GetY(); | |
630 | x[2] = trkl->GetLocalZ(); | |
631 | ||
632 | matrix->LocalToMaster(x, p); | |
633 | geo.RotateBack(trkl->GetDetector(), p, p2); | |
634 | SetPoint(0, p2[0], p2[1], p2[2]); | |
635 | ||
636 | x[0] -= length; | |
637 | x[1] -= length * trkl->GetdYdX(); | |
638 | matrix->LocalToMaster(x, p); | |
639 | p[2] *= p[0] / (p[0] + length); | |
640 | geo.RotateBack(trkl->GetDetector(), p, p2); | |
641 | SetPoint(1, p2[0], p2[1], p2[2]); | |
642 | } | |
643 | ||
644 | AliEveTRDTrackletOnline::AliEveTRDTrackletOnline(AliTRDtrackletWord *tracklet) : | |
645 | TEveLine(), | |
646 | fDetector(-1), | |
647 | fROB(-1), | |
648 | fMCM(-1) | |
649 | { | |
650 | AliTRDtrackletWord *trkl = new AliTRDtrackletWord(*tracklet); | |
651 | SetUserData(trkl); | |
652 | ||
653 | fDetector = trkl->GetDetector(); | |
654 | fROB = trkl->GetROB(); | |
655 | fMCM = trkl->GetMCM(); | |
656 | ||
657 | SetName("raw tracklet"); | |
658 | SetTitle(Form("Det: %i, ROB: %i, MCM: %i, Label: %i\n0x%08x", | |
659 | trkl->GetDetector(), fROB, fMCM, -1, | |
660 | trkl->GetTrackletWord())); | |
661 | SetLineColor(kRed); | |
662 | SetLineWidth(3); | |
663 | ||
664 | AliTRDgeometry geo; | |
665 | TGeoHMatrix *matrix = geo.GetClusterMatrix(trkl->GetDetector()); | |
666 | ||
667 | Float_t length = 3.; | |
668 | Double_t x[3]; | |
669 | Double_t p[3]; | |
670 | Double_t p2[3]; | |
671 | x[0] = AliTRDgeometry::AnodePos(); | |
672 | x[1] = trkl->GetY(); | |
673 | x[2] = trkl->GetLocalZ(); | |
674 | ||
675 | matrix->LocalToMaster(x, p); | |
676 | geo.RotateBack(trkl->GetDetector(), p, p2); | |
677 | SetPoint(0, p2[0], p2[1], p2[2]); | |
678 | ||
679 | x[0] -= length; | |
680 | x[1] -= length * trkl->GetdYdX(); | |
681 | matrix->LocalToMaster(x, p); | |
682 | p[2] *= p[0] / (p[0] + length); | |
683 | geo.RotateBack(trkl->GetDetector(), p, p2); | |
684 | SetPoint(1, p2[0], p2[1], p2[2]); | |
685 | } | |
686 | ||
687 | AliEveTRDTrackletOnline::~AliEveTRDTrackletOnline() | |
688 | { | |
689 | delete ((AliTRDtrackletBase*) GetUserData()); | |
690 | SetUserData(0x0); | |
691 | } | |
692 | ||
693 | void AliEveTRDTrackletOnline::ShowMCM(Option_t *opt) const | |
694 | { | |
695 | if (fDetector < 0 || fROB < 0 || fMCM < 0) | |
696 | return; | |
697 | ||
698 | AliEveTRDmcm *evemcm = new AliEveTRDmcm(); | |
699 | evemcm->Init(fDetector, fROB, fMCM); | |
700 | evemcm->LoadDigits(); | |
701 | evemcm->Draw(opt); | |
702 | ||
703 | TEveElementList *mcmlist = NULL; | |
704 | if (gEve->GetCurrentEvent()) | |
705 | mcmlist = (TEveElementList*) gEve->GetCurrentEvent()->FindChild("TRD MCMs"); | |
706 | if (!mcmlist) { | |
707 | mcmlist = new TEveElementList("TRD MCMs"); | |
708 | gEve->AddElement(mcmlist); | |
709 | } | |
710 | gEve->AddElement(evemcm, mcmlist); | |
711 | } | |
712 | ||
713 | ||
714 | AliEveTRDmcm::AliEveTRDmcm() : | |
715 | TEveElement(), | |
716 | TNamed(), | |
717 | fMCM(new AliTRDmcmSim()) | |
718 | { | |
719 | SetName("MCM"); | |
720 | SetTitle("Unknown MCM"); | |
721 | } | |
722 | ||
723 | AliEveTRDmcm::~AliEveTRDmcm() | |
724 | { | |
725 | delete fMCM; | |
726 | } | |
727 | ||
728 | Bool_t AliEveTRDmcm::Init(Int_t det, Int_t rob, Int_t mcm) | |
729 | { | |
730 | SetName(Form("MCM: Det. %i, ROB %i, MCM %i", det, rob, mcm)); | |
731 | SetTitle(Form("MCM: Det. %i, ROB %i, MCM %i", det, rob, mcm)); | |
732 | fMCM->Init(det, rob, mcm); | |
733 | fMCM->Reset(); | |
734 | return kTRUE; | |
735 | } | |
736 | ||
737 | Bool_t AliEveTRDmcm::LoadDigits() | |
738 | { | |
739 | AliRunLoader *rl = AliEveEventManager::AssertRunLoader(); | |
740 | return fMCM->LoadMCM(rl, fMCM->GetDetector(), | |
741 | fMCM->GetRobPos(), fMCM->GetMcmPos()); | |
742 | } | |
743 | ||
744 | Bool_t AliEveTRDmcm::Filter() | |
745 | { | |
746 | fMCM->Filter(); | |
747 | return kTRUE; | |
748 | } | |
749 | ||
750 | Bool_t AliEveTRDmcm::Tracklet() | |
751 | { | |
752 | fMCM->Tracklet(); | |
753 | return kTRUE; | |
754 | } | |
755 | ||
756 | void AliEveTRDmcm::Draw(Option_t* option) | |
757 | { | |
758 | const char *mcmname = Form("mcm_%i_%i_%i", fMCM->GetDetector(), | |
759 | fMCM->GetRobPos(), fMCM->GetMcmPos()); | |
760 | ||
761 | TCanvas *c = dynamic_cast<TCanvas*> (gROOT->FindObject(mcmname)); | |
762 | if (!c) | |
763 | c = gEve->AddCanvasTab("TRD MCM"); | |
764 | c->SetTitle(Form("MCM %i on ROB %i of det. %i", | |
765 | fMCM->GetMcmPos(), fMCM->GetRobPos(), fMCM->GetDetector())); | |
766 | c->SetName(mcmname); | |
767 | c->cd(); | |
768 | fMCM->Draw(option); | |
769 | } | |
770 | ||
771 | Bool_t AliEveTRDmcm::AssignPointer(const char* ptrname) | |
772 | { | |
773 | gROOT->ProcessLine(Form("AliTRDmcmSim* %s = (AliTRDmcmSim *)%p", ptrname, (void*)fMCM)); | |
774 | return kTRUE; | |
775 | } |