]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ITS/ShowITSHitsRecPoints.C
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Do not get main runloader from gAlice
[u/mrichter/AliRoot.git] / ITS / ShowITSHitsRecPoints.C
CommitLineData
3fd01708 1#if !defined(__CINT__) || defined(__MAKECINT__)
2#include <TCanvas.h>
3#include <TClassTable.h>
4#include <TGraph.h>
5#include <TGraph2D.h>
6#include <TGeoManager.h>
7#include <TH1.h>
8#include <TInterpreter.h>
9#include "AliCDBManager.h"
10#include "AliGeomManager.h"
11#include "AliHeader.h"
12#include "AliITS.h"
13#include "AliITSDetTypeRec.h"
14#include "AliITSgeom.h"
15#include "AliITSRecPoint.h"
16#include "AliRun.h"
17#endif
18
00fe805a 19Int_t ShowITSHitsRecPoints(Bool_t align=kFALSE,
3fd01708 20 TString alignfile="ITSfullv11Misalignment.root")
21{
22 ///////////////////////////////////////////////////////////////////////
23 // Macro to check clusters and hits in the 6 ITS layers //
24 ///////////////////////////////////////////////////////////////////////
25
26 if (gClassTable->GetID("AliRun") < 0) {
27 gInterpreter->ExecuteMacro("loadlibs.C");
28 }
29 else {
30 if(gAlice){
880d6abe 31 delete AliRunLoader::GetRunLoader();
3fd01708 32 delete gAlice;
33 gAlice=0;
34 }
35 }
36 // Set OCDB if needed
37 AliCDBManager* man = AliCDBManager::Instance();
38 if (!man->IsDefaultStorageSet()) {
39 printf("Setting a local default storage and run number 0\n");
40 man->SetDefaultStorage("local://$ALICE_ROOT");
41 man->SetRun(0);
42 }
43 else {
44 printf("Using deafult storage \n");
45 }
46
47 // retrives geometry
48 if(!gGeoManager){
49 AliGeomManager::LoadGeometry("geometry.root");
50 }
51 if(align) {
52 TFile f(alignfile.Data());
53 TClonesArray* ar = (TClonesArray*)f.Get("ITSAlignObjs");
54 AliGeomManager::ApplyAlignObjsToGeom(*ar);
55 f.Close();
56 }
57
58 AliRunLoader* rl = AliRunLoader::Open("galice.root");
59 if (rl == 0x0){
60 cerr<<"Can not open session RL=NULL"<< endl;
61 return -1;
62 }
63 Int_t retval = rl->LoadgAlice();
64 if (retval){
65 cerr<<"LoadgAlice returned error"<<endl;
66 return -1;
67 }
68 gAlice=rl->GetAliRun();
69
70 retval = rl->LoadHeader();
71 if (retval){
72 cerr<<"LoadHeader returned error"<<endl;
73 return -1;
74 }
75
76
77 AliITSLoader* ITSloader = (AliITSLoader*) rl->GetLoader("ITSLoader");
78 if(!ITSloader){
79 cerr<<"ITS loader not found"<<endl;
80 return -1;
81 }
82 ITSloader->LoadRecPoints("read");
83 ITSloader->LoadHits("read");
84
85
86 Float_t cluglo[3]={0.,0.,0.};
87 AliITS *ITS = (AliITS*)gAlice->GetModule("ITS");
88 ITS->SetTreeAddress();
89 AliITSgeom *geom = ITS->GetITSgeom();
90 AliITSDetTypeRec* detTypeRec = new AliITSDetTypeRec();
91 detTypeRec->SetITSgeom(ITSloader->GetITSgeom());
92 detTypeRec->SetDefaults();
93
94 Int_t modmin=geom->GetStartDet(0);
95 Int_t modmax=geom->GetLastDet(2);
96 Int_t totmod=modmax-modmin;
97 Float_t xlim[6]={4.5,7.5,16.,26.,40.,45.};
98 Float_t zlim[6]={15.,15.,22.,30.,45.,55.};
99
100 TH1F* hlayer=new TH1F("hlayer","",6,0.5,6.5);
101 TH1F** hmod=new TH1F*[6];
102 TH1F** hxl=new TH1F*[6];
103 TH1F** hzl=new TH1F*[6];
104 TH1F** hxg=new TH1F*[6];
105 TH1F** hyg=new TH1F*[6];
106 TH2F** hxyg=new TH2F*[6];
107 TH2F** hxygHits=new TH2F*[6];
108 TH1F** hzg=new TH1F*[6];
109 TH1F** hr=new TH1F*[6];
110 TH1F** hphi=new TH1F*[6];
111 TH1F** hq=new TH1F*[6];
112
113 Char_t name[10];
114 for(Int_t iLay=0;iLay<6;iLay++){
115 sprintf(name,"hmod%d",iLay+1);
116 hmod[iLay]=new TH1F(name,"",totmod,modmin-0.5,modmax+0.5);
117 hmod[iLay]->GetXaxis()->SetTitle("Module");
118 hmod[iLay]->GetXaxis()->CenterTitle();
119 sprintf(name,"hxloc%d",iLay+1);
120 hxl[iLay]=new TH1F(name,"",100,-4.,4.);
121 hxl[iLay]->GetXaxis()->SetTitle("Xloc (cm)");
122 hxl[iLay]->GetXaxis()->CenterTitle();
123 sprintf(name,"hzloc%d",iLay+1);
124 hzl[iLay]=new TH1F(name,"",100,-4.,4.);
125 hzl[iLay]->GetXaxis()->SetTitle("Zloc (cm)");
126 hzl[iLay]->GetXaxis()->CenterTitle();
127 sprintf(name,"hxgl%d",iLay+1);
128 hxg[iLay]=new TH1F(name,"",100,-xlim[iLay],xlim[iLay]);
129 hxg[iLay]->GetXaxis()->SetTitle("Xglob (cm)");
130 hxg[iLay]->GetXaxis()->CenterTitle();
131 sprintf(name,"hygl%d",iLay+1);
132 hyg[iLay]=new TH1F(name,"",100,-xlim[iLay],xlim[iLay]);
133 hyg[iLay]->GetXaxis()->SetTitle("Yglob (cm)");
134 hyg[iLay]->GetXaxis()->CenterTitle();
135 sprintf(name,"hxygl%d",iLay+1);
136 hxyg[iLay]=new TH2F(name,"",1000,-xlim[iLay],xlim[iLay],1000,-xlim[iLay],xlim[iLay]);
137 hxyg[iLay]->GetXaxis()->SetTitle("Xglob (cm)");
138 hxyg[iLay]->GetYaxis()->SetTitle("Yglob (cm)");
139 sprintf(name,"hxygHitsl%d",iLay+1);
140 hxygHits[iLay]=new TH2F(name,"",1000,-xlim[iLay],xlim[iLay],1000,-xlim[iLay],xlim[iLay]);
141 sprintf(name,"hzgl%d",iLay+1);
142 hzg[iLay]=new TH1F(name,"",100,-zlim[iLay],zlim[iLay]);
143 hzg[iLay]->GetXaxis()->SetTitle("Zglob (cm)");
144 hzg[iLay]->GetXaxis()->CenterTitle();
145 sprintf(name,"hr%d",iLay+1);
146 hr[iLay]=new TH1F(name,"",100,0.,50.);
147 hr[iLay]->GetXaxis()->SetTitle("r (cm)");
148 hr[iLay]->GetXaxis()->CenterTitle();
149 sprintf(name,"hphi%d",iLay+1);
150 hphi[iLay]=new TH1F(name,"",100,-TMath::Pi(),TMath::Pi());
151 hphi[iLay]->GetXaxis()->SetTitle("#varphi (rad)");
152 hphi[iLay]->GetXaxis()->CenterTitle();
153 sprintf(name,"hq%d",iLay+1);
154 hq[iLay]=new TH1F(name,"",100,0.,300.);
155 hq[iLay]->GetXaxis()->SetTitle("Charge (keV)");
156 hq[iLay]->GetXaxis()->CenterTitle();
157 }
158
159 TGraph **gpts=new TGraph*[3];
160 TGraph2D **gpts3d=new TGraph2D*[3];
161 TGraph **gRP=new TGraph*[6];
162 TGraph **gHend=new TGraph*[6];
163 TGraph **gHstart=new TGraph*[6];
164 for(Int_t i=0;i<3;i++){
165 gpts[i]=new TGraph(0);
166 gpts3d[i]=new TGraph2D(0);
167 }
168 for(Int_t i=0;i<6;i++){
169 gRP[i]=new TGraph(0);
170 gRP[i]->GetXaxis()->SetTitle("global x [cm]");
171 gRP[i]->GetYaxis()->SetTitle("global y [cm]");
172 gHend[i]=new TGraph(0);
173 gHstart[i]=new TGraph(0);
174 }
175 Int_t totev=rl->GetNumberOfEvents();
176 printf("Total Number of events = %d\n",totev);
177
178 Int_t iRP[6]={0},iH[6]={0};
179
180 for(Int_t iev=0;iev<totev;iev++){
181 rl->GetEvent(iev);
182 TTree *TR = ITSloader->TreeR();
183 TClonesArray *ITSrec = detTypeRec->RecPoints();
184 TBranch *branch = 0;
185 if(TR && ITSrec){
186 branch = ITSloader->TreeR()->GetBranch("ITSRecPoints");
187 if(branch)branch->SetAddress(&ITSrec);
188 }
189 TTree *TH = ITSloader->TreeH();
190 TClonesArray *hits=new TClonesArray("AliITShit",10000);
191 TH->SetBranchAddress("ITS",&hits);
192
193 Int_t nparticles = rl->GetHeader()->GetNtrack();
194 cout<<"Event #"<<iev<<" #Particles="<<nparticles<<endl;
195
196
197 Int_t ipt=0;
198 for(Int_t subd=0;subd<3;subd++){
199
200 Int_t first = geom->GetStartDet(subd);
201 Int_t last = geom->GetLastDet(subd);
202
203 for (Int_t mod=first; mod<=last; mod++){
204 detTypeRec->ResetRecPoints();
205 branch->GetEvent(mod);
206 Int_t nrecp = ITSrec->GetEntries();
207 if(nrecp>0){
208 for(Int_t irec=0;irec<nrecp;irec++) {
209 AliITSRecPoint *recp = (AliITSRecPoint*)ITSrec->At(irec);
210 Int_t lay=recp->GetLayer();
211 hlayer->Fill(lay);
212 recp->GetGlobalXYZ(cluglo);
213 Float_t rad=TMath::Sqrt(cluglo[0]*cluglo[0]+cluglo[1]*cluglo[1]);
214 Float_t phi=TMath::ATan2(cluglo[1],cluglo[0]);
215 if(iev<3){
216 gpts[iev]->SetPoint(ipt,cluglo[0],cluglo[1]);
217 gpts3d[iev]->SetPoint(ipt,cluglo[0],cluglo[1],cluglo[2]);
218 ipt++;
219 }
220 gRP[lay]->SetPoint(iRP[lay],cluglo[0],cluglo[1]); iRP[lay]++;
221 hmod[lay]->Fill(mod);
222 hzl[lay]->Fill(recp->GetDetLocalZ());
223 hxl[lay]->Fill(recp->GetDetLocalX());
224 hzg[lay]->Fill(cluglo[2]);
225 hyg[lay]->Fill(cluglo[1]);
226 hxyg[lay]->Fill(cluglo[0],cluglo[1]);
227 hxg[lay]->Fill(cluglo[0]);
228 hr[lay]->Fill(rad);
229 hphi[lay]->Fill(phi);
230 hq[lay]->Fill(recp->GetQ());
231 }
232 }
233 }
234 }
235
236
237 cout<<" Now read hits "<<endl;
238
239 Int_t nentrHits=(Int_t)TH->GetEntries();
240 for (Int_t i=0; i<nentrHits; i++) {
241 TH->GetEvent(i);
242 Int_t nhit=hits->GetEntriesFast();
243 for (Int_t ih=0; ih<nhit; ih++) {
244 AliITShit *h=(AliITShit*)hits->UncheckedAt(ih);
245 if(h->StatusExiting()) {
246 Double_t xl,yl,zl,tl,xl0,yl0,zl0,tl0;
247 h->GetPositionL(xl,yl,zl,tl);
248 h->GetPositionL0(xl0,yl0,zl0,tl0);
249 //if(TMath::Abs(yl-yl0)<0.0290) continue;
250 hxygHits[h->GetLayer()-1]->Fill(h->GetXG(),h->GetYG());
251 gHend[h->GetLayer()-1]->SetPoint(iH[h->GetLayer()-1],h->GetXG(),h->GetYG());
252 Double_t x0,y0,z0,t0;
253 h->GetPositionG0(x0,y0,z0,t0);
254 gHstart[h->GetLayer()-1]->SetPoint(iH[h->GetLayer()-1],x0,y0);
255 iH[h->GetLayer()-1]++;
256 //printf("layer %d hit length xy %f hit length yl %f\n",h->GetLayer()-1,TMath::Sqrt((x0-h->GetXG())*(x0-h->GetXG())+(y0-h->GetYG())*(y0-h->GetYG())),TMath::Abs(yl-yl0));
257 }
258 }
259 }
260 }
261
262
263 TCanvas **c=new TCanvas*[6];
264 Char_t ctit[30];
265 for(Int_t iLay=0;iLay<6;iLay++){
266 sprintf(name,"can%d",iLay+1);
267 sprintf(ctit,"Layer %d",iLay+1);
268 c[iLay]=new TCanvas(name,ctit,1200,900);
269 c[iLay]->Divide(3,3,0.001,0.001);
270 c[iLay]->cd(1);
271 hmod[iLay]->Draw();
272 c[iLay]->cd(2);
273 hxl[iLay]->Draw();
274 c[iLay]->cd(3);
275 hzl[iLay]->Draw();
276 c[iLay]->cd(4);
277 hxg[iLay]->Draw();
278 c[iLay]->cd(5);
279 hyg[iLay]->Draw();
280 c[iLay]->cd(6);
281 hzg[iLay]->Draw();
282 c[iLay]->cd(7);
283 hr[iLay]->Draw();
284 c[iLay]->cd(8);
285 hphi[iLay]->Draw();
286 c[iLay]->cd(9);
287 hxyg[iLay]->Draw();
288 }
289
290 TCanvas *cev0;
291 cev0=new TCanvas("cev0","Event 0",600,600);
292 gpts[0]->SetMarkerStyle(7);
293 gpts[0]->SetTitle(0);
294 gpts[0]->Draw("AP");
295
296 TCanvas *cev1;
297 cev1=new TCanvas("cev1","Event 1",600,600);
298 gpts[1]->SetMarkerStyle(7);
299 gpts[1]->SetTitle(0);
300 gpts[1]->Draw("AP");
301
302 TCanvas *cev2;
303 cev2=new TCanvas("cev2","Event 2",600,600);
304 gpts[2]->SetMarkerStyle(7);
305 gpts[2]->SetTitle(0);
306 gpts[2]->Draw("AP");
307
308 TCanvas *chr = new TCanvas("chr","chr");
309 chr->Divide(3,2);
310 for(Int_t i=0;i<6;i++) {
311 chr->cd(i+1);
312 gHend[i]->SetMarkerStyle(7);
313 gHend[i]->SetMarkerColor(1);
314 gHend[i]->Draw("A,P");
315 gHstart[i]->SetMarkerStyle(7);
316 gHstart[i]->SetMarkerColor(4);
317 gHstart[i]->Draw("P");
318 gRP[i]->SetMarkerStyle(7);
319 gRP[i]->SetMarkerColor(2);
320 gRP[i]->Draw("P");
321 }
322
323
324 return 0;
325}
326
327
00fe805a 328Int_t ShowITSHitsRecPointsNtuple(Bool_t align=kFALSE,
3fd01708 329 TString alignfile="ITSfullv11Misalignment.root")
330{
331 ///////////////////////////////////////////////////////////////////////
332 // Macro to check clusters and hits in the 6 ITS layers //
333 // Creates also ntuple //
334 ///////////////////////////////////////////////////////////////////////
335
336 if (gClassTable->GetID("AliRun") < 0) {
337 gInterpreter->ExecuteMacro("loadlibs.C");
338 }
339 else {
340 if(gAlice){
880d6abe 341 delete AliRunLoader::GetRunLoader();
3fd01708 342 delete gAlice;
343 gAlice=0;
344 }
345 }
346 // Set OCDB if needed
347 AliCDBManager* man = AliCDBManager::Instance();
348 if (!man->IsDefaultStorageSet()) {
349 printf("Setting a local default storage and run number 0\n");
350 man->SetDefaultStorage("local://$ALICE_ROOT");
351 man->SetRun(0);
352 }
353 else {
354 printf("Using deafult storage \n");
355 }
356
357 // retrives geometry
358 if(!gGeoManager){
359 AliGeomManager::LoadGeometry("geometry.root");
360 }
361 if(align) {
362 TFile f(alignfile.Data());
363 TClonesArray* ar = (TClonesArray*)f.Get("ITSAlignObjs");
364 AliGeomManager::ApplyAlignObjsToGeom(*ar);
365 f.Close();
366 }
367
368 AliRunLoader* rl = AliRunLoader::Open("galice.root");
369 if (rl == 0x0){
370 cerr<<"Can not open session RL=NULL"<< endl;
371 return -1;
372 }
373 Int_t retval = rl->LoadgAlice();
374 if (retval){
375 cerr<<"LoadgAlice returned error"<<endl;
376 return -1;
377 }
378 gAlice=rl->GetAliRun();
379
380 retval = rl->LoadHeader();
381 if (retval){
382 cerr<<"LoadHeader returned error"<<endl;
383 return -1;
384 }
385
386
387 AliITSLoader* ITSloader = (AliITSLoader*) rl->GetLoader("ITSLoader");
388 if(!ITSloader){
389 cerr<<"ITS loader not found"<<endl;
390 return -1;
391 }
392 ITSloader->LoadRecPoints("read");
393 ITSloader->LoadHits("read");
394
395
396 Float_t cluglo[3]={0.,0.,0.};
397 AliITS *ITS = (AliITS*)gAlice->GetModule("ITS");
398 ITS->SetTreeAddress();
399 AliITSgeom *geom = ITS->GetITSgeom();
400 AliITSDetTypeRec* detTypeRec = new AliITSDetTypeRec();
401 detTypeRec->SetITSgeom(ITSloader->GetITSgeom());
402 detTypeRec->SetDefaults();
403
404 Int_t modmin=geom->GetStartDet(0);
405 Int_t modmax=geom->GetLastDet(2);
406 Int_t totmod=modmax-modmin;
407 Float_t xlim[6]={4.5,7.5,16.,26.,40.,45.};
408 Float_t zlim[6]={15.,15.,22.,30.,45.,55.};
409
410 TH1F* hlayer=new TH1F("hlayer","",6,0.5,6.5);
411 TH1F** hmod=new TH1F*[6];
412 TH1F** hxl=new TH1F*[6];
413 TH1F** hzl=new TH1F*[6];
414 TH1F** hxg=new TH1F*[6];
415 TH1F** hyg=new TH1F*[6];
416 TH2F** hxyg=new TH2F*[6];
417 TH2F** hxygHits=new TH2F*[6];
418 TH1F** hzg=new TH1F*[6];
419 TH1F** hr=new TH1F*[6];
420 TH1F** hphi=new TH1F*[6];
421 TH1F** hq=new TH1F*[6];
422
423 Char_t name[10];
424 for(Int_t iLay=0;iLay<6;iLay++){
425 sprintf(name,"hmod%d",iLay+1);
426 hmod[iLay]=new TH1F(name,"",totmod,modmin-0.5,modmax+0.5);
427 hmod[iLay]->GetXaxis()->SetTitle("Module");
428 hmod[iLay]->GetXaxis()->CenterTitle();
429 sprintf(name,"hxloc%d",iLay+1);
430 hxl[iLay]=new TH1F(name,"",100,-4.,4.);
431 hxl[iLay]->GetXaxis()->SetTitle("Xloc (cm)");
432 hxl[iLay]->GetXaxis()->CenterTitle();
433 sprintf(name,"hzloc%d",iLay+1);
434 hzl[iLay]=new TH1F(name,"",100,-4.,4.);
435 hzl[iLay]->GetXaxis()->SetTitle("Zloc (cm)");
436 hzl[iLay]->GetXaxis()->CenterTitle();
437 sprintf(name,"hxgl%d",iLay+1);
438 hxg[iLay]=new TH1F(name,"",100,-xlim[iLay],xlim[iLay]);
439 hxg[iLay]->GetXaxis()->SetTitle("Xglob (cm)");
440 hxg[iLay]->GetXaxis()->CenterTitle();
441 sprintf(name,"hygl%d",iLay+1);
442 hyg[iLay]=new TH1F(name,"",100,-xlim[iLay],xlim[iLay]);
443 hyg[iLay]->GetXaxis()->SetTitle("Yglob (cm)");
444 hyg[iLay]->GetXaxis()->CenterTitle();
445 sprintf(name,"hxygl%d",iLay+1);
446 hxyg[iLay]=new TH2F(name,"",1000,-xlim[iLay],xlim[iLay],1000,-xlim[iLay],xlim[iLay]);
447 hxyg[iLay]->GetXaxis()->SetTitle("Xglob (cm)");
448 hxyg[iLay]->GetYaxis()->SetTitle("Yglob (cm)");
449 sprintf(name,"hxygHitsl%d",iLay+1);
450 hxygHits[iLay]=new TH2F(name,"",1000,-xlim[iLay],xlim[iLay],1000,-xlim[iLay],xlim[iLay]);
451 sprintf(name,"hzgl%d",iLay+1);
452 hzg[iLay]=new TH1F(name,"",100,-zlim[iLay],zlim[iLay]);
453 hzg[iLay]->GetXaxis()->SetTitle("Zglob (cm)");
454 hzg[iLay]->GetXaxis()->CenterTitle();
455 sprintf(name,"hr%d",iLay+1);
456 hr[iLay]=new TH1F(name,"",100,0.,50.);
457 hr[iLay]->GetXaxis()->SetTitle("r (cm)");
458 hr[iLay]->GetXaxis()->CenterTitle();
459 sprintf(name,"hphi%d",iLay+1);
460 hphi[iLay]=new TH1F(name,"",100,-TMath::Pi(),TMath::Pi());
461 hphi[iLay]->GetXaxis()->SetTitle("#varphi (rad)");
462 hphi[iLay]->GetXaxis()->CenterTitle();
463 sprintf(name,"hq%d",iLay+1);
464 hq[iLay]=new TH1F(name,"",100,0.,300.);
465 hq[iLay]->GetXaxis()->SetTitle("Charge (keV)");
466 hq[iLay]->GetXaxis()->CenterTitle();
467 }
468
469 TGraph **gpts=new TGraph*[3];
470 TGraph2D **gpts3d=new TGraph2D*[3];
471 TGraph **gRP=new TGraph*[6];
472 TGraph **gHend=new TGraph*[6];
473 TGraph **gHstart=new TGraph*[6];
474 for(Int_t i=0;i<3;i++){
475 gpts[i]=new TGraph(0);
476 gpts3d[i]=new TGraph2D(0);
477 }
478 for(Int_t i=0;i<6;i++){
479 gRP[i]=new TGraph(0);
480 gRP[i]->GetXaxis()->SetTitle("global x [cm]");
481 gRP[i]->GetYaxis()->SetTitle("global y [cm]");
482 gHend[i]=new TGraph(0);
483 gHstart[i]=new TGraph(0);
484 }
485 Int_t totev=rl->GetNumberOfEvents();
486 printf("Total Number of events = %d\n",totev);
487
488 Int_t iRP[6]={0},iH[6]={0};
489
490 TNtuple *nt = new TNtuple("nt","ntuple","lay:mod:deltaxl:deltax:deltay:deltaz:phi:z:xl:hlength");
491
492 for(Int_t iev=0;iev<totev;iev++){
493 rl->GetEvent(iev);
494 TTree *TR = ITSloader->TreeR();
495 TClonesArray *ITSrec = detTypeRec->RecPoints();
496 TBranch *branch = 0;
497 if(TR && ITSrec){
498 branch = ITSloader->TreeR()->GetBranch("ITSRecPoints");
499 if(branch)branch->SetAddress(&ITSrec);
500 }
501 TTree *TH = ITSloader->TreeH();
502 TClonesArray *hits=new TClonesArray("AliITShit",10000);
503 TH->SetBranchAddress("ITS",&hits);
504
505 Int_t nparticles = rl->GetHeader()->GetNtrack();
506 cout<<"Event #"<<iev<<" #Particles="<<nparticles<<endl;
507
508
509 Int_t ipt=0;
510 for(Int_t subd=0;subd<3;subd++){
511
512 Int_t first = geom->GetStartDet(subd);
513 Int_t last = geom->GetLastDet(subd);
514
515 for (Int_t mod=first; mod<=last; mod++){
516 detTypeRec->ResetRecPoints();
517 branch->GetEvent(mod);
518 Int_t nrecp = ITSrec->GetEntries();
519 if(nrecp>0){
520 for(Int_t irec=0;irec<nrecp;irec++) {
521 AliITSRecPoint *recp = (AliITSRecPoint*)ITSrec->At(irec);
522 Int_t lay=recp->GetLayer();
523 hlayer->Fill(lay);
524 recp->GetGlobalXYZ(cluglo);
525 Float_t rad=TMath::Sqrt(cluglo[0]*cluglo[0]+cluglo[1]*cluglo[1]);
526 Float_t phi=TMath::ATan2(cluglo[1],cluglo[0]);
527 if(cluglo[2]>0) {gRP[lay]->SetPoint(iRP[lay],cluglo[0],cluglo[1]); iRP[lay]++;}
528 hmod[lay]->Fill(mod);
529 hzl[lay]->Fill(recp->GetDetLocalZ());
530 hxl[lay]->Fill(recp->GetDetLocalX());
531 hzg[lay]->Fill(cluglo[2]);
532 hyg[lay]->Fill(cluglo[1]);
533 hxyg[lay]->Fill(cluglo[0],cluglo[1]);
534 hxg[lay]->Fill(cluglo[0]);
535 hr[lay]->Fill(rad);
536 hphi[lay]->Fill(phi);
537 hq[lay]->Fill(recp->GetQ());
538
539 Double_t hlength=0.;
540 for(Int_t jhits=0;jhits<TH->GetEntries();jhits++) {
541 TH->GetEvent(jhits);
542 Int_t nhit=hits->GetEntriesFast();
543 for (Int_t ih=0; ih<nhit; ih++) {
544 AliITShit *h=(AliITShit*)hits->UncheckedAt(ih);
545 if(h->GetTrack()!=recp->GetLabel(0)) continue;
546 if(h->GetLayer()-1!=lay) continue;
547 if(h->GetModule()!=mod) continue;
548 //hxygHits[h->GetLayer()-1]->Fill(h->GetXG(),h->GetYG());
549 //gHend[h->GetLayer()-1]->SetPoint(iH[h->GetLayer()-1],h->GetXG(),h->GetYG());
550 Double_t x0,y0,z0,t0;
551 Double_t x,y,z,t;
552 h->GetPositionG0(x0,y0,z0,t0);
553 h->GetPositionG(x,y,z,t);
554 //gHstart[h->GetLayer()-1]->SetPoint(iH[h->GetLayer()-1],x0,y0);
555 //iH[h->GetLayer()-1]++;
556 Double_t xl,yl,zl,tl,xl0,yl0,zl0,tl0;
557 h->GetPositionL(xl,yl,zl,tl);
558 h->GetPositionL0(xl0,yl0,zl0,tl0);
559 hlength += TMath::Abs(yl-yl0);
560 if(!h->StatusExiting()) continue;
561 nt->Fill(lay,
562 mod,
563 1.e4*(0.5*(xl+xl0)-recp->GetDetLocalX()),
564 1.e4*(0.5*(x+x0)-cluglo[0]),1.e4*(0.5*(y+y0)-cluglo[1]),
565 1.e4*(0.5*(z+z0)-cluglo[2]),
566 phi,
567 cluglo[2],
568 recp->GetDetLocalX(),
569 hlength);
570 hlength = 0.;
571 }
572 }
573
574
575 }
576 }
577 }
578 }
579
580
581 }
582
583
584 TCanvas **c=new TCanvas*[6];
585 Char_t ctit[30];
586 for(Int_t iLay=0;iLay<6;iLay++){
587 sprintf(name,"can%d",iLay+1);
588 sprintf(ctit,"Layer %d",iLay+1);
589 c[iLay]=new TCanvas(name,ctit,1200,900);
590 c[iLay]->Divide(3,3,0.001,0.001);
591 c[iLay]->cd(1);
592 hmod[iLay]->Draw();
593 c[iLay]->cd(2);
594 hxl[iLay]->Draw();
595 c[iLay]->cd(3);
596 hzl[iLay]->Draw();
597 c[iLay]->cd(4);
598 hxg[iLay]->Draw();
599 c[iLay]->cd(5);
600 hyg[iLay]->Draw();
601 c[iLay]->cd(6);
602 hzg[iLay]->Draw();
603 c[iLay]->cd(7);
604 hr[iLay]->Draw();
605 c[iLay]->cd(8);
606 hphi[iLay]->Draw();
607 c[iLay]->cd(9);
608 hxyg[iLay]->Draw();
609 }
610
611 TCanvas *cev0;
612 cev0=new TCanvas("cev0","Event 0",600,600);
613 gpts[0]->SetMarkerStyle(7);
614 gpts[0]->SetTitle(0);
615 gpts[0]->Draw("AP");
616
617 TCanvas *cev1;
618 cev1=new TCanvas("cev1","Event 1",600,600);
619 gpts[1]->SetMarkerStyle(7);
620 gpts[1]->SetTitle(0);
621 gpts[1]->Draw("AP");
622
623 TCanvas *cev2;
624 cev2=new TCanvas("cev2","Event 2",600,600);
625 gpts[2]->SetMarkerStyle(7);
626 gpts[2]->SetTitle(0);
627 gpts[2]->Draw("AP");
628
629 TFile *out= new TFile("ntHitsRecPoints.root","recreate");
630 nt->Write();
631 out->Close();
632
633 return 0;
634}