]>
Commit | Line | Data |
---|---|---|
1 | void MUONtestsim(Int_t ncham, Int_t evNumber1=0,Int_t evNumber2=0) | |
2 | { | |
3 | ///////////////////////////////////////////////////////////////////////// | |
4 | // This macro is a small example of a ROOT macro | |
5 | // illustrating how to read the output of GALICE | |
6 | // and do some analysis. | |
7 | // | |
8 | ///////////////////////////////////////////////////////////////////////// | |
9 | ||
10 | // Dynamically link some shared libs | |
11 | ||
12 | ||
13 | if (gClassTable->GetID("AliRun") < 0) { | |
14 | gROOT->LoadMacro("loadlibs.C"); | |
15 | loadlibs(); | |
16 | } | |
17 | ||
18 | // Connect the Root Galice file containing Geometry, Kine and Hits | |
19 | ||
20 | TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject("galice.root"); | |
21 | ||
22 | ||
23 | if (!file) { | |
24 | printf("\n Creating galice.root \n"); | |
25 | file = new TFile("galice.root"); | |
26 | } else { | |
27 | printf("\n galice.root found in file list"); | |
28 | } | |
29 | file->ls(); | |
30 | ||
31 | // Get AliRun object from file or create it if not on file | |
32 | if (!gAlice) { | |
33 | gAlice = (AliRun*)(file->Get("gAlice")); | |
34 | if (gAlice) printf("AliRun object found on file\n"); | |
35 | if (!gAlice) { | |
36 | printf("\n create new gAlice object"); | |
37 | gAlice = new AliRun("gAlice","Alice test program"); | |
38 | } | |
39 | } | |
40 | ||
41 | // Create some histograms | |
42 | ||
43 | ||
44 | TH1F *ccharge1 = new TH1F("ccharge1","Cluster Charge distribution" | |
45 | ,100,0.,5000.); | |
46 | TH1F *pcharge1 = new TH1F("pcharge1","Pad Charge distribution" | |
47 | ,100,0.,200.); | |
48 | TH1F *xresid1 = new TH1F("xresid1","x-Residuals" | |
49 | ,100,-10.,10); | |
50 | TH1F *yresid1 = new TH1F("yresid1","y-Residuals" | |
51 | ,100,-0.2,0.2); | |
52 | TH1F *npads1 = new TH1F("npads1" ,"Pads per Hit" | |
53 | ,20,-0.5,19.5); | |
54 | TH2F *xresys1 = new TH2F("xresys1","x-Residuals systematics" | |
55 | ,50,-2.0,2.0,100,-1.0,1.0); | |
56 | TH2F *yresys1 = new TH2F("yresys1","y-Residuals systematics" | |
57 | ,50,-0.4,0.4,100,-0.1,0.1); | |
58 | ||
59 | TH1F *ccharge2 = new TH1F("ccharge2","Cluster Charge distribution" | |
60 | ,100,0.,5000.); | |
61 | TH1F *pcharge2 = new TH1F("pcharge2","Pad Charge distribution" | |
62 | ,100,0.,200.); | |
63 | TH1F *xresid2 = new TH1F("xresid2","x-Residuals" | |
64 | ,100,-1,1); | |
65 | TH1F *yresid2 = new TH1F("yresid2","y-Residuals" | |
66 | ,100,-10, 10.); | |
67 | TH1F *npads2 = new TH1F("npads2" ,"Pads per Hit" | |
68 | ,20,-0.5,19.5); | |
69 | TH2F *xresys2 = new TH2F("xresys2","x-Residuals systematics" | |
70 | ,50,-2.0,2.0,100,-1.0,1.0); | |
71 | TH2F *yresys2 = new TH2F("yresys2","y-Residuals systematics" | |
72 | ,50,-0.4,0.4,100,-0.1,0.1); | |
73 | ||
74 | ||
75 | AliMUONChamber* iChamber; | |
76 | // | |
77 | // Loop over events | |
78 | // | |
79 | Int_t Nh=0; | |
80 | Int_t Nh1=0; | |
81 | for (Int_t nev=0; nev<= evNumber2; nev++) { | |
82 | cout << "nev " << nev <<endl; | |
83 | Int_t nparticles = gAlice->GetEvent(nev); | |
84 | cout << "nparticles " << nparticles <<endl; | |
85 | if (nev < evNumber1) continue; | |
86 | if (nparticles <= 0) return; | |
87 | ||
88 | AliMUON *MUON = (AliMUON*) gAlice->GetDetector("MUON"); | |
89 | printf("\n track %d %d \n ", nev, MUON); | |
90 | ||
91 | TTree *TH = gAlice->TreeH(); | |
92 | Int_t ntracks = TH->GetEntries(); | |
93 | Int_t Nc=0; | |
94 | Int_t npad[2]; | |
95 | Float_t Q[2], xmean[2],ymean[2],xres[2],yres[2], xonpad[2], yonpad[2]; | |
96 | // | |
97 | // Loop over events | |
98 | // | |
99 | for (Int_t track=0; track<ntracks;track++) { | |
100 | ||
101 | gAlice->ResetHits(); | |
102 | Int_t nbytes += TH->GetEvent(track); | |
103 | if (MUON) { | |
104 | for(AliMUONHit* mHit=(AliMUONHit*)MUON->FirstHit(-1); | |
105 | mHit; | |
106 | mHit=(AliMUONHit*)MUON->NextHit()) | |
107 | { | |
108 | Int_t nch = mHit->fChamber; // chamber number | |
109 | Float_t x = mHit->X(); // x-pos of hit | |
110 | Float_t y = mHit->Y(); // y-pos | |
111 | Float_t z = mHit->Z(); | |
112 | ||
113 | // | |
114 | // | |
115 | iChamber = & MUON->Chamber(nch-1); | |
116 | response=iChamber->ResponseModel(); | |
117 | // | |
118 | // | |
119 | if (nch == ncham) { | |
120 | for (Int_t i=0; i<2; i++) { | |
121 | xmean[i]=0; | |
122 | ymean[i]=0; | |
123 | Q[i]=0; | |
124 | npad[i]=0; | |
125 | } | |
126 | ||
127 | for (AliMUONPadHit* mPad= | |
128 | (AliMUONPadHit*)MUON->FirstPad(mHit,MUON->PadHits()); | |
129 | mPad; | |
130 | mPad=(AliMUONPadHit*)MUON->NextPad(MUON->PadHits())) | |
131 | { | |
132 | ||
133 | Int_t nseg = mPad->fCathode; // segmentation | |
134 | Int_t nhit = mPad->fHitNumber; // hit number | |
135 | Int_t qtot = mPad->fQ; // charge | |
136 | Int_t ipx = mPad->fPadX; // pad number on X | |
137 | Int_t ipy = mPad->fPadY; // pad number on Y | |
138 | Int_t iqpad = mPad->fQpad; // charge per pad | |
139 | Int_t secpad = mPad->fRSec; // r-pos of pad | |
140 | // | |
141 | // | |
142 | segmentation=iChamber->SegmentationModel(nseg); | |
143 | Int_t ind=nseg-1; | |
144 | Float_t xg,yg,zg; | |
145 | segmentation->GetPadC(ipx, ipy, xg, yg, zg); | |
146 | Int_t sec = segmentation->Sector(ipx, ipy); | |
147 | Float_t dpx=segmentation->Dpx(sec); | |
148 | Float_t dpy=segmentation->Dpy(sec); | |
149 | ||
150 | if (nseg == 1) { | |
151 | pcharge1->Fill(iqpad,(float) 1); | |
152 | ccharge1->Fill(qtot ,(float) 1); | |
153 | } else { | |
154 | pcharge2->Fill(iqpad,(float) 1); | |
155 | ccharge2->Fill(qtot ,(float) 1); | |
156 | } | |
157 | // Calculate centre of gravity | |
158 | // | |
159 | if (iqpad == 0 && ind==0) { | |
160 | printf("\n attention iqpad 0"); | |
161 | } | |
162 | ||
163 | if (iqpad > 0) { | |
164 | Float_t xc,yc,zc; | |
165 | npad[ind]++; | |
166 | segmentation->GetPadC(ipx,ipy,xc,yc,zc); | |
167 | xmean[ind]+=Float_t(iqpad*xc); | |
168 | ymean[ind]+=Float_t(iqpad*yc); | |
169 | Q[ind]+=iqpad; | |
170 | } | |
171 | ||
172 | } //pad hit loop | |
173 | for (Int_t i=0; i<2; i++) { | |
174 | segmentation = iChamber->SegmentationModel(i+1); | |
175 | if (Q[i] >0) { | |
176 | xmean[i] = xmean[i]/Q[i]; | |
177 | xres[i] = xmean[i]-x; | |
178 | ymean[i] = ymean[i]/Q[i]; | |
179 | yres[i] = ymean[i]-y; | |
180 | // Systematic Error | |
181 | // | |
182 | Int_t icx, icy; | |
183 | Float_t zonpad; | |
184 | ||
185 | segmentation-> | |
186 | GetPadI(x,y,z,icx,icy); | |
187 | segmentation-> | |
188 | GetPadC(icx,icy,xonpad[i],yonpad[i],zonpad); | |
189 | xonpad[i]-=x; | |
190 | yonpad[i]-=y; | |
191 | } // charge not 0 | |
192 | } // plane loop | |
193 | xresid1->Fill(xres[0],(float) 1); | |
194 | yresid1->Fill(yres[0],(float) 1); | |
195 | npads1->Fill(npad[0],(float) 1); | |
196 | if (npad[0] >=2 && Q[0] > 6 ) { | |
197 | xresys1->Fill(xonpad[0],xres[0],(float) 1); | |
198 | yresys1->Fill(yonpad[0],yres[0],(float) 1); | |
199 | } | |
200 | ||
201 | xresid2->Fill(xres[1],(float) 1); | |
202 | yresid2->Fill(yres[1],(float) 1); | |
203 | npads2->Fill(npad[1],(float) 1); | |
204 | if (npad[1] >=2 && Q[1] > 6) { | |
205 | xresys2->Fill(xonpad[1],xres[1],(float) 1); | |
206 | yresys2->Fill(yonpad[1],yres[1],(float) 1); | |
207 | } | |
208 | } // chamber 1 | |
209 | } // hit loop | |
210 | } // if MUON | |
211 | } // track loop | |
212 | } // event loop | |
213 | ||
214 | //Create a canvas, set the view range, show histograms | |
215 | TCanvas *c1 = new TCanvas("c1","Charge and Residuals (1)",400,10,600,700); | |
216 | pad11 = new TPad("pad11"," ",0.01,0.51,0.49,0.99); | |
217 | pad12 = new TPad("pad12"," ",0.51,0.51,0.99,0.99); | |
218 | pad13 = new TPad("pad13"," ",0.01,0.01,0.49,0.49); | |
219 | pad14 = new TPad("pad14"," ",0.51,0.01,0.99,0.49); | |
220 | pad11->SetFillColor(11); | |
221 | pad12->SetFillColor(11); | |
222 | pad13->SetFillColor(11); | |
223 | pad14->SetFillColor(11); | |
224 | pad11->Draw(); | |
225 | pad12->Draw(); | |
226 | pad13->Draw(); | |
227 | pad14->Draw(); | |
228 | ||
229 | pad11->cd(); | |
230 | ccharge1->SetFillColor(42); | |
231 | ccharge1->SetXTitle("ADC units"); | |
232 | ccharge1->Draw(); | |
233 | ||
234 | pad12->cd(); | |
235 | pcharge1->SetFillColor(42); | |
236 | pcharge1->SetXTitle("ADC units"); | |
237 | pcharge1->Draw(); | |
238 | ||
239 | pad13->cd(); | |
240 | xresid1->SetFillColor(42); | |
241 | xresid1->Draw(); | |
242 | ||
243 | pad14->cd(); | |
244 | yresid1->SetFillColor(42); | |
245 | yresid1->Draw(); | |
246 | ||
247 | TCanvas *c2 = new TCanvas("c2","Cluster Size (1)",400,10,600,700); | |
248 | pad21 = new TPad("pad21"," ",0.01,0.51,0.49,0.99); | |
249 | pad22 = new TPad("pad22"," ",0.51,0.51,0.99,0.99); | |
250 | pad23 = new TPad("pad23"," ",0.01,0.01,0.49,0.49); | |
251 | pad24 = new TPad("pad24"," ",0.51,0.01,0.99,0.49); | |
252 | pad21->SetFillColor(11); | |
253 | pad22->SetFillColor(11); | |
254 | pad23->SetFillColor(11); | |
255 | pad24->SetFillColor(11); | |
256 | pad21->Draw(); | |
257 | pad22->Draw(); | |
258 | pad23->Draw(); | |
259 | pad24->Draw(); | |
260 | ||
261 | pad21->cd(); | |
262 | npads1->SetFillColor(42); | |
263 | npads1->SetXTitle("Cluster Size"); | |
264 | npads1->Draw(); | |
265 | ||
266 | pad23->cd(); | |
267 | xresys1->SetXTitle("x on pad"); | |
268 | xresys1->SetYTitle("x-xcog"); | |
269 | xresys1->Draw(); | |
270 | ||
271 | pad24->cd(); | |
272 | yresys1->SetXTitle("y on pad"); | |
273 | yresys1->SetYTitle("y-ycog"); | |
274 | yresys1->Draw(); | |
275 | ||
276 | TCanvas *c3 = new TCanvas("c3","Charge and Residuals (2)",400,10,600,700); | |
277 | pad31 = new TPad("pad31"," ",0.01,0.51,0.49,0.99); | |
278 | pad32 = new TPad("pad32"," ",0.51,0.51,0.99,0.99); | |
279 | pad33 = new TPad("pad33"," ",0.01,0.01,0.49,0.49); | |
280 | pad34 = new TPad("pad34"," ",0.51,0.01,0.99,0.49); | |
281 | pad31->SetFillColor(11); | |
282 | pad32->SetFillColor(11); | |
283 | pad33->SetFillColor(11); | |
284 | pad34->SetFillColor(11); | |
285 | pad31->Draw(); | |
286 | pad32->Draw(); | |
287 | pad33->Draw(); | |
288 | pad34->Draw(); | |
289 | ||
290 | pad31->cd(); | |
291 | ccharge2->SetFillColor(42); | |
292 | ccharge2->SetXTitle("ADC units"); | |
293 | ccharge2->Draw(); | |
294 | ||
295 | pad32->cd(); | |
296 | pcharge2->SetFillColor(42); | |
297 | pcharge2->SetXTitle("ADC units"); | |
298 | pcharge2->Draw(); | |
299 | ||
300 | pad33->cd(); | |
301 | xresid2->SetFillColor(42); | |
302 | xresid2->Draw(); | |
303 | ||
304 | pad34->cd(); | |
305 | yresid2->SetFillColor(42); | |
306 | yresid2->Draw(); | |
307 | ||
308 | TCanvas *c4 = new TCanvas("c4","Cluster Size (2)",400,10,600,700); | |
309 | pad41 = new TPad("pad41"," ",0.01,0.51,0.49,0.99); | |
310 | pad42 = new TPad("pad42"," ",0.51,0.51,0.99,0.99); | |
311 | pad43 = new TPad("pad43"," ",0.01,0.01,0.49,0.49); | |
312 | pad44 = new TPad("pad44"," ",0.51,0.01,0.99,0.49); | |
313 | pad41->SetFillColor(11); | |
314 | pad42->SetFillColor(11); | |
315 | pad43->SetFillColor(11); | |
316 | pad44->SetFillColor(11); | |
317 | pad41->Draw(); | |
318 | pad42->Draw(); | |
319 | pad43->Draw(); | |
320 | pad44->Draw(); | |
321 | ||
322 | pad41->cd(); | |
323 | npads2->SetFillColor(42); | |
324 | npads2->SetXTitle("Cluster Size"); | |
325 | npads2->Draw(); | |
326 | ||
327 | pad43->cd(); | |
328 | xresys2->SetXTitle("x on pad"); | |
329 | xresys2->SetYTitle("x-xcog"); | |
330 | xresys2->Draw(); | |
331 | ||
332 | pad44->cd(); | |
333 | yresys2->SetXTitle("y on pad"); | |
334 | yresys2->SetYTitle("y-ycog"); | |
335 | yresys2->Draw(); | |
336 | ||
337 | } |