Added support to display the HLT ESD Tree. See the comment in visscan_local.C for How
[u/mrichter/AliRoot.git] / EVE / alice-macros / muon_trackRefs.C
1 // $Id$
2
3 /**************************************************************************
4  * Copyright(c) 1998-2008, ALICE Experiment at CERN, all rights reserved. *
5  * See http://aliceinfo.cern.ch/Offline/AliRoot/License.html for          *
6  * full copyright notice.                                                 *
7  **************************************************************************/
8
9 /// \ingroup evemacros
10 /// \file muon_trackRefs.C
11 /// \brief Macro to visualise trackRef in MUON spectrometer 
12 /// (both tracker and trigger).
13 ///
14 /// Use muon_trackRefs(Bool_t showSimClusters) in order to run it
15 ///
16 /// Needs that alieve_init() is already called
17 ///
18 /// \author P. Pillot, L. Aphecetche; Subatech
19
20 #if !defined(__CINT__) || defined(__MAKECINT__)
21 #include <TClonesArray.h>
22 #include <TTree.h>
23 #include <TParticle.h>
24 #include <TMath.h>
25 #include <TROOT.h>
26 #include <TEveManager.h>
27 #include <TEveUtil.h>
28 #include <TEveTrack.h>
29 #include <TEvePointSet.h>
30 #include <TEveVSDStructs.h>
31 #include <TEveTrackPropagator.h>
32
33 #include <AliMUONClusterStoreV2.h>
34 #include <AliMUONRawClusterV2.h>
35 #include <AliMUONVCluster.h>
36 #include <AliMUONConstants.h>
37 #include <AliMUONRecoParam.h>
38 #include <AliMUONCDB.h>
39 #include <AliStack.h>
40 #include <AliTrackReference.h>
41 #include <AliRunLoader.h>
42 #include <AliEveMagField.h>
43 #include <AliEveTrack.h>
44 #include <AliEveEventManager.h>
45 #endif
46
47 //______________________________________________________________________________
48 void muon_trackRef_propagator_setup(TEveTrackPropagator* trkProp, Bool_t showVertex)
49 {
50   // set magnetic field
51   trkProp->SetMagFieldObj(new AliEveMagField);
52   trkProp->SetStepper(TEveTrackPropagator::kRungeKutta);
53   
54   // set propagation range
55   trkProp->SetMaxR(500.);
56   trkProp->SetMaxZ(-AliMUONConstants::DefaultChamberZ(13)+3.5);
57   
58   // go through pathmarks
59   trkProp->SetFitDaughters(kFALSE);
60   trkProp->SetFitReferences(kTRUE);
61   trkProp->SetFitDecay(kTRUE);
62   trkProp->SetFitCluster2Ds(kFALSE);
63   
64   // Render first vertex if required
65   if (showVertex)
66   {
67     trkProp->SetRnrFV(kTRUE);
68     trkProp->RefFVAtt().SetMarkerSize(0.5);
69     trkProp->RefFVAtt().SetMarkerColor(kRed);
70   }
71 }
72
73 //______________________________________________________________________________
74 Bool_t isReconstructible(Bool_t* chHit)
75 {
76   // load recoParam from OCDB
77   static AliMUONRecoParam* gRecoParam = 0x0;
78   static UInt_t gRequestedStationMask = 0;
79   static Bool_t gRequest2ChInSameSt45 = kFALSE;
80   if (!gRecoParam)
81   {
82     gRecoParam = AliMUONCDB::LoadRecoParam();
83     if (!gRecoParam) exit(-1);
84     // compute the mask of requested stations
85     gRequestedStationMask = 0;
86     for (Int_t i = 0; i < 5; i++) if (gRecoParam->RequestStation(i)) gRequestedStationMask |= ( 1 << i );
87     // get whether a track need 2 chambers hit in the same station (4 or 5) or not to be reconstructible
88     gRequest2ChInSameSt45 = !gRecoParam->MakeMoreTrackCandidates();
89   }
90   
91   // check which chambers are hit
92   UInt_t presentStationMask(0);
93   Int_t nChHitInSt4 = 0, nChHitInSt5 = 0;
94   for (Int_t ich=0; ich<10; ich++)
95   {
96     if (!chHit[ich]) continue;
97     Int_t ist = ich/2;
98     presentStationMask |= ( 1 << ist );
99     if (ist == 3) nChHitInSt4++;
100     if (ist == 4) nChHitInSt5++;
101   }
102   
103   // at least one cluster per requested station
104   if ((gRequestedStationMask & presentStationMask) != gRequestedStationMask) return kFALSE;
105   
106   // 2 chambers hit in the same station (4 or 5)
107   if (gRequest2ChInSameSt45) return (nChHitInSt4 == 2 || nChHitInSt5 == 2);
108   // or 2 chambers hit in station 4 & 5 together
109   else return (nChHitInSt4+nChHitInSt5 >= 2);
110 }
111
112 //______________________________________________________________________________
113 void add_muon_trackRefs(AliStack* stack, TTree* treeTR, TEveTrackList* reco, TEveTrackList* other,
114                         TEvePointSet* RecoClusters, TEvePointSet* OtherClusters)
115 {
116   TClonesArray* trackRefs = 0;
117   treeTR->SetBranchAddress("TrackReferences", &trackRefs);
118   Int_t nSimTracks = stack->GetNtrack();
119   AliMUONClusterStoreV2 clusters;
120   TClonesArray clustersTrigger("AliMUONRawClusterV2", 10);
121   
122   // loop over simulated track
123   for (Int_t itr = 0; itr < nSimTracks; itr++)
124   {
125     treeTR->GetEntry(stack->TreeKEntry(itr));
126     Int_t nTrackRefs = trackRefs->GetEntriesFast();
127     if (nTrackRefs <= 0) continue;
128     
129     TEveTrack* track = 0x0;
130     Bool_t chHit[10];
131     for (Int_t ich=0; ich<10; ich++) chHit[ich] = kFALSE;
132     
133     // loop over simulated track hits
134     for (Int_t itrR = 0; itrR < nTrackRefs; ++itrR)
135     {
136       AliTrackReference* atr = static_cast<AliTrackReference*>(trackRefs->UncheckedAt(itrR));
137       
138       // skip trackRefs not in MUON
139       if (atr->DetectorId() != AliTrackReference::kMUON) continue;
140       
141       // record chamber hit
142       Int_t detElemId = atr->UserId();
143       Int_t chamberId = detElemId / 100 - 1;
144       if (chamberId < 0) continue;
145       if (chamberId < 10) chHit[chamberId] = kTRUE;
146       
147       // produce eve track if not already done
148       if (!track)
149       {
150         TParticle* p = stack->Particle(itr);
151         track = new AliEveTrack(p, itr, 0x0);
152         track->SetName(Form("%s [%d]", p->GetName(), itr));
153         track->SetStdTitle();
154         track->SetSourceObject(p);
155         
156         clusters.Clear();
157         clustersTrigger.Clear("C");
158       }
159       
160       // add path mark
161       track->AddPathMark(TEvePathMark(TEvePathMark::kReference,
162                                       TEveVector(atr->X(),  atr->Y(),  atr->Z()),
163                                       TEveVector(atr->Px(), atr->Py(), atr->Pz()),
164                                       atr->GetTime()));
165       
166       // produce clusters if required
167       if (RecoClusters || OtherClusters)
168       {
169         // from tracker
170         if (chamberId < 10)
171         {
172           // produce a new cluster on that DE or update existing one
173           AliMUONVCluster* cl = 0x0;
174           Int_t clNum = -1;
175           do cl = clusters.FindObject(AliMUONVCluster::BuildUniqueID(chamberId, detElemId, ++clNum));
176           while (cl && cl->GetNDigits() == 2);
177           if (cl) {
178             cl->SetXYZ((cl->GetX() + atr->X()) / 2., (cl->GetY() + atr->Y()) / 2., (cl->GetZ() + atr->Z()) / 2.);
179           }
180           else
181           {
182             cl = clusters.Add(chamberId, detElemId, clNum);
183             cl->SetXYZ(atr->X(), atr->Y(), atr->Z());
184           }
185           cl->AddDigitId((UInt_t)itrR);
186         }
187         else // from trigger
188         {
189           AliMUONVCluster* cl = new(clustersTrigger[clustersTrigger.GetLast()+1]) AliMUONRawClusterV2();
190           cl->SetXYZ(atr->X(), atr->Y(), atr->Z());
191         }
192       }
193     }
194     
195     if (track)
196     {
197       track->SortPathMarksByTime();
198       // stop track propagation at last path mark
199       track->RefPathMarks().back().fType = TEvePathMarkT<double>::EType_e(TEvePathMark::kDecay);
200       
201       // add the track and trackRefs to proper lists
202       if (isReconstructible(chHit)) {
203         track->SetPropagator(reco->GetPropagator());
204         track->SetAttLineAttMarker(reco);
205         reco->AddElement(track);
206         
207         // trackRefs
208         if (RecoClusters)
209         {
210           // tracker
211           TIter next(clusters.CreateIterator());
212           gROOT->ProcessLine(Form("add_muon_clusters((TIter*)%p, (TEvePointSet*)%p);",&next, RecoClusters));
213           // trigger
214           TIter next2(clustersTrigger.MakeIterator());
215           gROOT->ProcessLine(Form("add_muon_clusters((TIter*)%p, (TEvePointSet*)%p);",&next2, RecoClusters));
216         }
217       }
218       else {
219         track->SetPropagator(other->GetPropagator());
220         track->SetAttLineAttMarker(other);
221         other->AddElement(track);
222         
223         // trackRefs
224         if (OtherClusters)
225         {
226           // tracker
227           TIter next(clusters.CreateIterator());
228           gROOT->ProcessLine(Form("add_muon_clusters((TIter*)%p, (TEvePointSet*)%p);",&next, OtherClusters));
229           // trigger
230           TIter next2(clustersTrigger.MakeIterator());
231           gROOT->ProcessLine(Form("add_muon_clusters((TIter*)%p, (TEvePointSet*)%p);",&next2, OtherClusters));
232         }
233       }
234     }
235   }
236   delete trackRefs;
237 }
238
239 //______________________________________________________________________________
240 void muon_trackRefs(Bool_t showSimClusters)
241 {
242   // load kinematics and trackRefs
243   AliRunLoader* rl =  AliEveEventManager::AssertRunLoader();
244   rl->LoadKinematics();
245   AliStack* stack = rl->Stack();
246   if (!stack) return;
247   rl->LoadTrackRefs();
248   TTree* treeTR = rl->TreeTR();  
249   if (!treeTR) return;
250   
251   // track containers
252   TEveElementList* trackCont = new TEveElementList("Sim MUON Tracks");
253   
254   TEveTrackList* reco = new TEveTrackList("reconstructible");
255   reco->SetRnrPoints(kFALSE);
256   reco->SetRnrLine(kTRUE);
257   reco->SetLineColor(kRed);
258   reco->SetLineStyle(2);
259   muon_trackRef_propagator_setup(reco->GetPropagator(), kTRUE);
260   trackCont->AddElement(reco);
261   
262   TEveTrackList* other = new TEveTrackList("others");
263   other->SetRnrPoints(kFALSE);
264   other->SetRnrLine(kTRUE);
265   other->SetLineColor(kWhite);
266   other->SetLineStyle(3);
267   muon_trackRef_propagator_setup(other->GetPropagator(), kFALSE);
268   trackCont->AddElement(other);
269   
270   // cluster container
271   TEveElementList* clusterCont = 0x0;
272   TEvePointSet* RecoClusters = 0x0;
273   TEvePointSet* OtherClusters = 0x0;
274   if (showSimClusters)
275   {
276     clusterCont = new TEveElementList("Sim MUON Clusters");
277     
278     RecoClusters = new TEvePointSet(1000);
279     RecoClusters->SetName("Reconstructibles");
280     RecoClusters->SetPickable(kFALSE);
281     RecoClusters->SetMarkerStyle(2);
282     RecoClusters->SetMarkerColor(kRed);
283     RecoClusters->SetMarkerSize(0.5);
284     clusterCont->AddElement(RecoClusters);
285     
286     OtherClusters = new TEvePointSet(10000);
287     OtherClusters->SetName("Others");
288     OtherClusters->SetPickable(kFALSE);
289     OtherClusters->SetMarkerStyle(2);
290     OtherClusters->SetMarkerColor(kWhite);
291     OtherClusters->SetMarkerSize(0.5);
292     clusterCont->AddElement(OtherClusters);
293     
294     TEveUtil::LoadMacro("muon_clusters.C+");
295   }
296   
297   // add tracks to the proper list and propagate them. Add also clusters if required.
298   add_muon_trackRefs(stack, treeTR, reco, other, RecoClusters, OtherClusters);
299   reco->MakeTracks();
300   other->MakeTracks();
301   trackCont->SetTitle(Form("N=%d", reco->NumChildren()+other->NumChildren()));
302   reco->SetTitle(Form("N=%d",reco->NumChildren()));
303   other->SetTitle(Form("N=%d",other->NumChildren()));
304   if (showSimClusters)
305   {
306     clusterCont->SetTitle(Form("N=%d",RecoClusters->GetLastPoint()+OtherClusters->GetLastPoint()+2));
307     RecoClusters->SetTitle(Form("N=%d",RecoClusters->GetLastPoint()+1));
308     OtherClusters->SetTitle(Form("N=%d",OtherClusters->GetLastPoint()+1));
309   }
310   
311   // add graphic containers
312   gEve->DisableRedraw();
313   gEve->AddElement(trackCont);
314   if (clusterCont) gEve->AddElement(clusterCont);
315   gEve->EnableRedraw();
316   gEve->Redraw3D();
317 }