Replace GetTracks()[0] by GetTrack(0).
[u/mrichter/AliRoot.git] / ITS / ITSgeoplot.C
CommitLineData
79f97fd2 1#define __NOCOMPILED__
2#ifdef __COMPILED__
3#include<iostream.h>
4#include<TROOT.h>
5#include<TArrayI.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<TObject.h>
13#include<TObjArray.h>
14#include<TTree.h>
15#include <AliRun.h>
16#include <AliITS.h>
17#include <AliITSgeom.h>
18#include <AliITSDetType.h>
19#include <AliITSRecPoint.h>
20#include <AliITSdigit.h>
21#include <AliITShit.h>
22#include <AliITSmodule.h>
23#include <AliITSsegmentation.h>
24#include <AliITSsegmentationSPD.h>
25#include <AliITSsegmentationSDD.h>
26#include <AliITSsegmentationSSD.h>
27#endif
2f101c6a 28Int_t ITSgeoplot (char *opt="All+Rec", char *filename="galice.root") {
29 /*******************************************************************
30 * This macro displays geometrical information related to the
31 * hits, digits and rec points in ITS.
32 * There are histograms that are not displayed (i.e. energy
33 * deposition) but that are saved onto a file (see below)
34 *
35 * Options: Any combination of:
36 * 1) subdetector name: "SPD", "SDD", "SSD", "All" (default)
37 * 2) Printouts: "Verbose" Almost everything is printed out:
38 * it is wise to redirect the output onto a file
39 * e.g.: .x ITSgeoplot.C("All+Verbose") > out.log
40 * 3) Rec Points option: "Rec" ---> Uses Rec. Points (default)
41 * otherwise ---> uses hits and digits only
42 * Examples:
43 * .x ITSgeoplot(); (All subdetectors; no-verbose; no-recpoints)
44 * .x ITSgeoplot("SPD+SSD+Verbose+Rec");
45 *
46 * OUTPUT: It produces a root file with a list of histograms
47 * The list can be accessed either with the Root TBrowser
48 * or with the simple macro ITSgeoplotread.C
49 * WARNING: spatial information for SSD/DIGITS is obtained by pairing
50 * digits on p and n side originating from the same track, when
51 * possible. This (mis)use of DIGITS is tolerated for debugging
52 * purposes only !!!! The pairing in real life should be done
79f97fd2 53 * starting from info really available...
54 *
55 * COMPILATION: this macro can be compiled.
56 * 1) Change the first line above from its default
57 * value (#define __NOCOMPILED__) to
58 * #define __COMPILED__
59 * 2) You need to set your include path with
60 * gSystem->SetIncludePath("-I- -I$ALICE_ROOT/include -I$ALICE_ROOT/ITS -g");
61 * 3) If you are using root instead of aliroot you need to
62 * execute the macro loadlibs.C in advance
63 * 4) To compile this macro from root (or aliroot):
64 * --- .L ITSgeoplot.C++
65 * --- ITSgeoplot("ListOfParametersIfAny");
2f101c6a 66 *
faafc0eb 67 * M.Masera 14/05/2001 18:30
68 * Use DetToLocal instead of GetCxz for SPD and SDD
2f101c6a 69 ********************************************************************/
70
71 extern void GetHitsCoor(TObject *its, Int_t mod, TObjArray & histos, Int_t subd,Bool_t verb);
72 extern Int_t GetRecCoor(TObject *ge, TClonesArray *ITSrec, Int_t mod, TH2F *h2, TH1F *h1, Bool_t verb);
73 extern void GetDigits(TObject *tmps,TObject *ge,TClonesArray *ITSdigits, Int_t subd, Int_t mod, Bool_t verbose, TObjArray & histos);
2f101c6a 74 //Options
75 TString choice(opt);
76 Bool_t All = choice.Contains("All");
77 Bool_t verbose=choice.Contains("Verbose");
78 Bool_t userec=choice.Contains("Rec");
79 Int_t retcode=1; //return code
79f97fd2 80#ifdef __NOCOMPILED__
2f101c6a 81 if (gClassTable->GetID("AliRun") < 0) {
82 gROOT->LoadMacro("loadlibs.C");
83 loadlibs();
84 }
85 else {
79f97fd2 86#endif
2f101c6a 87 if(gAlice){
88 delete gAlice;
89 gAlice=0;
90 }
79f97fd2 91#ifdef __NOCOMPILED__
2f101c6a 92 }
79f97fd2 93#endif
2f101c6a 94 // Connect the Root input file containing Geometry, Kine and Hits
95 // galice.root file by default
96
97 TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(filename);
98 if (!file) file = new TFile(filename);
99 file->ls();
100
101 // Get AliRun object from file
102
103 if (!gAlice) {
104 gAlice = (AliRun*)file->Get("gAlice");
105 if (gAlice && verbose)cout<<"AliRun object found on file "<<filename<<endl;
106 if(!gAlice){
107 cout<<"Can't access AliRun object on file "<<filename<<endl;
108 cout<<"Macro execution stopped!!!"<<endl;
109 retcode=-1;
110 return retcode;
111 }
112 }
113 Int_t nparticles = gAlice->GetEvent(0);
114 if(verbose) {
115 cout<<" "<<endl<<" "<<endl;
116 cout<<"******* Event processing started *******"<<endl;
117 cout<<"In the following, hits with 'StatusEntering' flag active"<<endl;
118 cout<<"will not be processed"<<endl;
119 cout << "Number of particles= " << nparticles <<endl;
120 }
121
122 // HITS
123 TTree *TH = gAlice->TreeH();
124 Stat_t ntracks = TH->GetEntries();
125 if(verbose)cout<<"Number of primary tracks= "<<ntracks<<endl;
126
127 // ITS
128 AliITS *ITS = (AliITS*)gAlice->GetModule("ITS");
129 Int_t nmodules;
130 ITS->InitModules(-1,nmodules);
131 cout<<"Number of ITS modules= "<<nmodules<<endl;
132 cout<<"Filling modules... It takes a while, now. Please be patient"<<endl;
133 ITS->FillModules(0,0,nmodules," "," ");
134 cout<<"ITS modules .... DONE!"<<endl;
135
136 // DIGITS
137 TTree *TD = gAlice->TreeD();
138
139 //RECPOINTS
140 TTree *TR = gAlice->TreeR();
141 TClonesArray *ITSrec = ITS->RecPoints();
142 if(userec && (!TR || !ITSrec)){
143 userec = kFALSE;
144 cout<<" "<<endl<<"======================================================="<<endl;
145 cout<<"WARNING: there are no RECPOINTS on this file ! "<<endl;
146 cout<<"======================================================="<<endl<<" "<<endl;
147 }
148 //local variables
149 Int_t mod; //module number
150 Int_t nbytes = 0;
151 Double_t pos[3]; // Global position of the current module
152 Float_t ragdet; // Radius of detector (x y plane)
153 Int_t first; // first module
154 Int_t last; // last module
155 Int_t nrecp; //number of RecPoints for a given module
156
157 //List of histograms
158 TObjArray histos(26,0); // contains the pointers to the histograms
159 // Book histograms SPD
160 TH2F *rspd = new TH2F("rspd","Radii of digits - SPD",50,-10.,10.,50,-10.,10.);
161 TH2F *rhspd = new TH2F("rhspd","Radii of hits - SPD",50,-10.,10.,50,-10.,10.);
162 TH2F *rmspd = new TH2F("rmspd","Radii of SPD modules",50,-10.,10.,50,-10.,10.);
163 TH1F *zspd = new TH1F("zspd","Z of digits - SPD",100,-30.,30.);
164 TH1F *zhspd = new TH1F("zhspd","Z of hits - SPD",100,-30.,30.);
165 TH1F *zmspd = new TH1F("zmspd","Z of SPD modules",100,-30,30.);
166 TH2F *rrspd = new TH2F("rrspd","Radii of recpoints - SPD",50,-10.,10.,50,-10.,10.);
167 TH1F *zrspd = new TH1F("zrspd","Z of recpoints - SPD",100,-30.,30.);
168 TH1F *enespd = new TH1F("enespd","Energy deposition SPD (KeV)",100,0.,1000.);
169 histos.AddLast(rspd); // 0
170 histos.AddLast(rhspd); // 1
171 histos.AddLast(rmspd); // 2
172 histos.AddLast(zspd); // 3
173 histos.AddLast(zhspd); // 4
174 histos.AddLast(zmspd); // 5
175 histos.AddLast(rrspd); // 6
176 histos.AddLast(zrspd); // 7
177 histos.AddLast(enespd); // 8
178 // Book histograms SDD
179 TH2F *rsdd = new TH2F("rsdd","Radii of digits - SDD",50,-40.,40.,50,-40.,40.);
180 TH2F *rhsdd = new TH2F("rhsdd","Radii of hits - SDD",50,-40.,40.,50,-40.,40.);
181 TH2F *rmsdd = new TH2F("rmsdd","Radii of SDD modules",50,-40.,40.,50,-40.,40.);
182 TH1F *zsdd = new TH1F("zsdd","Z of digits - SDD",100,-40.,40.);
183 TH1F *zhsdd = new TH1F("zhsdd","Z of hits - SDD",100,-40.,40.);
184 TH1F *zmsdd = new TH1F("zmsdd","Z of SDD modules",100,-40,40.);
185 TH2F *rrsdd = new TH2F("rrsdd","Radii of recpoints - SDD",50,-40.,40.,50,-40.,40.);
186 TH1F *zrsdd = new TH1F("zrsdd","Z of recpoints - SDD",100,-40.,40.);
187 TH1F *enesdd = new TH1F("enesdd","Energy deposition SDD (KeV)",100,0.,1000.);
188 histos.AddLast(rsdd); // 9
189 histos.AddLast(rhsdd); // 10
190 histos.AddLast(rmsdd); // 11
191 histos.AddLast(zsdd); // 12
192 histos.AddLast(zhsdd); // 13
193 histos.AddLast(zmsdd); // 14
194 histos.AddLast(rrsdd); // 15
195 histos.AddLast(zrsdd); // 16
196 histos.AddLast(enesdd); // 17
197 // Book histogram SSD
198 TH2F *rssd = new TH2F("rssd","Radii of digits - SSD",50,-50.,50.,50,-50.,50.);
199 TH2F *rhssd = new TH2F("rhssd","Radii of hits - SSD",50,-50.,50.,50,-50.,50.);
200 TH2F *rmssd = new TH2F("rmssd","Radii of SSD modules",50,-50.,50.,50,-50.,50.);
201 TH1F *zssd = new TH1F("zssd","Z of digits - SSD",100,-70.,70.);
202 TH1F *zhssd = new TH1F("zhssd","Z of hits - SSD",100,-70.,70.);
203 TH1F *zmssd = new TH1F("zmssd","Z of SSD modules",100,-70,70.);
204 TH2F *rrssd = new TH2F("rrssd","Radii of recpoints - SSD",50,-50.,50.,50,-50.,50.);
205 TH1F *zrssd = new TH1F("zrssd","Z of recpoints - SSD",100,-70.,70.);
206 TH1F *enessd = new TH1F("enessd","Energy deposition SSD (KeV)",100,0.,1000.);
207 histos.AddLast(rssd); // 18
208 histos.AddLast(rhssd); // 19
209 histos.AddLast(rmssd); // 20
210 histos.AddLast(zssd); // 21
211 histos.AddLast(zhssd); // 22
212 histos.AddLast(zmssd); // 23
213 histos.AddLast(rrssd); // 24
214 histos.AddLast(zrssd); // 25
215 histos.AddLast(enessd); // 26
216 //
217 // Loop on subdetectors
218 //
219 AliITSgeom *geom = ITS->GetITSgeom();
79f97fd2 220 TString detna; // subdetector name
2f101c6a 221 for(Int_t subd=0;subd<3;subd++){
222 if(All || (choice.Contains("SPD") && subd==0) || (choice.Contains("SDD") && subd==1) || (choice.Contains("SSD") && subd==2)){
223 // Prepare array for the digits
224 TClonesArray *ITSdigits = ITS->DigitsAddress(subd);
225 Bool_t usedigits = kTRUE;
226 if(!ITSdigits){
227 usedigits = kFALSE;
228 cout<<" "<<endl<<"======================================================="<<endl;
229 cout<<"WARNING: there are no DIGITS on this file ! "<<endl;
230 cout<<"======================================================="<<endl<<" "<<endl;
231 }
232 // Get segmentation model
233 AliITSDetType *iDetType=ITS->DetType(subd);
2f101c6a 234 if(subd==0)detna="SPD";
235 if(subd==1)detna="SDD";
236 if(subd==2)detna="SSD";
237 AliITSsegmentation *seg=(AliITSsegmentation*)iDetType->GetSegmentationModel();
238 // Loop on modules
239 first = geom->GetStartDet(subd);
240 last = geom->GetLastDet(subd);
241 if(verbose){
242 cout<<" "<<endl<<"-------------------------------------"<<endl;
243 cout<<"Start processing subdetector "<<detna<<endl;
244 cout<<detna<<" first module "<<first<<endl;
245 cout<<detna<<" last module "<<last<<endl;
246 cout<<" "<<endl<<" "<<endl;
247 }
248 for (mod=first; mod<=last; mod++){
249 geom->GetTrans(mod,pos); // position of the module in the MRS
250 ragdet=sqrt(pos[0]*pos[0]+pos[1]*pos[1]);
251 // The following 2 histos are a check of the geometry
252 TH2F *bidi = (TH2F*)histos.At(2+subd*9);
253 TH1F *uni = (TH1F*)histos.At(5+subd*9);
254 bidi->Fill(pos[0],pos[1]);
255 uni->Fill(pos[2]);
256 if(verbose){
257 cout<<"==========================================================="<<endl;
258 cout<<detna<<" module="<<mod<<endl;
259 cout<<"Mod. coordinates: "<<pos[0]<<", "<<pos[1]<<", "<<pos[2]<<" Radius "<<ragdet<<endl;
260 }
261
262 // Hits
263 GetHitsCoor(ITS,mod,histos,subd,verbose);
264
265 //RecPoints
266 if(userec){
267 ITS->ResetRecPoints();
268 TR->GetEvent(mod);
269 TH2F *bidi=(TH2F*)histos.At(6+subd*9);
270 TH1F *uni=(TH1F*)histos.At(7+subd*9);
271 nrecp=GetRecCoor(geom,ITSrec,mod,bidi,uni,verbose);
272 }
273
274 // Digits
275 if(usedigits){
276 ITS->ResetDigits();
277 nbytes += TD->GetEvent(mod);
278 GetDigits(seg,geom,ITSdigits,subd,mod,verbose,histos);
279 }
280
281 } // End of loop on the modules
282 TH1F *h1tmp;
283 TH2F *h2tmp;
284 // Plot the histograms
285 TCanvas *current=0; // current Canvas (1--> SPD, 2---> SDD, 3---> SSD)
286 if(subd==0){
287 // Prepare canvas 1
288 TCanvas *c1 = new TCanvas("c1","SPD",10,10,600,900);
289 c1->Divide(2,3);
290 current=c1;
291 }
292 if(subd==1){
293 // Prepare canvas 2
294 TCanvas *c2 = new TCanvas("c2","SDD",40,40,600,900);
295 c2->Divide(2,3);
296 current=c2;
297 }
298 if(subd==2){
299 // Prepare canvas 3
300 TCanvas *c3 = new TCanvas("c3","SSD",70,70,600,900);
301 c3->Divide(2,3);
302 current=c3;
303 }
304 current->cd(1);
305 h2tmp = (TH2F*)histos.At(9*subd);
306 h2tmp->Draw();
307 current->cd(2);
308 h1tmp=(TH1F*)histos.At(3+subd*9);
309 h1tmp->Draw();
310 current->cd(3);
311 h2tmp=(TH2F*)histos.At(1+9*subd);
312 h2tmp->Draw();
313 current->cd(4);
314 h1tmp=(TH1F*)histos.At(4+subd*9);
315 h1tmp->Draw();
316
317 if(userec){
318 current->cd(5);
319 h2tmp=(TH2F*)histos.At(6+9*subd);
320 h2tmp->Draw();
321 current->cd(6);
322 h1tmp=(TH1F*)histos.At(7+subd*9);
323 h1tmp->Draw();
324 }
325
326 else {
327 current->cd(5);
328 h2tmp=(TH2F*)histos.At(2+9*subd);
329 h2tmp->Draw();
330 current->cd(6);
331 h2tmp=(TH2F*)histos.At(5+9*subd);
332 h2tmp->Draw();
333 }
334 } // if(All.....
335 } // end of loop on subdetectors
336 // Save the histograms
337 TFile *fh = new TFile("ITSgeoplot.root","recreate");
338 // The list is written to file as a single entry
339 TList *lihi = new TList();
340 // copy the pointers to the histograms to a TList object.
341 // The histograms concerning recpoints are not copied if
342 // 'userec' is false.
343 for(Int_t i=0;i<histos.GetEntriesFast();i++){
344 if(choice.Contains("All") || (choice.Contains("SPD") && i<8) || (choice.Contains("SDD") && i>7 && i<16) || (choice.Contains("SSD") && i>15)){
345 if(!(!userec && ((i+2)%9==0 || (i+1)%9==0)))lihi->Add(histos.At(i));
346 }
347 }
348 lihi->Write("Histograms ITS hits+digits+recpoints",TObject::kSingleKey);
349 fh->Close();
350
351 return retcode;
352}
353
354
355void GetHitsCoor(TObject *its, Int_t mod, TObjArray & histos, Int_t subd,Bool_t verb){
356 TH2F *h2=(TH2F*)histos.At(1+subd*9);
357 TH1F *h1=(TH1F*)histos.At(4+subd*9);
358 TH1F *ener=(TH1F*)histos.At(8+subd*9);
359 AliITS *ITS= (AliITS*)its;
360 AliITSmodule *modu = ITS->GetModule(mod);
361 TObjArray *fHits = modu->GetHits();
362 Int_t nhits = fHits->GetEntriesFast();
363 if(nhits>0){
364 if(verb){
365 cout<<"-------------------------------------------------------"<<endl;
366 cout<<"Number of HITS for module "<<mod<<": "<<nhits<<endl;
367 }
368 for (Int_t hit=0;hit<nhits;hit++) {
369 AliITShit *iHit = (AliITShit*) fHits->At(hit);
370 if(!iHit->StatusEntering()){
371 Float_t x=iHit->GetXG();
372 Float_t y=iHit->GetYG();
373 Float_t z=iHit->GetZG();
374 Float_t edep=iHit->GetIonization()*1000000;
375 h2->Fill(x,y);
376 h1->Fill(z);
377 ener->Fill(edep);
378 if(verb){
379 cout<<"hit # "<<hit<<" Global coordinates "<<x<<" "<<y<<" "<<z<<endl;
380 cout<<"track # "<<iHit->GetTrack()<<" energy deposition (KeV)= "<<edep<<endl;
381 }
382 }
383 }
384 }
385}
386
387
388Int_t GetRecCoor(TObject *ge, TClonesArray *ITSrec, Int_t mod, TH2F *h2, TH1F *h1, Bool_t verb){
389 AliITSgeom *geom = (AliITSgeom*)ge;
390 Int_t nrecp = ITSrec->GetEntries();
391 if(nrecp>0){
392 Float_t lc[3]; for(Int_t i=0; i<3; i++) lc[i]=0.;
393 Float_t gc[3]; for(Int_t i=0; i<3; i++) gc[i]=0.;
394 if(verb){
395 cout<<"-------------------------------------------------------"<<endl;
396 cout<<"Number of REC POINTS for module "<<mod<<": "<<nrecp<<endl;
397 }
398 for(Int_t irec=0;irec<nrecp;irec++) {
399 AliITSRecPoint *recp = (AliITSRecPoint*)ITSrec->At(irec);
400 lc[0]=recp->GetX();
401 lc[2]=recp->GetZ();
402 geom->LtoG(mod,lc,gc);
403 if(verb){
404 cout<<"recp # "<<irec<<" local coordinates. lx= "<<lc[0]<<" lz= "<<lc[2]<<endl;
405 Float_t r=sqrt(gc[0]*gc[0]+gc[1]*gc[1]);
406 cout<<"Global coor. + radius "<<gc[0]<<" "<<gc[1]<<" "<<gc[2]<<" "<<r<<endl;
407 }
408 h2->Fill(gc[0],gc[1]);
409 h1->Fill(gc[2]);
410 }
411 }
412 return nrecp;
413}
414
415void GetDigits(TObject *tmps,TObject *ge,TClonesArray *ITSdigits, Int_t subd, Int_t mod, Bool_t verbose, TObjArray & histos){
416 AliITSsegmentation *seg = (AliITSsegmentation*)tmps;
417 AliITSgeom *geom = (AliITSgeom*)ge;
79f97fd2 418 Int_t layer;
419 Int_t ladder;
420 Int_t detec;
421 if(subd==2){
422 geom->GetModuleId(mod,layer,ladder,detec);
423 seg->SetLayer(layer);
424 }
2f101c6a 425 Float_t lcoor[3]; for(Int_t j=0; j<3; j++) lcoor[j]=0.; //local coord dig.
426 Float_t gcoor[3]; for(Int_t j=0; j<3; j++) gcoor[j]=0.; // global coo. dig.
427 Float_t ragdig; // Radius digit
428 TArrayI ssdone(5000); // used to match p and n side digits of strips
429 TArrayI pair(5000); // as above
430 Int_t ndigits = ITSdigits->GetEntries();
431 AliITSdigit *digs;
432 if(ndigits){
433 if(verbose){
434 cout<<"-------------------------------------------------------"<<endl;
435 cout<<"Number of DIGITS for module "<<mod<<": "<<ndigits<<endl;
436 }
437 // Get the coordinates of the module
438 if(subd==2){
439 for (Int_t digit=0;digit<ndigits;digit++){
440 ssdone[digit]=0;
441 pair[digit]=0;
442 }
443 }
444 for (Int_t digit=0;digit<ndigits;digit++) {
445 digs = (AliITSdigit*)ITSdigits->UncheckedAt(digit);
446 Int_t iz=digs->fCoord1; // cell number z
447 Int_t ix=digs->fCoord2; // cell number x
faafc0eb 448 // Get local coordinates of the element
2f101c6a 449 if(subd<2){
faafc0eb 450 seg->DetToLocal(ix,iz,lcoor[0],lcoor[2]);
2f101c6a 451 }
452 else{
453 // SSD: if iz==0 ---> N side; if iz==1 P side
454 if(ssdone[digit]==0){
455 ssdone[digit]=1;
456 pair[digit]=-1;
457 Bool_t pside=(iz==1);
458 Bool_t impaired=kTRUE;
459 Int_t pstrip=0;
460 Int_t nstrip=0;
461 if(pside)pstrip=ix;
462 if(!pside)nstrip=ix;
463 for(Int_t digi2=0;digi2<ndigits;digi2++){
464 if(ssdone[digi2]==0 && impaired){
465 AliITSdigitSSD *dig2=(AliITSdigitSSD*)ITSdigits->UncheckedAt(digi2);
3a48961e 466 if(dig2->fCoord1 != iz && dig2->GetTrack(0)==digs->GetTrack(0) && dig2->GetTrack(0)>0){
2f101c6a 467 ssdone[digi2]=2;
468 pair[digit]=digi2;
469 if(pside)nstrip=dig2->fCoord2;
470 if(!pside)pstrip=dig2->fCoord2;
471 impaired=kFALSE;
472 }
473 }
474 }
475 if(!impaired)seg->GetPadCxz(pstrip,nstrip,lcoor[0],lcoor[2]);
476 }
477 }
2f101c6a 478 if(subd<2 || (subd==2 && ssdone[digit]==1)){
479 Int_t coor1=digs->fCoord1;
480 Int_t coor2=digs->fCoord2;
3a48961e 481 Int_t tra0=digs->GetTrack(0);
2f101c6a 482 if(verbose){
483 cout<<"digit # "<<digit<<" fCoord1= "<<coor1<<" fCoord2= "<<coor2<<" track "<<tra0<<" "<<endl;
484 if(subd<2)cout<<"local coordinates -- x="<<lcoor[0]<<", z="<<lcoor[2]<<endl;
485 if(subd==2){
486 if(pair[digit]==-1){
487 cout<<"(digit not paired)"<<endl;
488 }
489 else {
490 cout<<"local coordinates -- x="<<lcoor[0]<<", z="<<lcoor[2]<<endl;
491 Int_t dtmp=pair[digit];
492 AliITSdigitSSD *dig2=(AliITSdigitSSD*)ITSdigits->UncheckedAt(dtmp);
493 Int_t coor1b=dig2->fCoord1;
494 Int_t coor2b=dig2->fCoord2;
3a48961e 495 Int_t tra0b=dig2->GetTrack(0);
2f101c6a 496 cout<<"(digit paired with digit #"<<dtmp<<endl;
497 cout<<"with fCoord1= "<<coor1b<<" fCoord2= "<<coor2b<<" track "<<tra0b<<")"<<endl;
498 }
499 }
500 }
501 if(subd<2 || (subd==2 && pair[digit]!=-1)){
502 // Global coordinates of the element
faafc0eb 503 //SDD and SPD use cm, SSD microns (GetPadCxz)
504 if(subd==2)for(Int_t j=0;j<3;j++)lcoor[j]=lcoor[j]/10000.;
2f101c6a 505 lcoor[1]=0.;
506 geom->LtoG(mod,lcoor,gcoor); // global coord. in cm
507 ragdig=sqrt(gcoor[0]*gcoor[0]+gcoor[1]*gcoor[1]);
508 if(verbose)cout<<"global coordinates "<<gcoor[0]<<" "<<gcoor[1]<<" "<<gcoor[2]<<" Radius "<<ragdig<<endl;
509 //Fill histograms
510 TH2F *bidi = (TH2F*)histos.At(subd*9);
511 TH1F *uni = (TH1F*)histos.At(3+subd*9);
512 bidi->Fill(gcoor[0],gcoor[1]);
513 uni->Fill(gcoor[2]);
514 }
515 }
516 } // loop on digits for this module
517 } // if(ndigits>0....
518}
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538