For Pythia with tune don't switch off MI in ConfigHeavyFlavor
[u/mrichter/AliRoot.git] / ITS / AliITSCompareHitsRecPoints.C
CommitLineData
012f0f4c 1Bool_t AliITSCompareHitsRecPoints(Char_t *rfn="galice.root"){
2 // This macro compares the average location of the hits associated
3 // with a specific track and the best match RecPoint location.
4 // Inputs:
5 // none.
6 // Outputs:
7 // none.
8 // Return:
9 // kTRUE if no errors were encountered, otherwise kFALSE.
10 TProfile *pxg,*pyg,*pzg,*pxl,*pyl,*pzl;
11 Double_t hg[3],hl[3],tof,rg[3],rl[3];
12 Float_t rgf[3];
13
14 // Dynamically link some shared libs
15 if (gClassTable->GetID("AliRun") < 0) {
16 gROOT->LoadMacro("loadlibs.C");
17 loadlibs();
18 } else {
19 if(gAlice){
33c3c91a 20 delete AliRunLoader::Instance();
012f0f4c 21 delete gAlice;
22 gAlice=0;
23 } // end if gAlice
24 }// end if gClassTable->GetID()
25
26 gROOT->LoadMacro("$(ALICE_ROOT)/ITS/AliITSstandard.C");
27
28 AliRunLoader *runl = AccessFile(rfn); // Set up to read in Data
29 Int_t retval = runl->LoadHeader();
30 if (retval){
31 cerr<<"AliITSPrintRecPoints.C : LoadHeader returned error"<<endl;
32 return kFALSE;
33 } // end if retval
34
35 AliITSLoader* ITSloader = (AliITSLoader*) runl->GetLoader("ITSLoader");
36 if(!ITSloader){
37 cerr<<"AliITSPrintRecPoints.C : ITS loader not found"<<endl;
38 return kFALSE;
39 } // end if !ITSloader
40
41 if(!gGeoManager){
42 gGeoManger = new TGeoManager();
43 gGeoManager->Import("geometry.root","");
44 } // end if !gGeoManager
45 if(ITSloader->LoadHits("read")!=0){
46 cerr<<"Error Loading Hits"<<endl;
47 return kFALSE;
48 }// end if ITSloader->LoadHits
49 if(ITSloader->LoadRecPoints("read")!=0){
50 cerr<<"Error Loading RecPoints"<<endl;
51 return kFALSE;
52 }// end if ITSloader->LoadRecPoints
53 AliITS *ITS = 0;
54 ITS = (AliITS*)(gAlice->GetDetector("ITS"));
55 if(!ITS){
56 cout << "Error: no ITS found. Aborting"<<endl;
57 return kFALSE;
58 } // end if !ITS
59 if(!(ITS->GetDetTypeSim())){
60 cout <<"No AliITSDetTypeSim object found in ITS"<<endl;
61 return kFALSE;
62 } // end if
63 AliITSmodule *m = 0;
64 AliITShit *h = 0;
65 AliITSDetTypeSim *sim = ITS->GetDetTypeSim();
66 AliITSgeom *gm=0;
67 gm = ITS->GetITSgeom();
68 if(!gm){
69 cout <<"No AliITSgeom object found in ITS"<<endl;
70 if(!gGeoManager){
71 cout <<"No gGeoManger. Aborting"<<endl;
72 return kFALSE;
73 }else{
74 ITS->UpdateInternalGeometry(gGeoManager);
75 } // end if
76 } // end if !AliITSgeom
77 //
880d6abe 78 Int_t nMods= gm->GetIndexMax(),nEvents=AliRunLoader::GetNumberOfEvents();
012f0f4c 79 Int_t mod=0,evnt=0,size=-1,irp=0,ih=0,trkindexOld=-1;
80 Double_t xmod,nHitPerTrack;
81 TTree *rpt = 0;
82 TClonesArray *rpa = 0;
83 AliITSRecPoint *rp = 0;
84 AliITSDetTypeRec *rec = new AliITSDetTypeRec(ITSloader);
85 rec->SetDefaults();
86 // We are going to use the following
87 // <hit_x[i]> == { Sum_j hit_x[i][j] }/N_hits per track
88 // d_x[i] == <hit_x[i]>-recpoint_x[i]
89 // <d_x> == {sum_i d_x[i]}/N_RecPoints (excluding noise and no merging)
90 // = {sum_i (<hit_x[i]>-recpoint_x[i])}/N_RecPoints
91 // = {sum_i ((sum_j hit_x[i][j])/N_hits_per_track -recpoinnt_x[i])}/N_Recpoints
92 // = sum_i sum_j {hit_x[i][j]/N_hits_per_track}/N_Recpoints -
93 // sum_i {recpoint_x[i]}/N_Recpoints
94 pxg = new TProfile("XdiffGlobal","Mean displacement in gX direction",nMods,
95 -0.5,(Double_t)(nMods)+0.5," ");
96 pxg->SetXTitle("module");
97 pxg->SetYTitle("Global X [cm]");
98 pyg = new TProfile("YdiffGlobal","Mean displacement in gY direction",nMods,
99 -0.5,(Double_t)(nMods)+0.5," ");
100 pyg->SetXTitle("module");
101 pyg->SetYTitle("Global Y [cm]");
102 pzg = new TProfile("ZdiffGlobal","Mean displacement in gZ direction",nMods,
103 -0.5,(Double_t)(nMods)+0.5," ");
104 pzg->SetXTitle("module");
105 pzg->SetYTitle("Global Z [cm]");
106 pxl = new TProfile("XdiffLocal","Mean displacement in lX direction",nMods,
107 -0.5,(Double_t)(nMods)+0.5," ");
108 pxl->SetXTitle("module");
109 pxl->SetYTitle("local X [cm]");
110 pyl = new TProfile("YdiffLocal","Mean displacement in lY direction",nMods,
111 -0.5,(Double_t)(nMods)+0.5," ");
112 pyl->SetXTitle("module");
113 pyl->SetYTitle("local Y [cm]");
114 pzl = new TProfile("ZdiffLocal","Mean displacement in lZ direction",nMods,
115 -0.5,(Double_t)(nMods)+0.5," ");
116 pzl->SetXTitle("module");
117 pzl->SetYTitle("local Z [cm]");
118 //
119 for(evnt=0;evnt<nEvents;evnt++){
120 runl->GetEvent(evnt);
121 ITS->InitModules(size,nMods);
122 ITS->FillModules(evnt,0,-1," "," ");
123 rec->SetTreeAddress();
124 for(mod=0;mod<nMods;mod++){
125 xmod = (Double_t) mod;
126 m = ITS->GetModule(mod);
127 rec->ResetRecPoints();
128 rpt = ITSloader->TreeR();
129 rpt->GetEvent(mod);
130 rpa = rec->RecPoints();
131 for(irp=0;irp<rpa->GetEntriesFast();irp++){
132 rp = (AliITSRecPoint*)(rpa->At(irp));
133 rp->GetGlobalXYZ(rgf);rg[0]=rgf[0];rg[1]=rgf[1];rg[2]=rgf[2];
134 rl[0] = rp->GetDetLocalX();
135 rl[1] = 0.0;
136 rl[2] = rp->GetDetLocalZ();
137 pxg->Fill(xmod,-rg[0],0.5);
138 pyg->Fill(xmod,-rg[1],0.5);
139 pzg->Fill(xmod,-rg[2],0.5);
140 pxl->Fill(xmod,-rl[0],0.5);
141 //pyl->Fill(xmod,-rl[1],0.5); // assumed to be zero always.
142 pzl->Fill(xmod,-rl[2],0.5);
143 }// rnf got itp
144 trkindexOld=-1;
145 nHitPerTrack = 0.0;
146 for(ih=0;ih<m->GetNhits();ih++){ // We want the median hit location
147 // for each track.
148 h = m->GetHit(ih);
149 if(m->GetHitTrackIndex(ih)!=trkindexOld){// Enterence location
150 trkindexOld = m->GetHitTrackIndex(ih);
151 nHitPerTrack = 1.0;
152 do{
153 nHitPerTrack += 1.0;
154 }while (m->GetHitTrackIndex(ih+nHitPerTrack-1)==trkindexOld);
155 h->GetPositionG0(hg[0],hg[1],hg[2],tof);
156 h->GetPositionL0(hl[0],hl[1],hl[2],tof);
157 pxg->Fill(xmod,hg[0],0.5/nHitPerTrack);
158 pyg->Fill(xmod,hg[1],0.5/nHitPerTrack);
159 pzg->Fill(xmod,hg[2],0.5/nHitPerTrack);
160 pxl->Fill(xmod,hl[0],0.5/nHitPerTrack);
161 pyl->Fill(xmod,hl[1],1.0/nHitPerTrack); // rl[1]=0 always.
162 pzl->Fill(xmod,hl[2],0.5/nHitPerTrack);
163 } // end if
164 if(nHitPerTrack==0.0) continue;
165 h->GetPositionG(hg[0],hg[1],hg[2],tof);
166 h->GetPositionL(hl[0],hl[1],hl[2],tof);
167 pxg->Fill(xmod,hg[0],0.5/nHitPerTrack);
168 pyg->Fill(xmod,hg[1],0.5/nHitPerTrack);
169 pzg->Fill(xmod,hg[2],0.5/nHitPerTrack);
170 pxl->Fill(xmod,hl[0],0.5/nHitPerTrack);
171 pyl->Fill(xmod,hl[1],1.0/nHitPerTrack); // rl[1]=0 always.
172 pzl->Fill(xmod,hl[2],0.5/nHitPerTrack);
173 } // end for ih
174 } // end for mod
175 ITS->ClearModules();
176 } // end for evnt
177 //
178 Int_t wh=800;
179 TCanvas *c0 = new TCanvas("c0","Average displacements between hits and RecPoints"
180 ,1.5*wh,wh);
181 c0->Divide(3,2);
182 c0->cd(1);
183 pxg->Draw();
184 c0->cd(2);
185 pyg->Draw();
186 c0->cd(3);
187 pzg->Draw();
188 c0->cd(4);
189 pxl->Draw();
190 c0->cd(5);
191 pyl->Draw();
192 c0->cd(6);
193 pzl->Draw();
194 c0->Update();
195 return kTRUE;
196}