TrackRnrStyle: add pt-range memebers.
[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::ImportDaughters()
269 {
270   static const Exc_t eH("Track::ImportDaughters ");
271
272   if (fLabel < 0)
273     throw(eH + "label not set.");
274
275   Reve::LoadMacro("kine_daughter_tracks.C");
276   gROOT->ProcessLine(Form("kine_daughter_tracks(%d);", fLabel));
277 }
278
279 /**************************************************************************/
280
281 void Track::PrintKineStack()
282 {
283   Reve::LoadMacro("print_kine_from_label.C");
284   gROOT->ProcessLine(Form("print_kine_from_label(%d);", fLabel));
285 }
286
287
288 void Track::PrintPathMarks()
289 {
290   static const Exc_t eH("Track::PrintPathMarks ");
291
292   if (fLabel < 0)
293     throw(eH + "label not set.");
294
295   printf("Number of path marks %d label %d\n",
296          fPathMarks.size(), fLabel);
297
298   PathMark* pm;
299   for(vpPathMark_i i=fPathMarks.begin(); i!=fPathMarks.end(); i++) 
300   {
301     pm = *i;
302     printf("Reve::PathMark: %-9s  p: %8f %8f %8f Vertex: %8e %8e %8e %g \n",
303            pm->type_name(),
304            pm->P.x,  pm->P.y, pm->P.z,
305            pm->V.x,  pm->V.y, pm->V.z,
306            pm->time);
307   }
308 }
309
310 /**************************************************************************/
311
312 void Track::CtrlClicked(Reve::Track* track)
313 {
314   Emit("CtrlClicked(Reve::Track*)", (Long_t)track);
315 }
316
317
318 /**************************************************************************/
319 /**************************************************************************/
320
321 //______________________________________________________________________
322 // TrackRnrStyle
323 //
324
325 ClassImp(Reve::TrackRnrStyle)
326
327 Float_t       TrackRnrStyle::fgDefMagField = 5;
328 const Float_t TrackRnrStyle::fgkB2C        = 0.299792458e-3;
329 TrackRnrStyle TrackRnrStyle::fgDefStyle;
330
331 TrackRnrStyle::TrackRnrStyle() :
332   TObject(),
333
334   fColor(1),
335   fWidth(1),
336   fMagField(fgDefMagField),
337
338   fMaxR  (350),
339   fMaxZ  (450),
340
341   fMaxOrbs (0.5),
342   fMinAng  (45),
343   fDelta   (0.1),
344
345   fMinPt   (0.1),
346   fMaxPt   (10),
347
348   fFitDaughters  (kTRUE),
349   fFitReferences (kTRUE),
350   fFitDecay      (kTRUE),
351
352   fRnrDaughters  (kTRUE),
353   fRnrReferences (kTRUE),
354   fRnrDecay      (kTRUE)
355 {}
356 /**************************************************************************/
357 /**************************************************************************/
358
359 //______________________________________________________________________
360 // TrackList
361 //
362
363 ClassImp(Reve::TrackList)
364
365 void TrackList::Init()
366 {
367   fMarkerStyle = 5;
368   fMarkerColor = 5;
369   // fMarker->SetMarkerSize(0.05);
370
371   if (fRnrStyle== 0) fRnrStyle = new TrackRnrStyle;
372   SetMainColorPtr(&fRnrStyle->fColor);
373 }
374
375 TrackList::TrackList(Int_t n_tracks, TrackRnrStyle* rs) :
376   RenderElementListBase(),
377   TPolyMarker3D(n_tracks),
378
379   fTitle(),
380
381   fRnrStyle      (rs),
382   fRnrTracks     (kTRUE),
383   fEditPathMarks (kFALSE)
384 {
385   Init();
386 }
387
388 TrackList::TrackList(const Text_t* name, Int_t n_tracks, TrackRnrStyle* rs) :
389   RenderElementListBase(),
390   TPolyMarker3D(n_tracks),
391   
392   fTitle(),
393
394   fRnrStyle   (rs),
395   fRnrTracks  (kTRUE)
396 {
397   Init();
398   SetName(name);
399 }
400
401 void TrackList::Reset(Int_t n_tracks)
402 {
403   delete [] fP; fP = 0;
404   fN = n_tracks;
405   if(fN) fP = new Float_t [3*fN];
406   memset(fP, 0, 3*fN*sizeof(Float_t));
407   fLastPoint = -1;
408 }
409
410 /**************************************************************************/
411
412 void TrackList::Paint(Option_t* option)
413 {
414   if(fRnrElement) {
415     if(fRnrMarkers) {
416       TPolyMarker3D::Paint(option);
417     }
418     if(fRnrTracks) {
419       for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
420         if((*i)->GetRnrElement())
421           (*i)->GetObject()->Paint(option);
422       }
423     }
424   }
425 }
426
427 /**************************************************************************/
428
429 void TrackList::AddElement(RenderElement* el)
430 {
431   static const Exc_t eH("TrackList::AddElement ");
432   if (dynamic_cast<Track*>(el)  == 0)
433     throw(eH + "new element not a Track.");
434   RenderElementListBase::AddElement(el);
435 }
436
437 /**************************************************************************/
438
439 void TrackList::MakeTracks()
440 {
441   for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
442     ((Track*)(*i))->MakeTrack();
443   }
444   gReve->Redraw3D();
445 }
446
447
448 void TrackList::MakeMarkers()
449 {
450   Reset(fChildren.size());
451   for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
452     Track& t = *((Track*)(*i));
453     if(t.GetN() > 0)
454       SetNextPoint(t.fV.x, t.fV.y, t.fV.z);
455   }
456   gReve->Redraw3D();
457 }
458
459 /**************************************************************************/
460 /*************************************************************************/
461
462 void TrackList::SetWidth(Width_t w)
463 {
464   Width_t oldw = fRnrStyle->fWidth;
465   fRnrStyle->fWidth = w;
466   for (List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
467     Track& t = *((Track*)(*i));
468     if (t.GetLineWidth() == oldw)
469       t.SetLineWidth(w);
470   }
471 }
472
473 void TrackList::SetMaxR(Float_t x)
474 {
475   fRnrStyle->fMaxR = x;
476   MakeTracks();
477   MakeMarkers();
478 }
479
480 void TrackList::SetMaxZ(Float_t x)
481 {
482   fRnrStyle->fMaxZ = x;
483   MakeTracks();
484   MakeMarkers();
485 }
486
487 void TrackList::SetMaxOrbs(Float_t x)
488 {
489   fRnrStyle->fMaxOrbs = x;
490   MakeTracks();
491 }
492
493 void TrackList::SetMinAng(Float_t x)
494 {
495   fRnrStyle->fMinAng = x;
496   MakeTracks();
497 }
498
499 void TrackList::SetDelta(Float_t x)
500 {
501   fRnrStyle->fDelta = x;
502   MakeTracks();
503 }
504
505 void TrackList::SetFitDaughters(Bool_t x)
506 {
507   fRnrStyle->fFitDaughters = x;
508   MakeTracks();
509 }
510
511 void TrackList::SetFitReferences(Bool_t x)
512 {
513   fRnrStyle->fFitReferences = x;
514   MakeTracks();
515 }
516
517 void TrackList::SetFitDecay(Bool_t x)
518 {
519   fRnrStyle->fFitDecay = x;
520   MakeTracks();
521 }
522
523 void TrackList::SetRnrDecay(Bool_t rnr)
524 {
525   fRnrStyle->fRnrDecay = rnr;
526   MakeTracks();
527 }
528
529 void TrackList::SetRnrDaughters(Bool_t rnr)
530 {
531   fRnrStyle->fRnrDaughters = rnr;
532   MakeTracks();
533 }
534
535 void TrackList::SetRnrReferences(Bool_t rnr)
536 {
537   fRnrStyle->fRnrReferences = rnr;
538   MakeTracks();
539 }
540  
541 void TrackList::SetRnrMarkers(Bool_t rnr)
542 {
543   fRnrMarkers = rnr;
544   gReve->Redraw3D();
545 }
546
547 void TrackList::SetRnrTracks(Bool_t rnr)
548 {
549
550   fRnrTracks = rnr;
551   gReve->Redraw3D();
552 }
553
554 /**************************************************************************/
555 /**************************************************************************/
556
557 void TrackList::SelectByPt(Float_t min_pt, Float_t max_pt)
558 {
559   fRnrStyle->fMinPt = min_pt;
560   fRnrStyle->fMaxPt = max_pt;
561
562   Float_t minptsq = min_pt*min_pt;
563   Float_t maxptsq = max_pt*max_pt;
564   Float_t ptsq;
565
566   for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
567     ptsq = ((Track*)(*i))->fP.Perp2();
568     (*i)->SetRnrElement(ptsq >= minptsq && ptsq <= maxptsq);
569   }
570 }
571
572 /**************************************************************************/
573
574 void TrackList::ImportHits()
575 {
576   for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
577     ((Track*)(*i))->ImportHits();
578   }
579 }
580
581 void TrackList::ImportClusters()
582 {
583   for(List_i i=fChildren.begin(); i!=fChildren.end(); ++i) {
584     ((Track*)(*i))->ImportClusters();
585   }
586 }
587
588 /**************************************************************************/
589 /**************************************************************************/
590 /**************************************************************************/
591
592 #include "RGEditor.h"
593
594 //______________________________________________________________________
595 // TrackCounter
596 //
597
598 ClassImp(TrackCounter)
599
600 TrackCounter* TrackCounter::fgInstance = 0;
601
602 TrackCounter::TrackCounter(const Text_t* name, const Text_t* title) :
603   RenderElement(),
604   TNamed(name, title),
605
606   fBadLineStyle (6),
607   fClickAction  (CA_ToggleTrack),
608   fAllTracks    (0),
609   fGoodTracks   (0),
610   fTrackLists   ()
611 {
612   if (fgInstance == 0) fgInstance = this;
613   TQObject::Connect("Reve::Track", "CtrlClicked(Reve::Track*)",
614                     "Reve::TrackCounter", this, "DoTrackAction(Reve::Track*)");
615 }
616
617 TrackCounter::~TrackCounter()
618 {
619   TQObject::Disconnect("Reve::Track", "DoTrackAction(Reve::Track*)");
620   if (fgInstance == this) fgInstance = 0;
621 }
622
623 /**************************************************************************/
624
625 void TrackCounter::Reset()
626 {
627   printf("TrackCounter::Reset()\n");
628   fAllTracks  = 0;
629   fGoodTracks = 0;
630   TIter next(&fTrackLists);
631   TrackList* tlist;
632   while ((tlist = dynamic_cast<TrackList*>(next())))
633     tlist->RemoveParent(this);
634   fTrackLists.Clear();
635 }
636
637 void TrackCounter::RegisterTracks(TrackList* tlist, Bool_t goodTracks)
638 {
639   // printf("TrackCounter::RegisterTracks '%s', %s\n",
640   //   tlist->GetObject()->GetName(), goodTracks ? "good" : "bad");
641
642   tlist->AddParent(this);
643   fTrackLists.Add(tlist);
644
645   List_i i = tlist->BeginChildren();
646   while (i != tlist->EndChildren())
647   {
648     Track* t = dynamic_cast<Track*>(*i);
649     if (t != 0)
650     {
651       if (goodTracks)
652       {
653         ++fGoodTracks;
654       } else {
655         t->SetLineStyle(fBadLineStyle);
656       }
657       ++fAllTracks;
658     }
659     ++i;
660   }
661 }
662
663 void TrackCounter::DoTrackAction(Track* track)
664 {
665   // !!!! No check done if ok.
666   // !!!! Should also override RemoveElementLocal
667   // !!!! But then ... should also sotre local information if track is ok.
668
669   switch (fClickAction)
670   {
671
672     case CA_PrintTrackInfo:
673     {
674       printf("Track '%s'\n", track->GetObject()->GetName());
675       Vector &v = track->fV, &p = track->fP;
676       printf("  Vx=%f, Vy=%f, Vz=%f; Pt=%f, Pz=%f, phi=%f)\n",
677              v.x, v.y, v.z, p.Perp(), p.z, TMath::RadToDeg()*p.Phi());
678       printf("  <other information should be printed ... full AliESDtrack>\n");
679       break;
680     }
681
682     case CA_ToggleTrack:
683     {
684       if (track->GetLineStyle() == 1)
685       {
686         track->SetLineStyle(fBadLineStyle);
687         --fGoodTracks;
688       } else {
689         track->SetLineStyle(1);
690         ++fGoodTracks;
691       }
692       gReve->Redraw3D();
693
694       printf("TrackCounter::CountTrack All=%d, Good=%d, Bad=%d\n",
695              fAllTracks, fGoodTracks, fAllTracks-fGoodTracks);
696
697       if (gReve->GetEditor()->GetModel() == GetObject())
698         gReve->EditRenderElement(this);
699
700       break;
701     }
702
703   } // end switch fClickAction
704 }
705
706 /**************************************************************************/
707
708 void TrackCounter::OutputEventTracks(FILE* out)
709 {
710   if (out == 0)
711   {
712     out = stdout;
713     fprintf(out, "TrackCounter::FinalizeEvent()\n");
714   }
715
716   fprintf(out, "Event = %d  Ntracks = %d\n", fEventId, fGoodTracks);
717
718   TIter tlists(&fTrackLists);
719   TrackList* tlist;
720   Int_t cnt = 0;
721   while ((tlist = (TrackList*) tlists()) != 0)
722   {
723     List_i i = tlist->BeginChildren();
724     while (i != tlist->EndChildren())
725     {
726       Track* t = dynamic_cast<Track*>(*i);
727       if (t != 0 && t->GetLineStyle() == 1)
728       {
729         ++cnt;
730         fprintf(out, " %2d: chg=%+2d  pt=%8.5f  eta=%+8.5f\n",
731                cnt, t->fCharge, t->fP.Perp(), t->fP.Eta());
732       }
733       ++i;
734     }
735   }
736 }