Fixed bugs that made histograms being drawn on top of each other
[u/mrichter/AliRoot.git] / HLT / EVE / AliHLTEveHLT.cxx
1 /**************************************************************************
2  * This file is property of and copyright by the ALICE HLT Project        *
3  * ALICE Experiment at CERN, All rights reserved.                         *
4  *                                                                        *
5  * Primary Authors: Svein Lindal <slindal@fys.uio.no   >                  *
6  *                  for The ALICE HLT Project.                            *
7  *                                                                        *
8  * Permission to use, copy, modify and distribute this software and its   *
9  * documentation strictly for non-commercial purposes is hereby granted   *
10  * without fee, provided that the above copyright notice appears in all   *
11  * copies and that both the copyright notice and this permission notice   *
12  * appear in the supporting documentation. The authors make no claims     *
13  * about the suitability of this software for any purpose. It is          *
14  * provided "as is" without express or implied warranty.                  *
15  **************************************************************************/
16
17 /// @file   AliHLTEvePhos.cxx
18 /// @author Svein Lindal <slindal@fys.uio.no>
19 /// @brief  HLT class for the HLT EVE display
20
21 #include "AliHLTEveHLT.h"
22 #include "AliHLTHOMERBlockDesc.h"
23 #include "AliHLTEveBase.h"
24 #include "AliEveHOMERManager.h"
25 #include "TEveManager.h"
26 #include "TEvePointSet.h"
27 #include "TEveTrack.h"
28 #include "TCanvas.h"
29 #include "AliESDEvent.h"
30 #include "TEveTrackPropagator.h"
31 #include "AliEveTrack.h"
32 #include "TEveVSDStructs.h"
33 #include "TString.h"
34 #include "TPCLib/tracking-ca/AliHLTTPCCATrackParam.h"
35 #include "TPCLib/tracking-ca/AliHLTTPCCATrackConvertor.h"
36 #include "AliEveMagField.h"
37 #include "TH1F.h"
38 #include "TH2F.h"
39
40 ClassImp(AliHLTEveHLT)
41
42 AliHLTEveHLT::AliHLTEveHLT() : 
43   AliHLTEveBase(), 
44   fTrueField(kFALSE),
45   fUseIpOnFailedITS(kFALSE),
46   fUseRkStepper(kFALSE),
47   fTrackList(NULL)
48 {
49   // Constructor.
50 }
51
52 AliHLTEveHLT::~AliHLTEveHLT()
53 {
54   //Destructor, not implemented
55   if(fTrackList)
56     delete fTrackList;
57   fTrackList = NULL;
58 }
59
60
61 void AliHLTEveHLT::ProcessBlock(AliHLTHOMERBlockDesc * block) {
62   //See header file for documentation
63   if ( ! block->GetDataType().CompareTo("ALIESDV0") ) {
64     if(!fTrackList) CreateTrackList();
65     ProcessEsdBlock(block, fTrackList);
66   } 
67   
68   else if ( ! block->GetDataType().CompareTo("ROOTTOBJ") ) {
69     //processROOTTOBJ( block, gHLTText );
70   } 
71
72   else if ( ! block->GetDataType().CompareTo("HLTRDLST") ) {
73     //processHLTRDLST( block );
74   } 
75
76   else if ( !block->GetDataType().CompareTo("ROOTHIST") ) {      
77     if( !fCanvas ) { 
78       fCanvas = CreateCanvas("Primary Vertex", "Primary Vertex");
79       fCanvas->Divide(3, 2);
80     }
81     ProcessHistograms( block , fCanvas);
82   }
83   
84 }
85
86
87 void AliHLTEveHLT::UpdateElements() {
88   //See header file for documentation
89   if(fCanvas) fCanvas->Update();
90   if(fTrackList) fTrackList->ElementChanged();
91
92 }
93
94 void AliHLTEveHLT::ResetElements(){
95     //See header file for documentation
96   if(fTrackList) fTrackList->DestroyElements();
97   fHistoCount = 0;
98
99 }
100
101 void AliHLTEveHLT::ProcessHistograms(AliHLTHOMERBlockDesc * block, TCanvas * canvas) {
102   //See header file for documentation
103   if ( ! block->GetClassName().CompareTo("TH1F")) {
104     TH1F* histo = reinterpret_cast<TH1F*>(block->GetTObject());
105     if( histo ){
106       TString name(histo->GetName());
107       if( !name.CompareTo("primVertexZ") ){
108         canvas->cd(2);
109         histo->Draw();
110       }else if( !name.CompareTo("primVertexX") ){
111         canvas->cd(3);
112         histo->Draw();
113       }else if( !name.CompareTo("primVertexY") ){
114         canvas->cd(4);
115         histo->Draw();
116       }
117     }
118   }  else if ( ! block->GetClassName().CompareTo("TH2F")) {
119     TH2F *hista = reinterpret_cast<TH2F*>(block->GetTObject());
120     if (hista ){
121        TString name(hista->GetName());
122        if( !name.CompareTo("primVertexXY")) {      
123          canvas->cd(1);
124          hista->Draw();
125        }
126     }
127   }
128   canvas->cd();
129
130
131
132
133 }
134
135 void AliHLTEveHLT::CreateTrackList() {
136   //See header file for documentation
137   fTrackList = new TEveTrackList("ESD Tracks");
138   fTrackList->SetMainColor(6);
139   gEve->AddElement(fTrackList);
140 }
141
142
143 void AliHLTEveHLT::ProcessEsdBlock( AliHLTHOMERBlockDesc * block, TEveTrackList * cont ) {
144   //See header file for documentation
145
146   AliESDEvent* esd = (AliESDEvent *) (block->GetTObject());
147   esd->GetStdContent();
148   
149   SetUpTrackPropagator(cont->GetPropagator(),-0.1*esd->GetMagneticField(), 520);
150
151   for (Int_t iter = 0; iter < esd->GetNumberOfTracks(); ++iter) {
152     AliEveTrack* track = dynamic_cast<AliEveTrack*>(MakeEsdTrack(esd->GetTrack(iter), cont));
153     cont->AddElement(track);
154     
155 //     gTPCPt->Fill(esd->GetTrack(iter)->GetSignedPt()); // KK
156 //     gTPCEta->Fill(esd->GetTrack(iter)->GetSnp());
157 //     gTPCPsi->Fill(esd->GetTrack(iter)->GetTgl());
158 //     gTPCnClusters->Fill(esd->GetTrack(iter)->GetTPCNcls());  
159   }
160   
161 //   gTPCMult->Fill(esd->GetNumberOfTracks()); // KK
162   
163 //   Int_t icd = 0;
164 //   gTPCClustCanvas->Clear();
165 //   gTPCClustCanvas->Divide(2, 2);
166 //   gTPCClustCanvas->cd(icd++);
167 //   gTPCPt->Draw();
168 //   gTPCClustCanvas->cd(icd++);
169 //   gTPCEta->Draw();
170 //   gTPCClustCanvas->cd(icd++);
171 //   gTPCPsi->Draw();
172 //   gTPCClustCanvas->cd(icd++);
173 //   gTPCnClusters->Draw();
174 //   gTPCClustCanvas->cd(icd++);
175 //   gTPCMult->Draw();
176 //   gTPCClustCanvas->Update();
177
178   
179   cont->SetTitle(Form("N=%d", esd->GetNumberOfTracks()) );
180   cont->MakeTracks();
181
182 }
183
184 AliEveTrack* AliHLTEveHLT::MakeEsdTrack (AliESDtrack *at, TEveTrackList* cont) {
185   //See header file for documentation
186
187
188   const double kCLight = 0.000299792458;
189   double bz = - kCLight*10.*( cont->GetPropagator()->GetMagField(0,0,0).fZ);
190
191   Bool_t innerTaken = kFALSE;
192   if ( ! at->IsOn(AliESDtrack::kITSrefit) && fUseIpOnFailedITS)
193   {
194     //tp = at->GetInnerParam();
195     innerTaken = kTRUE;
196   }
197
198   // Add inner/outer track parameters as path-marks.
199
200   Double_t     pbuf[3], vbuf[3];
201
202   AliExternalTrackParam trackParam = *at;
203
204   // take parameters constrained to vertex (if they are)
205
206   if( at->GetConstrainedParam() ){
207     trackParam = *at->GetConstrainedParam();
208   }
209   else if( at->GetInnerParam() ){
210     trackParam = *(at->GetInnerParam());
211   }
212   if( at->GetStatus()&AliESDtrack::kTRDin ){
213     // transport to TRD in
214     trackParam = *at;
215     trackParam.PropagateTo( 290.45, -10.*( cont->GetPropagator()->GetMagField(0,0,0).fZ) );
216   }
217
218   TEveRecTrack rt;
219   {
220     rt.fLabel  = at->GetLabel();
221     rt.fIndex  = (Int_t) at->GetID();
222     rt.fStatus = (Int_t) at->GetStatus();
223     rt.fSign   = (Int_t) trackParam.GetSign();  
224     trackParam.GetXYZ(vbuf);
225     trackParam.GetPxPyPz(pbuf);    
226     rt.fV.Set(vbuf);
227     rt.fP.Set(pbuf);
228     Double_t ep = at->GetP(), mc = at->GetMass();
229     rt.fBeta = ep/TMath::Sqrt(ep*ep + mc*mc);
230   }
231
232   AliEveTrack* track = new AliEveTrack(&rt, cont->GetPropagator());
233   track->SetAttLineAttMarker(cont);
234   track->SetName(Form("AliEveTrack %d", at->GetID()));
235   track->SetElementTitle(CreateTrackTitle(at));
236   track->SetSourceObject(at);
237
238
239   // Set reference points along the trajectory
240   // and the last point
241
242   { 
243     TEvePathMark startPoint(TEvePathMark::kReference);
244     trackParam.GetXYZ(vbuf);
245     trackParam.GetPxPyPz(pbuf);    
246     startPoint.fV.Set(vbuf);
247     startPoint.fP.Set(pbuf);
248     rt.fV.Set(vbuf);
249     rt.fP.Set(pbuf);
250     Double_t ep = at->GetP(), mc = at->GetMass();
251     rt.fBeta = ep/TMath::Sqrt(ep*ep + mc*mc);
252
253     track->AddPathMark( startPoint );    
254   }
255
256
257   if( at->GetTPCPoints(2)>80 ){
258   
259     //
260     // use AliHLTTPCCATrackParam propagator 
261     // since AliExternalTrackParam:PropagateTo()
262     // has an offset at big distances
263     //
264     
265     AliHLTTPCCATrackParam t;
266     AliHLTTPCCATrackConvertor::SetExtParam( t, trackParam );
267     
268     Double_t x0 = trackParam.GetX();
269     Double_t dx = at->GetTPCPoints(2) - x0;
270     
271     //
272     // set a reference at the half of trajectory for better drawing
273     //
274     
275     for( double dxx=dx/2; TMath::Abs(dxx)>=1.; dxx*=.9 ){
276       if( !t.TransportToX(x0+dxx, bz, .99 ) ) continue;
277       AliHLTTPCCATrackConvertor::GetExtParam( t, trackParam, trackParam.GetAlpha() ); 
278       trackParam.GetXYZ(vbuf);
279       trackParam.GetPxPyPz(pbuf);
280       TEvePathMark midPoint(TEvePathMark::kReference);
281       midPoint.fV.Set(vbuf);
282       midPoint.fP.Set(pbuf);    
283       track->AddPathMark( midPoint );
284       break;
285     }
286     
287     //
288     // Set a reference at the end of the trajectory
289     // and a "decay point", to let the event display know where the track ends
290     //
291     
292     for( ; TMath::Abs(dx)>=1.; dx*=.9 ){
293       if( !t.TransportToX(x0+dx, bz, .99 ) ) continue;
294       AliHLTTPCCATrackConvertor::GetExtParam( t, trackParam, trackParam.GetAlpha() ); 
295       trackParam.GetXYZ(vbuf);
296       trackParam.GetPxPyPz(pbuf);
297       TEvePathMark endPoint(TEvePathMark::kReference);
298       TEvePathMark decPoint(TEvePathMark::kDecay);
299       endPoint.fV.Set(vbuf);
300       endPoint.fP.Set(pbuf);
301       decPoint.fV.Set(vbuf);
302       decPoint.fP.Set(pbuf);
303       track->AddPathMark( endPoint );
304       track->AddPathMark( decPoint );
305       break;
306     }  
307   }
308
309   if (at->IsOn(AliESDtrack::kTPCrefit))
310   {
311     if ( ! innerTaken)
312     {
313       AddTrackParamToTrack(track, at->GetInnerParam());
314     }
315     AddTrackParamToTrack(track, at->GetOuterParam());
316   }
317   return track;
318 }
319
320 void AliHLTEveHLT::SetUpTrackPropagator(TEveTrackPropagator* trkProp, Float_t magF, Float_t maxR) {
321   //See header file for documentation
322
323   if (fTrueField) {
324     trkProp->SetMagFieldObj(new AliEveMagField);
325   
326   } else {
327     trkProp->SetMagField(magF);
328   }
329  
330   if (fUseRkStepper) {
331     trkProp->SetStepper(TEveTrackPropagator::kRungeKutta);
332   }
333
334   trkProp->SetMaxR(maxR);
335 }
336
337
338 void AliHLTEveHLT::AddTrackParamToTrack(AliEveTrack* track, const AliExternalTrackParam* tp) {
339   //See header file for documentation
340
341   if (tp == 0)
342     return;
343
344   Double_t pbuf[3], vbuf[3];
345   tp->GetXYZ(vbuf);
346   tp->GetPxPyPz(pbuf);
347
348   TEvePathMark pm(TEvePathMark::kReference);
349   pm.fV.Set(vbuf);
350   pm.fP.Set(pbuf);
351   track->AddPathMark(pm);
352 }
353
354
355
356 TString AliHLTEveHLT::CreateTrackTitle(AliESDtrack* t) {
357   // Add additional track parameters as a path-mark to track.
358
359   TString s;
360
361   Int_t label = t->GetLabel(), index = t->GetID();
362   TString idx(index == kMinInt ? "<undef>" : Form("%d", index));
363   TString lbl(label == kMinInt ? "<undef>" : Form("%d", label));
364
365   Double_t p[3], v[3];
366   t->GetXYZ(v);
367   t->GetPxPyPz(p);
368   Double_t pt    = t->Pt();
369   Double_t ptsig = TMath::Sqrt(t->GetSigma1Pt2());
370   Double_t ptsq  = pt*pt;
371   Double_t ptm   = pt / (1.0 + pt*ptsig);
372   Double_t ptM   = pt / (1.0 - pt*ptsig);
373
374   s = Form("Index=%s, Label=%s\nChg=%d, Pdg=%d\n"
375            "pT = %.3f + %.3f - %.3f [%.3f]\n"
376            "P  = (%.3f, %.3f, %.3f)\n"
377            "V  = (%.3f, %.3f, %.3f)\n",
378            idx.Data(), lbl.Data(), t->Charge(), 0,
379            pt, ptM - pt, pt - ptm, ptsig*ptsq,
380            p[0], p[1], p[2],
381            v[0], v[1], v[2]);
382
383   Int_t   o;
384   s += "Det (in,out,refit,pid):\n";
385   o  = AliESDtrack::kITSin;
386   s += Form("ITS (%d,%d,%d,%d)  ",  t->IsOn(o), t->IsOn(o<<1), t->IsOn(o<<2), t->IsOn(o<<3));
387   o  = AliESDtrack::kTPCin;
388   s += Form("TPC(%d,%d,%d,%d)\n",   t->IsOn(o), t->IsOn(o<<1), t->IsOn(o<<2), t->IsOn(o<<3));
389   o  = AliESDtrack::kTRDin;
390   s += Form("TRD(%d,%d,%d,%d) ",    t->IsOn(o), t->IsOn(o<<1), t->IsOn(o<<2), t->IsOn(o<<3));
391   o  = AliESDtrack::kTOFin;
392   s += Form("TOF(%d,%d,%d,%d)\n",   t->IsOn(o), t->IsOn(o<<1), t->IsOn(o<<2), t->IsOn(o<<3));
393   o  = AliESDtrack::kHMPIDout;
394   s += Form("HMPID(out=%d,pid=%d)\n", t->IsOn(o), t->IsOn(o<<1));
395   s += Form("ESD pid=%d", t->IsOn(AliESDtrack::kESDpid));
396
397   if (t->IsOn(AliESDtrack::kESDpid))
398   {
399     Double_t pid[5];
400     t->GetESDpid(pid);
401     s += Form("\n[%.2f %.2f %.2f %.2f %.2f]", pid[0], pid[1], pid[2], pid[3], pid[4]);
402   }
403
404   return s;
405 }
406