a12f71ff2562237b24e41965d9269c873262e569
[u/mrichter/AliRoot.git] / TRD / AliTRD.cxx
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.26  2001/05/07 08:03:22  cblume
19 Generate also hits in the amplification region
20
21 Revision 1.25  2001/03/13 09:30:35  cblume
22 Update of digitization. Moved digit branch definition to AliTRD
23
24 Revision 1.24  2001/01/26 19:56:49  hristov
25 Major upgrade of AliRoot code
26
27 Revision 1.23  2000/11/01 14:53:20  cblume
28 Merge with TRD-develop
29
30
31 Revision 1.17.2.6  2000/10/15 23:29:08  cblume
32 Introduced more detailed geometry for the display
33
34 Revision 1.17.2.5  2000/10/06 16:49:46  cblume
35 Made Getters const
36
37 Revision 1.17.2.4  2000/10/04 16:34:57  cblume
38 Replace include files by forward declarations
39
40 Revision 1.17.2.3  2000/09/22 14:45:17  cblume
41 Included changes for the tracking
42
43 Revision 1.17.2.2  2000/09/18 13:25:13  cblume
44 Included LoadPoints() method to display the TR photons
45
46 Revision 1.22  2000/10/02 21:28:19  fca
47 Removal of useless dependecies via forward declarations
48
49 Revision 1.21  2000/06/09 11:10:07  cblume
50 Compiler warnings and coding conventions, next round
51
52 Revision 1.20  2000/06/08 18:32:57  cblume
53 Make code compliant to coding conventions
54
55 Revision 1.19  2000/06/07 16:25:37  cblume
56 Try to remove compiler warnings on Sun and HP
57
58 Revision 1.18  2000/05/08 16:17:27  cblume
59 Merge TRD-develop
60
61 Revision 1.21  2000/06/09 11:10:07  cblume
62 Compiler warnings and coding conventions, next round
63
64 Revision 1.20  2000/06/08 18:32:57  cblume
65 Make code compliant to coding conventions
66
67 Revision 1.19  2000/06/07 16:25:37  cblume
68 Try to remove compiler warnings on Sun and HP
69
70 Revision 1.18  2000/05/08 16:17:27  cblume
71 Merge TRD-develop
72
73 Revision 1.17.2.1  2000/05/08 14:28:59  cblume
74 Introduced SetPHOShole() and SetRICHhole(). AliTRDrecPoint container is now a TObjArray
75
76 Revision 1.17  2000/02/28 19:10:26  cblume
77 Include the new TRD classes
78
79 Revision 1.16.2.2  2000/02/28 17:53:24  cblume
80 Introduce TRD geometry classes
81
82 Revision 1.16.2.1  2000/02/28 17:04:19  cblume
83 Include functions and data members for AliTRDrecPoint
84
85 Revision 1.16  2000/01/19 17:17:35  fca
86 Introducing a list of lists of hits -- more hits allowed for detector now
87
88 Revision 1.15  1999/11/02 17:04:25  fca
89 Small syntax change for HP compiler
90
91 Revision 1.14  1999/11/02 16:57:02  fca
92 Avoid non ansi warnings on HP compilers
93
94 Revision 1.13  1999/11/02 16:35:56  fca
95 New version of TRD introduced
96
97 Revision 1.12  1999/11/01 20:41:51  fca
98 Added protections against using the wrong version of FRAME
99
100 Revision 1.11  1999/09/29 09:24:34  fca
101 Introduction of the Copyright and cvs Log
102
103 */
104
105 ///////////////////////////////////////////////////////////////////////////////
106 //                                                                           //
107 //  Transition Radiation Detector                                            //
108 //  This class contains the basic functions for the Transition Radiation     //
109 //  Detector.                                                                //
110 //                                                                           //
111 ///////////////////////////////////////////////////////////////////////////////
112
113 #include <stdlib.h>
114 #include <iostream.h>
115
116 #include <TMath.h>
117 #include <TNode.h>
118 #include <TGeometry.h>
119 #include <TTree.h>                                                              
120 #include <TPGON.h> 
121
122 #include "AliRun.h"
123 #include "AliConst.h"
124 #include "AliDigit.h"
125 #include "AliMagF.h"
126 #include "AliMC.h"                                                              
127  
128 #include "AliTRD.h"
129 #include "AliTRDhit.h"
130 #include "AliTRDpoints.h"
131 #include "AliTRDdigit.h"
132 #include "AliTRDdigitizer.h"
133 #include "AliTRDclusterizer.h"
134 #include "AliTRDgeometryHole.h"
135 #include "AliTRDgeometryFull.h"
136 #include "AliTRDrecPoint.h"
137 #include "AliTRDcluster.h"
138 #include "AliTRDdigitsManager.h"
139 #include "AliTRDdataArrayI.h"
140 #include "AliTRDsegmentArray.h"
141
142 ClassImp(AliTRD)
143  
144 //_____________________________________________________________________________
145 AliTRD::AliTRD()
146 {
147   //
148   // Default constructor
149   //
150
151   fIshunt        = 0;
152   fGasMix        = 0;
153   fHits          = 0;
154   fDigits        = 0;
155
156   fRecPoints     = 0;
157   fNRecPoints    = 0;
158
159   fGeometry      = 0;
160
161   fGasDensity    = 0;
162   fFoilDensity   = 0;
163
164   fDrawTR        = 0;
165   fDisplayType   = 0; 
166
167   fDigitsArray   = 0; 
168   for (Int_t iDict = 0; iDict < AliTRDdigitsManager::NDict(); iDict++) {
169     fDictionaryArray[iDict] = 0; 
170   }
171
172 }
173  
174 //_____________________________________________________________________________
175 AliTRD::AliTRD(const char *name, const char *title)
176        : AliDetector(name,title)
177 {
178   //
179   // Standard constructor for the TRD
180   //
181
182   // Check that FRAME is there otherwise we have no place where to
183   // put TRD
184   AliModule* frame = gAlice->GetModule("FRAME");
185   if (!frame) {
186     Error("Ctor","TRD needs FRAME to be present\n");
187     exit(1);
188   } 
189
190   // Define the TRD geometry according to the FRAME geometry
191   if      (frame->IsVersion() == 0) {
192     // Geometry with hole
193     fGeometry = new AliTRDgeometryHole();
194   }
195   else if (frame->IsVersion() == 1) {
196     // Geometry without hole
197     fGeometry = new AliTRDgeometryFull();
198   }
199   else {
200     Error("Ctor","Could not find valid FRAME version\n");
201     exit(1);
202   }
203
204   // Allocate the hit array
205   fHits           = new TClonesArray("AliTRDhit"     ,405);
206   gAlice->AddHitList(fHits);
207
208   // Allocate the digits array
209   fDigits         = 0;
210
211   // Allocate the rec point array
212   fRecPoints     = new TObjArray(400);
213   fNRecPoints    = 0;
214    
215   fIshunt        = 0;
216   fGasMix        = 1;
217
218   fGasDensity    = 0;
219   fFoilDensity   = 0;
220
221   fDrawTR        = 0;
222   fDisplayType   = 0;
223
224   fDigitsArray   = 0; 
225   for (Int_t iDict = 0; iDict < AliTRDdigitsManager::NDict(); iDict++) {
226     fDictionaryArray[iDict] = 0; 
227   }
228
229   SetMarkerColor(kWhite);   
230
231 }
232
233 //_____________________________________________________________________________
234 AliTRD::AliTRD(const AliTRD &trd)
235 {
236   //
237   // Copy constructor
238   //
239
240   ((AliTRD &) trd).Copy(*this);
241
242 }
243
244 //_____________________________________________________________________________
245 AliTRD::~AliTRD()
246 {
247   //
248   // TRD destructor
249   //
250
251   fIshunt = 0;
252
253   delete fGeometry;
254   delete fHits;
255   delete fRecPoints;
256   if (fDigitsArray) delete fDigitsArray;
257   for (Int_t iDict = 0; iDict < AliTRDdigitsManager::NDict(); iDict++) {
258     if (fDictionaryArray[iDict]) delete fDictionaryArray[iDict];
259   }
260
261 }
262
263 //_____________________________________________________________________________
264 void AliTRD::AddCluster(Float_t *pos, Int_t *digits, Int_t det, Float_t amp
265                        , Int_t *tracks, Int_t iType)
266 {
267   //
268   // Add a cluster for the TRD
269   //
270
271   Int_t   pla = fGeometry->GetPlane(det);
272   Int_t   cha = fGeometry->GetChamber(det);
273   Int_t   sec = fGeometry->GetSector(det);
274
275   Float_t padRow  = pos[0]; 
276   Float_t padCol  = pos[1]; 
277   Float_t row0    = fGeometry->GetRow0(pla,cha,sec);
278   Float_t col0    = fGeometry->GetCol0(pla);
279   Float_t rowSize = fGeometry->GetRowPadSize(pla,cha,sec);
280   Float_t colSize = fGeometry->GetColPadSize(pla);
281
282   AliTRDcluster *c = new AliTRDcluster();
283
284   c->SetDetector(det);
285   c->AddTrackIndex(tracks);
286
287   c->SetQ(amp);
288
289   c->SetLocalTimeBin(((Int_t) pos[2]));
290   c->SetY(- (col0 + padCol * colSize));
291   c->SetZ(   row0 + padRow * rowSize);
292   
293   c->SetSigmaY2(0.05 * 0.05);
294   c->SetSigmaZ2(rowSize * rowSize / 12.);
295
296   if (iType == 1) {
297     c->SetUnfolding();
298   }
299
300   fRecPoints->Add(c);
301
302 }
303
304 //_____________________________________________________________________________
305 void AliTRD::Hits2Digits()
306 {
307   //
308   // Create digits
309   //
310
311   AliTRDdigitizer *digitizer = new AliTRDdigitizer("TRDdigitizer"
312                                                   ,"TRD digitizer class");
313
314   // Set the parameter
315   digitizer->SetDiffusion();
316   digitizer->SetExB();
317   digitizer->SetEvent(gAlice->GetEvNumber());
318   // Initialization
319   //digitizer->InitDetector();
320     
321   // Create the digits
322   digitizer->MakeDigits();
323   
324   // Write the digits into the input file
325   if (digitizer->Digits()->MakeBranch(fDigitsFile)) {
326
327     digitizer->WriteDigits();
328
329     // Save the digitizer class in the AliROOT 
330     digitizer->Write();
331
332   }
333
334 }
335
336 //_____________________________________________________________________________
337 void AliTRD::Hits2SDigits()
338 {
339   //
340   // Create summable digits
341   //
342
343   AliTRDdigitizer *digitizer = new AliTRDdigitizer("TRDdigitizer"
344                                                   ,"TRD digitizer class");
345
346   // For the summable digits
347   digitizer->SetSDigits(kTRUE);
348
349   // Set the parameter
350   digitizer->SetDiffusion();
351   digitizer->SetExB();
352   digitizer->SetEvent(gAlice->GetEvNumber());
353
354   // Initialization
355   //digitizer->InitDetector();
356     
357   // Create the digits
358   digitizer->MakeDigits();
359   
360   // Write the digits into the input file
361   if (digitizer->Digits()->MakeBranch(fDigitsFile)) {
362
363     digitizer->WriteDigits();
364
365     // Save the digitizer class in the AliROOT 
366     digitizer->Write();
367
368   }
369
370 }
371
372 //_____________________________________________________________________________
373 void AliTRD::SDigits2Digits()
374 {
375   //
376   // Create final digits from summable digits
377   //
378
379 }
380
381 //_____________________________________________________________________________
382 void AliTRD::AddDigit(Int_t *digits, Int_t *amp)
383 {
384   //
385   // Add a digit for the TRD
386   //
387
388   TClonesArray &ldigits = *fDigits;
389   new(ldigits[fNdigits++]) AliTRDdigit(kFALSE,digits,amp);
390
391 }
392
393 //_____________________________________________________________________________
394 void AliTRD::AddHit(Int_t track, Int_t det, Float_t *hits, Int_t q)
395 {
396   //
397   // Add a hit for the TRD
398   //
399
400   TClonesArray &lhits = *fHits;
401   new(lhits[fNhits++]) AliTRDhit(fIshunt,track,det,hits,q);
402
403 }
404
405 //_____________________________________________________________________________
406 void AliTRD::BuildGeometry()
407 {
408   //
409   // Create the ROOT TNode geometry for the TRD
410   //
411
412   TNode *node, *top;
413   TPGON *pgon;
414
415   Float_t rmin, rmax;
416   Float_t zmax1, zmax2;
417  
418   const Int_t kColorTRD = 46;
419   
420   // Find the top node alice
421   top = gAlice->GetGeometry()->GetNode("alice");
422   
423   switch (fDisplayType) {
424
425   case 0:
426
427     pgon = new TPGON("S_TRD","TRD","void",0,360,AliTRDgeometry::Nsect(),4);
428     rmin = AliTRDgeometry::Rmin();
429     rmax = AliTRDgeometry::Rmax();
430     pgon->DefineSection(0,-AliTRDgeometry::Zmax1(),rmax,rmax);
431     pgon->DefineSection(1,-AliTRDgeometry::Zmax2(),rmin,rmax);
432     pgon->DefineSection(2, AliTRDgeometry::Zmax2(),rmin,rmax);
433     pgon->DefineSection(3, AliTRDgeometry::Zmax1(),rmax,rmax);
434     top->cd();
435     node = new TNode("TRD","TRD","S_TRD",0,0,0,"");
436     node->SetLineColor(kColorTRD);
437     fNodes->Add(node);
438
439     break;
440
441   case 1:
442
443     Char_t name[7];
444
445     Float_t slope = (AliTRDgeometry::Zmax1() - AliTRDgeometry::Zmax2())
446                   / (AliTRDgeometry::Rmax()  - AliTRDgeometry::Rmin());
447
448     rmin  = AliTRDgeometry::Rmin() + AliTRDgeometry::Ccframe()/2.
449                                    + AliTRDgeometry::DrZpos()
450                                    - AliTRDgeometry::DrThick()/2.;
451     rmax  = rmin + AliTRDgeometry::DrThick();
452
453     Float_t thickness = rmin - AliTRDgeometry::Rmin();
454     zmax2 = AliTRDgeometry::Zmax2() + slope * thickness;
455     zmax1 = zmax2 + slope * AliTRDgeometry::DrThick();
456
457     for (Int_t iPlan = 0; iPlan < AliTRDgeometry::Nplan(); iPlan++) {
458
459       sprintf(name,"S_TR1%d",iPlan);
460       pgon  = new TPGON(name,"TRD","void",0,360,AliTRDgeometry::Nsect(),4);
461       pgon->DefineSection(0,-zmax1,rmax,rmax);
462       pgon->DefineSection(1,-zmax2,rmin,rmax);
463       pgon->DefineSection(2, zmax2,rmin,rmax);
464       pgon->DefineSection(3, zmax1,rmax,rmax);
465       top->cd();
466       node = new TNode("TRD","TRD",name,0,0,0,"");
467       node->SetLineColor(kColorTRD);
468       fNodes->Add(node);
469
470       Float_t height = AliTRDgeometry::Cheight() + AliTRDgeometry::Cspace(); 
471       rmin  = rmin  + height;
472       rmax  = rmax  + height;
473       zmax1 = zmax1 + slope * height;
474       zmax2 = zmax2 + slope * height;
475
476     }
477
478     thickness += AliTRDgeometry::DrThick();
479     rmin  = AliTRDgeometry::Rmin() + thickness;
480     rmax  = rmin + AliTRDgeometry::AmThick();
481     zmax2 = AliTRDgeometry::Zmax2() + slope * thickness;
482     zmax1 = zmax2 + slope * AliTRDgeometry::AmThick();
483
484     for (Int_t iPlan = 0; iPlan < AliTRDgeometry::Nplan(); iPlan++) {
485
486       sprintf(name,"S_TR2%d",iPlan);
487       pgon  = new TPGON(name,"TRD","void",0,360,AliTRDgeometry::Nsect(),4);
488       pgon->DefineSection(0,-zmax1,rmax,rmax);
489       pgon->DefineSection(1,-zmax2,rmin,rmax);
490       pgon->DefineSection(2, zmax2,rmin,rmax);
491       pgon->DefineSection(3, zmax1,rmax,rmax);
492       top->cd();
493       node = new TNode("TRD","TRD",name,0,0,0,"");
494       node->SetLineColor(kColorTRD);
495       fNodes->Add(node);
496
497       Float_t height = AliTRDgeometry::Cheight() + AliTRDgeometry::Cspace(); 
498       rmin  = rmin  + height;
499       rmax  = rmax  + height;
500       zmax1 = zmax1 + slope * height;
501       zmax2 = zmax2 + slope * height;
502
503     }
504
505     break;
506
507   };
508
509 }
510  
511 //_____________________________________________________________________________
512 void AliTRD::Copy(TObject &trd)
513 {
514   //
515   // Copy function
516   //
517
518   ((AliTRD &) trd).fGasMix      = fGasMix;
519   ((AliTRD &) trd).fGeometry    = fGeometry;       
520   ((AliTRD &) trd).fRecPoints   = fRecPoints;
521   ((AliTRD &) trd).fNRecPoints  = fNRecPoints;
522   ((AliTRD &) trd).fGasDensity  = fGasDensity;
523   ((AliTRD &) trd).fFoilDensity = fFoilDensity;
524   ((AliTRD &) trd).fDrawTR      = fDrawTR;
525   ((AliTRD &) trd).fDisplayType = fDisplayType;
526
527   //AliDetector::Copy(trd);
528
529 }
530
531 //_____________________________________________________________________________
532 void AliTRD::CreateGeometry()
533 {
534   //
535   // Creates the volumes for the TRD chambers
536   //
537
538   // Check that FRAME is there otherwise we have no place where to put the TRD
539   AliModule* frame = gAlice->GetModule("FRAME");
540   if (!frame) {
541     printf(" The TRD needs the FRAME to be defined first\n");
542     return;
543   }
544
545   fGeometry->CreateGeometry(fIdtmed->GetArray() - 1299);
546
547 }
548  
549 //_____________________________________________________________________________
550 void AliTRD::CreateMaterials()
551 {
552   //
553   // Create the materials for the TRD
554   // Origin Y.Foka
555   //
556
557   Int_t   isxfld = gAlice->Field()->Integ();
558   Float_t sxmgmx = gAlice->Field()->Max();
559   
560   // For polyethilene (CH2) 
561   Float_t ape[2] = { 12., 1. };
562   Float_t zpe[2] = {  6., 1. };
563   Float_t wpe[2] = {  1., 2. };
564   Float_t dpe    = 0.95;
565
566   // For mylar (C5H4O2) 
567   Float_t amy[3] = { 12., 1., 16. };
568   Float_t zmy[3] = {  6., 1.,  8. };
569   Float_t wmy[3] = {  5., 4.,  2. };
570   Float_t dmy    = 1.39;
571
572   // For CO2 
573   Float_t aco[2] = { 12., 16. };
574   Float_t zco[2] = {  6.,  8. };
575   Float_t wco[2] = {  1.,  2. };
576   Float_t dco    = 0.001977;
577
578   // For water
579   Float_t awa[2] = {  1., 16. };
580   Float_t zwa[2] = {  1.,  8. };
581   Float_t wwa[2] = {  2.,  1. };
582   Float_t dwa    = 1.0;
583
584   // For isobutane (C4H10)
585   Float_t ais[2] = { 12.,  1. };
586   Float_t zis[2] = {  6.,  1. };
587   Float_t wis[2] = {  4., 10. };
588   Float_t dis    = 0.00267;
589
590   // For Xe/CO2-gas-mixture 
591   // Xe-content of the Xe/CO2-mixture (90% / 10%) 
592   Float_t fxc    = .90;
593   // Xe-content of the Xe/Isobutane-mixture (97% / 3%) 
594   Float_t fxi    = .97;
595   Float_t dxe    = .005858;
596   
597   // General tracking parameter
598   Float_t tmaxfd = -10.;
599   Float_t stemax = -1e10;
600   Float_t deemax = -0.1;
601   Float_t epsil  =  1e-4;
602   Float_t stmin  = -0.001;
603   
604   Float_t absl, radl, d, buf[1];
605   Float_t agm[2], zgm[2], wgm[2];
606   Float_t dgm1, dgm2;
607   Int_t   nbuf;
608   
609   //////////////////////////////////////////////////////////////////////////
610   //     Define Materials 
611   //////////////////////////////////////////////////////////////////////////
612
613   AliMaterial( 1, "Al $",  26.98, 13.0, 2.7     ,     8.9 ,    37.2);
614   AliMaterial( 2, "Air$",  14.61,  7.3, 0.001205, 30420.0 , 67500.0);
615   AliMaterial( 4, "Xe $", 131.29, 54.0, dxe     ,  1447.59,     0.0);
616   AliMaterial( 5, "Cu $",  63.54, 29.0, 8.96    ,     1.43,    14.8);
617   AliMaterial( 6, "C  $",  12.01,  6.0, 2.265   ,    18.8 ,    74.4);
618   AliMaterial(12, "G10$",  20.00, 10.0, 1.7     ,    19.4 ,   999.0);
619
620   // Mixtures 
621   AliMixture(3, "Polyethilene$",   ape, zpe, dpe, -2, wpe);
622   AliMixture(7, "Mylar$",          amy, zmy, dmy, -3, wmy);
623   AliMixture(8, "CO2$",            aco, zco, dco, -2, wco);
624   AliMixture(9, "Isobutane$",      ais, zis, dis, -2, wis);
625   AliMixture(13,"Water$",          awa, zwa, dwa, -2, wwa);
626
627   // Gas mixtures
628   Char_t namate[21];
629   // Xe/CO2-mixture
630   // Get properties of Xe 
631   gMC->Gfmate((*fIdmate)[4], namate, agm[0], zgm[0], d, radl, absl, buf, nbuf);
632   // Get properties of CO2 
633   gMC->Gfmate((*fIdmate)[8], namate, agm[1], zgm[1], d, radl, absl, buf, nbuf);
634   // Create gas mixture 
635   wgm[0] = fxc;
636   wgm[1] = 1. - fxc;
637   dgm1   = wgm[0] * dxe + wgm[1] * dco;
638   AliMixture(10, "Gas mixture 1$", agm, zgm, dgm1,  2, wgm);
639   // Xe/Isobutane-mixture
640   // Get properties of Xe 
641   gMC->Gfmate((*fIdmate)[4], namate, agm[0], zgm[0], d, radl, absl, buf, nbuf);
642   // Get properties of Isobutane
643   gMC->Gfmate((*fIdmate)[9], namate, agm[1], zgm[1], d, radl, absl, buf, nbuf);
644   // Create gas mixture 
645   wgm[0] = fxi;
646   wgm[1] = 1. - fxi;
647   dgm2   = wgm[0] * dxe + wgm[1] * dis;
648   AliMixture(11, "Gas mixture 2$", agm, zgm, dgm2,  2, wgm);
649  
650   //////////////////////////////////////////////////////////////////////////
651   //     Tracking Media Parameters 
652   //////////////////////////////////////////////////////////////////////////
653
654   // Al Frame 
655   AliMedium(1, "Al Frame$",   1, 0, isxfld, sxmgmx
656                 , tmaxfd, stemax, deemax, epsil, stmin);
657   // Air 
658   AliMedium(2, "Air$",        2, 0, isxfld, sxmgmx
659                 , tmaxfd, stemax, deemax, epsil, stmin);
660   // Polyethilene 
661   AliMedium(3, "Radiator$",   3, 0, isxfld, sxmgmx
662                 , tmaxfd, stemax, deemax, epsil, stmin);
663   // Xe 
664   AliMedium(4, "Xe$",         4, 1, isxfld, sxmgmx
665                 , tmaxfd, stemax, deemax, epsil, stmin);
666   // Cu pads 
667   AliMedium(5, "Padplane$",   5, 1, isxfld, sxmgmx
668                 , tmaxfd, stemax, deemax, epsil, stmin);
669   // Fee + cables 
670   AliMedium(6, "Readout$",    1, 0, isxfld, sxmgmx
671                 , tmaxfd, stemax, deemax, epsil, stmin);
672   // C frame 
673   AliMedium(7, "C Frame$",    6, 0, isxfld, sxmgmx
674                 , tmaxfd, stemax, deemax, epsil, stmin);
675   // Mylar foils 
676   AliMedium(8, "Mylar$",      7, 0, isxfld, sxmgmx
677                 , tmaxfd, stemax, deemax, epsil, stmin);
678   if (fGasMix == 1) {
679     // Gas-mixture (Xe/CO2) 
680     AliMedium(9, "Gas-mix$",   10, 1, isxfld, sxmgmx
681                   , tmaxfd, stemax, deemax, epsil, stmin);
682   }
683   else {
684     // Gas-mixture (Xe/Isobutane) 
685     AliMedium(9, "Gas-mix$",   11, 1, isxfld, sxmgmx
686                   , tmaxfd, stemax, deemax, epsil, stmin);
687   }
688   // Nomex-honeycomb (use carbon for the time being) 
689   AliMedium(10, "Nomex$",      6, 0, isxfld, sxmgmx
690                 , tmaxfd, stemax, deemax, epsil, stmin);
691   // Kapton foils (use Mylar for the time being) 
692   AliMedium(11, "Kapton$",     7, 0, isxfld, sxmgmx
693                 , tmaxfd, stemax, deemax, epsil, stmin);
694   // Gas-filling of the radiator 
695   AliMedium(12, "CO2$",        8, 0, isxfld, sxmgmx
696                 , tmaxfd, stemax, deemax, epsil, stmin);
697   // G10-plates
698   AliMedium(13, "G10-plates$",12, 0, isxfld, sxmgmx
699                 , tmaxfd, stemax, deemax, epsil, stmin);
700   // Cooling water
701   AliMedium(14, "Water$",     13, 0, isxfld, sxmgmx
702                 , tmaxfd, stemax, deemax, epsil, stmin);
703
704   // Save the density values for the TRD absorbtion
705   fFoilDensity = dmy;
706   if (fGasMix == 1)
707     fGasDensity = dgm1;
708   else
709     fGasDensity = dgm2;
710
711 }
712
713 //_____________________________________________________________________________
714 void AliTRD::DrawModule()
715 {
716   //
717   // Draw a shaded view of the Transition Radiation Detector version 0
718   //
719
720   // Set everything unseen
721   gMC->Gsatt("*"   ,"SEEN",-1);
722   
723   // Set ALIC mother transparent
724   gMC->Gsatt("ALIC","SEEN", 0);
725   
726   // Set the volumes visible
727   if (fGeometry->IsVersion() == 0) {
728     gMC->Gsatt("B071","SEEN", 0);
729     gMC->Gsatt("B074","SEEN", 0);
730     gMC->Gsatt("B075","SEEN", 0);
731     gMC->Gsatt("B077","SEEN", 0);
732     gMC->Gsatt("BTR1","SEEN", 0);
733     gMC->Gsatt("BTR2","SEEN", 0);
734     gMC->Gsatt("BTR3","SEEN", 0);
735     gMC->Gsatt("TRD1","SEEN", 0);
736     gMC->Gsatt("TRD2","SEEN", 0);
737     gMC->Gsatt("TRD3","SEEN", 0);
738   }
739   else {
740     gMC->Gsatt("B071","SEEN", 0);
741     gMC->Gsatt("B074","SEEN", 0);
742     gMC->Gsatt("B075","SEEN", 0);
743     gMC->Gsatt("B077","SEEN", 0);
744     gMC->Gsatt("BTR1","SEEN", 0);
745     gMC->Gsatt("BTR2","SEEN", 0);
746     gMC->Gsatt("BTR3","SEEN", 0);
747     gMC->Gsatt("TRD1","SEEN", 0);
748     if (fGeometry->GetPHOShole())
749       gMC->Gsatt("TRD2","SEEN", 0);
750     if (fGeometry->GetRICHhole())
751       gMC->Gsatt("TRD3","SEEN", 0);
752   }
753   gMC->Gsatt("UCII","SEEN", 0);
754   gMC->Gsatt("UCIM","SEEN", 0);
755   gMC->Gsatt("UCIO","SEEN", 0);
756   gMC->Gsatt("UL02","SEEN", 1);
757   gMC->Gsatt("UL05","SEEN", 1);
758   gMC->Gsatt("UL06","SEEN", 1);
759   
760   gMC->Gdopt("hide", "on");
761   gMC->Gdopt("shad", "on");
762   gMC->Gsatt("*", "fill", 7);
763   gMC->SetClipBox(".");
764   gMC->SetClipBox("*", 0, 2000, -2000, 2000, -2000, 2000);
765   gMC->DefaultRange();
766   gMC->Gdraw("alic", 40, 30, 0, 12, 9.4, .021, .021);
767   gMC->Gdhead(1111, "Transition Radiation Detector");
768   gMC->Gdman(18, 4, "MAN");
769
770 }
771
772 //_____________________________________________________________________________
773 Int_t AliTRD::DistancetoPrimitive(Int_t , Int_t )
774 {
775   //
776   // Distance between the mouse and the TRD detector on the screen
777   // Dummy routine
778   
779   return 9999;
780
781 }
782  
783 //_____________________________________________________________________________
784 void AliTRD::Init()
785 {
786   //
787   // Initialize the TRD detector after the geometry has been created
788   //
789
790   Int_t i;
791
792   if(fDebug) {
793     printf("\n%s: ",ClassName());
794     for (i = 0; i < 35; i++) printf("*");
795     printf(" TRD_INIT ");
796     for (i = 0; i < 35; i++) printf("*");
797     printf("\n");
798   }
799
800   if      (fGeometry->IsVersion() == 0) {
801     printf("%s: Geometry for spaceframe with holes initialized\n",ClassName());
802   }
803   else if (fGeometry->IsVersion() == 1) {
804     printf("%s: Geometry for spaceframe without holes initialized\n",ClassName());
805     if (fGeometry->GetPHOShole())
806       printf("%s: Leave space in front of PHOS free\n",ClassName());
807     if (fGeometry->GetRICHhole())
808       printf("%s: Leave space in front of RICH free\n",ClassName());
809   }
810   
811   if (fGasMix == 1)
812     printf("%s: Gas Mixture: 90%% Xe + 10%% CO2\n",ClassName());
813   else
814     printf("%s: Gas Mixture: 97%% Xe + 3%% Isobutane\n",ClassName());
815
816 }
817
818 //_____________________________________________________________________________
819 void AliTRD::LoadPoints(Int_t track)
820 {
821   //
822   // Store x, y, z of all hits in memory.
823   // Hit originating from TR photons are given a different color
824   //
825
826   if (!fDrawTR) {
827     AliDetector::LoadPoints(track);
828     return;
829   }
830
831   if (fHits == 0) return;
832
833   Int_t nhits  = fHits->GetEntriesFast();
834   if (nhits == 0) return;
835
836   Int_t tracks = gAlice->GetNtrack();
837   if (fPoints == 0) fPoints = new TObjArray(tracks);
838
839   AliTRDhit *ahit;
840   
841   Int_t    *ntrkE = new Int_t[tracks];
842   Int_t    *ntrkT = new Int_t[tracks];
843   Int_t    *limiE = new Int_t[tracks];
844   Int_t    *limiT = new Int_t[tracks];
845   Float_t **coorE = new Float_t*[tracks];
846   Float_t **coorT = new Float_t*[tracks];
847   for(Int_t i = 0; i < tracks; i++) {
848     ntrkE[i] = 0;
849     ntrkT[i] = 0;
850     coorE[i] = 0;
851     coorT[i] = 0;
852     limiE[i] = 0;
853     limiT[i] = 0;
854   }
855   
856   AliTRDpoints  *points = 0;
857   Float_t       *fp     = 0;
858   Int_t          trk;
859   Int_t          chunk  = nhits / 4 + 1;
860
861   // Loop over all the hits and store their position
862   for (Int_t hit = 0; hit < nhits; hit++) {
863
864     ahit = (AliTRDhit *) fHits->UncheckedAt(hit);
865
866     // dEdx hits
867     if (ahit->FromDrift() || ahit->FromAmplification()) {
868
869       trk = ahit->GetTrack();
870       if (ntrkE[trk] == limiE[trk]) {
871         // Initialise a new track
872         fp = new Float_t[3*(limiE[trk]+chunk)];
873         if (coorE[trk]) {
874           memcpy(fp,coorE[trk],sizeof(Float_t)*3*limiE[trk]);
875           delete [] coorE[trk];
876         }
877         limiE[trk] += chunk;
878         coorE[trk]  = fp;
879       } 
880       else {
881         fp = coorE[trk];
882       }
883       fp[3*ntrkE[trk]  ] = ahit->X();
884       fp[3*ntrkE[trk]+1] = ahit->Y();
885       fp[3*ntrkE[trk]+2] = ahit->Z();
886       ntrkE[trk]++;
887
888     }
889     // TR photon hits
890     else if (ahit->FromTRphoton()) {
891
892       trk = ahit->GetTrack();
893       if (ntrkT[trk] == limiT[trk]) {
894         // Initialise a new track
895         fp = new Float_t[3*(limiT[trk]+chunk)];
896         if (coorT[trk]) {
897           memcpy(fp,coorT[trk],sizeof(Float_t)*3*limiT[trk]);
898           delete [] coorT[trk];
899         }
900         limiT[trk] += chunk;
901         coorT[trk]  = fp;
902       } 
903       else {
904         fp = coorT[trk];
905       }
906       fp[3*ntrkT[trk]  ] = ahit->X();
907       fp[3*ntrkT[trk]+1] = ahit->Y();
908       fp[3*ntrkT[trk]+2] = ahit->Z();
909       ntrkT[trk]++;
910
911     }
912
913   }
914
915   for (trk = 0; trk < tracks; ++trk) {
916
917     if (ntrkE[trk] || ntrkT[trk]) {
918
919       points = new AliTRDpoints();
920       points->SetDetector(this);
921       points->SetParticle(trk);
922
923       // Set the dEdx points
924       if (ntrkE[trk]) {
925         points->SetMarkerColor(GetMarkerColor());
926         points->SetMarkerSize(GetMarkerSize());
927         points->SetPolyMarker(ntrkE[trk],coorE[trk],GetMarkerStyle());
928         delete [] coorE[trk];
929         coorE[trk] = 0;
930       }
931
932       // Set the TR photon points
933       if (ntrkT[trk]) {
934         points->SetTRpoints(ntrkT[trk],coorT[trk]);
935         delete [] coorT[trk];
936         coorT[trk] = 0;
937       }
938
939       fPoints->AddAt(points,trk);
940
941     }
942
943   }
944
945   delete [] coorE;
946   delete [] coorT;
947   delete [] ntrkE;
948   delete [] ntrkT;
949   delete [] limiE;
950   delete [] limiT;
951
952 }
953
954 //_____________________________________________________________________________
955 void AliTRD::MakeBranch(Option_t* option, const char *file)
956 {
957   //
958   // Create Tree branches for the TRD digits and cluster.
959   //
960
961   //Int_t  buffersize = 4000;
962   //Char_t branchname[15];
963
964   AliDetector::MakeBranch(option,file);
965
966   Int_t buffersize = 64000;
967
968   fDigitsArray = new AliTRDdataArrayI();
969   MakeBranchInTree(gAlice->TreeD() 
970                    ,"TRDdigits", fDigitsArray->IsA()->GetName()
971                    ,&fDigitsArray,buffersize,1,file);
972
973   for (Int_t iDict = 0; iDict < AliTRDdigitsManager::NDict(); iDict++) {
974     Char_t branchname[15];
975     sprintf(branchname,"TRDdictionary%d",iDict);
976     fDictionaryArray[iDict] = new AliTRDdataArrayI();
977     MakeBranchInTree(gAlice->TreeD() 
978                      ,branchname,fDictionaryArray[iDict]->IsA()->GetName()
979                      ,&fDictionaryArray[iDict],buffersize,1,file);
980   }
981 }
982
983 //_____________________________________________________________________________
984 void AliTRD::ResetDigits()
985 {
986   //
987   // Resets the digits
988   //
989
990   if (gAlice->TreeD()) {
991     TBranch *branch;
992     branch = gAlice->TreeD()->GetBranch("TRDdigits");
993     if (branch) {
994       branch->Reset();
995       for (Int_t iDict = 0; iDict < AliTRDdigitsManager::NDict(); iDict++) {
996         Char_t branchname[15];
997         sprintf(branchname,"TRDdictionary%d",iDict);
998         branch = gAlice->TreeD()->GetBranch(branchname);
999         branch->Reset();
1000       }
1001     }
1002   }
1003
1004 }
1005
1006 //_____________________________________________________________________________
1007 void AliTRD::ResetRecPoints()
1008 {
1009   //
1010   // Reset number of reconstructed points and the point array
1011   //
1012
1013   if (fRecPoints) {
1014     fNRecPoints = 0;
1015     Int_t nentr = fRecPoints->GetEntriesFast();
1016     for (Int_t i = 0; i < nentr; i++) delete fRecPoints->RemoveAt(i);
1017   }
1018
1019 }
1020
1021 //_____________________________________________________________________________
1022 void AliTRD::SetTreeAddress()
1023 {
1024   //
1025   // Set the branch addresses for the trees.
1026   //
1027
1028   Char_t branchname[15];
1029
1030   AliDetector::SetTreeAddress();
1031
1032   TBranch *branch;
1033   TTree   *treeR = gAlice->TreeR();
1034
1035   if (treeR) {
1036     sprintf(branchname,"%scluster",GetName());
1037     if (fRecPoints) {
1038       branch = treeR->GetBranch(branchname);
1039       if (branch) {
1040         branch->SetAddress(&fRecPoints);
1041       }
1042     }
1043   }
1044
1045 }
1046
1047 //_____________________________________________________________________________
1048 void AliTRD::SetGasMix(Int_t imix)
1049 {
1050   //
1051   // Defines the gas mixture (imix=0:  Xe/Isobutane imix=1: Xe/CO2)
1052   //
1053   
1054   if ((imix < 0) || (imix > 1)) {
1055     printf("Wrong input value: %d\n",imix);
1056     printf("Use standard setting\n");
1057     fGasMix = 1;
1058     return;
1059   }
1060
1061   fGasMix = imix;
1062
1063 }
1064
1065 //_____________________________________________________________________________
1066 void AliTRD::SetPHOShole()
1067 {
1068   //
1069   // Selects a geometry with a hole in front of the PHOS
1070   //
1071
1072   fGeometry->SetPHOShole();
1073
1074 }
1075
1076 //_____________________________________________________________________________
1077 void AliTRD::SetRICHhole()
1078 {
1079   //
1080   // Selects a geometry with a hole in front of the RICH
1081   //
1082
1083   fGeometry->SetRICHhole();
1084
1085 }
1086
1087 //_____________________________________________________________________________
1088 AliTRD &AliTRD::operator=(const AliTRD &trd)
1089 {
1090   //
1091   // Assignment operator
1092   //
1093
1094   if (this != &trd) ((AliTRD &) trd).Copy(*this);
1095   return *this;
1096
1097