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