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