]>
Commit | Line | Data |
---|---|---|
9123a941 | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: The ALICE Off-line Project. * | |
5 | * Contributors are mentioned in the code where appropriate. * | |
6 | * * | |
7 | * Permission to use, copy, modify and distribute this software and its * | |
8 | * documentation strictly for non-commercial purposes is hereby granted * | |
9 | * without fee, provided that the above copyright notice appears in all * | |
10 | * copies and that both the copyright notice and this permission notice * | |
11 | * appear in the supporting documentation. The authors make no claims * | |
12 | * about the suitability of this software for any purpose. It is * | |
13 | * provided "as is" without express or implied warranty. * | |
14 | **************************************************************************/ | |
15 | ||
9d6c0962 | 16 | |
7118aef0 | 17 | #include <Riostream.h> |
9123a941 | 18 | |
88cb7938 | 19 | #include <TBRIK.h> |
20 | #include <TGeometry.h> | |
9d6c0962 | 21 | #include <TLorentzVector.h> |
88cb7938 | 22 | #include <TNode.h> |
9d6c0962 | 23 | #include <TParticle.h> |
88cb7938 | 24 | #include <TVector3.h> |
25 | #include <TVirtualMC.h> | |
cd1cf354 | 26 | #include <TPDGCode.h> //for kNuetron |
c021cb15 | 27 | #include <TCanvas.h> |
28 | #include <TF1.h> | |
29 | #include <TH1.h> | |
30 | #include <TH2.h> | |
31 | #include <TStyle.h> | |
9d6c0962 | 32 | |
88cb7938 | 33 | #include "AliConst.h" |
34 | #include "AliMagF.h" | |
35 | #include "AliPDG.h" | |
ae714751 | 36 | #include "AliRICHGeometry.h" |
88cb7938 | 37 | #include "AliRICHResponseV0.h" |
38 | #include "AliRICHSegmentationV1.h" | |
39 | #include "AliRICHv3.h" | |
40 | #include "AliRun.h" | |
c021cb15 | 41 | #include "AliRICHRawCluster.h" |
42 | #include "AliRICHDigit.h" | |
43 | #include "AliRICHRecHit1D.h" | |
44 | ||
ae714751 | 45 | |
9d6c0962 | 46 | ClassImp(AliRICHv3) |
ae714751 | 47 | |
9d6c0962 | 48 | //______________________________________________________________ |
49 | // Implementation of the RICH version 3 with azimuthal rotation | |
9123a941 | 50 | |
9123a941 | 51 | |
9d6c0962 | 52 | AliRICHv3::AliRICHv3(const char *sName, const char *sTitle) |
53 | :AliRICH(sName,sTitle) | |
9123a941 | 54 | { |
9d6c0962 | 55 | // The named ctor currently creates a single copy of |
56 | // AliRICHGeometry AliRICHSegmentationV1 AliRICHResponseV0 | |
57 | // and initialises the corresponding models of all 7 chambers with these stuctures. | |
58 | // Note: all chambers share the single copy of models. MUST be changed later (???). | |
cd1cf354 | 59 | if(GetDebug())Info("named ctor","Start."); |
9123a941 | 60 | |
cd1cf354 | 61 | fCkovNumber=fFreonProd=0; |
ae714751 | 62 | |
9d6c0962 | 63 | AliRICHGeometry *pRICHGeometry =new AliRICHGeometry; // ??? to be moved to AlRICHChamber::named ctor |
64 | AliRICHSegmentationV1 *pRICHSegmentation=new AliRICHSegmentationV1; // ??? to be moved to AlRICHChamber::named ctor | |
65 | AliRICHResponseV0 *pRICHResponse =new AliRICHResponseV0; // ??? to be moved to AlRICHChamber::named ctor | |
ae714751 | 66 | |
c021cb15 | 67 | for (Int_t i=1; i<=kNCH; i++){ |
ae714751 | 68 | SetGeometryModel(i,pRICHGeometry); |
69 | SetSegmentationModel(i,pRICHSegmentation); | |
70 | SetResponseModel(i,pRICHResponse); | |
c021cb15 | 71 | C(i)->Init(i); // ??? to be removed |
9123a941 | 72 | } |
cd1cf354 | 73 | if(GetDebug())Info("named ctor","Stop."); |
9123a941 | 74 | }//AliRICHv3::ctor(const char *pcName, const char *pcTitle) |
75 | ||
9d6c0962 | 76 | AliRICHv3::~AliRICHv3() |
77 | { | |
78 | // Dtor deletes RICH models. In future (???) AliRICHChamber will be responsible for that. | |
cd1cf354 | 79 | if(GetDebug()) cout<<ClassName()<<"::dtor()>\n"; |
9d6c0962 | 80 | |
00ec1fc7 | 81 | if(fChambers) { |
c021cb15 | 82 | AliRICHChamber *ch =C(1); |
00ec1fc7 | 83 | if(ch) { |
84 | delete ch->GetGeometryModel(); | |
85 | delete ch->GetResponseModel(); | |
86 | delete ch->GetSegmentationModel(); | |
87 | } | |
88 | } | |
9d6c0962 | 89 | }//AliRICHv3::dtor() |
90 | ||
9123a941 | 91 | |
92 | void AliRICHv3::CreateGeometry() | |
93 | { | |
9d6c0962 | 94 | // Provides geometry structure for simulation (currently GEANT volumes tree) |
cd1cf354 | 95 | if(GetDebug()) cout<<ClassName()<<"::CreateGeometry()>\n"; |
9123a941 | 96 | |
97 | AliRICH *pRICH = (AliRICH *) gAlice->GetDetector("RICH"); | |
98 | AliRICHSegmentationV0* segmentation; | |
99 | AliRICHGeometry* geometry; | |
100 | AliRICHChamber* iChamber; | |
101 | ||
102 | iChamber = &(pRICH->Chamber(0)); | |
cc23c5c6 | 103 | segmentation=(AliRICHSegmentationV0*) iChamber->GetSegmentationModel(); |
9123a941 | 104 | geometry=iChamber->GetGeometryModel(); |
105 | ||
106 | Float_t distance; | |
107 | distance = geometry->GetFreonThickness()/2 + geometry->GetQuartzThickness() + geometry->GetGapThickness(); | |
108 | geometry->SetRadiatorToPads(distance); | |
109 | ||
110 | //Opaque quartz thickness | |
111 | Float_t oqua_thickness = .5; | |
112 | //CsI dimensions | |
113 | ||
114 | ||
115 | Float_t csi_width = segmentation->Npx()*segmentation->Dpx() + segmentation->DeadZone(); | |
116 | Float_t csi_length = segmentation->Npy()*segmentation->Dpy() + 2*segmentation->DeadZone(); | |
117 | ||
118 | ||
119 | Int_t *idtmed = fIdtmed->GetArray()-999; | |
120 | ||
121 | Int_t i; | |
122 | Float_t zs; | |
123 | Int_t idrotm[1099]; | |
124 | Float_t par[3]; | |
125 | ||
126 | // --- Define the RICH detector | |
127 | // External aluminium box | |
128 | par[0] = 68.8; | |
129 | par[1] = 13; //Original Settings | |
130 | par[2] = 70.86; | |
131 | gMC->Gsvolu("RICH", "BOX ", idtmed[1009], par, 3); | |
132 | ||
133 | // Air | |
134 | par[0] = 66.3; | |
135 | par[1] = 13; //Original Settings | |
136 | par[2] = 68.35; | |
137 | gMC->Gsvolu("SRIC", "BOX ", idtmed[1000], par, 3); | |
138 | ||
139 | // Air 2 (cutting the lower part of the box) | |
140 | ||
141 | par[0] = 1.25; | |
142 | par[1] = 3; //Original Settings | |
143 | par[2] = 70.86; | |
144 | gMC->Gsvolu("AIR2", "BOX ", idtmed[1000], par, 3); | |
145 | ||
146 | // Air 3 (cutting the lower part of the box) | |
147 | ||
148 | par[0] = 66.3; | |
149 | par[1] = 3; //Original Settings | |
150 | par[2] = 1.2505; | |
151 | gMC->Gsvolu("AIR3", "BOX ", idtmed[1000], par, 3); | |
152 | ||
153 | // Honeycomb | |
154 | par[0] = 66.3; | |
155 | par[1] = .188; //Original Settings | |
156 | par[2] = 68.35; | |
157 | gMC->Gsvolu("HONE", "BOX ", idtmed[1001], par, 3); | |
158 | ||
159 | // Aluminium sheet | |
160 | par[0] = 66.3; | |
161 | par[1] = .025; //Original Settings | |
162 | par[2] = 68.35; | |
163 | /*par[0] = 66.5; | |
164 | par[1] = .025; | |
165 | par[2] = 63.1;*/ | |
166 | gMC->Gsvolu("ALUM", "BOX ", idtmed[1009], par, 3); | |
167 | ||
168 | // Quartz | |
169 | par[0] = geometry->GetQuartzWidth()/2; | |
170 | par[1] = geometry->GetQuartzThickness()/2; | |
171 | par[2] = geometry->GetQuartzLength()/2; | |
172 | gMC->Gsvolu("QUAR", "BOX ", idtmed[1002], par, 3); | |
173 | ||
174 | // Spacers (cylinders) | |
175 | par[0] = 0.; | |
176 | par[1] = .5; | |
177 | par[2] = geometry->GetFreonThickness()/2; | |
178 | gMC->Gsvolu("SPAC", "TUBE", idtmed[1002], par, 3); | |
179 | ||
180 | // Feet (freon slabs supports) | |
181 | ||
182 | par[0] = .7; | |
183 | par[1] = .3; | |
184 | par[2] = 1.9; | |
185 | gMC->Gsvolu("FOOT", "BOX", idtmed[1009], par, 3); | |
186 | ||
187 | // Opaque quartz | |
188 | par[0] = geometry->GetQuartzWidth()/2; | |
189 | par[1] = .2; | |
190 | par[2] = geometry->GetQuartzLength()/2; | |
191 | gMC->Gsvolu("OQUA", "BOX ", idtmed[1007], par, 3); | |
192 | ||
193 | // Frame of opaque quartz | |
194 | par[0] = geometry->GetOuterFreonWidth()/2; | |
195 | par[1] = geometry->GetFreonThickness()/2; | |
196 | par[2] = geometry->GetOuterFreonLength()/2; | |
197 | gMC->Gsvolu("OQF1", "BOX ", idtmed[1007], par, 3); | |
198 | ||
199 | par[0] = geometry->GetInnerFreonWidth()/2; | |
200 | par[1] = geometry->GetFreonThickness()/2; | |
201 | par[2] = geometry->GetInnerFreonLength()/2; | |
202 | gMC->Gsvolu("OQF2", "BOX ", idtmed[1007], par, 3); | |
203 | ||
204 | ||
205 | // Freon | |
206 | par[0] = geometry->GetOuterFreonWidth()/2 - oqua_thickness; | |
207 | par[1] = geometry->GetFreonThickness()/2; | |
208 | par[2] = geometry->GetOuterFreonLength()/2 - 2*oqua_thickness; | |
209 | gMC->Gsvolu("FRE1", "BOX ", idtmed[1003], par, 3); | |
210 | ||
211 | par[0] = geometry->GetInnerFreonWidth()/2 - oqua_thickness; | |
212 | par[1] = geometry->GetFreonThickness()/2; | |
213 | par[2] = geometry->GetInnerFreonLength()/2 - 2*oqua_thickness; | |
214 | gMC->Gsvolu("FRE2", "BOX ", idtmed[1003], par, 3); | |
215 | ||
216 | // Methane | |
217 | par[0] = csi_width/2; | |
218 | par[1] = geometry->GetGapThickness()/2; | |
219 | par[2] = csi_length/2; | |
220 | gMC->Gsvolu("META", "BOX ", idtmed[1004], par, 3); | |
221 | ||
222 | // Methane gap | |
223 | par[0] = csi_width/2; | |
224 | par[1] = geometry->GetProximityGapThickness()/2; | |
225 | par[2] = csi_length/2; | |
226 | gMC->Gsvolu("GAP ", "BOX ", idtmed[1008], par, 3); | |
227 | ||
228 | // CsI photocathode | |
229 | par[0] = csi_width/2; | |
230 | par[1] = .25; | |
231 | par[2] = csi_length/2; | |
232 | gMC->Gsvolu("CSI ", "BOX ", idtmed[1005], par, 3); | |
233 | ||
234 | // Anode grid | |
235 | par[0] = 0.; | |
236 | par[1] = .001; | |
237 | par[2] = 20.; | |
238 | gMC->Gsvolu("GRID", "TUBE", idtmed[1006], par, 3); | |
239 | ||
240 | // Wire supports | |
241 | // Bar of metal | |
242 | ||
243 | par[0] = csi_width/2; | |
244 | par[1] = 1.05; | |
245 | par[2] = 1.05; | |
246 | gMC->Gsvolu("WSMe", "BOX ", idtmed[1009], par, 3); | |
247 | ||
248 | // Ceramic pick up (base) | |
249 | ||
250 | par[0] = csi_width/2; | |
251 | par[1] = .25; | |
252 | par[2] = 1.05; | |
253 | gMC->Gsvolu("WSG1", "BOX ", idtmed[1010], par, 3); | |
254 | ||
255 | // Ceramic pick up (head) | |
256 | ||
257 | par[0] = csi_width/2; | |
258 | par[1] = .1; | |
259 | par[2] = .1; | |
260 | gMC->Gsvolu("WSG2", "BOX ", idtmed[1010], par, 3); | |
261 | ||
262 | // Aluminium supports for methane and CsI | |
263 | // Short bar | |
264 | ||
265 | par[0] = csi_width/2; | |
266 | par[1] = geometry->GetGapThickness()/2 + .25; | |
267 | par[2] = (68.35 - csi_length/2)/2; | |
268 | gMC->Gsvolu("SMSH", "BOX", idtmed[1009], par, 3); | |
269 | ||
270 | // Long bar | |
271 | ||
272 | par[0] = (66.3 - csi_width/2)/2; | |
273 | par[1] = geometry->GetGapThickness()/2 + .25; | |
274 | par[2] = csi_length/2 + 68.35 - csi_length/2; | |
275 | gMC->Gsvolu("SMLG", "BOX", idtmed[1009], par, 3); | |
276 | ||
277 | // Aluminium supports for freon | |
278 | // Short bar | |
279 | ||
280 | par[0] = geometry->GetQuartzWidth()/2; | |
281 | par[1] = .3; | |
282 | par[2] = (68.35 - geometry->GetQuartzLength()/2)/2; | |
283 | gMC->Gsvolu("SFSH", "BOX", idtmed[1009], par, 3); | |
284 | ||
285 | // Long bar | |
286 | ||
287 | par[0] = (66.3 - geometry->GetQuartzWidth()/2)/2; | |
288 | par[1] = .3; | |
289 | par[2] = geometry->GetQuartzLength()/2 + 68.35 - geometry->GetQuartzLength()/2; | |
290 | gMC->Gsvolu("SFLG", "BOX", idtmed[1009], par, 3); | |
291 | ||
292 | // PCB backplane | |
293 | ||
294 | par[0] = csi_width/2; | |
295 | par[1] = .25; | |
296 | par[2] = csi_length/4 -.5025; | |
297 | gMC->Gsvolu("PCB ", "BOX", idtmed[1011], par, 3); | |
298 | ||
299 | ||
300 | // Backplane supports | |
301 | ||
302 | // Aluminium slab | |
303 | ||
304 | par[0] = 33.15; | |
305 | par[1] = 2; | |
306 | par[2] = 21.65; | |
307 | gMC->Gsvolu("BACK", "BOX", idtmed[1009], par, 3); | |
308 | ||
309 | // Big hole | |
310 | ||
311 | par[0] = 9.05; | |
312 | par[1] = 2; | |
313 | par[2] = 4.4625; | |
314 | gMC->Gsvolu("BKHL", "BOX", idtmed[1000], par, 3); | |
315 | ||
316 | // Small hole | |
317 | ||
318 | par[0] = 5.7; | |
319 | par[1] = 2; | |
320 | par[2] = 4.4625; | |
321 | gMC->Gsvolu("BKHS", "BOX", idtmed[1000], par, 3); | |
322 | ||
323 | // Place holes inside backplane support | |
324 | ||
325 | gMC->Gspos("BKHS", 1, "BACK", .8 + 5.7,0., .6 + 4.4625, 0, "ONLY"); | |
326 | gMC->Gspos("BKHS", 2, "BACK", -.8 - 5.7,0., .6 + 4.4625, 0, "ONLY"); | |
327 | gMC->Gspos("BKHS", 3, "BACK", .8 + 5.7,0., -.6 - 4.4625, 0, "ONLY"); | |
328 | gMC->Gspos("BKHS", 4, "BACK", -.8 - 5.7,0., -.6 - 4.4625, 0, "ONLY"); | |
329 | gMC->Gspos("BKHS", 5, "BACK", .8 + 5.7,0., .6 + 8.925 + 1.2 + 4.4625, 0, "ONLY"); | |
330 | gMC->Gspos("BKHS", 6, "BACK", -.8 - 5.7,0., .6 + 8.925 + 1.2 + 4.4625, 0, "ONLY"); | |
331 | gMC->Gspos("BKHS", 7, "BACK", .8 + 5.7,0., -.6 - 8.925 - 1.2 - 4.4625, 0, "ONLY"); | |
332 | gMC->Gspos("BKHS", 8, "BACK", -.8 - 5.7,0., -.6 - 8.925 - 1.2 - 4.4625, 0, "ONLY"); | |
333 | gMC->Gspos("BKHL", 1, "BACK", .8 + 11.4 + 1.6 + 9.05, 0., .6 + 4.4625, 0, "ONLY"); | |
334 | gMC->Gspos("BKHL", 2, "BACK", -.8 - 11.4 - 1.6 - 9.05, 0., .6 + 4.4625, 0, "ONLY"); | |
335 | gMC->Gspos("BKHL", 3, "BACK", .8 + 11.4 + 1.6 + 9.05, 0., -.6 - 4.4625, 0, "ONLY"); | |
336 | gMC->Gspos("BKHL", 4, "BACK", -.8 - 11.4 - 1.6 - 9.05, 0., -.6 - 4.4625, 0, "ONLY"); | |
337 | gMC->Gspos("BKHL", 5, "BACK", .8 + 11.4+ 1.6 + 9.05, 0., .6 + 8.925 + 1.2 + 4.4625, 0, "ONLY"); | |
338 | gMC->Gspos("BKHL", 6, "BACK", -.8 - 11.4 - 1.6 - 9.05, 0., .6 + 8.925 + 1.2 + 4.4625, 0, "ONLY"); | |
339 | gMC->Gspos("BKHL", 7, "BACK", .8 + 11.4 + 1.6 + 9.05, 0., -.6 - 8.925 - 1.2 - 4.4625, 0, "ONLY"); | |
340 | gMC->Gspos("BKHL", 8, "BACK", -.8 - 11.4 - 1.6 - 9.05, 0., -.6 - 8.925 - 1.2 - 4.4625, 0, "ONLY"); | |
341 | ||
342 | ||
343 | ||
344 | // --- Places the detectors defined with GSVOLU | |
345 | // Place material inside RICH | |
346 | gMC->Gspos("SRIC", 1, "RICH", 0.,0., 0., 0, "ONLY"); | |
347 | gMC->Gspos("AIR2", 1, "RICH", 66.3 + 1.2505, 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .6 - .05 - .376 -.5 - 3.35, 0., 0, "ONLY"); | |
348 | gMC->Gspos("AIR2", 2, "RICH", -66.3 - 1.2505, 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .6 - .05 - .376 -.5 - 3.35, 0., 0, "ONLY"); | |
349 | gMC->Gspos("AIR3", 1, "RICH", 0., 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .6 - .05 - .376 -.5 - 3.35, -68.35 - 1.25, 0, "ONLY"); | |
350 | gMC->Gspos("AIR3", 2, "RICH", 0., 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .6 - .05 - .376 -.5 - 3.35, 68.35 + 1.25, 0, "ONLY"); | |
351 | ||
352 | ||
353 | gMC->Gspos("ALUM", 1, "SRIC", 0., 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .6 - .05 - .376 -.025, 0., 0, "ONLY"); | |
354 | gMC->Gspos("HONE", 1, "SRIC", 0., 1.276- geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .6 - .05 - .188, 0., 0, "ONLY"); | |
355 | gMC->Gspos("ALUM", 2, "SRIC", 0., 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .6 - .025, 0., 0, "ONLY"); | |
356 | gMC->Gspos("FOOT", 1, "SRIC", 64.95, 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .3, 36.9, 0, "ONLY"); | |
357 | gMC->Gspos("FOOT", 2, "SRIC", 21.65, 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .3 , 36.9, 0, "ONLY"); | |
358 | gMC->Gspos("FOOT", 3, "SRIC", -21.65, 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .3, 36.9, 0, "ONLY"); | |
359 | gMC->Gspos("FOOT", 4, "SRIC", -64.95, 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .3, 36.9, 0, "ONLY"); | |
360 | gMC->Gspos("FOOT", 5, "SRIC", 64.95, 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .3, -36.9, 0, "ONLY"); | |
361 | gMC->Gspos("FOOT", 6, "SRIC", 21.65, 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .3, -36.9, 0, "ONLY"); | |
362 | gMC->Gspos("FOOT", 7, "SRIC", -21.65, 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .3, -36.9, 0, "ONLY"); | |
363 | gMC->Gspos("FOOT", 8, "SRIC", -64.95, 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .4 - .3, -36.9, 0, "ONLY"); | |
364 | gMC->Gspos("OQUA", 1, "SRIC", 0., 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()- .2, 0., 0, "ONLY"); | |
365 | ||
366 | // Supports placing | |
367 | ||
368 | // Methane supports | |
369 | gMC->Gspos("SMLG", 1, "SRIC", csi_width/2 + (66.3 - csi_width/2)/2, 1.276 + .25, 0., 0, "ONLY"); | |
370 | gMC->Gspos("SMLG", 2, "SRIC", - csi_width/2 - (66.3 - csi_width/2)/2, 1.276 + .25, 0., 0, "ONLY"); | |
371 | gMC->Gspos("SMSH", 1, "SRIC", 0., 1.276 + .25, csi_length/2 + (68.35 - csi_length/2)/2, 0, "ONLY"); | |
372 | gMC->Gspos("SMSH", 2, "SRIC", 0., 1.276 + .25, - csi_length/2 - (68.35 - csi_length/2)/2, 0, "ONLY"); | |
373 | ||
374 | //Freon supports | |
375 | ||
376 | Float_t supp_y = 1.276 - geometry->GetGapThickness()/2- geometry->GetQuartzThickness() -geometry->GetFreonThickness() - .2 + .3; //y position of freon supports | |
377 | ||
378 | gMC->Gspos("SFLG", 1, "SRIC", geometry->GetQuartzWidth()/2 + (66.3 - geometry->GetQuartzWidth()/2)/2, supp_y, 0., 0, "ONLY"); | |
379 | gMC->Gspos("SFLG", 2, "SRIC", - geometry->GetQuartzWidth()/2 - (66.3 - geometry->GetQuartzWidth()/2)/2, supp_y, 0., 0, "ONLY"); | |
380 | gMC->Gspos("SFSH", 1, "SRIC", 0., supp_y, geometry->GetQuartzLength()/2 + (68.35 - geometry->GetQuartzLength()/2)/2, 0, "ONLY"); | |
381 | gMC->Gspos("SFSH", 2, "SRIC", 0., supp_y, - geometry->GetQuartzLength()/2 - (68.35 - geometry->GetQuartzLength()/2)/2, 0, "ONLY"); | |
382 | ||
383 | AliMatrix(idrotm[1019], 0., 0., 90., 0., 90., 90.); | |
384 | ||
385 | ||
386 | Int_t nspacers = 30; | |
387 | ||
388 | for (i = 0; i < nspacers/3; i++) { | |
389 | zs = -11.6/2 + (TMath::Abs(nspacers/6) - i) * 12.2; | |
390 | gMC->Gspos("SPAC", i, "FRE1", 10.5, 0., zs, idrotm[1019], "ONLY"); //Original settings | |
391 | } | |
392 | ||
393 | for (i = nspacers/3; i < (nspacers*2)/3; i++) { | |
394 | zs = -11.6/2 + (nspacers/3 + TMath::Abs(nspacers/6) - i) * 12.2; | |
395 | gMC->Gspos("SPAC", i, "FRE1", 0, 0., zs, idrotm[1019], "ONLY"); //Original settings | |
396 | } | |
397 | ||
398 | for (i = (nspacers*2)/3; i < nspacers; ++i) { | |
399 | zs = -11.6/2 + ((nspacers*2)/3 + TMath::Abs(nspacers/6) - i) * 12.2; | |
400 | gMC->Gspos("SPAC", i, "FRE1", -10.5, 0., zs, idrotm[1019], "ONLY"); //Original settings | |
401 | } | |
402 | ||
403 | for (i = 0; i < nspacers/3; i++) { | |
404 | zs = -11.6/2 + (TMath::Abs(nspacers/6) - i) * 12.2; | |
405 | gMC->Gspos("SPAC", i, "FRE2", 10.5, 0., zs, idrotm[1019], "ONLY"); //Original settings | |
406 | } | |
407 | ||
408 | for (i = nspacers/3; i < (nspacers*2)/3; i++) { | |
409 | zs = -11.6/2 + (nspacers/3 + TMath::Abs(nspacers/6) - i) * 12.2; | |
410 | gMC->Gspos("SPAC", i, "FRE2", 0, 0., zs, idrotm[1019], "ONLY"); //Original settings | |
411 | } | |
412 | ||
413 | for (i = (nspacers*2)/3; i < nspacers; ++i) { | |
414 | zs = -11.6/2 + ((nspacers*2)/3 + TMath::Abs(nspacers/6) - i) * 12.2; | |
415 | gMC->Gspos("SPAC", i, "FRE2", -10.5, 0., zs, idrotm[1019], "ONLY"); //Original settings | |
416 | } | |
417 | ||
418 | ||
419 | gMC->Gspos("FRE1", 1, "OQF1", 0., 0., 0., 0, "ONLY"); | |
420 | gMC->Gspos("FRE2", 1, "OQF2", 0., 0., 0., 0, "ONLY"); | |
421 | gMC->Gspos("OQF1", 1, "SRIC", geometry->GetOuterFreonWidth()/2 + geometry->GetInnerFreonWidth()/2 + 2, 1.276 - geometry->GetGapThickness()/2- geometry->GetQuartzThickness() -geometry->GetFreonThickness()/2, 0., 0, "ONLY"); //Original settings (31.3) | |
422 | gMC->Gspos("OQF2", 2, "SRIC", 0., 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()/2, 0., 0, "ONLY"); //Original settings | |
423 | gMC->Gspos("OQF1", 3, "SRIC", - (geometry->GetOuterFreonWidth()/2 + geometry->GetInnerFreonWidth()/2) - 2, 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness() - geometry->GetFreonThickness()/2, 0., 0, "ONLY"); //Original settings (-31.3) | |
424 | gMC->Gspos("QUAR", 1, "SRIC", 0., 1.276 - geometry->GetGapThickness()/2 - geometry->GetQuartzThickness()/2, 0., 0, "ONLY"); | |
425 | gMC->Gspos("GAP ", 1, "META", 0., geometry->GetGapThickness()/2 - geometry->GetProximityGapThickness()/2 - 0.0001, 0., 0, "ONLY"); | |
426 | gMC->Gspos("META", 1, "SRIC", 0., 1.276, 0., 0, "ONLY"); | |
427 | gMC->Gspos("CSI ", 1, "SRIC", 0., 1.276 + geometry->GetGapThickness()/2 + .25, 0., 0, "ONLY"); | |
428 | printf("CSI pos: %f\n",1.276 + geometry->GetGapThickness()/2 + .25); | |
429 | ||
430 | // Wire support placing | |
431 | ||
432 | gMC->Gspos("WSG2", 1, "GAP ", 0., geometry->GetProximityGapThickness()/2 - .1, 0., 0, "ONLY"); | |
433 | gMC->Gspos("WSG1", 1, "CSI ", 0., 0., 0., 0, "ONLY"); | |
434 | gMC->Gspos("WSMe", 1, "SRIC ", 0., 1.276 + geometry->GetGapThickness()/2 + .5 + 1.05, 0., 0, "ONLY"); | |
435 | ||
436 | // Backplane placing | |
437 | ||
438 | gMC->Gspos("BACK", 1, "SRIC ", -33.15, 1.276 + geometry->GetGapThickness()/2 + .5 + 2.1 + 2, 43.3, 0, "ONLY"); | |
439 | gMC->Gspos("BACK", 2, "SRIC ", 33.15, 1.276 + geometry->GetGapThickness()/2 + .5 + 2.1 + 2 , 43.3, 0, "ONLY"); | |
440 | gMC->Gspos("BACK", 3, "SRIC ", -33.15, 1.276 + geometry->GetGapThickness()/2 + .5 + 2.1 + 2, 0., 0, "ONLY"); | |
441 | gMC->Gspos("BACK", 4, "SRIC ", 33.15, 1.276 + geometry->GetGapThickness()/2 + .5 + 2.1 + 2, 0., 0, "ONLY"); | |
442 | gMC->Gspos("BACK", 5, "SRIC ", 33.15, 1.276 + geometry->GetGapThickness()/2 + .5 + 2.1 + 2, -43.3, 0, "ONLY"); | |
443 | gMC->Gspos("BACK", 6, "SRIC ", -33.15, 1.276 + geometry->GetGapThickness()/2 + .5 + 2.1 + 2, -43.3, 0, "ONLY"); | |
444 | ||
445 | // PCB placing | |
446 | ||
447 | gMC->Gspos("PCB ", 1, "SRIC ", 0., 1.276 + geometry->GetGapThickness()/2 + .5 + 1.05, csi_width/4 + .5025 + 2.5, 0, "ONLY"); | |
448 | gMC->Gspos("PCB ", 2, "SRIC ", 0., 1.276 + geometry->GetGapThickness()/2 + .5 + 1.05, -csi_width/4 - .5025 - 2.5, 0, "ONLY"); | |
449 | ||
450 | // Place chambers into mother volume ALIC | |
451 | ||
9d6c0962 | 452 | Double_t dOffset = geometry->GetOffset() - geometry->GetGapThickness()/2; // distance from center of mother volume ALIC to methane |
453 | ||
454 | Double_t dAlpha = geometry->GetAlphaAngle(); // angle between centers of chambers - y-z plane | |
455 | Double_t dAlphaRad = dAlpha*kDegrad; | |
456 | ||
457 | Double_t dBeta = geometry->GetBetaAngle(); // angle between center of chambers - y-x plane | |
458 | Double_t dBetaRad = dBeta*kDegrad; | |
459 | ||
460 | Double_t dRotAngle = geometry->GetRotationAngle(); // the whole RICH is to be rotated in x-y plane + means clockwise rotation | |
461 | Double_t dRotAngleRad = dRotAngle*kDegrad; | |
9d6c0962 | 462 | |
4b7382d1 | 463 | |
9d6c0962 | 464 | TRotMatrix *pRotMatrix; // tmp pointer |
465 | ||
466 | TVector3 vector(0,dOffset,0); // Position of chamber 2 without rotation | |
467 | ||
9123a941 | 468 | // Chamber 0 standalone (no other chambers in this row) |
4b7382d1 | 469 | pRotMatrix = new TRotMatrix("rot993","rot993", 0., 0., 0.,0.,0.,0.); |
470 | const Double_t* r = pRotMatrix->SetAngles(90., 0., 90.-dAlpha , 90., dAlpha, -90.); | |
471 | Double_t* rr = RotateXY(r, -dRotAngleRad); | |
472 | AliMatrix(idrotm[1000], rr[0], rr[1], rr[2], rr[3], rr[4], rr[5]); | |
473 | pRotMatrix->SetAngles(rr[0], rr[1], rr[2], rr[3], rr[4], rr[5]); | |
9123a941 | 474 | |
9d6c0962 | 475 | vector.SetXYZ(0,dOffset,0); vector.RotateX(dAlphaRad); |
476 | vector.RotateZ(-dRotAngleRad); | |
9123a941 | 477 | |
9d6c0962 | 478 | gMC->Gspos("RICH",1,"ALIC",vector.X(),vector.Y(),vector.Z(),idrotm[1000], "ONLY"); |
479 | Chamber(0).SetChamberTransform(vector.X(),vector.Y(),vector.Z(),pRotMatrix); | |
e1826cc2 | 480 | if(GetDebug()) Info("CreateGeometry 0","%8.3f %8.3f %8.3f %8.3f %8.3f %8.3f",rr[0],rr[1],rr[2],rr[3],rr[4],rr[5]); |
481 | if(GetDebug()) Info("CreateGeometry 0","x=%8.3f y=%8.3f z=%8.3f",vector.X(),vector.Y(),vector.Z()); | |
9123a941 | 482 | // Chamber 1 |
4b7382d1 | 483 | pRotMatrix = new TRotMatrix("rot994","rot994", 0., 0., 0.,0.,0.,0.); |
484 | r = pRotMatrix->SetAngles(90., -dBeta, 90., 90.-dBeta, 0., 0.); | |
485 | rr = RotateXY(r, -dRotAngleRad); | |
486 | AliMatrix(idrotm[1001], rr[0], rr[1], rr[2], rr[3], rr[4], rr[5]); | |
487 | pRotMatrix->SetAngles(rr[0], rr[1], rr[2], rr[3], rr[4], rr[5]); | |
9d6c0962 | 488 | vector.SetXYZ(0,dOffset,0); vector.RotateZ(-dBetaRad); |
489 | vector.RotateZ(-dRotAngleRad); | |
9123a941 | 490 | |
9d6c0962 | 491 | gMC->Gspos("RICH",2,"ALIC",vector.X(),vector.Y(),vector.Z(),idrotm[1001], "ONLY"); |
492 | Chamber(1).SetChamberTransform(vector.X(),vector.Y(),vector.Z(),pRotMatrix); | |
e1826cc2 | 493 | if(GetDebug()) Info("CreateGeometry 1","%8.3f %8.3f %8.3f %8.3f %8.3f %8.3f",rr[0],rr[1],rr[2],rr[3],rr[4],rr[5]); |
494 | if(GetDebug()) Info("CreateGeometry 1","x=%8.3f y=%8.3f z=%8.3f",vector.X(),vector.Y(),vector.Z()); | |
9123a941 | 495 | // Chamber 2 the top one with no Alpha-Beta rotation |
4b7382d1 | 496 | pRotMatrix = new TRotMatrix("rot995","rot995", 0., 0., 0.,0.,0.,0.); |
497 | r = pRotMatrix->SetAngles(90., 0., 90., 90., 0., 0.); | |
498 | rr = RotateXY(r, -dRotAngleRad); | |
499 | AliMatrix(idrotm[1002], rr[0], rr[1], rr[2], rr[3], rr[4], rr[5]); | |
500 | pRotMatrix->SetAngles(rr[0], rr[1], rr[2], rr[3], rr[4], rr[5]); | |
9d6c0962 | 501 | vector.SetXYZ(0,dOffset,0); |
502 | vector.RotateZ(-dRotAngleRad); | |
9d6c0962 | 503 | gMC->Gspos("RICH",3,"ALIC",vector.X(),vector.Y(),vector.Z(),idrotm[1002], "ONLY"); |
504 | Chamber(2).SetChamberTransform(vector.X(),vector.Y(),vector.Z(),pRotMatrix); | |
e1826cc2 | 505 | if(GetDebug()) Info("CreateGeometry 2","%8.3f %8.3f %8.3f %8.3f %8.3f %8.3f",rr[0],rr[1],rr[2],rr[3],rr[4],rr[5]); |
506 | if(GetDebug()) Info("CreateGeometry 2","x=%8.3f y=%8.3f z=%8.3f",vector.X(),vector.Y(),vector.Z()); | |
9123a941 | 507 | // Chamber 3 |
4b7382d1 | 508 | pRotMatrix = new TRotMatrix("rot996","rot996", 0., 0., 0.,0.,0.,0.); |
509 | r = pRotMatrix->SetAngles(90., dBeta, 90., 90.+dBeta, 0., 0.); | |
510 | rr = RotateXY(r, -dRotAngleRad); | |
511 | AliMatrix(idrotm[1003], rr[0], rr[1], rr[2], rr[3], rr[4], rr[5]); | |
512 | pRotMatrix->SetAngles(rr[0], rr[1], rr[2], rr[3], rr[4], rr[5]); | |
9d6c0962 | 513 | vector.SetXYZ(0,dOffset,0); vector.RotateZ(dBetaRad); |
514 | vector.RotateZ(-dRotAngleRad); | |
9123a941 | 515 | |
9d6c0962 | 516 | gMC->Gspos("RICH",4,"ALIC",vector.X(),vector.Y(),vector.Z(),idrotm[1003], "ONLY"); |
517 | Chamber(3).SetChamberTransform(vector.X(),vector.Y(),vector.Z(),pRotMatrix); | |
e1826cc2 | 518 | if(GetDebug()) Info("CreateGeometry 3","%8.3f %8.3f %8.3f %8.3f %8.3f %8.3f",rr[0],rr[1],rr[2],rr[3],rr[4],rr[5]); |
519 | if(GetDebug()) Info("CreateGeometry 3","x=%8.3f y=%8.3f z=%8.3f",vector.X(),vector.Y(),vector.Z()); | |
9123a941 | 520 | // Chamber 4 |
4b7382d1 | 521 | pRotMatrix = new TRotMatrix("rot997","rot997", 0., 0., 0.,0.,0.,0.); |
522 | r = pRotMatrix->SetAngles(90., 360.-dBeta, 108.2, 90.-dBeta, 18.2, 90.-dBeta); | |
523 | rr = RotateXY(r, -dRotAngleRad); | |
524 | AliMatrix(idrotm[1004], rr[0], rr[1], rr[2], rr[3], rr[4], rr[5]); | |
525 | pRotMatrix->SetAngles(rr[0], rr[1], rr[2], rr[3], rr[4], rr[5]); | |
9d6c0962 | 526 | vector.SetXYZ(0,dOffset,0); vector.RotateZ(-dBetaRad); vector.RotateX(-dAlphaRad); |
527 | vector.RotateZ(-dRotAngleRad); | |
9123a941 | 528 | |
9d6c0962 | 529 | gMC->Gspos("RICH",5,"ALIC",vector.X(),vector.Y(),vector.Z(),idrotm[1004], "ONLY"); |
530 | Chamber(4).SetChamberTransform(vector.X(),vector.Y(),vector.Z(),pRotMatrix); | |
e1826cc2 | 531 | if(GetDebug()) Info("CreateGeometry 4","%8.3f %8.3f %8.3f %8.3f %8.3f %8.3f",rr[0],rr[1],rr[2],rr[3],rr[4],rr[5]); |
532 | if(GetDebug()) Info("CreateGeometry 4","x=%8.3f y=%8.3f z=%8.3f",vector.X(),vector.Y(),vector.Z()); | |
9123a941 | 533 | // Chamber 5 |
4b7382d1 | 534 | pRotMatrix = new TRotMatrix("rot998","rot998", 0., 0., 0.,0.,0.,0.); |
535 | r = pRotMatrix->SetAngles(90., 0., 90.+dAlpha, 90., dAlpha, 90.); | |
536 | rr = RotateXY(r, -dRotAngleRad); | |
537 | AliMatrix(idrotm[1005], rr[0], rr[1], rr[2], rr[3], rr[4], rr[5]); | |
538 | pRotMatrix->SetAngles(rr[0], rr[1], rr[2], rr[3], rr[4], rr[5]); | |
9d6c0962 | 539 | vector.SetXYZ(0,dOffset,0); vector.RotateX(-dAlphaRad); |
540 | vector.RotateZ(-dRotAngleRad); | |
541 | ||
542 | gMC->Gspos("RICH",6,"ALIC",vector.X(),vector.Y(),vector.Z(),idrotm[1005], "ONLY"); | |
543 | Chamber(5).SetChamberTransform(vector.X(),vector.Y(),vector.Z(),pRotMatrix); | |
e1826cc2 | 544 | if(GetDebug()) Info("CreateGeometry 5","%8.3f %8.3f %8.3f %8.3f %8.3f %8.3f",rr[0],rr[1],rr[2],rr[3],rr[4],rr[5]); |
545 | if(GetDebug()) Info("CreateGeometry 5","x=%8.3f y=%8.3f z=%8.3f",vector.X(),vector.Y(),vector.Z()); | |
4b7382d1 | 546 | // Chamber 6 |
547 | pRotMatrix = new TRotMatrix("rot999","rot999", 0., 0., 0.,0.,0.,0.); | |
548 | r = pRotMatrix->SetAngles(90., dBeta, 108.2, 90.+dBeta, 18.2, 90.+dBeta); | |
549 | rr = RotateXY(r, -dRotAngleRad); | |
550 | AliMatrix(idrotm[1006], rr[0], rr[1], rr[2], rr[3], rr[4], rr[5]); | |
551 | pRotMatrix->SetAngles(rr[0], rr[1], rr[2], rr[3], rr[4], rr[5]); | |
698656f6 | 552 | vector.SetXYZ(0,dOffset,0); vector.RotateZ(dBetaRad); vector.RotateX(-dAlphaRad); |
9d6c0962 | 553 | vector.RotateZ(-dRotAngleRad); |
554 | ||
555 | gMC->Gspos("RICH",7,"ALIC",vector.X(),vector.Y(),vector.Z(),idrotm[1006], "ONLY"); | |
556 | Chamber(6).SetChamberTransform(vector.X(),vector.Y(),vector.Z(),pRotMatrix); | |
e1826cc2 | 557 | if(GetDebug()) Info("CreateGeometry 6","%8.3f %8.3f %8.3f %8.3f %8.3f %8.3f",rr[0],rr[1],rr[2],rr[3],rr[4],rr[5]); |
558 | if(GetDebug()) Info("CreateGeometry 6","x=%8.3f y=%8.3f z=%8.3f",vector.X(),vector.Y(),vector.Z()); | |
9123a941 | 559 | |
560 | }//void AliRICHv3::CreateGeometry() | |
e1826cc2 | 561 | //______________________________________________________________________________ |
9123a941 | 562 | void AliRICHv3::Init() |
cd1cf354 | 563 | {//Makes nothing for a while |
564 | if(GetDebug())Info("Init","Start."); | |
565 | if(GetDebug())Info("Init","Stop."); | |
9123a941 | 566 | } |
e1826cc2 | 567 | //______________________________________________________________________________ |
4b7382d1 | 568 | Double_t* AliRICHv3::RotateXY(const Double_t* r, Double_t a) |
569 | { | |
570 | // Rotatation in xy-plane | |
571 | // by angle a | |
572 | // The resulting rotation matrix is given back in the G3 notation. | |
573 | Double_t* rr = new Double_t[6]; | |
574 | Double_t m[9]; | |
575 | Int_t i,j,k; | |
576 | ||
577 | for (i = 0; i < 3; i++) { | |
578 | j = 3*i; | |
579 | m[j] = r[j] * TMath::Cos(a) - r[j+1] * TMath::Sin(a); | |
580 | m[j+1] = r[j] * TMath::Sin(a) + r[j+1] * TMath::Cos(a); | |
581 | m[j+2] = r[j+2]; | |
582 | } | |
583 | ||
584 | for (i = 0; i < 3; i++) { | |
585 | j = 3*i; | |
586 | k = 2*i; | |
587 | rr[k] = TMath::ACos(m[j+2]) * kRaddeg; | |
588 | rr[k+1] = TMath::ATan2(m[j+1], m[j]) * kRaddeg; | |
589 | } | |
590 | return rr; | |
e1826cc2 | 591 | }//Double_t* AliRICHv3::RotateXY(const Double_t* r, Double_t a) |
cd1cf354 | 592 | //______________________________________________________________________________ |
593 | void AliRICHv3::StepManager() | |
594 | {//Full Step Manager | |
595 | ||
596 | Int_t copy, id; | |
597 | static Int_t idvol; | |
598 | static Int_t vol[2]; | |
599 | Int_t ipart; | |
600 | static Float_t hits[22]; | |
601 | static Float_t ckovData[19]; | |
602 | TLorentzVector position; | |
603 | TLorentzVector momentum; | |
604 | Float_t pos[3]; | |
605 | Float_t mom[4]; | |
606 | Float_t localPos[3]; | |
607 | Float_t localMom[4]; | |
608 | Float_t localTheta,localPhi; | |
609 | Float_t theta,phi; | |
610 | Float_t destep, step; | |
611 | Double_t ranf[2]; | |
cd1cf354 | 612 | Float_t coscerenkov; |
613 | static Float_t eloss, xhit, yhit, tlength; | |
614 | const Float_t kBig=1.e10; | |
615 | ||
616 | TClonesArray &lhits = *fHits; | |
642f15cf | 617 | TParticle *current = (TParticle*)(*gAlice->Particles())[gAlice->GetCurrentTrackNumber()]; |
cd1cf354 | 618 | |
619 | //if (current->Energy()>1) | |
620 | //{ | |
621 | ||
622 | // Only gas gap inside chamber | |
623 | // Tag chambers and record hits when track enters | |
624 | ||
625 | ||
626 | id=gMC->CurrentVolID(copy); | |
627 | idvol = copy-1; | |
628 | Float_t cherenkovLoss=0; | |
642f15cf | 629 | //gAlice->KeepTrack(gAlice->GetCurrentTrackNumber()); |
cd1cf354 | 630 | |
631 | gMC->TrackPosition(position); | |
632 | pos[0]=position(0); | |
633 | pos[1]=position(1); | |
634 | pos[2]=position(2); | |
635 | //bzero((char *)ckovData,sizeof(ckovData)*19); | |
636 | ckovData[1] = pos[0]; // X-position for hit | |
637 | ckovData[2] = pos[1]; // Y-position for hit | |
638 | ckovData[3] = pos[2]; // Z-position for hit | |
639 | ckovData[6] = 0; // dummy track length | |
642f15cf | 640 | //ckovData[11] = gAlice->GetCurrentTrackNumber(); |
cd1cf354 | 641 | |
642f15cf | 642 | //printf("\n+++++++++++\nTrack: %d\n++++++++++++\n",gAlice->GetCurrentTrackNumber()); |
cd1cf354 | 643 | |
644 | //AliRICH *RICH = (AliRICH *) gAlice->GetDetector("RICH"); | |
645 | ||
646 | /********************Store production parameters for Cerenkov photons************************/ | |
647 | //is it a Cerenkov photon? | |
648 | if (gMC->TrackPid() == 50000050) { | |
649 | ||
650 | //if (gMC->VolId("GAP ")==gMC->CurrentVolID(copy)) | |
651 | //{ | |
652 | Float_t ckovEnergy = current->Energy(); | |
653 | //energy interval for tracking | |
654 | if (ckovEnergy > 5.6e-09 && ckovEnergy < 7.8e-09 ) | |
655 | //if (ckovEnergy > 0) | |
656 | { | |
657 | if (gMC->IsTrackEntering()){ //is track entering? | |
658 | //printf("Track entered (1)\n"); | |
659 | if (gMC->VolId("FRE1")==gMC->CurrentVolID(copy) || gMC->VolId("FRE2")==gMC->CurrentVolID(copy)) | |
660 | { //is it in freo? | |
661 | if (gMC->IsNewTrack()){ //is it the first step? | |
662 | //printf("I'm in!\n"); | |
663 | Int_t mother = current->GetFirstMother(); | |
664 | ||
665 | //printf("Second Mother:%d\n",current->GetSecondMother()); | |
666 | ||
667 | ckovData[10] = mother; | |
642f15cf | 668 | ckovData[11] = gAlice->GetCurrentTrackNumber(); |
cd1cf354 | 669 | ckovData[12] = 1; //Media where photon was produced 1->Freon, 2->Quarz |
670 | //printf("Produced in FREO\n"); | |
671 | fCkovNumber++; | |
672 | fFreonProd=1; | |
673 | //printf("Index: %d\n",fCkovNumber); | |
674 | } //first step question | |
675 | } //freo question | |
676 | ||
677 | if (gMC->IsNewTrack()){ //is it first step? | |
678 | if (gMC->VolId("QUAR")==gMC->CurrentVolID(copy)) //is it in quarz? | |
679 | { | |
680 | ckovData[12] = 2; | |
681 | //printf("Produced in QUAR\n"); | |
682 | } //quarz question | |
683 | } //first step question | |
684 | ||
685 | //printf("Before %d\n",fFreonProd); | |
686 | } //track entering question | |
687 | ||
688 | if (ckovData[12] == 1) //was it produced in Freon? | |
689 | //if (fFreonProd == 1) | |
690 | { | |
691 | if (gMC->IsTrackEntering()){ //is track entering? | |
692 | //printf("Track entered (2)\n"); | |
693 | //printf("Current volume (should be META): %s\n",gMC->CurrentVolName()); | |
694 | //printf("VolId: %d, CurrentVolID: %d\n",gMC->VolId("META"),gMC->CurrentVolID(copy)); | |
695 | if (gMC->VolId("META")==gMC->CurrentVolID(copy)) //is it in gap? | |
696 | { | |
697 | //printf("Got in META\n"); | |
698 | gMC->TrackMomentum(momentum); | |
699 | mom[0]=momentum(0); | |
700 | mom[1]=momentum(1); | |
701 | mom[2]=momentum(2); | |
702 | mom[3]=momentum(3); | |
703 | ||
704 | gMC->Gmtod(mom,localMom,2); | |
705 | Float_t cophi = TMath::Cos(TMath::ATan2(localMom[0], localMom[1])); | |
706 | Float_t t = (1. - .025 / cophi) * (1. - .05 / cophi); | |
707 | /**************** Photons lost in second grid have to be calculated by hand************/ | |
708 | gMC->GetRandom()->RndmArray(1,ranf); | |
709 | if (ranf[0] > t) { | |
710 | gMC->StopTrack(); | |
711 | ckovData[13] = 5; | |
642f15cf | 712 | AddCerenkov(gAlice->GetCurrentTrackNumber(),vol,ckovData); |
cd1cf354 | 713 | //printf("Added One (1)!\n"); |
714 | //printf("Lost one in grid\n"); | |
715 | } | |
716 | /**********************************************************************************/ | |
717 | } //gap | |
718 | ||
719 | //printf("Current volume (should be CSI) (1): %s\n",gMC->CurrentVolName()); | |
720 | //printf("VolId: %d, CurrentVolID: %d\n",gMC->VolId("CSI "),gMC->CurrentVolID(copy)); | |
721 | if (gMC->VolId("CSI ")==gMC->CurrentVolID(copy)) //is it in csi? | |
722 | { | |
723 | //printf("Got in CSI\n"); | |
724 | gMC->TrackMomentum(momentum); | |
725 | mom[0]=momentum(0); | |
726 | mom[1]=momentum(1); | |
727 | mom[2]=momentum(2); | |
728 | mom[3]=momentum(3); | |
729 | ||
730 | gMC->Gmtod(mom,localMom,2); | |
731 | /********* Photons lost by Fresnel reflection have to be calculated by hand********/ | |
732 | /***********************Cerenkov phtons (always polarised)*************************/ | |
733 | Double_t localTc = localMom[0]*localMom[0]+localMom[2]*localMom[2]; | |
734 | Double_t localRt = TMath::Sqrt(localTc); | |
735 | localTheta = Float_t(TMath::ATan2(localRt,Double_t(localMom[1]))); | |
736 | Double_t cotheta = TMath::Abs(cos(localTheta)); | |
737 | Float_t t = Fresnel(ckovEnergy*1e9,cotheta,1); | |
738 | gMC->GetRandom()->RndmArray(1,ranf); | |
739 | if (ranf[0] < t) { | |
740 | gMC->StopTrack(); | |
741 | ckovData[13] = 6; | |
642f15cf | 742 | AddCerenkov(gAlice->GetCurrentTrackNumber(),vol,ckovData); |
cd1cf354 | 743 | |
744 | //printf("Added One (2)!\n"); | |
745 | //printf("Lost by Fresnel\n"); | |
746 | } | |
747 | /**********************************************************************************/ | |
748 | } | |
749 | } //track entering? | |
750 | ||
751 | ||
752 | /********************Evaluation of losses************************/ | |
753 | /******************still in the old fashion**********************/ | |
754 | ||
755 | TArrayI procs; | |
756 | Int_t i1 = gMC->StepProcesses(procs); //number of physics mechanisms acting on the particle | |
757 | for (Int_t i = 0; i < i1; ++i) { | |
758 | // Reflection loss | |
759 | if (procs[i] == kPLightReflection) { //was it reflected | |
760 | ckovData[13]=10; | |
761 | if (gMC->VolId("FRE1")==gMC->CurrentVolID(copy) || gMC->VolId("FRE2")==gMC->CurrentVolID(copy)) | |
762 | ckovData[13]=1; | |
763 | if (gMC->CurrentVolID(copy) == gMC->VolId("QUAR")) | |
764 | ckovData[13]=2; | |
765 | //gMC->StopTrack(); | |
642f15cf | 766 | //AddCerenkov(gAlice->GetCurrentTrackNumber(),vol,ckovData); |
cd1cf354 | 767 | } //reflection question |
768 | ||
769 | // Absorption loss | |
770 | else if (procs[i] == kPLightAbsorption) { //was it absorbed? | |
771 | //printf("Got in absorption\n"); | |
772 | ckovData[13]=20; | |
773 | if (gMC->VolId("FRE1")==gMC->CurrentVolID(copy) || gMC->VolId("FRE2")==gMC->CurrentVolID(copy)) | |
774 | ckovData[13]=11; | |
775 | if (gMC->CurrentVolID(copy) == gMC->VolId("QUAR")) | |
776 | ckovData[13]=12; | |
777 | if (gMC->CurrentVolID(copy) == gMC->VolId("META")) | |
778 | ckovData[13]=13; | |
779 | if (gMC->CurrentVolID(copy) == gMC->VolId("GAP ")) | |
780 | ckovData[13]=13; | |
781 | ||
782 | if (gMC->CurrentVolID(copy) == gMC->VolId("SRIC")) | |
783 | ckovData[13]=15; | |
784 | ||
785 | // CsI inefficiency | |
786 | if (gMC->CurrentVolID(copy) == gMC->VolId("CSI ")) { | |
787 | ckovData[13]=16; | |
788 | } | |
789 | gMC->StopTrack(); | |
642f15cf | 790 | AddCerenkov(gAlice->GetCurrentTrackNumber(),vol,ckovData); |
cd1cf354 | 791 | //printf("Added One (3)!\n"); |
792 | //printf("Added cerenkov %d\n",fCkovNumber); | |
793 | } //absorption question | |
794 | ||
795 | ||
796 | // Photon goes out of tracking scope | |
797 | else if (procs[i] == kPStop) { //is it below energy treshold? | |
798 | ckovData[13]=21; | |
799 | gMC->StopTrack(); | |
642f15cf | 800 | AddCerenkov(gAlice->GetCurrentTrackNumber(),vol,ckovData); |
cd1cf354 | 801 | //printf("Added One (4)!\n"); |
802 | } // energy treshold question | |
803 | } //number of mechanisms cycle | |
804 | /**********************End of evaluation************************/ | |
805 | } //freon production question | |
806 | } //energy interval question | |
807 | //}//inside the proximity gap question | |
808 | } //cerenkov photon question | |
809 | ||
810 | /**************************************End of Production Parameters Storing*********************/ | |
811 | ||
812 | ||
813 | /*******************************Treat photons that hit the CsI (Ckovs and Feedbacks)************/ | |
814 | ||
815 | if (gMC->TrackPid() == 50000050 || gMC->TrackPid() == 50000051) { | |
816 | //printf("Cerenkov\n"); | |
817 | ||
818 | //if (gMC->TrackPid() == 50000051) | |
819 | //printf("Tracking a feedback\n"); | |
820 | ||
821 | if (gMC->VolId("CSI ")==gMC->CurrentVolID(copy)) | |
822 | { | |
823 | //printf("Current volume (should be CSI) (2): %s\n",gMC->CurrentVolName()); | |
824 | //printf("VolId: %d, CurrentVolID: %d\n",gMC->VolId("CSI "),gMC->CurrentVolID(copy)); | |
825 | //printf("Got in CSI\n"); | |
826 | //printf("Tracking a %d\n",gMC->TrackPid()); | |
827 | if (gMC->Edep() > 0.){ | |
828 | gMC->TrackPosition(position); | |
829 | gMC->TrackMomentum(momentum); | |
830 | pos[0]=position(0); | |
831 | pos[1]=position(1); | |
832 | pos[2]=position(2); | |
833 | mom[0]=momentum(0); | |
834 | mom[1]=momentum(1); | |
835 | mom[2]=momentum(2); | |
836 | mom[3]=momentum(3); | |
837 | Double_t tc = mom[0]*mom[0]+mom[1]*mom[1]; | |
838 | Double_t rt = TMath::Sqrt(tc); | |
839 | theta = Float_t(TMath::ATan2(rt,Double_t(mom[2])))*kRaddeg; | |
840 | phi = Float_t(TMath::ATan2(Double_t(mom[1]),Double_t(mom[0])))*kRaddeg; | |
841 | ||
842 | gMC->CurrentVolOffID(2,copy); | |
843 | vol[0]=copy; | |
844 | idvol=vol[0]-1; | |
845 | ||
846 | ||
847 | gMC->Gmtod(pos,localPos,1); | |
848 | ||
849 | //Chamber(idvol).GlobaltoLocal(pos,localPos); | |
850 | ||
851 | gMC->Gmtod(mom,localMom,2); | |
852 | ||
853 | //Chamber(idvol).GlobaltoLocal(mom,localMom); | |
854 | ||
855 | gMC->CurrentVolOffID(2,copy); | |
856 | vol[0]=copy; | |
857 | idvol=vol[0]-1; | |
858 | ||
859 | //Int_t sector=((AliRICHChamber*) (*fChambers)[idvol]) | |
860 | //->Sector(localPos[0], localPos[2]); | |
861 | //printf("Sector:%d\n",sector); | |
862 | ||
863 | /*if (gMC->TrackPid() == 50000051){ | |
864 | fFeedbacks++; | |
865 | printf("Feedbacks:%d\n",fFeedbacks); | |
866 | }*/ | |
867 | ||
868 | //PH ((AliRICHChamber*) (*fChambers)[idvol]) | |
869 | ((AliRICHChamber*)fChambers->At(idvol)) | |
870 | ->SigGenInit(localPos[0], localPos[2], localPos[1]); | |
871 | if(idvol<kNCH) { | |
872 | ckovData[0] = gMC->TrackPid(); // particle type | |
873 | ckovData[1] = pos[0]; // X-position for hit | |
874 | ckovData[2] = pos[1]; // Y-position for hit | |
875 | ckovData[3] = pos[2]; // Z-position for hit | |
876 | ckovData[4] = theta; // theta angle of incidence | |
877 | ckovData[5] = phi; // phi angle of incidence | |
853634d3 | 878 | ckovData[8] = (Float_t) fNsdigits; // first sdigit |
cd1cf354 | 879 | ckovData[9] = -1; // last pad hit |
880 | ckovData[13] = 4; // photon was detected | |
881 | ckovData[14] = mom[0]; | |
882 | ckovData[15] = mom[1]; | |
883 | ckovData[16] = mom[2]; | |
884 | ||
885 | destep = gMC->Edep(); | |
886 | gMC->SetMaxStep(kBig); | |
887 | cherenkovLoss += destep; | |
888 | ckovData[7]=cherenkovLoss; | |
889 | ||
c021cb15 | 890 | ckovData[17] = Hits2SDigits(localPos[0],localPos[2],cherenkovLoss,idvol,kPhoton);//for photons in CsI |
cd1cf354 | 891 | |
853634d3 | 892 | if (fNsdigits > (Int_t)ckovData[8]) { |
cd1cf354 | 893 | ckovData[8]= ckovData[8]+1; |
853634d3 | 894 | ckovData[9]= (Float_t) fNsdigits; |
cd1cf354 | 895 | } |
896 | ||
cd1cf354 | 897 | |
898 | //TClonesArray *Hits = RICH->Hits(); | |
853634d3 | 899 | AliRICHhit *mipHit = (AliRICHhit*) (fHits->UncheckedAt(0)); |
cd1cf354 | 900 | if (mipHit) |
901 | { | |
902 | mom[0] = current->Px(); | |
903 | mom[1] = current->Py(); | |
904 | mom[2] = current->Pz(); | |
905 | Float_t mipPx = mipHit->MomX(); | |
906 | Float_t mipPy = mipHit->MomY(); | |
907 | Float_t mipPz = mipHit->MomZ(); | |
908 | ||
909 | Float_t r = mom[0]*mom[0] + mom[1]*mom[1] + mom[2]*mom[2]; | |
910 | Float_t rt = TMath::Sqrt(r); | |
911 | Float_t mipR = mipPx*mipPx + mipPy*mipPy + mipPz*mipPz; | |
912 | Float_t mipRt = TMath::Sqrt(mipR); | |
913 | if ((rt*mipRt) > 0) | |
914 | { | |
915 | coscerenkov = (mom[0]*mipPx + mom[1]*mipPy + mom[2]*mipPz)/(rt*mipRt); | |
916 | } | |
917 | else | |
918 | { | |
919 | coscerenkov = 0; | |
920 | } | |
921 | Float_t cherenkov = TMath::ACos(coscerenkov); | |
922 | ckovData[18]=cherenkov; | |
923 | } | |
924 | //if (sector != -1) | |
925 | //{ | |
642f15cf | 926 | AddHit(gAlice->GetCurrentTrackNumber(),vol,ckovData); |
927 | AddCerenkov(gAlice->GetCurrentTrackNumber(),vol,ckovData); | |
cd1cf354 | 928 | //printf("Added One (5)!\n"); |
929 | //} | |
930 | } | |
931 | } | |
932 | } | |
933 | } | |
934 | ||
935 | /***********************************************End of photon hits*********************************************/ | |
936 | ||
937 | ||
938 | /**********************************************Charged particles treatment*************************************/ | |
939 | ||
dfa42ff8 | 940 | else if (gMC->TrackCharge()){ |
cd1cf354 | 941 | //If MIP |
942 | /*if (gMC->IsTrackEntering()) | |
943 | { | |
944 | hits[13]=20;//is track entering? | |
945 | }*/ | |
946 | if (gMC->VolId("FRE1")==gMC->CurrentVolID(copy) || gMC->VolId("FRE2")==gMC->CurrentVolID(copy)) | |
947 | { | |
948 | gMC->TrackMomentum(momentum); | |
949 | mom[0]=momentum(0); | |
950 | mom[1]=momentum(1); | |
951 | mom[2]=momentum(2); | |
952 | mom[3]=momentum(3); | |
953 | hits [19] = mom[0]; | |
954 | hits [20] = mom[1]; | |
955 | hits [21] = mom[2]; | |
956 | fFreonProd=1; | |
957 | } | |
958 | ||
dfa42ff8 | 959 | if (gMC->VolId("GAP ")== gMC->CurrentVolID(copy)) {//is in GAP? |
cd1cf354 | 960 | // Get current particle id (ipart), track position (pos) and momentum (mom) |
961 | ||
962 | gMC->CurrentVolOffID(3,copy); | |
963 | vol[0]=copy; | |
964 | idvol=vol[0]-1; | |
965 | ||
966 | //Int_t sector=((AliRICHChamber*) (*fChambers)[idvol]) | |
967 | //->Sector(localPos[0], localPos[2]); | |
968 | //printf("Sector:%d\n",sector); | |
969 | ||
970 | gMC->TrackPosition(position); | |
971 | gMC->TrackMomentum(momentum); | |
972 | pos[0]=position(0); | |
973 | pos[1]=position(1); | |
974 | pos[2]=position(2); | |
975 | mom[0]=momentum(0); | |
976 | mom[1]=momentum(1); | |
977 | mom[2]=momentum(2); | |
978 | mom[3]=momentum(3); | |
979 | ||
980 | gMC->Gmtod(pos,localPos,1); | |
981 | ||
982 | //Chamber(idvol).GlobaltoLocal(pos,localPos); | |
983 | ||
984 | gMC->Gmtod(mom,localMom,2); | |
985 | ||
986 | //Chamber(idvol).GlobaltoLocal(mom,localMom); | |
987 | ||
988 | ipart = gMC->TrackPid(); | |
989 | // | |
990 | // momentum loss and steplength in last step | |
991 | destep = gMC->Edep(); | |
992 | step = gMC->TrackStep(); | |
993 | ||
994 | // | |
995 | // record hits when track enters ... | |
996 | if( gMC->IsTrackEntering()) { | |
997 | // gMC->SetMaxStep(fMaxStepGas); | |
998 | Double_t tc = mom[0]*mom[0]+mom[1]*mom[1]; | |
999 | Double_t rt = TMath::Sqrt(tc); | |
1000 | theta = Float_t(TMath::ATan2(rt,Double_t(mom[2])))*kRaddeg; | |
1001 | phi = Float_t(TMath::ATan2(Double_t(mom[1]),Double_t(mom[0])))*kRaddeg; | |
1002 | ||
1003 | ||
1004 | Double_t localTc = localMom[0]*localMom[0]+localMom[2]*localMom[2]; | |
1005 | Double_t localRt = TMath::Sqrt(localTc); | |
1006 | localTheta = Float_t(TMath::ATan2(localRt,Double_t(localMom[1])))*kRaddeg; | |
1007 | localPhi = Float_t(TMath::ATan2(Double_t(localMom[2]),Double_t(localMom[0])))*kRaddeg; | |
1008 | ||
1009 | hits[0] = Float_t(ipart); // particle type | |
1010 | hits[1] = localPos[0]; // X-position for hit | |
1011 | hits[2] = localPos[1]; // Y-position for hit | |
1012 | hits[3] = localPos[2]; // Z-position for hit | |
1013 | hits[4] = localTheta; // theta angle of incidence | |
1014 | hits[5] = localPhi; // phi angle of incidence | |
853634d3 | 1015 | hits[8] = (Float_t) fNsdigits; // first sdigit |
cd1cf354 | 1016 | hits[9] = -1; // last pad hit |
1017 | hits[13] = fFreonProd; // did id hit the freon? | |
1018 | hits[14] = mom[0]; | |
1019 | hits[15] = mom[1]; | |
1020 | hits[16] = mom[2]; | |
1021 | hits[18] = 0; // dummy cerenkov angle | |
1022 | ||
1023 | tlength = 0; | |
1024 | eloss = 0; | |
1025 | fFreonProd = 0; | |
1026 | ||
1027 | Chamber(idvol).LocaltoGlobal(localPos,hits+1); | |
1028 | ||
1029 | ||
1030 | //To make chamber coordinates x-y had to pass localPos[0], localPos[2] | |
1031 | xhit = localPos[0]; | |
1032 | yhit = localPos[2]; | |
1033 | // Only if not trigger chamber | |
1034 | if(idvol<kNCH) { | |
1035 | // | |
1036 | // Initialize hit position (cursor) in the segmentation model | |
1037 | //PH ((AliRICHChamber*) (*fChambers)[idvol]) | |
1038 | ((AliRICHChamber*)fChambers->At(idvol)) | |
1039 | ->SigGenInit(localPos[0], localPos[2], localPos[1]); | |
1040 | } | |
1041 | } | |
1042 | ||
1043 | // | |
1044 | // Calculate the charge induced on a pad (disintegration) in case | |
1045 | // | |
1046 | // Mip left chamber ... | |
1047 | if( gMC->IsTrackExiting() || gMC->IsTrackStop() || gMC->IsTrackDisappeared()){ | |
1048 | gMC->SetMaxStep(kBig); | |
1049 | eloss += destep; | |
1050 | tlength += step; | |
1051 | ||
1052 | ||
1053 | // Only if not trigger chamber | |
1054 | if(idvol<kNCH) { | |
1055 | if (eloss > 0) | |
1056 | { | |
1057 | if(gMC->TrackPid() == kNeutron) | |
1058 | printf("\n\n\n\n\n Neutron Making Pad Hit!!! \n\n\n\n"); | |
dfa42ff8 | 1059 | hits[17] = Hits2SDigits(xhit,yhit,eloss,idvol,kMip); //for MIP |
cd1cf354 | 1060 | } |
1061 | } | |
1062 | ||
1063 | hits[6]=tlength; | |
1064 | hits[7]=eloss; | |
853634d3 | 1065 | if (fNsdigits > (Int_t)hits[8]) { |
cd1cf354 | 1066 | hits[8]= hits[8]+1; |
853634d3 | 1067 | hits[9]= (Float_t) fNsdigits; |
cd1cf354 | 1068 | } |
1069 | ||
1070 | //if(sector !=-1) | |
853634d3 | 1071 | new(lhits[fNhits++]) AliRICHhit(fIshunt,gAlice->GetCurrentTrackNumber(),vol,hits); |
cd1cf354 | 1072 | eloss = 0; |
1073 | // | |
1074 | // Check additional signal generation conditions | |
1075 | // defined by the segmentation | |
1076 | // model (boundary crossing conditions) | |
dfa42ff8 | 1077 | }else if(((AliRICHChamber*)fChambers->At(idvol))->SigGenCond(localPos[0], localPos[2], localPos[1])){ |
1078 | ((AliRICHChamber*)fChambers->At(idvol))->SigGenInit(localPos[0], localPos[2], localPos[1]); | |
cd1cf354 | 1079 | if (eloss > 0) |
1080 | { | |
1081 | if(gMC->TrackPid() == kNeutron) | |
1082 | printf("\n\n\n\n\n Neutron Making Pad Hit!!! \n\n\n\n"); | |
dfa42ff8 | 1083 | hits[17] = Hits2SDigits(xhit,yhit,eloss,idvol,kMip);//for n |
cd1cf354 | 1084 | } |
1085 | xhit = localPos[0]; | |
1086 | yhit = localPos[2]; | |
1087 | eloss = destep; | |
1088 | tlength += step ; | |
1089 | // | |
1090 | // nothing special happened, add up energy loss | |
1091 | } else { | |
1092 | eloss += destep; | |
1093 | tlength += step ; | |
1094 | } | |
dfa42ff8 | 1095 | }//is in GAP? |
1096 | }//is MIP? | |
cd1cf354 | 1097 | /*************************************************End of MIP treatment**************************************/ |
1098 | }//void AliRICHv3::StepManager() | |
c021cb15 | 1099 | //__________________________________________________________________________________________________ |
1100 | Int_t AliRICHv3::Hits2SDigits(Float_t xhit,Float_t yhit,Float_t eloss, Int_t idvol, ResponseType res) | |
1101 | {//calls the charge disintegration method of the current chamber and adds all generated sdigits to the list of digits | |
1102 | ||
c60862bf | 1103 | Float_t newclust[4][500]; |
1104 | Int_t clhits[5]; | |
c021cb15 | 1105 | Int_t iNdigits; |
c60862bf | 1106 | clhits[0]=fNhits+1; |
1107 | ||
1108 | ((AliRICHChamber*)fChambers->At(idvol))->DisIntegration(eloss, xhit, yhit, iNdigits,newclust, res); | |
1109 | ||
1110 | for (Int_t i=0; i<iNdigits; i++) { | |
1111 | if (Int_t(newclust[0][i]) > 0) { | |
1112 | clhits[1] = Int_t(newclust[0][i]);// Cluster Charge | |
1113 | clhits[2] = Int_t(newclust[1][i]);// Pad: ix | |
1114 | clhits[3] = Int_t(newclust[2][i]);// Pad: iy | |
1115 | clhits[4] = Int_t(newclust[3][i]);// Pad: chamber sector | |
543d5224 | 1116 | AddSpecialOld(clhits); |
c60862bf | 1117 | } |
1118 | } | |
c60862bf | 1119 | return iNdigits; |
c021cb15 | 1120 | }//Int_t AliRICHv3::Hits2SDigits(Float_t xhit,Float_t yhit,Float_t eloss, Int_t idvol, ResponseType res) |
1121 | //__________________________________________________________________________________________________ | |
1122 | void AliRICHv3::DiagnosticsFE(Int_t evNumber1,Int_t evNumber2) | |
1123 | { | |
1124 | ||
1125 | Int_t NpadX = 162; // number of pads on X | |
1126 | Int_t NpadY = 162; // number of pads on Y | |
1127 | ||
1128 | Int_t Pad[162][162]; | |
1129 | for (Int_t i=0;i<NpadX;i++) { | |
1130 | for (Int_t j=0;j<NpadY;j++) { | |
1131 | Pad[i][j]=0; | |
1132 | } | |
1133 | } | |
1134 | ||
1135 | // Create some histograms | |
1136 | ||
1137 | TH1F *pionspectra1 = new TH1F("pionspectra1","Pion Spectra",200,-4,2); | |
1138 | TH1F *pionspectra2 = new TH1F("pionspectra2","Pion Spectra",200,-4,2); | |
1139 | TH1F *pionspectra3 = new TH1F("pionspectra3","Pion Spectra",200,-4,2); | |
1140 | TH1F *protonspectra1 = new TH1F("protonspectra1","Proton Spectra",200,-4,2); | |
1141 | TH1F *protonspectra2 = new TH1F("protonspectra2","Proton Spectra",200,-4,2); | |
1142 | TH1F *protonspectra3 = new TH1F("protonspectra3","Proton Spectra",200,-4,2); | |
1143 | TH1F *kaonspectra1 = new TH1F("kaonspectra1","Kaon Spectra",100,-4,2); | |
1144 | TH1F *kaonspectra2 = new TH1F("kaonspectra2","Kaon Spectra",100,-4,2); | |
1145 | TH1F *kaonspectra3 = new TH1F("kaonspectra3","Kaon Spectra",100,-4,2); | |
1146 | TH1F *electronspectra1 = new TH1F("electronspectra1","Electron Spectra",100,-4,2); | |
1147 | TH1F *electronspectra2 = new TH1F("electronspectra2","Electron Spectra",100,-4,2); | |
1148 | TH1F *electronspectra3 = new TH1F("electronspectra3","Electron Spectra",100,-4,2); | |
1149 | TH1F *muonspectra1 = new TH1F("muonspectra1","Muon Spectra",100,-4,2); | |
1150 | TH1F *muonspectra2 = new TH1F("muonspectra2","Muon Spectra",100,-4,2); | |
1151 | TH1F *muonspectra3 = new TH1F("muonspectra3","Muon Spectra",100,-4,2); | |
1152 | TH1F *neutronspectra1 = new TH1F("neutronspectra1","Neutron Spectra",100,-4,2); | |
1153 | TH1F *neutronspectra2 = new TH1F("neutronspectra2","Neutron Spectra",100,-4,2); | |
1154 | TH1F *neutronspectra3 = new TH1F("neutronspectra2","Neutron Spectra",100,-4,2); | |
1155 | TH1F *chargedspectra1 = new TH1F("chargedspectra1","Charged particles above 1 GeV Spectra",100,-1,3); | |
1156 | TH1F *chargedspectra2 = new TH1F("chargedspectra2","Charged particles above 1 GeV Spectra",100,-1,3); | |
1157 | TH1F *chargedspectra3 = new TH1F("chargedspectra2","Charged particles above 1 GeV Spectra",100,-1,3); | |
1158 | TH1F *pionptspectrafinal = new TH1F("pionptspectrafinal","Primary Pions Transverse Momenta at HMPID",20,0,5); | |
1159 | TH1F *pionptspectravertex = new TH1F("pionptspectravertex","Primary Pions Transverse Momenta at vertex",20,0,5); | |
1160 | TH1F *kaonptspectrafinal = new TH1F("kaonptspectrafinal","Primary Kaons Transverse Momenta at HMPID",20,0,5); | |
1161 | TH1F *kaonptspectravertex = new TH1F("kaonptspectravertex","Primary Kaons Transverse Momenta at vertex",20,0,5); | |
1162 | //TH1F *hitsPhi = new TH1F("hitsPhi","Distribution of phi angle of incidence",100,-180,180); | |
1163 | TH1F *hitsTheta = new TH1F("hitsTheta","Distribution of Theta angle of incidence, all tracks",100,0,50); | |
1164 | TH1F *hitsTheta500MeV = new TH1F("hitsTheta500MeV","Distribution of Theta angle of incidence, 0.5-1 GeV primary tracks",100,0,50); | |
1165 | TH1F *hitsTheta1GeV = new TH1F("hitsTheta1GeV","Distribution of Theta angle of incidence, 1-2 GeV primary tracks",100,0,50); | |
1166 | TH1F *hitsTheta2GeV = new TH1F("hitsTheta2GeV","Distribution of Theta angle of incidence, 2-3 GeV primary tracks",100,0,50); | |
1167 | TH1F *hitsTheta3GeV = new TH1F("hitsTheta3GeV","Distribution of Theta angle of incidence, >3 GeV primary tracks",100,0,50); | |
1168 | TH2F *production = new TH2F("production","Mother production vertices",100,-300,300,100,0,600); | |
1169 | ||
1170 | ||
1171 | ||
1172 | ||
1173 | // Start loop over events | |
1174 | ||
1175 | Int_t pion=0, kaon=0, proton=0, electron=0, positron=0, neutron=0, highneutrons=0, muon=0; | |
1176 | Int_t chargedpions=0,primarypions=0,highprimarypions=0,chargedkaons=0,primarykaons=0,highprimarykaons=0; | |
1177 | Int_t photons=0, primaryphotons=0, highprimaryphotons=0; | |
1178 | TRandom* random=0; | |
1179 | ||
1180 | for (int nev=0; nev<= evNumber2; nev++) { | |
1181 | Int_t nparticles = gAlice->GetEvent(nev); | |
1182 | ||
1183 | ||
1184 | if (nev < evNumber1) continue; | |
1185 | if (nparticles <= 0) return; | |
1186 | ||
1187 | // Get pointers to RICH detector and Hits containers | |
1188 | ||
1189 | AliRICH *pRICH = (AliRICH *) gAlice->GetDetector("RICH"); | |
1190 | ||
1191 | TTree *treeH = TreeH(); | |
1192 | Int_t ntracks =(Int_t) treeH->GetEntries(); | |
1193 | ||
1194 | // Start loop on tracks in the hits containers | |
1195 | ||
1196 | for (Int_t track=0; track<ntracks;track++) { | |
1197 | printf ("Processing Track: %d\n",track); | |
1198 | gAlice->ResetHits(); | |
1199 | treeH->GetEvent(track); | |
1200 | ||
1201 | for(AliRICHhit* mHit=(AliRICHhit*)pRICH->FirstHit(-1); | |
1202 | mHit; | |
1203 | mHit=(AliRICHhit*)pRICH->NextHit()) | |
1204 | { | |
1205 | //Int_t nch = mHit->fChamber; // chamber number | |
1206 | //Float_t x = mHit->X(); // x-pos of hit | |
1207 | //Float_t y = mHit->Z(); // y-pos | |
1208 | //Float_t z = mHit->Y(); | |
1209 | //Float_t phi = mHit->Phi(); //Phi angle of incidence | |
1210 | Float_t theta = mHit->Theta(); //Theta angle of incidence | |
1211 | Float_t px = mHit->MomX(); | |
1212 | Float_t py = mHit->MomY(); | |
1213 | Int_t index = mHit->Track(); | |
1214 | Int_t particle = (Int_t)(mHit->Particle()); | |
1215 | Float_t R; | |
1216 | Float_t PTfinal; | |
1217 | Float_t PTvertex; | |
1218 | ||
1219 | TParticle *current = gAlice->Particle(index); | |
1220 | ||
1221 | //Float_t energy=current->Energy(); | |
1222 | ||
1223 | R=TMath::Sqrt(current->Vx()*current->Vx() + current->Vy()*current->Vy()); | |
1224 | PTfinal=TMath::Sqrt(px*px + py*py); | |
1225 | PTvertex=TMath::Sqrt(current->Px()*current->Px() + current->Py()*current->Py()); | |
1226 | ||
1227 | ||
1228 | ||
1229 | if (TMath::Abs(particle) < 10000000) | |
1230 | { | |
1231 | hitsTheta->Fill(theta,(float) 1); | |
1232 | if (R<5) | |
1233 | { | |
1234 | if (PTvertex>.5 && PTvertex<=1) | |
1235 | { | |
1236 | hitsTheta500MeV->Fill(theta,(float) 1); | |
1237 | } | |
1238 | if (PTvertex>1 && PTvertex<=2) | |
1239 | { | |
1240 | hitsTheta1GeV->Fill(theta,(float) 1); | |
1241 | } | |
1242 | if (PTvertex>2 && PTvertex<=3) | |
1243 | { | |
1244 | hitsTheta2GeV->Fill(theta,(float) 1); | |
1245 | } | |
1246 | if (PTvertex>3) | |
1247 | { | |
1248 | hitsTheta3GeV->Fill(theta,(float) 1); | |
1249 | } | |
1250 | } | |
1251 | ||
1252 | } | |
1253 | ||
1254 | //if (nch == 3) | |
1255 | //{ | |
1256 | ||
1257 | if (TMath::Abs(particle) < 50000051) | |
1258 | { | |
1259 | //if (TMath::Abs(particle) == 50000050 || TMath::Abs(particle) == 2112) | |
1260 | if (TMath::Abs(particle) == 2112 || TMath::Abs(particle) == 50000050) | |
1261 | { | |
1262 | //gMC->Rndm(&random, 1); | |
1263 | if (random->Rndm() < .1) | |
1264 | production->Fill(current->Vz(),R,(float) 1); | |
1265 | if (TMath::Abs(particle) == 50000050) | |
1266 | //if (TMath::Abs(particle) > 50000000) | |
1267 | { | |
1268 | photons +=1; | |
1269 | if (R<5) | |
1270 | { | |
1271 | primaryphotons +=1; | |
1272 | if (current->Energy()>0.001) | |
1273 | highprimaryphotons +=1; | |
1274 | } | |
1275 | } | |
1276 | if (TMath::Abs(particle) == 2112) | |
1277 | { | |
1278 | neutron +=1; | |
1279 | if (current->Energy()>0.0001) | |
1280 | highneutrons +=1; | |
1281 | } | |
1282 | } | |
1283 | if (TMath::Abs(particle) < 50000000) | |
1284 | { | |
1285 | production->Fill(current->Vz(),R,(float) 1); | |
1286 | } | |
1287 | //mip->Fill(x,y,(float) 1); | |
1288 | } | |
1289 | ||
1290 | if (TMath::Abs(particle)==211 || TMath::Abs(particle)==111) | |
1291 | { | |
1292 | if (R<5) | |
1293 | { | |
1294 | pionptspectravertex->Fill(PTvertex,(float) 1); | |
1295 | pionptspectrafinal->Fill(PTfinal,(float) 1); | |
1296 | } | |
1297 | } | |
1298 | ||
1299 | if (TMath::Abs(particle)==321 || TMath::Abs(particle)==130 || TMath::Abs(particle)==310 | |
1300 | || TMath::Abs(particle)==311) | |
1301 | { | |
1302 | if (R<5) | |
1303 | { | |
1304 | kaonptspectravertex->Fill(PTvertex,(float) 1); | |
1305 | kaonptspectrafinal->Fill(PTfinal,(float) 1); | |
1306 | } | |
1307 | } | |
1308 | ||
1309 | ||
1310 | if (TMath::Abs(particle)==211 || TMath::Abs(particle)==111) | |
1311 | { | |
1312 | pionspectra1->Fill(TMath::Log10(current->Energy() - current->GetCalcMass()),(float) 1); | |
1313 | if (current->Vx()>5 && current->Vy()>5 && current->Vz()>5) | |
1314 | pionspectra2->Fill(TMath::Log10(current->Energy() - current->GetCalcMass()),(float) 1); | |
1315 | if (R>250 && R<450) | |
1316 | { | |
1317 | pionspectra3->Fill(TMath::Log10(current->Energy() - current->GetCalcMass()),(float) 1); | |
1318 | } | |
1319 | pion +=1; | |
1320 | if (TMath::Abs(particle)==211) | |
1321 | { | |
1322 | chargedpions +=1; | |
1323 | if (R<5) | |
1324 | { | |
1325 | primarypions +=1; | |
1326 | if (current->Energy()>1) | |
1327 | highprimarypions +=1; | |
1328 | } | |
1329 | } | |
1330 | } | |
1331 | if (TMath::Abs(particle)==2212) | |
1332 | { | |
1333 | protonspectra1->Fill(TMath::Log10(current->Energy() - current->GetCalcMass()),(float) 1); | |
1334 | //ptspectra->Fill(Pt,(float) 1); | |
1335 | if (current->Vx()>5 && current->Vy()>5 && current->Vz()>5) | |
1336 | protonspectra2->Fill(TMath::Log10(current->Energy() - current->GetCalcMass()),(float) 1); | |
1337 | if (R>250 && R<450) | |
1338 | protonspectra3->Fill(TMath::Log10(current->Energy() - current->GetCalcMass()),(float) 1); | |
1339 | proton +=1; | |
1340 | } | |
1341 | if (TMath::Abs(particle)==321 || TMath::Abs(particle)==130 || TMath::Abs(particle)==310 | |
1342 | || TMath::Abs(particle)==311) | |
1343 | { | |
1344 | kaonspectra1->Fill(TMath::Log10(current->Energy() - current->GetCalcMass()),(float) 1); | |
1345 | //ptspectra->Fill(Pt,(float) 1); | |
1346 | if (current->Vx()>5 && current->Vy()>5 && current->Vz()>5) | |
1347 | kaonspectra2->Fill(TMath::Log10(current->Energy() - current->GetCalcMass()),(float) 1); | |
1348 | if (R>250 && R<450) | |
1349 | kaonspectra3->Fill(TMath::Log10(current->Energy() - current->GetCalcMass()),(float) 1); | |
1350 | kaon +=1; | |
1351 | if (TMath::Abs(particle)==321) | |
1352 | { | |
1353 | chargedkaons +=1; | |
1354 | if (R<5) | |
1355 | { | |
1356 | primarykaons +=1; | |
1357 | if (current->Energy()>1) | |
1358 | highprimarykaons +=1; | |
1359 | } | |
1360 | } | |
1361 | } | |
1362 | if (TMath::Abs(particle)==11) | |
1363 | { | |
1364 | electronspectra1->Fill(TMath::Log10(current->Energy() - current->GetCalcMass()),(float) 1); | |
1365 | //ptspectra->Fill(Pt,(float) 1); | |
1366 | if (current->Vx()>5 && current->Vy()>5 && current->Vz()>5) | |
1367 | electronspectra2->Fill(TMath::Log10(current->Energy() - current->GetCalcMass()),(float) 1); | |
1368 | if (R>250 && R<450) | |
1369 | electronspectra3->Fill(TMath::Log10(current->Energy() - current->GetCalcMass()),(float) 1); | |
1370 | if (particle == 11) | |
1371 | electron +=1; | |
1372 | if (particle == -11) | |
1373 | positron +=1; | |
1374 | } | |
1375 | if (TMath::Abs(particle)==13) | |
1376 | { | |
1377 | muonspectra1->Fill(TMath::Log10(current->Energy() - current->GetCalcMass()),(float) 1); | |
1378 | //ptspectra->Fill(Pt,(float) 1); | |
1379 | if (current->Vx()>5 && current->Vy()>5 && current->Vz()>5) | |
1380 | muonspectra2->Fill(TMath::Log10(current->Energy() - current->GetCalcMass()),(float) 1); | |
1381 | if (R>250 && R<450) | |
1382 | muonspectra3->Fill(TMath::Log10(current->Energy() - current->GetCalcMass()),(float) 1); | |
1383 | muon +=1; | |
1384 | } | |
1385 | if (TMath::Abs(particle)==2112) | |
1386 | { | |
1387 | neutronspectra1->Fill(TMath::Log10(current->Energy() - current->GetCalcMass()),(float) 1); | |
1388 | //ptspectra->Fill(Pt,(float) 1); | |
1389 | if (current->Vx()>5 && current->Vy()>5 && current->Vz()>5) | |
1390 | neutronspectra2->Fill(TMath::Log10(current->Energy() - current->GetCalcMass()),(float) 1); | |
1391 | if (R>250 && R<450) | |
1392 | { | |
1393 | neutronspectra3->Fill(TMath::Log10(current->Energy() - current->GetCalcMass()),(float) 1); | |
1394 | } | |
1395 | neutron +=1; | |
1396 | } | |
1397 | if(TMath::Abs(particle)==211 || TMath::Abs(particle)==2212 || TMath::Abs(particle)==321) | |
1398 | { | |
1399 | if (current->Energy()-current->GetCalcMass()>1) | |
1400 | { | |
1401 | chargedspectra1->Fill(TMath::Log10(current->Energy() - current->GetCalcMass()),(float) 1); | |
1402 | if (current->Vx()>5 && current->Vy()>5 && current->Vz()>5) | |
1403 | chargedspectra2->Fill(TMath::Log10(current->Energy() - current->GetCalcMass()),(float) 1); | |
1404 | if (R>250 && R<450) | |
1405 | chargedspectra3->Fill(TMath::Log10(current->Energy() - current->GetCalcMass()),(float) 1); | |
1406 | } | |
1407 | } | |
1408 | // Fill the histograms | |
1409 | //Nh1+=nhits; | |
1410 | //h->Fill(x,y,(float) 1); | |
1411 | //} | |
1412 | //} | |
1413 | } | |
1414 | ||
1415 | } | |
1416 | ||
1417 | } | |
1418 | // } | |
1419 | ||
1420 | TStyle *mystyle=new TStyle("Plain","mystyle"); | |
1421 | mystyle->SetPalette(1,0); | |
1422 | mystyle->cd(); | |
1423 | ||
1424 | //Create canvases, set the view range, show histograms | |
1425 | ||
1426 | TCanvas *c2 = new TCanvas("c2","Angles of incidence",150,150,100,150); | |
1427 | c2->Divide(2,2); | |
1428 | //c2->SetFillColor(42); | |
1429 | ||
1430 | c2->cd(1); | |
1431 | hitsTheta500MeV->SetFillColor(5); | |
1432 | hitsTheta500MeV->Draw(); | |
1433 | c2->cd(2); | |
1434 | hitsTheta1GeV->SetFillColor(5); | |
1435 | hitsTheta1GeV->Draw(); | |
1436 | c2->cd(3); | |
1437 | hitsTheta2GeV->SetFillColor(5); | |
1438 | hitsTheta2GeV->Draw(); | |
1439 | c2->cd(4); | |
1440 | hitsTheta3GeV->SetFillColor(5); | |
1441 | hitsTheta3GeV->Draw(); | |
1442 | ||
1443 | ||
1444 | ||
1445 | TCanvas *c15 = new TCanvas("c15","Mothers Production Vertices",50,50,600,600); | |
1446 | c15->cd(); | |
1447 | production->SetFillColor(42); | |
1448 | production->SetXTitle("z (m)"); | |
1449 | production->SetYTitle("R (m)"); | |
1450 | production->Draw(); | |
1451 | ||
1452 | TCanvas *c10 = new TCanvas("c10","Pt Spectra",50,50,600,700); | |
1453 | c10->Divide(2,2); | |
1454 | c10->cd(1); | |
1455 | pionptspectravertex->SetFillColor(5); | |
1456 | pionptspectravertex->SetXTitle("Pt (GeV)"); | |
1457 | pionptspectravertex->Draw(); | |
1458 | c10->cd(2); | |
1459 | pionptspectrafinal->SetFillColor(5); | |
1460 | pionptspectrafinal->SetXTitle("Pt (GeV)"); | |
1461 | pionptspectrafinal->Draw(); | |
1462 | c10->cd(3); | |
1463 | kaonptspectravertex->SetFillColor(5); | |
1464 | kaonptspectravertex->SetXTitle("Pt (GeV)"); | |
1465 | kaonptspectravertex->Draw(); | |
1466 | c10->cd(4); | |
1467 | kaonptspectrafinal->SetFillColor(5); | |
1468 | kaonptspectrafinal->SetXTitle("Pt (GeV)"); | |
1469 | kaonptspectrafinal->Draw(); | |
1470 | ||
1471 | ||
1472 | TCanvas *c16 = new TCanvas("c16","Particles Spectra II",150,150,600,350); | |
1473 | c16->Divide(2,1); | |
1474 | ||
1475 | c16->cd(1); | |
1476 | //TCanvas *c13 = new TCanvas("c13","Electron Spectra",400,10,600,700); | |
1477 | electronspectra1->SetFillColor(5); | |
1478 | electronspectra1->SetXTitle("log(GeV)"); | |
1479 | electronspectra2->SetFillColor(46); | |
1480 | electronspectra2->SetXTitle("log(GeV)"); | |
1481 | electronspectra3->SetFillColor(10); | |
1482 | electronspectra3->SetXTitle("log(GeV)"); | |
1483 | //c13->SetLogx(); | |
1484 | electronspectra1->Draw(); | |
1485 | electronspectra2->Draw("same"); | |
1486 | electronspectra3->Draw("same"); | |
1487 | ||
1488 | c16->cd(2); | |
1489 | //TCanvas *c14 = new TCanvas("c14","Muon Spectra",400,10,600,700); | |
1490 | muonspectra1->SetFillColor(5); | |
1491 | muonspectra1->SetXTitle("log(GeV)"); | |
1492 | muonspectra2->SetFillColor(46); | |
1493 | muonspectra2->SetXTitle("log(GeV)"); | |
1494 | muonspectra3->SetFillColor(10); | |
1495 | muonspectra3->SetXTitle("log(GeV)"); | |
1496 | //c14->SetLogx(); | |
1497 | muonspectra1->Draw(); | |
1498 | muonspectra2->Draw("same"); | |
1499 | muonspectra3->Draw("same"); | |
1500 | ||
1501 | //c16->cd(3); | |
1502 | //TCanvas *c16 = new TCanvas("c16","Neutron Spectra",400,10,600,700); | |
1503 | //neutronspectra1->SetFillColor(42); | |
1504 | //neutronspectra1->SetXTitle("log(GeV)"); | |
1505 | //neutronspectra2->SetFillColor(46); | |
1506 | //neutronspectra2->SetXTitle("log(GeV)"); | |
1507 | //neutronspectra3->SetFillColor(10); | |
1508 | //neutronspectra3->SetXTitle("log(GeV)"); | |
1509 | //c16->SetLogx(); | |
1510 | //neutronspectra1->Draw(); | |
1511 | //neutronspectra2->Draw("same"); | |
1512 | //neutronspectra3->Draw("same"); | |
1513 | ||
1514 | TCanvas *c9 = new TCanvas("c9","Particles Spectra",150,150,600,700); | |
1515 | //TCanvas *c9 = new TCanvas("c9","Pion Spectra",400,10,600,700); | |
1516 | c9->Divide(2,2); | |
1517 | ||
1518 | c9->cd(1); | |
1519 | pionspectra1->SetFillColor(5); | |
1520 | pionspectra1->SetXTitle("log(GeV)"); | |
1521 | pionspectra2->SetFillColor(46); | |
1522 | pionspectra2->SetXTitle("log(GeV)"); | |
1523 | pionspectra3->SetFillColor(10); | |
1524 | pionspectra3->SetXTitle("log(GeV)"); | |
1525 | //c9->SetLogx(); | |
1526 | pionspectra1->Draw(); | |
1527 | pionspectra2->Draw("same"); | |
1528 | pionspectra3->Draw("same"); | |
1529 | ||
1530 | c9->cd(2); | |
1531 | //TCanvas *c10 = new TCanvas("c10","Proton Spectra",400,10,600,700); | |
1532 | protonspectra1->SetFillColor(5); | |
1533 | protonspectra1->SetXTitle("log(GeV)"); | |
1534 | protonspectra2->SetFillColor(46); | |
1535 | protonspectra2->SetXTitle("log(GeV)"); | |
1536 | protonspectra3->SetFillColor(10); | |
1537 | protonspectra3->SetXTitle("log(GeV)"); | |
1538 | //c10->SetLogx(); | |
1539 | protonspectra1->Draw(); | |
1540 | protonspectra2->Draw("same"); | |
1541 | protonspectra3->Draw("same"); | |
1542 | ||
1543 | c9->cd(3); | |
1544 | //TCanvas *c11 = new TCanvas("c11","Kaon Spectra",400,10,600,700); | |
1545 | kaonspectra1->SetFillColor(5); | |
1546 | kaonspectra1->SetXTitle("log(GeV)"); | |
1547 | kaonspectra2->SetFillColor(46); | |
1548 | kaonspectra2->SetXTitle("log(GeV)"); | |
1549 | kaonspectra3->SetFillColor(10); | |
1550 | kaonspectra3->SetXTitle("log(GeV)"); | |
1551 | //c11->SetLogx(); | |
1552 | kaonspectra1->Draw(); | |
1553 | kaonspectra2->Draw("same"); | |
1554 | kaonspectra3->Draw("same"); | |
1555 | ||
1556 | c9->cd(4); | |
1557 | //TCanvas *c12 = new TCanvas("c12","Charged Particles Spectra",400,10,600,700); | |
1558 | chargedspectra1->SetFillColor(5); | |
1559 | chargedspectra1->SetXTitle("log(GeV)"); | |
1560 | chargedspectra2->SetFillColor(46); | |
1561 | chargedspectra2->SetXTitle("log(GeV)"); | |
1562 | chargedspectra3->SetFillColor(10); | |
1563 | chargedspectra3->SetXTitle("log(GeV)"); | |
1564 | //c12->SetLogx(); | |
1565 | chargedspectra1->Draw(); | |
1566 | chargedspectra2->Draw("same"); | |
1567 | chargedspectra3->Draw("same"); | |
1568 | ||
1569 | ||
1570 | ||
1571 | printf("*****************************************\n"); | |
1572 | printf("* Particle * Counts *\n"); | |
1573 | printf("*****************************************\n"); | |
1574 | ||
1575 | printf("* Pions: * %4d *\n",pion); | |
1576 | printf("* Charged Pions: * %4d *\n",chargedpions); | |
1577 | printf("* Primary Pions: * %4d *\n",primarypions); | |
1578 | printf("* Primary Pions (p>1GeV/c): * %4d *\n",highprimarypions); | |
1579 | printf("* Kaons: * %4d *\n",kaon); | |
1580 | printf("* Charged Kaons: * %4d *\n",chargedkaons); | |
1581 | printf("* Primary Kaons: * %4d *\n",primarykaons); | |
1582 | printf("* Primary Kaons (p>1GeV/c): * %4d *\n",highprimarykaons); | |
1583 | printf("* Muons: * %4d *\n",muon); | |
1584 | printf("* Electrons: * %4d *\n",electron); | |
1585 | printf("* Positrons: * %4d *\n",positron); | |
1586 | printf("* Protons: * %4d *\n",proton); | |
1587 | printf("* All Charged: * %4d *\n",(chargedpions+chargedkaons+muon+electron+positron+proton)); | |
1588 | printf("*****************************************\n"); | |
1589 | //printf("* Photons: * %3.1f *\n",photons); | |
1590 | //printf("* Primary Photons: * %3.1f *\n",primaryphotons); | |
1591 | //printf("* Primary Photons (p>1MeV/c):* %3.1f *\n",highprimaryphotons); | |
1592 | //printf("*****************************************\n"); | |
1593 | //printf("* Neutrons: * %3.1f *\n",neutron); | |
1594 | //printf("* Neutrons (p>100keV/c): * %3.1f *\n",highneutrons); | |
1595 | //printf("*****************************************\n"); | |
1596 | ||
1597 | if (gAlice->TreeD()) | |
1598 | { | |
1599 | gAlice->TreeD()->GetEvent(0); | |
1600 | ||
1601 | Float_t occ[7]; | |
1602 | Float_t sum=0; | |
1603 | Float_t mean=0; | |
1604 | printf("\n*****************************************\n"); | |
1605 | printf("* Chamber * Digits * Occupancy *\n"); | |
1606 | printf("*****************************************\n"); | |
1607 | ||
1608 | for (Int_t ich=0;ich<7;ich++) | |
1609 | { | |
1610 | TClonesArray *Digits = DigitsAddress(ich); // Raw clusters branch | |
1611 | Int_t ndigits = Digits->GetEntriesFast(); | |
1612 | occ[ich] = Float_t(ndigits)/(160*144); | |
1613 | sum += Float_t(ndigits)/(160*144); | |
1614 | printf("* %d * %d * %3.1f%% *\n",ich,ndigits,occ[ich]*100); | |
1615 | } | |
1616 | mean = sum/7; | |
1617 | printf("*****************************************\n"); | |
1618 | printf("* Mean occupancy * %3.1f%% *\n",mean*100); | |
1619 | printf("*****************************************\n"); | |
1620 | } | |
1621 | ||
1622 | printf("\nEnd of analysis\n"); | |
1623 | ||
1624 | }//void AliRICHv3::DiagnosticsFE(Int_t evNumber1,Int_t evNumber2) | |
1625 | //__________________________________________________________________________________________________ | |
c021cb15 | 1626 |