]>
Commit | Line | Data |
---|---|---|
1 | #if !defined(__CINT__) || defined(__MAKECINT__) | |
2 | #include<iostream.h> | |
3 | #include<TROOT.h> | |
4 | #include<TArrayI.h> | |
5 | #include<TBranch.h> | |
6 | #include<TCanvas.h> | |
7 | #include<TClassTable.h> | |
8 | #include<TClonesArray.h> | |
9 | #include<TFile.h> | |
10 | #include<TH1.h> | |
11 | #include<TH2.h> | |
12 | #include<TLegend.h> | |
13 | #include<TObject.h> | |
14 | #include<TObjArray.h> | |
15 | #include<TTree.h> | |
16 | #include <AliRun.h> | |
17 | #include <AliITS.h> | |
18 | #include <AliITSgeom.h> | |
19 | #include <AliITSDetType.h> | |
20 | #include <AliITSRecPoint.h> | |
21 | #include <AliITSdigit.h> | |
22 | #include <AliITShit.h> | |
23 | #include <AliITSmodule.h> | |
24 | #include <AliITSsegmentation.h> | |
25 | #include <AliITSsegmentationSPD.h> | |
26 | #include <AliITSsegmentationSDD.h> | |
27 | #include <AliITSsegmentationSSD.h> | |
28 | #endif | |
29 | ||
30 | Int_t GetModuleHits (TObject* its, Int_t ID, Float_t*& X, Float_t*& Y, Float_t*& Z, Bool_t*& St); | |
31 | Int_t GetModuleRecPoints (TObject *its, Int_t ID, Float_t*& X, Float_t*& Z, Int_t isfastpoints); | |
32 | Int_t GetModuleDigits(TObject *its, Int_t ID, Int_t dtype, Float_t*& X, Float_t*& Z); | |
33 | ||
34 | Int_t AliITSReadPlotData( Int_t evNum = 0, char *filename = "galice.root", TString FileDigits="galice.root", TString FileRec="galice.root", Int_t isfastpoints = 0) { | |
35 | ||
36 | /********************************************************************* | |
37 | * * | |
38 | * Macro used to read hits, digits and recpoints for a module * | |
39 | * It draws 3 TH2Fs where stores the 2-D coordinates of these * | |
40 | * data for a specified module (for which the macro asks when * | |
41 | * it starts, with some checks to avoid wrong detector coords. * | |
42 | * * | |
43 | * Only a little 'experimental' warning: * | |
44 | * with all the tests I have made, I wasn't able to plot the * | |
45 | * digits fo the 5th layer... * | |
46 | * I skipped this problem with an alert to beware th user, while * | |
47 | * in this situation the macro plots only recpoints and hits * | |
48 | * * | |
49 | * Author: Alberto Pulvirenti * | |
50 | * * | |
51 | *********************************************************************/ | |
52 | ||
53 | ||
54 | // First of all, here are put some variable declarations | |
55 | // that are useful in the following part: | |
56 | Int_t nparticles; // number of particles | |
57 | // ITS module coordinates [layer = 1, ladder = 2, det = 3] and absolute ID[0] of module [0] | |
58 | TArrayI ID(4); | |
59 | Int_t nmodules, dtype; // Total number of modules and module type (SSD, SPD, SDD) | |
60 | Float_t *x = 0, *y = 0, *z = 0; // Arrays where to store read coords | |
61 | Bool_t *St = 0; // Status of the track (hits only) | |
62 | /* | |
63 | // It's necessary to load the gAlice shared libs | |
64 | // if they aren't already stored in memory... | |
65 | if (gClassTable->GetID("AliRun") < 0) { | |
66 | gROOT->LoadMacro("loadlibs.C"); | |
67 | loadlibs(); | |
68 | } | |
69 | // Anyway, this macro needs to read a gAlice file, so it | |
70 | // clears the gAlice object if there is already one in memory... | |
71 | else { | |
72 | */ | |
73 | if(gAlice){ | |
74 | delete gAlice; | |
75 | gAlice = 0; | |
76 | } | |
77 | // } | |
78 | ||
79 | // Now is opened the Root input file containing Geometry, Kine and Hits | |
80 | // by default its name must be "galice.root". | |
81 | // When the file is opened, its contens are shown. | |
82 | TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(filename); | |
83 | if (!file) file = new TFile(filename); | |
84 | file->ls(); | |
85 | ||
86 | // Then, the macro gets the AliRun object from file. | |
87 | // If this object is not present, an error occurs | |
88 | // and the execution is stopped. | |
89 | // Anyway, this operation needs some time, | |
90 | // don't worry about an apparent absence of output and actions... | |
91 | cout << "\nSearching in '" << filename << "' for an AliRun object ... " << flush; | |
92 | gAlice = (AliRun*)file->Get("gAlice"); | |
93 | if (gAlice) | |
94 | cout << "FOUND!" << endl; | |
95 | else { | |
96 | cout<<"NOT FOUND! The Macro can't continue!" << endl; | |
97 | return 0; | |
98 | } | |
99 | if(!(FileDigits.Data() == filename)){ | |
100 | gAlice->SetTreeDFileName(FileDigits); | |
101 | } | |
102 | if(!(FileRec.Data() == filename)){ | |
103 | gAlice->SetTreeRFileName(FileRec); | |
104 | } | |
105 | ||
106 | // Then, the macro selects the event number specified. Default is 0. | |
107 | nparticles = gAlice->GetEvent(evNum); | |
108 | cout << "\nNumber of particles = " << nparticles << endl; | |
109 | if (!nparticles) { | |
110 | cout << "With no particles I can't do much... Goodbye!" << endl; | |
111 | return 0; | |
112 | } | |
113 | ||
114 | // The next step is filling the ITS from the AliRun object. | |
115 | AliITS *ITS = (AliITS*)gAlice->GetModule("ITS"); | |
116 | ITS->InitModules(-1, nmodules); | |
117 | cout << "Number of ITS modules = " << nmodules << endl; | |
118 | cout << "\nFilling modules (it takes a while, now)..." << flush; | |
119 | ITS->FillModules(0, 0, nmodules, " ", " "); | |
120 | cout << "DONE!" << endl; | |
121 | AliITSgeom *gm = ITS->GetITSgeom(); | |
122 | // AliITSDetType *det = ITS->DetType(dtype); | |
123 | // AliITSsegmentation *seg = det->GetSegmentationModel(); | |
124 | ||
125 | for(;;) { | |
126 | ||
127 | // Input phase. | |
128 | // The macro asks if the user wants to put a detector ID[0] | |
129 | // or prefers to input layer, ladder and detector. | |
130 | for (Int_t i = 0; i < 4; i++) ID[i] = -1; | |
131 | Int_t answ; | |
132 | do { | |
133 | cout << "\nSelection modes:" << endl; | |
134 | cout << "1) by layer - ladder - detector numbers" << endl; | |
135 | cout << "2) by unique detector ID" << endl; | |
136 | cout << "0) exit macro" << endl; | |
137 | cout << "\nEnter your choice: "; | |
138 | cin >> answ; | |
139 | } while (answ < 0 || answ > 2); | |
140 | switch (answ) { | |
141 | case 0: | |
142 | // input = 0 ---> EXIT | |
143 | return 0; | |
144 | break; | |
145 | case 1: | |
146 | // input = 1 ---> SELECTION BY COORD | |
147 | do { | |
148 | cout << "\nLayer number [1-6, 0 to exit]: "; | |
149 | cin >> ID[1]; | |
150 | if (!ID[1]) return 0; | |
151 | } while (ID[1] < 0 || ID[1] > 6); | |
152 | ||
153 | // Detector type: 0 --> SPD, 1 --> SDD, 2 --> SSD. | |
154 | // Layer 1,2 --> 0 / Layer 3,4 --> 1 / Layer 5,6 --> 2 | |
155 | dtype = (ID[1] - 1) / 2; | |
156 | ||
157 | // Once fixed the layer number, the macro calculates the max number | |
158 | // for ladder and detector from geometry, and accepts only suitable values. | |
159 | do { | |
160 | ID[2] = gm->GetNladders(ID[1]); | |
161 | cout << "Ladder number [1-" << ID[2] << ", 0 to exit]: "; | |
162 | cin >> ID[2]; | |
163 | if (!ID[2]) return 0; | |
164 | } while (ID[2] < 0 || ID[2] > gm->GetNladders(ID[1])); | |
165 | do { | |
166 | ID[3] = gm->GetNdetectors(ID[1]); | |
167 | cout << "Detector number [1-" << ID[3] << ", 0 to exit]: "; | |
168 | cin >> ID[3]; | |
169 | if (!ID[3]) return 0; | |
170 | } while (ID[3] < 0 || ID[3] > gm->GetNdetectors(ID[1])); | |
171 | break; | |
172 | case 2: | |
173 | // input = 2 ---> SELECTION BY ID[0] | |
174 | do { | |
175 | ID[0] = gm->GetIndexMax(); | |
176 | cout << "\n Detector ID number [0-" << ID[0]-1 << ", -1 to exit]: "; | |
177 | cin >> ID[0]; | |
178 | ID[0]++; | |
179 | if (ID[0] == -1) return 0; | |
180 | } while (ID[0] < 0 || ID[0] > gm->GetIndexMax()); | |
181 | break; | |
182 | }; | |
183 | ||
184 | if (ID[0] == -1) | |
185 | // If ID[0] = -1 the chioce was by coords, so it's necessary to assign the ID: | |
186 | ID[0] = gm->GetModuleIndex(ID[1], ID[2], ID[3]); | |
187 | else { | |
188 | // Else we must get the layer, ladder and detector by the ID; | |
189 | // ...but first we must remember that the ID goest from 0 to NModules - 1! | |
190 | ID[0]--; | |
191 | ID[1] = ITS->GetModule(ID[0])->GetLayer(); | |
192 | ID[2] = ITS->GetModule(ID[0])->GetLadder(); | |
193 | ID[3] = ITS->GetModule(ID[0])->GetDet(); | |
194 | } | |
195 | ||
196 | // Defines the histograms inside the `for' cycle, so they are destroyed at the end | |
197 | // of every read sequqnce, in order to mek another withour segmentation faults | |
198 | Text_t msg[250]; | |
199 | Float_t xm = 0.0, zm = 0.0; | |
200 | dtype = (ID[1] - 1) / 2; | |
201 | switch (dtype) { | |
202 | case 0: xm = 1.5; zm = 7.0; break; | |
203 | case 1: xm = 7.5; zm = 8.0; break; | |
204 | case 2: xm = 7.5; zm = 4.5; break; | |
205 | } | |
206 | sprintf(msg, "Module index=%d lay=%d lad=%d det=%d", ID[0], ID[1], ID[2], ID[3]); | |
207 | TH2F *hhits = new TH2F("hhits", msg, 500, -xm, xm, 500, -zm, zm); // Histogram of hits | |
208 | TH2F *hrecs = new TH2F("hrecs", msg, 500, -xm, xm, 500, -zm, zm); // Histogram of recpoints | |
209 | TH2F *hdigits = new TH2F("hdigits", msg, 500, -xm, xm, 500, -zm, zm); // Histogram of digits | |
210 | ||
211 | cout << endl; | |
212 | ||
213 | // Reads hits... | |
214 | Int_t tmpint = ID[0]; | |
215 | Int_t hits = GetModuleHits(ITS, tmpint, x, y, z, St); | |
216 | if (!hits) { | |
217 | cout << "No hits in module!" << endl; | |
218 | continue; | |
219 | } | |
220 | else | |
221 | cout << "Hits scanned..." << endl; | |
222 | for (Int_t i = 0; i < hits; i++) if (!St[i]) hhits->Fill(x[i], z[i]); | |
223 | ||
224 | // Reads recpoints... | |
225 | Int_t recs = GetModuleRecPoints(ITS, ID[0], x, z, isfastpoints); | |
226 | if (!recs) { | |
227 | cout << "No recpoints in module!" << endl; | |
228 | continue; | |
229 | } | |
230 | else | |
231 | cout << "Recpoints scanned..." << endl; | |
232 | for (Int_t i = 0; i < recs; i++) hrecs->Fill(x[i], z[i]); | |
233 | ||
234 | // Reads digits... | |
235 | Int_t digits = GetModuleDigits(ITS, ID[0], dtype, x, z); | |
236 | if (!digits) { | |
237 | cout << "No digits in module!" << endl; | |
238 | //continue; | |
239 | } | |
240 | else | |
241 | cout << "Digits scanned..." << endl; | |
242 | for (Int_t i = 0; i < digits; i++) hdigits->Fill(x[i], z[i]); | |
243 | ||
244 | // Draws read data... | |
245 | // HITS -------> red (2) crosses. | |
246 | // DIGITS -----> green (8) boxes. | |
247 | // REC-POINTS -> blue (4) St. Andrew's crosses. | |
248 | ||
249 | TCanvas *viewer = new TCanvas("viewer", "Module viewer canvas", 0, 0, 800, 800); | |
250 | viewer->cd(); | |
251 | ||
252 | hdigits->SetMarkerStyle(25); | |
253 | hdigits->SetMarkerColor(8); | |
254 | hdigits->SetMarkerSize(2); | |
255 | hdigits->SetStats(kFALSE); | |
256 | hdigits->SetXTitle("Local X (cm)"); | |
257 | hdigits->SetYTitle("Local Z (cm)"); | |
258 | hdigits->Draw(); | |
259 | ||
260 | hhits->SetMarkerStyle(5); | |
261 | hhits->SetMarkerColor(2); | |
262 | hhits->SetMarkerSize(3); | |
263 | hhits->SetStats(kFALSE); | |
264 | hhits->Draw("same"); | |
265 | ||
266 | hrecs->SetMarkerStyle(2); | |
267 | hrecs->SetMarkerColor(4); | |
268 | hrecs->SetMarkerSize(3); | |
269 | hrecs->SetStats(kFALSE); | |
270 | hrecs->Draw("same"); | |
271 | ||
272 | TLegend *legend = new TLegend(0.7, 0.8, 0.9, 0.9); | |
273 | legend->SetMargin(0.2); | |
274 | legend->AddEntry(hhits, "hits","P"); | |
275 | legend->AddEntry(hrecs, "recpoints","P"); | |
276 | legend->AddEntry(hdigits, "digits","P"); | |
277 | legend->SetTextSizePixels(14); | |
278 | legend->Draw(); | |
279 | ||
280 | viewer->Update(); | |
281 | ||
282 | ||
283 | Text_t fname[250],ans; | |
284 | cout << "Do you want to save the current canvas on a file (y/n) ? "; | |
285 | cin >> ans; | |
286 | if(ans == 'y' || ans == 'Y') { | |
287 | cout << "Enter filename: "; | |
288 | cin >> fname; | |
289 | TString *control = new TString(fname); | |
290 | Bool_t ok=control->Contains(".C") || control->Contains(".root") || control->Contains(".ps") || control->Contains(".eps") || control->Contains(".gif"); | |
291 | if(!ok){ | |
292 | cout << "File extension is not recognized. The canvas will be saved as Postscript file"; | |
293 | strcat(fname,".ps"); | |
294 | } | |
295 | viewer->SaveAs(fname); | |
296 | } | |
297 | } | |
298 | ||
299 | cout << "Done. Goodbye" << endl; | |
300 | return 0; | |
301 | } | |
302 | ||
303 | Int_t GetModuleHits(TObject* its, Int_t ID, Float_t*& X, Float_t*& Y, Float_t*& Z, Bool_t*& St) { | |
304 | // First of all, the macro selects the specified module, | |
305 | // then gets the array of hits in it and their number. | |
306 | AliITS *ITS = (AliITS*) its; | |
307 | AliITSmodule *module = ITS->GetModule(ID); | |
308 | TObjArray *hits_array = module->GetHits(); | |
309 | Int_t hits_num = hits_array->GetEntriesFast(); | |
310 | ||
311 | // Now, if this count returns 0, there's nothing to do, | |
312 | // while, if it doesn't, the first thing to do is dimensioning | |
313 | // the coordinate arrays, and then the loop can start. | |
314 | if (!hits_num) | |
315 | return 0; | |
316 | else { | |
317 | if (X) delete [] X; | |
318 | if (Y) delete [] Y; | |
319 | if (Z) delete [] Z; | |
320 | if (St) delete [] St; | |
321 | X = new Float_t[hits_num]; | |
322 | Y = new Float_t[hits_num]; | |
323 | Z = new Float_t[hits_num]; | |
324 | St = new Bool_t[hits_num]; | |
325 | } | |
326 | ||
327 | for (Int_t j = 0; j < hits_num; j++) { | |
328 | AliITShit *hit = (AliITShit*) hits_array->At(j); | |
329 | X[j] = hit->GetXL(); | |
330 | Y[j] = hit->GetYL(); | |
331 | Z[j] = hit->GetZL(); | |
332 | St[j] = hit->StatusEntering(); | |
333 | } | |
334 | return hits_num; | |
335 | } | |
336 | ||
337 | Int_t GetModuleRecPoints (TObject *its, Int_t ID, Float_t*& X, Float_t*& Z, Int_t isfastpoints) { | |
338 | ||
339 | // First of all, the macro selects the specified module, | |
340 | // then gets the array of recpoints in it and their number. | |
341 | AliITS *ITS = (AliITS*) its; | |
342 | TTree *TR = gAlice->TreeR(); | |
343 | ITS->ResetRecPoints(); | |
344 | ||
345 | TClonesArray *recs_array = ITS->RecPoints(); | |
346 | ||
347 | TBranch *branch = 0; | |
348 | if(TR && recs_array){ | |
349 | if(isfastpoints==1){ | |
350 | branch = gAlice->TreeR()->GetBranch("ITSRecPointsF"); | |
351 | cout<<"using fast points\n"; | |
352 | } | |
353 | else { | |
354 | branch = gAlice->TreeR()->GetBranch("ITSRecPoints"); | |
355 | } | |
356 | if(branch)branch->SetAddress(&recs_array); | |
357 | } | |
358 | ||
359 | branch->GetEvent(ID); | |
360 | ||
361 | Int_t recs_num = recs_array->GetEntries(); | |
362 | ||
363 | ||
364 | // Now, if this count returns 0, there's nothing to do, | |
365 | // while, if it doesn't, the first thing to do is dimensioning | |
366 | // the coordinate and energy loss arrays, and then the loop can start. | |
367 | if (!recs_num) | |
368 | return 0; | |
369 | else { | |
370 | if (X) delete [] X; | |
371 | if (Z) delete [] Z; | |
372 | X = new Float_t[recs_num]; | |
373 | Z = new Float_t[recs_num]; | |
374 | } | |
375 | for(Int_t j = 0; j < recs_num; j++) { | |
376 | AliITSRecPoint *recp = (AliITSRecPoint*)recs_array->At(j); | |
377 | X[j] = recp->GetX(); | |
378 | Z[j] = recp->GetZ(); | |
379 | } | |
380 | return recs_num; | |
381 | } | |
382 | ||
383 | Int_t GetModuleDigits(TObject *its, Int_t ID, Int_t dtype, Float_t*& X, Float_t*& Z) { | |
384 | ||
385 | // First of all, the macro selects the specified module, | |
386 | // then gets the array of recpoints in it and their number. | |
387 | AliITSdigit *digit; | |
388 | AliITS *ITS = (AliITS*)its; | |
389 | TTree *TD = gAlice->TreeD(); | |
390 | ITS->ResetDigits(); | |
391 | TD->GetEvent(ID); | |
392 | TClonesArray *digits_array = ITS->DigitsAddress(dtype); | |
393 | AliITSgeom *gm = ITS->GetITSgeom(); | |
394 | AliITSDetType *det = ITS->DetType(dtype); | |
395 | AliITSsegmentation *seg = det->GetSegmentationModel(); | |
396 | TArrayI ssdone(5000); // used to match p and n side digits of strips | |
397 | TArrayI pair(5000); // as above | |
398 | Int_t digits_num = digits_array->GetEntries(); | |
399 | // Now, if this count returns 0, there's nothing to do, | |
400 | // while, if it doesn't, the first thing to do is dimensioning | |
401 | // the coordinate and energy loss arrays, and then the loop can start. | |
402 | ||
403 | if(dtype==2){ | |
404 | Int_t layer, ladder, detec; | |
405 | gm->GetModuleId(ID,layer,ladder,detec); | |
406 | seg->SetLayer(layer); | |
407 | } | |
408 | ||
409 | if (!digits_num) | |
410 | return 0; | |
411 | else { | |
412 | cout << "Digits to scan: " << digits_num << endl; | |
413 | if (X) delete [] X; | |
414 | if (Z) delete [] Z; | |
415 | X = new Float_t[digits_num]; | |
416 | Z = new Float_t[digits_num]; | |
417 | } | |
418 | ||
419 | // Get the coordinates of the module | |
420 | if (dtype == 2) { | |
421 | for (Int_t j = 0; j < digits_num; j++) { | |
422 | ssdone[j] = 0; | |
423 | pair[j] = 0; | |
424 | } | |
425 | } | |
426 | // AliITSdigit *digit; | |
427 | for (Int_t j = 0; j < digits_num; j++) { | |
428 | digit = (AliITSdigit*)digits_array->UncheckedAt(j); | |
429 | Int_t iz=digit->fCoord1; // cell number z | |
430 | Int_t ix=digit->fCoord2; // cell number x | |
431 | // Get local coordinates of the element (microns) | |
432 | // ******************************* PARTE CORRETTA *************************************** | |
433 | if(dtype < 2) { | |
434 | Float_t xx, zz; // aggiunta | |
435 | seg->DetToLocal(ix, iz, xx, zz); | |
436 | X[j] = xx; // aggiunta | |
437 | Z[j] = zz; // aggiunta | |
438 | } | |
439 | // ******************************* FINE PARTE CORRETTA *************************************** | |
440 | else { | |
441 | // SSD: if iz==0 ---> N side; if iz==1 P side | |
442 | if (ssdone[j] == 0) { | |
443 | ssdone[j]=1; | |
444 | pair[j]=-1; | |
445 | Bool_t pside = (iz == 1); | |
446 | Bool_t impaired = kTRUE; | |
447 | Int_t pstrip = 0; | |
448 | Int_t nstrip = 0; | |
449 | if (pside) pstrip = ix; else nstrip = ix; | |
450 | for (Int_t k = 0; k < digits_num; k++) { | |
451 | if (ssdone[k] == 0 && impaired) { | |
452 | AliITSdigitSSD *sdigit=(AliITSdigitSSD*)digits_array->UncheckedAt(k); | |
453 | if (sdigit->fCoord1 != iz && sdigit->GetTracks()[0] == digit->GetTracks()[0]) { | |
454 | ssdone[k]=2; | |
455 | pair[j]=k; | |
456 | if (pside) nstrip = sdigit->fCoord2; else pstrip = sdigit->fCoord2; | |
457 | impaired=kFALSE; | |
458 | } | |
459 | } | |
460 | } | |
461 | if (!impaired) seg->GetPadCxz(pstrip, nstrip, X[j], Z[j]); | |
462 | X[j] /= 10000.0; // convert microns to cm | |
463 | Z[j] /= 10000.0; // convert microns to cm | |
464 | } | |
465 | } | |
466 | } | |
467 | return digits_num; | |
468 | } | |
469 | ||
470 | ||
471 | ||
472 | ||
473 | ||
474 | ||
475 | ||
476 | ||
477 | ||
478 |