07af35322d01b73ff3fe60f38318ae014e0f2726
[u/mrichter/AliRoot.git] / EVE / Reve / Track.cxx
1 // $Header$
2
3 #include "Track.h"
4 #include "MCHelixLine.hi"
5
6 #include <TPolyLine3D.h>
7 #include <TPolyMarker3D.h>
8 #include <TColor.h>
9
10 // Updates
11 #include <Reve/RGTopFrame.h>
12 #include <TCanvas.h>
13
14 #include <vector>
15
16 using namespace Reve;
17
18 //______________________________________________________________________
19 // Track
20 //
21
22 ClassImp(Reve::Track)
23
24 Track::Track() :
25   Line(),
26
27   fV(),
28   fP(),
29   fBeta(0),
30   fCharge(0),
31   fLabel(0),
32   fPathMarks(),
33
34   fRnrStyle(0)
35 {}
36
37 Track::Track(TParticle* t, Int_t label, TrackRnrStyle* rs):
38   Line(),
39
40   fV(t->Vx(), t->Vy(), t->Vz()),
41   fP(t->Px(), t->Py(), t->Pz()),
42   fBeta(t->P()/t->Energy()),
43   fCharge(0),
44   fLabel(label),
45   fPathMarks(),
46
47   fRnrStyle(rs)
48 {
49   fLineColor = fRnrStyle->GetColor();
50   fMainColorPtr = &fLineColor;
51
52   TParticlePDG* pdgp = t->GetPDG();
53   if (pdgp)
54     fCharge = (Int_t) TMath::Nint(pdgp->Charge()/3);
55
56   SetName(t->GetName());
57 }
58
59 Track::Track(Reve::MCTrack* t, TrackRnrStyle* rs):
60   Line(),
61
62   fV(t->Vx(), t->Vy(), t->Vz()),
63   fP(t->Px(), t->Py(), t->Pz()),
64   fBeta(t->P()/t->Energy()),
65   fCharge(0),
66   fLabel(t->label),
67   fIndex(t->index),
68   fPathMarks(),
69
70   fRnrStyle(rs)
71 {
72   fLineColor = fRnrStyle->GetColor();
73   fMainColorPtr = &fLineColor;
74
75   TParticlePDG* pdgp = t->GetPDG();
76   if(pdgp == 0) {
77     t->ResetPdgCode(); pdgp = t->GetPDG();
78   }
79   fCharge = (Int_t) TMath::Nint(pdgp->Charge()/3);
80
81   SetName(t->GetName());
82 }
83
84 Track::Track(Reve::RecTrack* t, TrackRnrStyle* rs) :
85   Line(),
86
87   fV(t->V),
88   fP(t->P),
89   fBeta(t->beta),
90   fCharge(t->sign),
91   fLabel(t->label),
92   fIndex(t->index),
93   fPathMarks(),
94
95   fRnrStyle(rs)
96 {
97   fLineColor = fRnrStyle->GetColor();
98   fMainColorPtr = &fLineColor;
99
100   SetName(t->GetName());
101 }
102
103 Track::~Track()
104 {
105   for (vpPathMark_i i=fPathMarks.begin(); i!=fPathMarks.end(); ++i)
106     delete *i;
107 }
108
109 /*
110 void Track::Reset(Int_t n_points)
111 {
112   delete [] TPolyLine3D::fP; TPolyLine3D::fP = 0;
113   fN = n_points;
114   if(fN) TPolyLine3D::fP = new Float_t [3*fN];
115   memset(TPolyLine3D::fP, 0, 3*fN*sizeof(Float_t));
116   fLastPoint = -1;
117 }
118 */
119
120 /**************************************************************************/
121
122 void Track::MakeTrack()
123 {
124   
125   TrackRnrStyle& RS((fRnrStyle != 0) ? *fRnrStyle : TrackRnrStyle::fgDefStyle);
126
127   Float_t px = fP.x, py = fP.y, pz = fP.z;  
128
129   MCVertex  mc_v0;
130   mc_v0.x = fV.x;
131   mc_v0.y = fV.y; 
132   mc_v0.z = fV.z; 
133   mc_v0.t = 0;
134
135   std::vector<MCVertex> track_points;
136   Bool_t decay = kFALSE;
137
138   if ((TMath::Abs(fV.z) > RS.fMaxZ) || (fV.x*fV.x + fV.y*fV.y > RS.fMaxR*RS.fMaxR)) 
139     goto make_polyline;
140   
141   if (fCharge != 0 && TMath::Abs(RS.fMagField) > 1e-5) {
142
143     // Charged particle in magnetic field
144
145     Float_t a = RS.fgkB2C * RS.fMagField * fCharge;
146    
147     MCHelix helix(fRnrStyle, &mc_v0, TMath::C()*fBeta, &track_points, a); //m->cm
148     helix.Init(TMath::Sqrt(px*px+py*py), pz);
149    
150     if(!fPathMarks.empty()){
151       for(std::vector<Reve::PathMark*>::iterator i=fPathMarks.begin(); i!=fPathMarks.end(); ++i) {
152         Reve::PathMark* pm = *i;
153         
154         if(RS.fFitDaughters &&  pm->type == Reve::PathMark::Daughter){
155           if(TMath::Abs(pm->V.z) > RS.fMaxZ 
156              || TMath::Sqrt(pm->V.x*pm->V.x + pm->V.y*pm->V.y) > RS.fMaxR )
157             goto helix_bounds;
158
159           //printf("%s fit daughter  \n", fName.Data()); 
160           helix.LoopToVertex(fP.x, fP.y, fP.z, pm->V.x, pm->V.y, pm->V.z);
161           fP.x -=  pm->P.x;
162           fP.y -=  pm->P.y;
163           fP.z -=  pm->P.z;
164         }
165         if(RS.fFitDecay &&  pm->type == Reve::PathMark::Decay){
166           
167           if(TMath::Abs(pm->V.z) > RS.fMaxZ 
168              || TMath::Sqrt(pm->V.x*pm->V.x + pm->V.y*pm->V.y) > RS.fMaxR )
169             goto helix_bounds;
170           helix.LoopToVertex(fP.x, fP.y, fP.z, pm->V.x, pm->V.y, pm->V.z);
171           decay = true;
172           break;
173         }
174       }
175     }
176   helix_bounds:
177     //go to bounds
178     if(!decay || RS.fFitDecay == kFALSE){
179       helix.LoopToBounds(px,py,pz);
180       // printf("%s loop to bounds  \n",fName.Data() );
181     }
182
183   } else {
184
185     // Neutral particle or no field
186
187     MCLine line(fRnrStyle, &mc_v0, TMath::C()*fBeta, &track_points);
188    
189     if(!fPathMarks.empty()){
190       for(std::vector<Reve::PathMark*>::iterator i=fPathMarks.begin(); i!=fPathMarks.end(); ++i) {
191         Reve::PathMark* pm = *i;
192
193         if(RS.fFitDaughters &&  pm->type == Reve::PathMark::Daughter){
194           if(TMath::Abs(pm->V.z) > RS.fMaxZ 
195              || TMath::Sqrt(pm->V.x*pm->V.x + pm->V.y*pm->V.y) > RS.fMaxR )
196             goto line_bounds;
197           line.GotoVertex(pm->V.x, pm->V.y, pm->V.z);
198           fP.x -=  pm->P.x;
199           fP.y -=  pm->P.y;
200           fP.z -=  pm->P.z;
201         }
202
203         if(RS.fFitDecay &&  pm->type == Reve::PathMark::Decay){
204           if(TMath::Abs(pm->V.z) > RS.fMaxZ 
205              || TMath::Sqrt(pm->V.x*pm->V.x + pm->V.y*pm->V.y) > RS.fMaxR )
206             goto line_bounds;
207           line.GotoVertex(pm->V.x, pm->V.y, pm->V.z);
208           decay = true;
209           break;
210         }
211       }
212     }
213
214   line_bounds:
215     if(!decay || RS.fFitDecay == kFALSE)
216       line.GotoBounds(px,py,pz);
217
218   }
219 make_polyline:
220   Reset(track_points.size());
221   for(std::vector<MCVertex>::iterator i=track_points.begin(); i!=track_points.end(); ++i)
222     SetNextPoint(i->x, i->y, i->z);
223 }
224
225 /**************************************************************************/
226
227 void Track::ImportHits()
228 {
229   Reve::LoadMacro("hits_from_label.C");
230   gROOT->ProcessLine(Form("hits_from_label(%d);", fLabel));
231 }
232
233 void Track::ImportClusters()
234 {
235   Reve::LoadMacro("clusters_from_label.C");
236   gROOT->ProcessLine(Form("clusters_from_label(%d);", fLabel));
237 }
238
239 void Track::ImportClustersFromIndex()
240 {
241   Reve::LoadMacro("clusters_from_index.C");
242   gROOT->ProcessLine(Form("clusters_from_index(%d);", fIndex));
243 }
244
245 /**************************************************************************/
246
247 void Track::CtrlClicked(Reve::Track* track)
248 {
249   Emit("CtrlClicked(Reve::Track*)", (Long_t)track);
250 }
251
252
253 /**************************************************************************/
254 /**************************************************************************/
255
256 //______________________________________________________________________
257 // TrackRnrStyle
258 //
259
260 ClassImp(Reve::TrackRnrStyle)
261
262 Float_t       TrackRnrStyle::fgDefMagField = 5;
263 const Float_t TrackRnrStyle::fgkB2C        = 0.299792458e-3;
264 TrackRnrStyle TrackRnrStyle::fgDefStyle;
265
266 TrackRnrStyle::TrackRnrStyle() :
267   TObject(),
268
269   fColor(1),
270   fWidth(1),
271   fMagField(fgDefMagField),
272
273   fMaxR  (350),
274   fMaxZ  (450),
275
276   fMaxOrbs (0.5),
277   fMinAng  (45),
278   fDelta   (0.1),
279
280   fFitDaughters(kTRUE),
281   fFitDecay    (kTRUE)
282 {}
283
284 /**************************************************************************/
285 /**************************************************************************/
286
287 //______________________________________________________________________
288 // TrackList
289 //
290
291 ClassImp(Reve::TrackList)
292
293 void TrackList::Init()
294 {
295   fMarkerStyle = 6;
296   fMarkerColor = 5;
297   // fMarker->SetMarkerSize(0.05);
298
299   if (fRnrStyle== 0) fRnrStyle = new TrackRnrStyle;
300   SetMainColorPtr(&fRnrStyle->fColor);
301 }
302
303 TrackList::TrackList(Int_t n_tracks, TrackRnrStyle* rs) :
304   RenderElementListBase(),
305   TPolyMarker3D(n_tracks),
306
307   fTitle(),
308
309   fRnrStyle   (rs),
310   fRnrMarkers (kTRUE),
311   fRnrTracks  (kTRUE)
312 {
313   Init();
314 }
315
316 TrackList::TrackList(const Text_t* name, Int_t n_tracks, TrackRnrStyle* rs) :
317   RenderElementListBase(),
318   TPolyMarker3D(n_tracks),
319   
320   fTitle(),
321
322   fRnrStyle   (rs),
323   fRnrMarkers (kTRUE),
324   fRnrTracks  (kTRUE)
325 {
326   Init();
327   SetName(name);
328 }
329
330 void TrackList::Reset(Int_t n_tracks)
331 {
332   delete [] fP; fP = 0;
333   fN = n_tracks;
334   if(fN) fP = new Float_t [3*fN];
335   memset(fP, 0, 3*fN*sizeof(Float_t));
336   fLastPoint = -1;
337 }
338
339 /**************************************************************************/
340
341 void TrackList::Paint(Option_t* option)
342 {
343   if(fRnrElement) {
344     if(fRnrMarkers) {
345       TPolyMarker3D::Paint(option);
346     }
347     if(fRnrTracks) {
348       for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
349         if((*i)->GetRnrElement())
350           (*i)->GetObject()->Paint(option);
351       }
352     }
353   }
354 }
355
356 /**************************************************************************/
357
358 void TrackList::AddElement(RenderElement* el)
359 {
360   static const Exc_t eH("TrackList::AddElement ");
361   if (dynamic_cast<Track*>(el)  == 0)
362     throw(eH + "new element not a Track.");
363   RenderElementListBase::AddElement(el);
364 }
365
366 /**************************************************************************/
367
368 void TrackList::SetRnrMarkers(Bool_t rnr)
369 {
370   fRnrMarkers = rnr;
371   gReve->Redraw3D();
372 }
373
374 void TrackList::SetRnrTracks(Bool_t rnr)
375 {
376
377   fRnrTracks = rnr;
378   gReve->Redraw3D();
379 }
380
381 /**************************************************************************/
382
383 void TrackList::MakeTracks()
384 {
385   for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
386     ((Track*)(*i))->MakeTrack();
387   }
388   gReve->Redraw3D();
389 }
390
391
392 void TrackList::MakeMarkers()
393 {
394   Reset(fChildren.size());
395   for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
396     Track& t = *((Track*)(*i));
397     if(t.GetN() > 0)
398       SetNextPoint(t.fV.x, t.fV.y, t.fV.z);
399   }
400   gReve->Redraw3D();
401 }
402
403 /**************************************************************************/
404 /*************************************************************************/
405
406 void TrackList::SetWidth(Width_t w)
407 {
408   Width_t oldw = fRnrStyle->fWidth;
409   fRnrStyle->fWidth = w;
410   for (List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
411     Track& t = *((Track*)(*i));
412     if (t.GetLineWidth() == oldw)
413       t.SetLineWidth(w);
414   }
415 }
416
417 void TrackList::SetMaxR(Float_t x)
418 {
419   fRnrStyle->fMaxR = x;
420   MakeTracks();
421   MakeMarkers();
422 }
423
424 void TrackList::SetMaxZ(Float_t x)
425 {
426   fRnrStyle->fMaxZ = x;
427   MakeTracks();
428   MakeMarkers();
429 }
430
431 void TrackList::SetMaxOrbs(Float_t x)
432 {
433   fRnrStyle->fMaxOrbs = x;
434   MakeTracks();
435 }
436
437 void TrackList::SetMinAng(Float_t x)
438 {
439   fRnrStyle->fMinAng = x;
440   MakeTracks();
441 }
442
443 void TrackList::SetDelta(Float_t x)
444 {
445   fRnrStyle->fDelta = x;
446   MakeTracks();
447 }
448
449 void TrackList::SetFitDaughters(Bool_t x)
450 {
451   fRnrStyle->fFitDaughters = x;
452   MakeTracks();
453 }
454
455 void TrackList::SetFitDecay(Bool_t x)
456 {
457   fRnrStyle->fFitDecay = x;
458   MakeTracks();
459 }
460
461 /**************************************************************************/
462 /**************************************************************************/
463
464 void TrackList::SelectByPt(Float_t min_pt, Float_t max_pt)
465 {
466   Float_t minptsq = min_pt*min_pt;
467   Float_t maxptsq = max_pt*max_pt;
468   Float_t ptsq;
469
470   for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
471     ptsq = ((Track*)(*i))->fP.Perp2();
472     (*i)->SetRnrElement(ptsq >= minptsq && ptsq <= maxptsq);
473   }
474 }
475
476 /**************************************************************************/
477
478 void TrackList::ImportHits()
479 {
480   for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
481     ((Track*)(*i))->ImportHits();
482   }
483 }
484
485 void TrackList::ImportClusters()
486 {
487   for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
488     ((Track*)(*i))->ImportClusters();
489   }
490 }
491
492 /**************************************************************************/
493 /**************************************************************************/
494 /**************************************************************************/
495
496 #include "RGEditor.h"
497
498 //______________________________________________________________________
499 // TrackCounter
500 //
501
502 ClassImp(TrackCounter)
503
504 TrackCounter::TrackCounter(const Text_t* name, const Text_t* title) :
505   RenderElement(),
506   TNamed(name, title),
507
508   fBadLineStyle (6),
509   fClickAction  (CA_ToggleTrack),
510   fAllTracks    (0),
511   fGoodTracks   (0),
512   fTrackLists   ()
513 {
514   TQObject::Connect("Reve::Track", "CtrlClicked(Reve::Track*)",
515                     "Reve::TrackCounter", this, "DoTrackAction(Reve::Track*)");
516 }
517
518 TrackCounter::~TrackCounter()
519 {
520   TQObject::Disconnect("Reve::Track", "DoTrackAction(Reve::Track*)");
521 }
522
523 /**************************************************************************/
524
525 void TrackCounter::Reset()
526 {
527   printf("TrackCounter::Reset()\n");
528   fAllTracks  = 0;
529   fGoodTracks = 0;
530   TIter next(&fTrackLists);
531   TrackList* tlist;
532   while ((tlist = dynamic_cast<TrackList*>(next())))
533     tlist->RemoveParent(this);
534   fTrackLists.Clear();
535 }
536
537 void TrackCounter::RegisterTracks(TrackList* tlist, Bool_t goodTracks)
538 {
539   // printf("TrackCounter::RegisterTracks '%s', %s\n",
540   //   tlist->GetObject()->GetName(), goodTracks ? "good" : "bad");
541
542   tlist->AddParent(this);
543   fTrackLists.Add(tlist);
544
545   List_i i = tlist->BeginChildren();
546   while (i != tlist->EndChildren())
547   {
548     Track* t = dynamic_cast<Track*>(*i);
549     if (t != 0)
550     {
551       if (goodTracks)
552       {
553         ++fGoodTracks;
554       } else {
555         t->SetLineStyle(fBadLineStyle);
556       }
557       ++fAllTracks;
558     }
559     ++i;
560   }
561 }
562
563 void TrackCounter::DoTrackAction(Track* track)
564 {
565   // !!!! No check done if ok.
566   // !!!! Should also override RemoveElementLocal
567   // !!!! But then ... should also sore local information if track is ok.
568
569   switch (fClickAction)
570   {
571
572     case CA_PrintTrackInfo:
573     {
574       printf("Track '%s'\n", track->GetObject()->GetName());
575       Vector &v = track->fV, &p = track->fP;
576       printf("  Vx=%f, Vy=%f, Vz=%f; Pt=%f, Pz=%f, phi=%f)\n",
577              v.x, v.y, v.z, p.Perp(), p.z, TMath::RadToDeg()*p.Phi());
578       printf("  <other information should be printed ... full AliESDtrack>\n");
579       break;
580     }
581
582     case CA_ToggleTrack:
583     {
584       if (track->GetLineStyle() == 1)
585       {
586         track->SetLineStyle(fBadLineStyle);
587         --fGoodTracks;
588       } else {
589         track->SetLineStyle(1);
590         ++fGoodTracks;
591       }
592       gReve->Redraw3D();
593
594       printf("TrackCounter::CountTrack All=%d, Good=%d, Bad=%d\n",
595              fAllTracks, fGoodTracks, fAllTracks-fGoodTracks);
596
597       if (gReve->GetEditor()->GetModel() == GetObject())
598         gReve->EditRenderElement(this);
599
600       break;
601     }
602
603   } // end switch fClickAction
604 }