From Bogdan: new macro files for MUON display.
[u/mrichter/AliRoot.git] / EVE / alice-macros / MUON_display.C
CommitLineData
f6ee486c 1#include "TGLViewer.h"
2
3namespace Alieve {
4class MUONData;
5class Event;
6}
7
8Alieve::MUONData* g_muon_data = 0;
9Alieve::Event* g_muon_last_event = 0;
10
11Int_t g_currentEvent = -1;
12Bool_t g_fromRaw = kFALSE;
13
14void MUON_display(Bool_t fromRaw = kFALSE, Bool_t showTracks = kTRUE)
15{
16
17 TTree* dt = 0;
18 TTree* ct = 0;
19
20 if (Alieve::gEvent == 0) {
21 printf("No alieve event: use alieve_init(...) \n");
22 return;
23 }
24
25 if (g_currentEvent == Alieve::gEvent->GetEventId()) {
26 if (g_fromRaw == fromRaw) {
27 printf("Same event... \n");
28 return;
29 } else {
30 if (g_fromRaw) {
31 printf("Same event with digits.\n");
32 Alieve::gEvent->GotoEvent(g_currentEvent);
33 } else {
34 printf("Same event with raw.\n");
35 Alieve::gEvent->GotoEvent(g_currentEvent);
36 }
37 }
38 }
39
40 g_fromRaw = fromRaw;
41
42 TString dataPath = TString(Alieve::gEvent->GetTitle());
43 dataPath.Append("/raw.root");
44
45 AliRunLoader* rl = Alieve::Event::AssertRunLoader();
46 g_muon_data = new Alieve::MUONData;
47
48 if (!fromRaw) {
49 rl->LoadDigits("MUON");
50 dt = rl->GetTreeD("MUON", false);
51 if (dt == 0) {
52 cout << "No digits produced!" << endl;
53 } else {
54 cout << "Display aliroot digits!" << endl;
55 g_muon_data->LoadDigits(dt);
56 }
57 } else {
58 if (gSystem->AccessPathName(dataPath.Data(),kFileExists)) {
59 cout << "No raw data produced!" << endl;
60 } else {
61 cout << "Display raw digits!" << endl;
62 g_muon_data->LoadRaw(dataPath.Data());
63 }
64 }
65
66 rl->LoadRecPoints("MUON");
67 ct = rl->GetTreeR("MUON", false);
68 g_muon_data->LoadRecPoints(ct);
69
70 g_muon_last_event = Alieve::gEvent;
71
72 g_currentEvent = g_muon_last_event->GetEventId();
73
74 gStyle->SetPalette(1, 0);
75
76 gReve->DisableRedraw();
77
78 Reve::RenderElementList* l = new Reve::RenderElementList("MUONChambers");
79 l->SetTitle("MUON chambers");
80 l->SetMainColor(Color_t(2));
81 gReve->AddRenderElement(l);
82
83 for (Int_t ic = 0; ic < 14; ic++) {
84
85 Alieve::MUONChamber* mucha = new Alieve::MUONChamber();
86
87 mucha->SetFrameColor(2);
88 mucha->SetChamberID(ic);
89 mucha->SetDataSource(g_muon_data);
90
91 gReve->AddRenderElement(l,mucha);
92
93 }
94
95 if (showTracks) MUON_tracks();
96
97 gReve->EnableRedraw();
98 gReve->Redraw3D(kTRUE);
99 /*
100 TGLViewer* view = dynamic_cast<TGLViewer*>(gReve->GetGLCanvas()->GetViewer3D());
101 view->ResetCamerasAfterNextUpdate();
102 gReve->GetGLCanvas()->Modified();
103 gReve->GetGLCanvas()->Update();
104 */
105}
106
107//_____________________________________________________________________________
108void MUON_tracks() {
109
110 AliRunLoader* rl = Alieve::Event::AssertRunLoader();
111 rl->LoadTracks("MUON");
112 TTree* tt = rl->GetTreeT("MUON", false);
113
114 AliMagF* mag_field = Alieve::Event::AssertMagField();
115
116 TClonesArray *tracks = 0;
117 tt->SetBranchAddress("MUONTrack",&tracks);
118 tt->GetEntry(0);
119
120 Int_t ntracks = tracks->GetEntriesFast();
121 //printf("Found %d tracks. \n",ntracks);
122
123 Reve::TrackList* lt = new Reve::TrackList("M-Tracks");
124 lt->SetMainColor(Color_t(6));
125
126 gReve->AddRenderElement(lt);
127
128 TMatrixD smatrix(2,2);
129 TMatrixD sums(2,1);
130 TMatrixD res(2,1);
131
132 Float_t xRec, xRec0;
133 Float_t yRec, yRec0;
134 Float_t zRec, zRec0;
135
136 Float_t zg[4] = { -1603.5, -1620.5, -1703.5, -1720.5 };
137
138 AliMUONTrack *mt;
139 Reve::RecTrack rt;
140 Int_t count;
141 for (Int_t n = 0; n < ntracks; n++) {
142
143 count = 0;
144
145 mt = (AliMUONTrack*) tracks->At(n);
146
147 //printf("Match trigger %d \n",mt->GetMatchTrigger());
148
149 rt.label = n;
150
151 Alieve::MUONTrack* track = new Alieve::MUONTrack(&rt, lt->GetRnrStyle());
152
153 track->MakeTrack(mt, mag_field);
154
155 gReve->AddRenderElement(lt, track);
156
157 continue;
158
159 /* old tracks display */
160
161 Reve::Track* rtrack = new Reve::Track(&rt, lt->GetRnrStyle());
162
163 //PH The line below is replaced waiting for a fix in Root
164 //PH which permits to use variable siza arguments in CINT
165 //PH on some platforms (alphalinuxgcc, solariscc5, etc.)
166 char form[1000];
167 if (mt->GetMatchTrigger()) {
168 //PH track->SetName(Form("MUONTrack %2d (MT)", rt.label));
169 sprintf(form,"MUONTrack %2d (MT)", rt.label);
170 rtrack->SetName(form);
171 rtrack->SetLineStyle(1);
172 rtrack->SetLineColor(4);
173 } else {
174 //PH track->SetName(Form("MUONTrack %2d ", rt.label));
175 sprintf(form,"MUONTrack %2d ", rt.label);
176 rtrack->SetName(form);
177 rtrack->SetLineStyle(2);
178 rtrack->SetLineColor(4);
179 }
180
181 AliMUONTrackParam *trackParam = mt->GetTrackParamAtVertex();
182 xRec0 = trackParam->GetNonBendingCoor();
183 yRec0 = trackParam->GetBendingCoor();
184 zRec0 = trackParam->GetZ();
185
186 rtrack->SetPoint(count,xRec0,yRec0,zRec0);
187 count++;
188
189 Float_t xr[20], yr[20], zr[20];
190 for (Int_t i = 0; i < 10; i++) xr[i]=yr[i]=zr[i]=0.0;
191
192 Int_t nTrackHits = mt->GetNTrackHits();
193 //printf("Nhits = %d \n",nTrackHits);
194 TClonesArray* trackParamAtHit;
195 for (Int_t iHit = 0; iHit < nTrackHits; iHit++){
196 trackParamAtHit = mt->GetTrackParamAtHit();
197 trackParam = (AliMUONTrackParam*) trackParamAtHit->At(iHit);
198 xRec = trackParam->GetNonBendingCoor();
199 yRec = trackParam->GetBendingCoor();
200 zRec = trackParam->GetZ();
201
202 //printf("Hit %d x %f y %f z %f \n",iHit,xRec,yRec,zRec);
203
204 xr[iHit] = xRec;
205 yr[iHit] = yRec;
206 zr[iHit] = zRec;
207
208 rtrack->SetPoint(count,xRec,yRec,zRec);
209 count++;
210
211 }
212
213 Float_t xrc[20], yrc[20], zrc[20];
214 Int_t nrc = 0;
215 if (mt->GetMatchTrigger() && 1) {
216
217 for (Int_t i = 0; i < nTrackHits; i++) {
218 if (TMath::Abs(zr[i]) > 1000.0) {
219 //printf("Hit %d x %f y %f z %f \n",iHit,xr[i],yr[i],zr[i]);
220 xrc[nrc] = xr[i];
221 yrc[nrc] = yr[i];
222 zrc[nrc] = zr[i];
223 nrc++;
224 }
225 }
226
227 Double_t xv, yv;
228 Float_t ax, bx, ay, by;
229
230 if (nrc < 2) continue;
231
232 // fit x-z
233 smatrix.Zero();
234 sums.Zero();
235 for (Int_t i = 0; i < nrc; i++) {
236 xv = (Double_t)zrc[i];
237 yv = (Double_t)xrc[i];
238 //printf("x-z: xv %f yv %f \n",xv,yv);
239 smatrix(0,0) += 1.0;
240 smatrix(1,1) += xv*xv;
241 smatrix(0,1) += xv;
242 smatrix(1,0) += xv;
243 sums(0,0) += yv;
244 sums(1,0) += xv*yv;
245 }
246 res = smatrix.Invert() * sums;
247 ax = res(0,0);
248 bx = res(1,0);
249
250 // fit y-z
251 smatrix.Zero();
252 sums.Zero();
253 for (Int_t i = 0; i < nrc; i++) {
254 xv = (Double_t)zrc[i];
255 yv = (Double_t)yrc[i];
256 //printf("y-z: xv %f yv %f \n",xv,yv);
257 smatrix(0,0) += 1.0;
258 smatrix(1,1) += xv*xv;
259 smatrix(0,1) += xv;
260 smatrix(1,0) += xv;
261 sums(0,0) += yv;
262 sums(1,0) += xv*yv;
263 }
264 res = smatrix.Invert() * sums;
265 ay = res(0,0);
266 by = res(1,0);
267
268 Float_t xtc, ytc, ztc;
269 for (Int_t ii = 0; ii < 4; ii++) {
270
271 ztc = zg[ii];
272 ytc = ay+by*zg[ii];
273 xtc = ax+bx*zg[ii];
274
275 //printf("tc: x %f y %f z %f \n",xtc,ytc,ztc);
276
277 rtrack->SetPoint(count,xtc,ytc,ztc);
278 count++;
279
280 }
281
282 } // end match trigger
283
284 gReve->AddRenderElement(lt, rtrack);
285
286 } // end track loop
287
288}
289