]>
Commit | Line | Data |
---|---|---|
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 | ||
16 | /* | |
17 | $Log$ | |
18 | Revision 1.18 2000/05/08 16:17:27 cblume | |
19 | Merge TRD-develop | |
20 | ||
21 | Revision 1.17.2.1 2000/05/08 14:28:59 cblume | |
22 | Introduced SetPHOShole() and SetRICHhole(). AliTRDrecPoint container is now a TObjArray | |
23 | ||
24 | Revision 1.17 2000/02/28 19:10:26 cblume | |
25 | Include the new TRD classes | |
26 | ||
27 | Revision 1.16.2.2 2000/02/28 17:53:24 cblume | |
28 | Introduce TRD geometry classes | |
29 | ||
30 | Revision 1.16.2.1 2000/02/28 17:04:19 cblume | |
31 | Include functions and data members for AliTRDrecPoint | |
32 | ||
33 | Revision 1.16 2000/01/19 17:17:35 fca | |
34 | Introducing a list of lists of hits -- more hits allowed for detector now | |
35 | ||
36 | Revision 1.15 1999/11/02 17:04:25 fca | |
37 | Small syntax change for HP compiler | |
38 | ||
39 | Revision 1.14 1999/11/02 16:57:02 fca | |
40 | Avoid non ansi warnings on HP compilers | |
41 | ||
42 | Revision 1.13 1999/11/02 16:35:56 fca | |
43 | New version of TRD introduced | |
44 | ||
45 | Revision 1.12 1999/11/01 20:41:51 fca | |
46 | Added protections against using the wrong version of FRAME | |
47 | ||
48 | Revision 1.11 1999/09/29 09:24:34 fca | |
49 | Introduction of the Copyright and cvs Log | |
50 | ||
51 | */ | |
52 | ||
53 | /////////////////////////////////////////////////////////////////////////////// | |
54 | // // | |
55 | // Transition Radiation Detector // | |
56 | // This class contains the basic functions for the Transition Radiation // | |
57 | // Detector. // | |
58 | // // | |
59 | /////////////////////////////////////////////////////////////////////////////// | |
60 | ||
61 | #include <stdlib.h> | |
62 | ||
63 | #include <TMath.h> | |
64 | #include <TNode.h> | |
65 | #include <TPGON.h> | |
66 | ||
67 | #include "AliTRD.h" | |
68 | #include "AliRun.h" | |
69 | #include "AliConst.h" | |
70 | #include "AliTRDdigitizer.h" | |
71 | #include "AliTRDclusterizer.h" | |
72 | #include "AliTRDgeometryHole.h" | |
73 | #include "AliTRDgeometryFull.h" | |
74 | #include "AliTRDrecPoint.h" | |
75 | ||
76 | ClassImp(AliTRD) | |
77 | ||
78 | //_____________________________________________________________________________ | |
79 | AliTRD::AliTRD() | |
80 | { | |
81 | // | |
82 | // Default constructor | |
83 | // | |
84 | ||
85 | fIshunt = 0; | |
86 | fGasMix = 0; | |
87 | fHits = 0; | |
88 | fDigits = 0; | |
89 | ||
90 | fRecPoints = 0; | |
91 | fNRecPoints = 0; | |
92 | ||
93 | fGeometry = 0; | |
94 | ||
95 | } | |
96 | ||
97 | //_____________________________________________________________________________ | |
98 | AliTRD::AliTRD(const char *name, const char *title) | |
99 | : AliDetector(name,title) | |
100 | { | |
101 | // | |
102 | // Standard constructor for the TRD | |
103 | // | |
104 | ||
105 | // Check that FRAME is there otherwise we have no place where to | |
106 | // put TRD | |
107 | AliModule* FRAME=gAlice->GetModule("FRAME"); | |
108 | if (!FRAME) { | |
109 | Error("Ctor","TRD needs FRAME to be present\n"); | |
110 | exit(1); | |
111 | } | |
112 | ||
113 | // Define the TRD geometry according to the FRAME geometry | |
114 | if (FRAME->IsVersion() == 0) { | |
115 | // Geometry with hole | |
116 | fGeometry = new AliTRDgeometryHole(); | |
117 | } | |
118 | else if (FRAME->IsVersion() == 1) { | |
119 | // Geometry without hole | |
120 | fGeometry = new AliTRDgeometryFull(); | |
121 | } | |
122 | else { | |
123 | Error("Ctor","Could not find valid FRAME version\n"); | |
124 | exit(1); | |
125 | } | |
126 | ||
127 | // Allocate the hit array | |
128 | fHits = new TClonesArray("AliTRDhit" ,405); | |
129 | gAlice->AddHitList(fHits); | |
130 | ||
131 | // Allocate the digits array | |
132 | fDigits = 0; | |
133 | ||
134 | // Allocate the rec point array | |
135 | fRecPoints = new TObjArray(400); | |
136 | fNRecPoints = 0; | |
137 | ||
138 | fIshunt = 0; | |
139 | fGasMix = 0; | |
140 | ||
141 | SetMarkerColor(kWhite); | |
142 | ||
143 | } | |
144 | ||
145 | //_____________________________________________________________________________ | |
146 | AliTRD::~AliTRD() | |
147 | { | |
148 | // | |
149 | // TRD destructor | |
150 | // | |
151 | ||
152 | fIshunt = 0; | |
153 | ||
154 | delete fGeometry; | |
155 | delete fHits; | |
156 | delete fRecPoints; | |
157 | ||
158 | } | |
159 | ||
160 | //_____________________________________________________________________________ | |
161 | void AliTRD::AddRecPoint(Float_t *pos, Int_t *digits, Int_t det, Float_t amp) | |
162 | { | |
163 | // | |
164 | // Add a reconstructed point for the TRD | |
165 | // | |
166 | ||
167 | AliTRDrecPoint *RecPoint = new AliTRDrecPoint(); | |
168 | TVector3 posVec(pos[0],pos[1],pos[2]); | |
169 | RecPoint->SetLocalPosition(posVec); | |
170 | RecPoint->SetDetector(det); | |
171 | RecPoint->SetEnergy(amp); | |
172 | for (Int_t iDigit = 0; iDigit < 3; iDigit++) { | |
173 | RecPoint->AddDigit(digits[iDigit]); | |
174 | } | |
175 | ||
176 | fRecPoints->Add(RecPoint); | |
177 | ||
178 | } | |
179 | ||
180 | //_____________________________________________________________________________ | |
181 | void AliTRD::AddDigit(Int_t *digits, Int_t *amp) | |
182 | { | |
183 | // | |
184 | // Add a digit for the TRD | |
185 | // | |
186 | ||
187 | TClonesArray &ldigits = *fDigits; | |
188 | new(ldigits[fNdigits++]) AliTRDdigit(kFALSE,digits,amp); | |
189 | ||
190 | } | |
191 | ||
192 | //_____________________________________________________________________________ | |
193 | void AliTRD::AddHit(Int_t track, Int_t* det, Float_t *hits) | |
194 | { | |
195 | // | |
196 | // Add a hit for the TRD | |
197 | // | |
198 | ||
199 | TClonesArray &lhits = *fHits; | |
200 | new(lhits[fNhits++]) AliTRDhit(fIshunt,track,det,hits); | |
201 | ||
202 | } | |
203 | ||
204 | //_____________________________________________________________________________ | |
205 | void AliTRD::BuildGeometry() | |
206 | { | |
207 | // | |
208 | // Create the ROOT TNode geometry for the TRD | |
209 | // | |
210 | ||
211 | TNode *Node, *Top; | |
212 | TPGON *pgon; | |
213 | const Int_t kColorTRD = 46; | |
214 | ||
215 | // Find the top node alice | |
216 | Top = gAlice->GetGeometry()->GetNode("alice"); | |
217 | ||
218 | pgon = new TPGON("S_TRD","TRD","void",0,360,kNsect,4); | |
219 | Float_t ff = TMath::Cos(kDegrad * 180 / kNsect); | |
220 | Float_t rrmin = kRmin / ff; | |
221 | Float_t rrmax = kRmax / ff; | |
222 | pgon->DefineSection(0,-kZmax1,rrmax,rrmax); | |
223 | pgon->DefineSection(1,-kZmax2,rrmin,rrmax); | |
224 | pgon->DefineSection(2, kZmax2,rrmin,rrmax); | |
225 | pgon->DefineSection(3, kZmax1,rrmax,rrmax); | |
226 | Top->cd(); | |
227 | Node = new TNode("TRD","TRD","S_TRD",0,0,0,""); | |
228 | Node->SetLineColor(kColorTRD); | |
229 | fNodes->Add(Node); | |
230 | ||
231 | } | |
232 | ||
233 | //_____________________________________________________________________________ | |
234 | void AliTRD::CreateGeometry() | |
235 | { | |
236 | // | |
237 | // Creates the volumes for the TRD chambers | |
238 | // | |
239 | ||
240 | // Check that FRAME is there otherwise we have no place where to put the TRD | |
241 | AliModule* FRAME = gAlice->GetModule("FRAME"); | |
242 | if (!FRAME) { | |
243 | printf(" The TRD needs the FRAME to be defined first\n"); | |
244 | return; | |
245 | } | |
246 | ||
247 | fGeometry->CreateGeometry(fIdtmed->GetArray() - 1299); | |
248 | ||
249 | } | |
250 | ||
251 | //_____________________________________________________________________________ | |
252 | void AliTRD::CreateMaterials() | |
253 | { | |
254 | // | |
255 | // Create the materials for the TRD | |
256 | // Origin Y.Foka | |
257 | // | |
258 | ||
259 | Int_t ISXFLD = gAlice->Field()->Integ(); | |
260 | Float_t SXMGMX = gAlice->Field()->Max(); | |
261 | ||
262 | // For polyethilene (CH2) | |
263 | Float_t ape[2] = { 12., 1. }; | |
264 | Float_t zpe[2] = { 6., 1. }; | |
265 | Float_t wpe[2] = { 1., 2. }; | |
266 | Float_t dpe = 0.95; | |
267 | ||
268 | // For mylar (C5H4O2) | |
269 | Float_t amy[3] = { 12., 1., 16. }; | |
270 | Float_t zmy[3] = { 6., 1., 8. }; | |
271 | Float_t wmy[3] = { 5., 4., 2. }; | |
272 | Float_t dmy = 1.39; | |
273 | ||
274 | // For CO2 | |
275 | Float_t aco[2] = { 12., 16. }; | |
276 | Float_t zco[2] = { 6., 8. }; | |
277 | Float_t wco[2] = { 1., 2. }; | |
278 | Float_t dco = 0.001977; | |
279 | ||
280 | // For water | |
281 | Float_t awa[2] = { 1., 16. }; | |
282 | Float_t zwa[2] = { 1., 8. }; | |
283 | Float_t wwa[2] = { 2., 1. }; | |
284 | Float_t dwa = 1.0; | |
285 | ||
286 | // For isobutane (C4H10) | |
287 | Float_t ais[2] = { 12., 1. }; | |
288 | Float_t zis[2] = { 6., 1. }; | |
289 | Float_t wis[2] = { 4., 10. }; | |
290 | Float_t dis = 0.00267; | |
291 | ||
292 | // For Xe/CO2-gas-mixture | |
293 | // Xe-content of the Xe/CO2-mixture (90% / 10%) | |
294 | Float_t fxc = .90; | |
295 | // Xe-content of the Xe/Isobutane-mixture (97% / 3%) | |
296 | Float_t fxi = .97; | |
297 | Float_t dxe = .005858; | |
298 | ||
299 | // General tracking parameter | |
300 | Float_t tmaxfd = -10.; | |
301 | Float_t stemax = -1e10; | |
302 | Float_t deemax = -0.1; | |
303 | Float_t epsil = 1e-4; | |
304 | Float_t stmin = -0.001; | |
305 | ||
306 | Float_t absl, radl, d, buf[1]; | |
307 | Float_t agm[2], dgm, zgm[2], wgm[2]; | |
308 | Int_t nbuf; | |
309 | ||
310 | ////////////////////////////////////////////////////////////////////////// | |
311 | // Define Materials | |
312 | ////////////////////////////////////////////////////////////////////////// | |
313 | ||
314 | AliMaterial( 1, "Al $", 26.98, 13.0, 2.7 , 8.9 , 37.2); | |
315 | AliMaterial( 2, "Air$", 14.61, 7.3, 0.001205, 30420.0 , 67500.0); | |
316 | AliMaterial( 4, "Xe $", 131.29, 54.0, dxe , 1447.59, 0.0); | |
317 | AliMaterial( 5, "Cu $", 63.54, 29.0, 8.96 , 1.43, 14.8); | |
318 | AliMaterial( 6, "C $", 12.01, 6.0, 2.265 , 18.8 , 74.4); | |
319 | AliMaterial(12, "G10$", 20.00, 10.0, 1.7 , 19.4 , 999.0); | |
320 | ||
321 | // Mixtures | |
322 | AliMixture(3, "Polyethilene$", ape, zpe, dpe, -2, wpe); | |
323 | AliMixture(7, "Mylar$", amy, zmy, dmy, -3, wmy); | |
324 | AliMixture(8, "CO2$", aco, zco, dco, -2, wco); | |
325 | AliMixture(9, "Isobutane$", ais, zis, dis, -2, wis); | |
326 | AliMixture(13,"Water$", awa, zwa, dwa, -2, wwa); | |
327 | ||
328 | // Gas mixtures | |
329 | Char_t namate[21]; | |
330 | // Xe/CO2-mixture | |
331 | // Get properties of Xe | |
332 | gMC->Gfmate((*fIdmate)[4], namate, agm[0], zgm[0], d, radl, absl, buf, nbuf); | |
333 | // Get properties of CO2 | |
334 | gMC->Gfmate((*fIdmate)[8], namate, agm[1], zgm[1], d, radl, absl, buf, nbuf); | |
335 | // Create gas mixture | |
336 | wgm[0] = fxc; | |
337 | wgm[1] = 1. - fxc; | |
338 | dgm = wgm[0] * dxe + wgm[1] * dco; | |
339 | AliMixture(10, "Gas mixture 1$", agm, zgm, dgm, 2, wgm); | |
340 | // Xe/Isobutane-mixture | |
341 | // Get properties of Xe | |
342 | gMC->Gfmate((*fIdmate)[4], namate, agm[0], zgm[0], d, radl, absl, buf, nbuf); | |
343 | // Get properties of Isobutane | |
344 | gMC->Gfmate((*fIdmate)[9], namate, agm[1], zgm[1], d, radl, absl, buf, nbuf); | |
345 | // Create gas mixture | |
346 | wgm[0] = fxi; | |
347 | wgm[1] = 1. - fxi; | |
348 | dgm = wgm[0] * dxe + wgm[1] * dis; | |
349 | AliMixture(11, "Gas mixture 2$", agm, zgm, dgm, 2, wgm); | |
350 | ||
351 | ////////////////////////////////////////////////////////////////////////// | |
352 | // Tracking Media Parameters | |
353 | ////////////////////////////////////////////////////////////////////////// | |
354 | ||
355 | // Al Frame | |
356 | AliMedium(1, "Al Frame$", 1, 0, ISXFLD, SXMGMX | |
357 | , tmaxfd, stemax, deemax, epsil, stmin); | |
358 | // Air | |
359 | AliMedium(2, "Air$", 2, 0, ISXFLD, SXMGMX | |
360 | , tmaxfd, stemax, deemax, epsil, stmin); | |
361 | // Polyethilene | |
362 | AliMedium(3, "Radiator$", 3, 0, ISXFLD, SXMGMX | |
363 | , tmaxfd, stemax, deemax, epsil, stmin); | |
364 | // Xe | |
365 | AliMedium(4, "Xe$", 4, 1, ISXFLD, SXMGMX | |
366 | , tmaxfd, stemax, deemax, epsil, stmin); | |
367 | // Cu pads | |
368 | AliMedium(5, "Padplane$", 5, 1, ISXFLD, SXMGMX | |
369 | , tmaxfd, stemax, deemax, epsil, stmin); | |
370 | // Fee + cables | |
371 | AliMedium(6, "Readout$", 1, 0, ISXFLD, SXMGMX | |
372 | , tmaxfd, stemax, deemax, epsil, stmin); | |
373 | // C frame | |
374 | AliMedium(7, "C Frame$", 6, 0, ISXFLD, SXMGMX | |
375 | , tmaxfd, stemax, deemax, epsil, stmin); | |
376 | // Mylar foils | |
377 | AliMedium(8, "Mylar$", 7, 0, ISXFLD, SXMGMX | |
378 | , tmaxfd, stemax, deemax, epsil, stmin); | |
379 | if (fGasMix == 1) { | |
380 | // Gas-mixture (Xe/CO2) | |
381 | AliMedium(9, "Gas-mix$", 10, 1, ISXFLD, SXMGMX | |
382 | , tmaxfd, stemax, deemax, epsil, stmin); | |
383 | } | |
384 | else { | |
385 | // Gas-mixture (Xe/Isobutane) | |
386 | AliMedium(9, "Gas-mix$", 11, 1, ISXFLD, SXMGMX | |
387 | , tmaxfd, stemax, deemax, epsil, stmin); | |
388 | } | |
389 | // Nomex-honeycomb (use carbon for the time being) | |
390 | AliMedium(10, "Nomex$", 6, 0, ISXFLD, SXMGMX | |
391 | , tmaxfd, stemax, deemax, epsil, stmin); | |
392 | // Kapton foils (use Mylar for the time being) | |
393 | AliMedium(11, "Kapton$", 7, 0, ISXFLD, SXMGMX | |
394 | , tmaxfd, stemax, deemax, epsil, stmin); | |
395 | // Gas-filling of the radiator | |
396 | AliMedium(12, "CO2$", 8, 0, ISXFLD, SXMGMX | |
397 | , tmaxfd, stemax, deemax, epsil, stmin); | |
398 | // G10-plates | |
399 | AliMedium(13, "G10-plates$",12, 0, ISXFLD, SXMGMX | |
400 | , tmaxfd, stemax, deemax, epsil, stmin); | |
401 | // Cooling water | |
402 | AliMedium(14, "Water$", 13, 0, ISXFLD, SXMGMX | |
403 | , tmaxfd, stemax, deemax, epsil, stmin); | |
404 | ||
405 | } | |
406 | ||
407 | //_____________________________________________________________________________ | |
408 | void AliTRD::DrawModule() | |
409 | { | |
410 | // | |
411 | // Draw a shaded view of the Transition Radiation Detector version 0 | |
412 | // | |
413 | ||
414 | // Set everything unseen | |
415 | gMC->Gsatt("*" ,"SEEN",-1); | |
416 | ||
417 | // Set ALIC mother transparent | |
418 | gMC->Gsatt("ALIC","SEEN", 0); | |
419 | ||
420 | // Set the volumes visible | |
421 | if (fGeometry->IsVersion() == 0) { | |
422 | gMC->Gsatt("B071","SEEN", 0); | |
423 | gMC->Gsatt("B074","SEEN", 0); | |
424 | gMC->Gsatt("B075","SEEN", 0); | |
425 | gMC->Gsatt("B077","SEEN", 0); | |
426 | gMC->Gsatt("BTR1","SEEN", 0); | |
427 | gMC->Gsatt("BTR2","SEEN", 0); | |
428 | gMC->Gsatt("BTR3","SEEN", 0); | |
429 | gMC->Gsatt("TRD1","SEEN", 0); | |
430 | gMC->Gsatt("TRD2","SEEN", 0); | |
431 | gMC->Gsatt("TRD3","SEEN", 0); | |
432 | } | |
433 | else { | |
434 | gMC->Gsatt("B071","SEEN", 0); | |
435 | gMC->Gsatt("B074","SEEN", 0); | |
436 | gMC->Gsatt("B075","SEEN", 0); | |
437 | gMC->Gsatt("B077","SEEN", 0); | |
438 | gMC->Gsatt("BTR1","SEEN", 0); | |
439 | gMC->Gsatt("BTR2","SEEN", 0); | |
440 | gMC->Gsatt("BTR3","SEEN", 0); | |
441 | gMC->Gsatt("TRD1","SEEN", 0); | |
442 | if (fGeometry->GetPHOShole()) | |
443 | gMC->Gsatt("TRD2","SEEN", 0); | |
444 | if (fGeometry->GetRICHhole()) | |
445 | gMC->Gsatt("TRD3","SEEN", 0); | |
446 | } | |
447 | gMC->Gsatt("UCII","SEEN", 0); | |
448 | gMC->Gsatt("UCIM","SEEN", 0); | |
449 | gMC->Gsatt("UCIO","SEEN", 0); | |
450 | gMC->Gsatt("UL02","SEEN", 1); | |
451 | gMC->Gsatt("UL05","SEEN", 1); | |
452 | gMC->Gsatt("UL06","SEEN", 1); | |
453 | ||
454 | gMC->Gdopt("hide", "on"); | |
455 | gMC->Gdopt("shad", "on"); | |
456 | gMC->Gsatt("*", "fill", 7); | |
457 | gMC->SetClipBox("."); | |
458 | gMC->SetClipBox("*", 0, 2000, -2000, 2000, -2000, 2000); | |
459 | gMC->DefaultRange(); | |
460 | gMC->Gdraw("alic", 40, 30, 0, 12, 9.4, .021, .021); | |
461 | gMC->Gdhead(1111, "Transition Radiation Detector"); | |
462 | gMC->Gdman(18, 4, "MAN"); | |
463 | ||
464 | } | |
465 | ||
466 | //_____________________________________________________________________________ | |
467 | Int_t AliTRD::DistancetoPrimitive(Int_t , Int_t ) | |
468 | { | |
469 | // | |
470 | // Distance between the mouse and the TRD detector on the screen | |
471 | // Dummy routine | |
472 | ||
473 | return 9999; | |
474 | ||
475 | } | |
476 | ||
477 | //_____________________________________________________________________________ | |
478 | void AliTRD::Init() | |
479 | { | |
480 | // | |
481 | // Initialize the TRD detector after the geometry has been created | |
482 | // | |
483 | ||
484 | Int_t i; | |
485 | ||
486 | printf("\n"); | |
487 | for (i = 0; i < 35; i++) printf("*"); | |
488 | printf(" TRD_INIT "); | |
489 | for (i = 0; i < 35; i++) printf("*"); | |
490 | printf("\n"); | |
491 | printf("\n"); | |
492 | ||
493 | if (fGeometry->IsVersion() == 0) { | |
494 | printf(" Geometry for spaceframe with holes initialized.\n\n"); | |
495 | } | |
496 | else if (fGeometry->IsVersion() == 1) { | |
497 | printf(" Geometry for spaceframe without holes initialized.\n"); | |
498 | if (fGeometry->GetPHOShole()) | |
499 | printf(" Leave space in front of PHOS free.\n"); | |
500 | if (fGeometry->GetRICHhole()) | |
501 | printf(" Leave space in front of RICH free.\n"); | |
502 | printf("\n"); | |
503 | } | |
504 | ||
505 | if (fGasMix == 1) | |
506 | printf(" Gas Mixture: 90%% Xe + 10%% CO2\n\n"); | |
507 | else | |
508 | printf(" Gas Mixture: 97%% Xe + 3%% Isobutane\n\n"); | |
509 | ||
510 | } | |
511 | ||
512 | //_____________________________________________________________________________ | |
513 | void AliTRD::MakeBranch(Option_t* option) | |
514 | { | |
515 | // | |
516 | // Create Tree branches for the TRD digits and cluster. | |
517 | // | |
518 | ||
519 | Int_t buffersize = 4000; | |
520 | Char_t branchname[15]; | |
521 | ||
522 | AliDetector::MakeBranch(option); | |
523 | ||
524 | Char_t *R = strstr(option,"R"); | |
525 | sprintf(branchname,"%srecPoints",GetName()); | |
526 | if (fRecPoints && gAlice->TreeR() && R) { | |
527 | gAlice->TreeR()->Branch(branchname,fRecPoints->IsA()->GetName() | |
528 | ,&fRecPoints,buffersize,0); | |
529 | printf("* AliTRD::MakeBranch * Making Branch %s for points in TreeR\n",branchname); | |
530 | } | |
531 | ||
532 | } | |
533 | ||
534 | //_____________________________________________________________________________ | |
535 | void AliTRD::ResetRecPoints() | |
536 | { | |
537 | // | |
538 | // Reset number of reconstructed points and the point array | |
539 | // | |
540 | ||
541 | fNRecPoints = 0; | |
542 | if (fRecPoints) fRecPoints->Delete(); | |
543 | ||
544 | } | |
545 | ||
546 | //_____________________________________________________________________________ | |
547 | void AliTRD::SetTreeAddress() | |
548 | { | |
549 | // | |
550 | // Set the branch addresses for the trees. | |
551 | // | |
552 | ||
553 | Char_t branchname[15]; | |
554 | ||
555 | AliDetector::SetTreeAddress(); | |
556 | ||
557 | TBranch *branch; | |
558 | TTree *treeR = gAlice->TreeR(); | |
559 | ||
560 | if (treeR) { | |
561 | sprintf(branchname,"%srecPoints",GetName()); | |
562 | if (fRecPoints) { | |
563 | branch = treeR->GetBranch(branchname); | |
564 | if (branch) { | |
565 | branch->SetAddress(&fRecPoints); | |
566 | } | |
567 | } | |
568 | } | |
569 | ||
570 | } | |
571 | ||
572 | //_____________________________________________________________________________ | |
573 | void AliTRD::SetGasMix(Int_t imix) | |
574 | { | |
575 | // | |
576 | // Defines the gas mixture (imix=0: Xe/Isobutane imix=1: Xe/CO2) | |
577 | // | |
578 | ||
579 | if ((imix < 0) || (imix > 1)) { | |
580 | printf("Wrong input value: %d\n",imix); | |
581 | printf("Use standard setting\n"); | |
582 | fGasMix = 0; | |
583 | return; | |
584 | } | |
585 | ||
586 | fGasMix = imix; | |
587 | ||
588 | } | |
589 |