]>
Commit | Line | Data |
---|---|---|
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 | 28 | Int_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 | ||
355 | void 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 | ||
388 | Int_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 | ||
415 | void 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 |