Fix compiler warning
[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 /* $Id$ */
17
18 ///////////////////////////////////////////////////////////////////////////////
19 //                                                                           //
20 //  Transition Radiation Detector                                            //
21 //  This class contains the basic functions for the Transition Radiation     //
22 //  Detector.                                                                //
23 //                                                                           //
24 ///////////////////////////////////////////////////////////////////////////////
25
26 #include <stdlib.h>
27 #include <Riostream.h>
28
29 #include <TClonesArray.h>
30 #include <TFile.h>
31 #include <TGeometry.h>
32 #include <TLorentzVector.h>
33 #include <TMath.h>
34 #include <TNode.h>
35 #include <TPGON.h> 
36 #include <TParticle.h>
37 #include <TROOT.h>
38 #include <TTree.h>
39 #include <TVirtualMC.h>
40  
41 #include "AliConst.h"
42 #include "AliDigit.h"
43 #include "AliLoader.h"
44 #include "AliLog.h"
45 #include "AliMC.h"
46 #include "AliMagF.h"
47 #include "AliRun.h"
48 #include "AliTrackReference.h"
49 #include "AliRawReader.h"
50
51 #include "AliTRD.h"
52 #include "AliTRDdigit.h"
53 #include "AliTRDdigitizer.h"
54 #include "AliTRDdigitsManager.h"
55 #include "AliTRDgeometry.h"
56 #include "AliTRDhit.h"
57 #include "AliTRDpoints.h"
58 #include "AliTRDrawData.h"
59 #include "AliTRDSimParam.h"
60 #include "AliTRDCommonParam.h"
61 #include "AliTRDcalibDB.h"
62
63 ClassImp(AliTRD)
64  
65 //_____________________________________________________________________________
66 AliTRD::AliTRD()
67   :AliDetector()
68   ,fGeometry(0)
69   ,fGasDensity(0)
70   ,fFoilDensity(0)
71   ,fGasNobleFraction(0)
72   ,fDrawTR(0)
73   ,fDisplayType(0)
74 {
75   //
76   // Default constructor
77   //
78  
79 }
80  
81 //_____________________________________________________________________________
82 AliTRD::AliTRD(const char *name, const char *title)
83   :AliDetector(name,title)
84   ,fGeometry(0)
85   ,fGasDensity(0)
86   ,fFoilDensity(0)
87   ,fGasNobleFraction(0)
88   ,fDrawTR(0)
89   ,fDisplayType(0)
90 {
91   //
92   // Standard constructor for the TRD
93   //
94
95   // Check that FRAME is there otherwise we have no place where to put TRD
96   AliModule *frame = gAlice->GetModule("FRAME");
97   if (!frame) {
98     AliError("TRD needs FRAME to be present\n");
99     exit(1);
100   } 
101
102   // Define the TRD geometry
103   if ((frame->IsVersion() == 0) ||
104       (frame->IsVersion() == 1)) {
105     fGeometry = new AliTRDgeometry();
106   }
107   else {
108     AliError("Could not find valid FRAME version\n");
109     exit(1);
110   }
111
112   // Allocate the hit array
113   fHits = new TClonesArray("AliTRDhit",405);
114   gAlice->GetMCApp()->AddHitList(fHits);
115
116 }
117
118 //_____________________________________________________________________________
119 AliTRD::~AliTRD()
120 {
121   //
122   // TRD destructor
123   //
124
125   if (fGeometry) {
126     delete fGeometry;
127     fGeometry = 0;
128   }
129
130   if (fHits) {
131     delete fHits;
132     fHits     = 0;
133   }
134
135 }
136
137 //_____________________________________________________________________________
138 void AliTRD::Hits2Digits()
139 {
140   //
141   // Create digits
142   //
143
144   AliTRDdigitizer digitizer("TRDdigitizer","TRD digitizer class");
145   AliLog::SetClassDebugLevel("TRDdigitizer",AliDebugLevel());
146
147   // Initialization
148   digitizer.InitDetector();
149     
150   if (!fLoader->TreeH()) {
151     fLoader->LoadHits("read");
152   }
153   fLoader->LoadDigits("recreate");
154
155   AliRunLoader *runLoader = fLoader->GetRunLoader(); 
156
157   for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
158     runLoader->GetEvent(iEvent);
159     digitizer.Open(runLoader,iEvent);
160     digitizer.MakeDigits();
161     digitizer.WriteDigits();
162   }
163
164   fLoader->UnloadHits();
165   fLoader->UnloadDigits();
166
167 }
168
169 //_____________________________________________________________________________
170 void AliTRD::Hits2SDigits()
171 {
172   //
173   // Create summable digits
174   //
175
176   AliTRDdigitizer digitizer("TRDdigitizer","TRD digitizer class");
177   // For the summable digits
178   digitizer.SetSDigits(kTRUE);
179   AliLog::SetClassDebugLevel("TRDdigitizer",AliDebugLevel());
180
181   // Initialization
182   digitizer.InitDetector();
183     
184   if (!fLoader->TreeH()) {
185     fLoader->LoadHits("read");
186   }
187   fLoader->LoadSDigits("recreate");
188
189   AliRunLoader *runLoader = fLoader->GetRunLoader(); 
190
191   for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
192     runLoader->GetEvent(iEvent);
193     digitizer.Open(runLoader,iEvent);
194     digitizer.MakeDigits();
195     digitizer.WriteDigits();
196   }
197
198   fLoader->UnloadHits();
199   fLoader->UnloadSDigits();
200   
201 }
202
203 //_____________________________________________________________________________
204 AliDigitizer *AliTRD::CreateDigitizer(AliRunDigitizer *manager) const
205 {
206   //
207   // Creates a new digitizer object
208   //
209
210   return new AliTRDdigitizer(manager);
211
212 }
213
214 //_____________________________________________________________________________
215 void AliTRD::SDigits2Digits()
216 {
217   //
218   // Create final digits from summable digits
219   //
220
221   // Create the TRD digitizer
222   AliTRDdigitizer digitizer("TRDdigitizer","TRD digitizer class");  
223   AliLog::SetClassDebugLevel("TRDdigitizer",AliDebugLevel());
224
225   // Set the parameter
226   digitizer.SetEvent(gAlice->GetEvNumber());
227
228   // Initialization
229   digitizer.InitDetector();
230
231   // Read the s-digits via digits manager
232   AliTRDdigitsManager sdigitsManager;
233  
234   AliLog::SetClassDebugLevel("TRDdigitisManager",AliDebugLevel());
235   sdigitsManager.SetSDigits(kTRUE);
236   sdigitsManager.CreateArrays();
237   
238   if (!fLoader->TreeS()) { 
239     if (fLoader->LoadSDigits("read")) {
240       return;
241     }
242   }
243   if (!fLoader->TreeS()) {
244     AliError(Form("Error while reading SDigits for event %d",gAlice->GetEvNumber()));
245     return;
246   }
247   
248   sdigitsManager.ReadDigits(fLoader->TreeS());
249
250   // Add the s-digits to the input list 
251   digitizer.AddSDigitsManager(&sdigitsManager);
252
253   // Convert the s-digits to normal digits
254   digitizer.SDigits2Digits();
255
256   // Store the digits
257   if (!fLoader->TreeD()) {
258     fLoader->MakeTree("D");
259   }
260   if (digitizer.MakeBranch(fLoader->TreeD())){
261     digitizer.WriteDigits();
262   }
263
264 }
265
266 //_____________________________________________________________________________
267 void AliTRD::Digits2Raw() 
268 {
269   //
270   // Convert digits of the current event to raw data
271   //
272
273   fLoader->LoadDigits();
274   TTree *digits = fLoader->TreeD();
275   if (!digits) {
276     AliError("No digits tree");
277     return;
278   }
279
280   AliTRDrawData rawWriter;
281   if (!rawWriter.Digits2Raw(digits)) {
282     AliError("The raw writer could not load the digits tree");
283   }
284
285   fLoader->UnloadDigits();
286
287 }
288
289 //_____________________________________________________________________________
290 void AliTRD::AddHit(Int_t track, Int_t det, Float_t *hits, Int_t q
291                   , Float_t time, Bool_t inDrift)
292 {
293   //
294   // Add a hit for the TRD
295   // 
296
297   TClonesArray &lhits = *fHits;
298   AliTRDhit *hit = new(lhits[fNhits++]) AliTRDhit(fIshunt
299                                                  ,track
300                                                  ,det
301                                                  ,hits
302                                                  ,q
303                                                  ,time);
304
305   if (inDrift) {
306     hit->SetDrift();
307   }
308   else {
309     hit->SetAmplification();
310   }
311
312   if (q < 0) {
313     hit->SetTRphoton();
314   }
315
316 }
317
318 //_____________________________________________________________________________
319 void AliTRD::BuildGeometry()
320 {
321   //
322   // Create the ROOT TNode geometry for the TRD
323   //
324
325   TNode *node;
326   TNode *top;
327   TPGON *pgon;
328
329   // The dimensions of the TRD super module
330   const Float_t kRmin  = 291.0;
331   const Float_t kRmax  = 369.0;
332   const Float_t kZmax1 = 378.35;
333   const Float_t kZmax2 = 302.0;
334
335   Float_t rmin;
336   Float_t rmax;
337   Float_t zmax1;
338   Float_t zmax2;
339
340   Int_t   iLayer;
341  
342   const Int_t kColorTRD = 46;
343   
344   // Find the top node alice
345   top = gAlice->GetGeometry()->GetNode("alice");
346   
347   if      (fDisplayType == 0) {
348
349     pgon = new TPGON("S_TRD","TRD","void",0,360,AliTRDgeometry::Nsector(),4);
350     rmin = kRmin;
351     rmax = kRmax;
352     pgon->DefineSection(0,-kZmax1,rmax,rmax);
353     pgon->DefineSection(1,-kZmax2,rmin,rmax);
354     pgon->DefineSection(2, kZmax2,rmin,rmax);
355     pgon->DefineSection(3, kZmax1,rmax,rmax);
356     top->cd();
357     node = new TNode("TRD","TRD","S_TRD",0,0,0,"");
358     node->SetLineColor(kColorTRD);
359     fNodes->Add(node);
360
361   }
362   else if (fDisplayType == 1) {
363
364     Char_t name[7];
365
366     Float_t slope = (kZmax1 - kZmax2) / (kRmax  - kRmin);
367
368     rmin  = kRmin + AliTRDgeometry::CraHght();
369     rmax  = rmin  + AliTRDgeometry::CdrHght();
370
371     Float_t thickness = rmin - kRmin;
372     zmax2 = kZmax2 + slope * thickness;
373     zmax1 = zmax2 + slope * AliTRDgeometry::DrThick();
374
375     for (iLayer = 0; iLayer < AliTRDgeometry::Nlayer(); iLayer++) {
376
377       sprintf(name,"S_TR1%d",iLayer);
378       pgon  = new TPGON(name,"TRD","void",0,360,AliTRDgeometry::Nsector(),4);
379       pgon->DefineSection(0,-zmax1,rmax,rmax);
380       pgon->DefineSection(1,-zmax2,rmin,rmax);
381       pgon->DefineSection(2, zmax2,rmin,rmax);
382       pgon->DefineSection(3, zmax1,rmax,rmax);
383       top->cd();
384       node = new TNode("TRD","TRD",name,0,0,0,"");
385       node->SetLineColor(kColorTRD);
386       fNodes->Add(node);
387
388       Float_t height = AliTRDgeometry::Cheight() + AliTRDgeometry::Cspace(); 
389       rmin  = rmin  + height;
390       rmax  = rmax  + height;
391       zmax1 = zmax1 + slope * height;
392       zmax2 = zmax2 + slope * height;
393
394     }
395
396     thickness += AliTRDgeometry::DrThick();
397     rmin       = kRmin  + thickness;
398     rmax       = rmin   + AliTRDgeometry::AmThick();
399     zmax2      = kZmax2 + slope * thickness;
400     zmax1      = zmax2  + slope * AliTRDgeometry::AmThick();
401
402     for (iLayer = 0; iLayer < AliTRDgeometry::Nlayer(); iLayer++) {
403
404       sprintf(name,"S_TR2%d",iLayer);
405       pgon  = new TPGON(name,"TRD","void",0,360,AliTRDgeometry::Nsector(),4);
406       pgon->DefineSection(0,-zmax1,rmax,rmax);
407       pgon->DefineSection(1,-zmax2,rmin,rmax);
408       pgon->DefineSection(2, zmax2,rmin,rmax);
409       pgon->DefineSection(3, zmax1,rmax,rmax);
410       top->cd();
411       node = new TNode("TRD","TRD",name,0,0,0,"");
412       node->SetLineColor(kColorTRD);
413       fNodes->Add(node);
414
415       Float_t height = AliTRDgeometry::Cheight() + AliTRDgeometry::Cspace(); 
416       rmin  = rmin  + height;
417       rmax  = rmax  + height;
418       zmax1 = zmax1 + slope * height;
419       zmax2 = zmax2 + slope * height;
420
421     }
422
423   }
424
425 }
426  
427 //_____________________________________________________________________________
428 void AliTRD::CreateGeometry()
429 {
430   //
431   // Creates the volumes for the TRD chambers
432   //
433
434   // Check that FRAME is there otherwise we have no place where to put the TRD
435   AliModule *frame = gAlice->GetModule("FRAME");
436   if (!frame) {
437     AliFatal("The TRD needs the FRAME to be defined first");
438   }
439
440   fGeometry->CreateGeometry(fIdtmed->GetArray() - 1299);
441
442 }
443  
444 //_____________________________________________________________________________
445 void AliTRD::CreateMaterials()
446 {
447   //
448   // Create the materials for the TRD
449   //
450
451   Int_t   isxfld = gAlice->Field()->Integ();
452   Float_t sxmgmx = gAlice->Field()->Max();
453   
454   // For polyethilene (CH2) 
455   Float_t ape[2]    = { 12.011 ,  1.0079 };
456   Float_t zpe[2]    = {  6.0   ,  1.0    };
457   Float_t wpe[2]    = {  1.0   ,  2.0    };
458   Float_t dpe       = 0.95;
459
460   // For CO2 
461   Float_t aco[2]    = { 12.011 , 15.9994 };
462   Float_t zco[2]    = {  6.0   ,  8.0    };
463   Float_t wco[2]    = {  1.0   ,  2.0    };
464   Float_t dco       = 0.00186; // at 20C
465
466   // For water
467   Float_t awa[2]    = {  1.0079, 15.9994 };
468   Float_t zwa[2]    = {  1.0   ,  8.0    };
469   Float_t wwa[2]    = {  2.0   ,  1.0    };
470   Float_t dwa       = 1.0;
471
472   // For plexiglas (C5H8O2)
473   Float_t apg[3]    = { 12.011 ,  1.0079, 15.9994 };
474   Float_t zpg[3]    = {  6.0   ,  1.0   ,  8.0    };
475   Float_t wpg[3]    = {  5.0   ,  8.0   ,  2.0    };
476   Float_t dpg       = 1.18; 
477   
478   // For epoxy (C18H19O3)
479   Float_t aEpoxy[3] = { 15.9994,  1.0079, 12.011  }; 
480   Float_t zEpoxy[3] = {  8.0   ,  1.0   ,  6.0    }; 
481   Float_t wEpoxy[3] = {  3.0   , 19.0   , 18.0    }; 
482   Float_t dEpoxy    = 1.8 ; 
483
484   // For Araldite, low density epoxy (C18H19O3)
485   Float_t aAral[3]  = { 15.9994,  1.0079, 12.011  }; 
486   Float_t zAral[3]  = {  8.0   ,  1.0   ,  6.0    }; 
487   Float_t wAral[3]  = {  3.0   , 19.0   , 18.0    }; 
488   Float_t dAral     = 1.05; 
489
490   // For air  
491   Float_t aAir[4]   = { 12.011   , 14.0     , 15.9994  , 36.0      };
492   Float_t zAir[4]   = {  6.0     ,  7.0     ,  8.0     , 18.0      };
493   Float_t wAir[4]   = {  0.000124,  0.755267,  0.231781,  0.012827 };
494   Float_t dAir      = 1.20479e-03;
495
496   // For G10
497   Float_t aG10[4]   = {  1.0079  , 12.011   , 15.9994  , 28.086    };
498   Float_t zG10[4]   = {  1.0     ,  6.0     ,  8.0     , 14.0      };
499   Float_t wG10[4]   = {  0.15201 ,  0.10641 ,  0.49444 ,  0.24714  };
500   Float_t dG10      = 1.7;
501
502   // For Xe/CO2-gas-mixture 
503   Float_t aXeCO2[3] = { 131.29   ,  12.0107 ,  15.9994  };
504   Float_t zXeCO2[3] = {  54.0    ,   6.0    ,   8.0     };
505   Float_t wXeCO2[3] = {   8.5    ,   1.5    ,   3.0     }; 
506   // Xe-content of the Xe/CO2-mixture (85% / 15%) 
507   Float_t fxc       = 0.85;
508   Float_t dxe       = 0.00549; // at 20C
509   Float_t dgmXe     = fxc * dxe + (1.0 - fxc) * dco;
510
511   // For Ar/CO2-gas-mixture
512   Float_t aArCO2[3] = {  39.948  ,  12.0107 ,  15.9994  };
513   Float_t zArCO2[3] = {  18.0    ,   6.0    ,   8.0     };
514   Float_t wArCO2[3] = {   8.2    ,   1.8    ,   3.6     }; 
515   // Ar-content of the Ar/CO2-mixture (82% / 18%) 
516   Float_t fac       = 0.82;
517   Float_t dar       = 0.00166; // at 20C
518   Float_t dgmAr     = fac * dar + (1.0 - fac) * dco;
519   
520   // General tracking parameter
521   Float_t tmaxfd    = -10.0;
522   Float_t stemax    = -1.0e10;
523   Float_t deemax    = -0.1;
524   Float_t epsil     =  1.0e-4;
525   Float_t stmin     = -0.001;
526   
527   //////////////////////////////////////////////////////////////////////////
528   //     Define Materials 
529   //////////////////////////////////////////////////////////////////////////
530
531   AliMaterial( 1, "Al"   ,  26.98, 13.0, 2.7     ,     8.9 ,    37.2);
532   AliMaterial( 4, "Xe"   , 131.29, 54.0, dxe     ,  1546.16,     0.0);
533   AliMaterial( 5, "Cu"   ,  63.54, 29.0, 8.96    ,     1.43,    14.8);
534   AliMaterial( 6, "C"    ,  12.01,  6.0, 2.265   ,    18.8 ,    74.4);
535   AliMaterial(15, "Sn"   , 118.71, 50.0, 7.31    ,     1.21,    14.8);
536   AliMaterial(16, "Si"   ,  28.09, 14.0, 2.33    ,     9.36,    37.2);
537   AliMaterial(18, "Fe"   ,  55.85, 26.0, 7.87    ,     1.76,    14.8);
538
539   // Mixtures 
540   AliMixture(2, "Air"         , aAir,   zAir,   dAir,    4, wAir  );
541   AliMixture(3, "Polyethilene", ape,    zpe,    dpe,    -2, wpe   );
542   AliMixture(8, "CO2",          aco,    zco,    dco,    -2, wco   );
543   if      (AliTRDSimParam::Instance()->IsXenon()) {
544     AliMixture(10,"XeCO2",        aXeCO2, zXeCO2, dgmXe,  -3, wXeCO2);
545   }
546   else if (AliTRDSimParam::Instance()->IsArgon()) {
547     AliInfo("Gas mixture: Ar C02 (80/20)");
548     AliMixture(10,"ArCO2",        aArCO2, zArCO2, dgmAr,  -3, wArCO2);
549   }
550   else {
551     AliFatal("Wrong gas mixture");
552     exit(1);
553   }
554   AliMixture(12,"G10",          aG10,   zG10,   dG10,    4, wG10  );
555   AliMixture(13,"Water",        awa,    zwa,    dwa,    -2, wwa   );
556   AliMixture(14,"Plexiglas",    apg,    zpg,    dpg,    -3, wpg   );
557   AliMixture(17,"Epoxy",        aEpoxy, zEpoxy, dEpoxy, -3, wEpoxy);
558   AliMixture(19,"Araldite",     aAral,  zAral,  dAral,  -3, wAral );
559
560   //////////////////////////////////////////////////////////////////////////
561   //     Tracking Media Parameters 
562   //////////////////////////////////////////////////////////////////////////
563
564   // Al Frame 
565   AliMedium( 1,"Al Frame"   , 1,0,isxfld,sxmgmx
566               ,tmaxfd,stemax,deemax,epsil,stmin);
567   // Air 
568   AliMedium( 2,"Air"        , 2,0,isxfld,sxmgmx
569               ,tmaxfd,stemax,deemax,epsil,stmin);
570   // Wires
571   AliMedium( 3,"Wires"      , 5,0,isxfld,sxmgmx
572               ,tmaxfd,stemax,deemax,epsil,stmin);
573   // All other ROB materials (caps, etc.)
574   AliMedium( 4,"ROB Other"  , 5,0,isxfld,sxmgmx
575               ,tmaxfd,stemax,deemax,epsil,stmin);
576   // Cu pads 
577   AliMedium( 5,"Padplane"   , 5,1,isxfld,sxmgmx
578               ,tmaxfd,stemax,deemax,epsil,stmin);
579   // Fee + cables 
580   AliMedium( 6,"Readout"    , 5,0,isxfld,sxmgmx
581               ,tmaxfd,stemax,deemax,epsil,stmin);
582   // C frame 
583   AliMedium( 7,"C Frame"    , 6,0,isxfld,sxmgmx
584               ,tmaxfd,stemax,deemax,epsil,stmin);
585   // INOX of cooling bus bars
586   AliMedium( 8,"Cooling bus",18,0,isxfld,sxmgmx
587               ,tmaxfd,stemax,deemax,epsil,stmin);
588   // Gas-mixture (Xe/CO2) 
589   AliMedium( 9,"Gas-mix"    ,10,1,isxfld,sxmgmx
590               ,tmaxfd,stemax,deemax,epsil,stmin);
591   // Nomex-honeycomb
592   AliMedium(10,"Nomex"      ,12,0,isxfld,sxmgmx
593               ,tmaxfd,stemax,deemax,epsil,stmin);
594   // Araldite glue
595   AliMedium(11,"Glue"       ,19,0,isxfld,sxmgmx
596               ,tmaxfd,stemax,deemax,epsil,stmin);
597   // G10-plates
598   AliMedium(13,"G10-plates" ,12,0,isxfld,sxmgmx
599               ,tmaxfd,stemax,deemax,epsil,stmin);
600   // Cooling water
601   AliMedium(14,"Water"      ,13,0,isxfld,sxmgmx
602               ,tmaxfd,stemax,deemax,epsil,stmin);
603   // Rohacell (plexiglas) for the radiator
604   AliMedium(15,"Rohacell"   ,14,0,isxfld,sxmgmx
605               ,tmaxfd,stemax,deemax,epsil,stmin);
606   // Al layer in MCMs
607   AliMedium(16,"MCM-Al"     , 1,0,isxfld,sxmgmx
608               ,tmaxfd,stemax,deemax,epsil,stmin);
609   // Sn layer in MCMs
610   AliMedium(17,"MCM-Sn"     ,15,0,isxfld,sxmgmx
611               ,tmaxfd,stemax,deemax,epsil,stmin);
612   // Cu layer in MCMs
613   AliMedium(18,"MCM-Cu"     , 5,0,isxfld,sxmgmx
614               ,tmaxfd,stemax,deemax,epsil,stmin);
615   // G10 layer in MCMs
616   AliMedium(19,"MCM-G10"    ,12,0,isxfld,sxmgmx
617               ,tmaxfd,stemax,deemax,epsil,stmin);
618   // Si in readout chips
619   AliMedium(20,"Chip-Si"    ,16,0,isxfld,sxmgmx
620               ,tmaxfd,stemax,deemax,epsil,stmin);
621   // Epoxy in readout chips
622   AliMedium(21,"Chip-Ep"    ,17,0,isxfld,sxmgmx
623               ,tmaxfd,stemax,deemax,epsil,stmin);
624   // PE in connectors
625   AliMedium(22,"Conn-PE"    , 3,0,isxfld,sxmgmx
626               ,tmaxfd,stemax,deemax,epsil,stmin);
627   // Cu in connectors
628   AliMedium(23,"Chip-Cu"    , 5,0,isxfld,sxmgmx
629               ,tmaxfd,stemax,deemax,epsil,stmin);
630   // Al of cooling pipes
631   AliMedium(24,"Cooling"    , 1,0,isxfld,sxmgmx
632               ,tmaxfd,stemax,deemax,epsil,stmin);
633   // Cu in services
634   AliMedium(25,"Serv-Cu"    , 5,0,isxfld,sxmgmx
635               ,tmaxfd,stemax,deemax,epsil,stmin);
636
637   // Save the density values for the TRD absorbtion
638   Float_t dmy  = 1.39;
639   fFoilDensity = dmy;
640   if      (AliTRDSimParam::Instance()->IsXenon()) {
641     fGasDensity       = dgmXe;
642     fGasNobleFraction = fxc;
643   }
644   else if (AliTRDSimParam::Instance()->IsArgon()) {
645     fGasDensity       = dgmAr;
646     fGasNobleFraction = fac;
647   }
648
649 }
650
651 //_____________________________________________________________________________
652 void AliTRD::DrawModule() const
653 {
654   //
655   // Draw a shaded view of the Transition Radiation Detector version 0
656   //
657
658   // Set everything unseen
659   gMC->Gsatt("*"   ,"SEEN",-1);
660   
661   // Set ALIC mother transparent
662   gMC->Gsatt("ALIC","SEEN", 0);
663   
664   // Set the volumes visible
665   if (fGeometry->IsVersion() == 0) {
666     gMC->Gsatt("B071","SEEN", 0);
667     gMC->Gsatt("B074","SEEN", 0);
668     gMC->Gsatt("B075","SEEN", 0);
669     gMC->Gsatt("B077","SEEN", 0);
670     gMC->Gsatt("BTR1","SEEN", 0);
671     gMC->Gsatt("BTR2","SEEN", 0);
672     gMC->Gsatt("BTR3","SEEN", 0);
673     gMC->Gsatt("UTR1","SEEN", 0);
674     gMC->Gsatt("UTR2","SEEN", 0);
675     gMC->Gsatt("UTR3","SEEN", 0);
676   }
677   else {
678     gMC->Gsatt("B071","SEEN", 0);
679     gMC->Gsatt("B074","SEEN", 0);
680     gMC->Gsatt("B075","SEEN", 0);
681     gMC->Gsatt("B077","SEEN", 0);
682     gMC->Gsatt("BTR1","SEEN", 0);
683     gMC->Gsatt("BTR2","SEEN", 0);
684     gMC->Gsatt("BTR3","SEEN", 0);
685     gMC->Gsatt("UTR1","SEEN", 0);
686   }
687   
688   gMC->Gdopt("hide", "on");
689   gMC->Gdopt("shad", "on");
690   gMC->Gsatt("*", "fill", 7);
691   gMC->SetClipBox(".");
692   gMC->SetClipBox("*", 0, 2000, -2000, 2000, -2000, 2000);
693   gMC->DefaultRange();
694   gMC->Gdraw("alic", 40, 30, 0, 12, 9.4, .021, .021);
695   gMC->Gdhead(1111, "Transition Radiation Detector");
696   gMC->Gdman(18, 4, "MAN");
697
698 }
699  
700 //_____________________________________________________________________________
701 void AliTRD::Init()
702 {
703   //
704   // Initialize the TRD detector after the geometry has been created
705   //
706
707   AliDebug(1,"++++++++++++++++++++++++++++++++++++++++++++++");
708
709   if (fGeometry->IsVersion() != 1) {
710     AliError("Not a valid geometry");
711   }
712
713   // Special tracking options for charged particles for XeCO2
714   gMC->Gstpar((* fIdtmed)[9],"DRAY"    , 1.0);
715   gMC->Gstpar((* fIdtmed)[9],"STRA"    , 1.0); 
716   gMC->Gstpar((* fIdtmed)[9],"LOSS"    ,13.0);      // Specific energy loss
717   gMC->Gstpar((* fIdtmed)[9],"PRIMIO_E",23.53);     // 1st ionisation potential
718   gMC->Gstpar((* fIdtmed)[9],"PRIMIO_N",19.344431); // Number of primaries
719
720 }
721
722 //_____________________________________________________________________________
723 void AliTRD::LoadPoints(Int_t )
724 {
725   //
726   // Store x, y, z of all hits in memory.
727   // Hit originating from TR photons are given a different color
728   //
729
730   if (fHits == 0) {
731     return;
732   }
733
734   Int_t nhits  = fHits->GetEntriesFast();
735   if (nhits == 0) {
736     return;
737   }
738
739   Int_t tracks = gAlice->GetMCApp()->GetNtrack();
740   if (fPoints == 0) {
741     fPoints = new TObjArray(tracks);
742   }
743
744   AliTRDhit *ahit;
745   
746   Int_t    *ntrkE = new Int_t[tracks];
747   Int_t    *ntrkT = new Int_t[tracks];
748   Int_t    *limiE = new Int_t[tracks];
749   Int_t    *limiT = new Int_t[tracks];
750   Float_t **coorE = new Float_t*[tracks];
751   Float_t **coorT = new Float_t*[tracks];
752   for(Int_t i = 0; i < tracks; i++) {
753     ntrkE[i] = 0;
754     ntrkT[i] = 0;
755     coorE[i] = 0;
756     coorT[i] = 0;
757     limiE[i] = 0;
758     limiT[i] = 0;
759   }
760   
761   AliTRDpoints *points = 0;
762   Float_t      *fp     = 0;
763   Int_t         trk;
764   Int_t         chunk  = nhits / 4 + 1;
765
766   // Loop over all the hits and store their position
767   ahit = (AliTRDhit *) FirstHit(-1);
768   while (ahit) {
769
770     // dEdx hits
771     if (ahit->GetCharge() >= 0) {
772
773       trk = ahit->GetTrack();
774       if (ntrkE[trk] == limiE[trk]) {
775         // Initialise a new track
776         fp = new Float_t[3*(limiE[trk]+chunk)];
777         if (coorE[trk]) {
778           memcpy(fp,coorE[trk],sizeof(Float_t)*3*limiE[trk]);
779           delete [] coorE[trk];
780         }
781         limiE[trk] += chunk;
782         coorE[trk]  = fp;
783       } 
784       else {
785         fp = coorE[trk];
786       }
787       fp[3*ntrkE[trk]  ] = ahit->X();
788       fp[3*ntrkE[trk]+1] = ahit->Y();
789       fp[3*ntrkE[trk]+2] = ahit->Z();
790       ntrkE[trk]++;
791
792     }
793     // TR photon hits
794     else if ((ahit->GetCharge() < 0) && 
795              (fDrawTR)) {
796
797       trk = ahit->GetTrack();
798       if (ntrkT[trk] == limiT[trk]) {
799         // Initialise a new track
800         fp = new Float_t[3*(limiT[trk]+chunk)];
801         if (coorT[trk]) {
802           memcpy(fp,coorT[trk],sizeof(Float_t)*3*limiT[trk]);
803           delete [] coorT[trk];
804         }
805         limiT[trk] += chunk;
806         coorT[trk]  = fp;
807       } 
808       else {
809         fp = coorT[trk];
810       }
811       fp[3*ntrkT[trk]  ] = ahit->X();
812       fp[3*ntrkT[trk]+1] = ahit->Y();
813       fp[3*ntrkT[trk]+2] = ahit->Z();
814       ntrkT[trk]++;
815
816     }
817
818     ahit = (AliTRDhit *) NextHit();
819
820   }
821
822   for (trk = 0; trk < tracks; ++trk) {
823
824     if (ntrkE[trk] || ntrkT[trk]) {
825
826       points = new AliTRDpoints();
827       points->SetDetector(this);
828       points->SetParticle(trk);
829
830       // Set the dEdx points
831       if (ntrkE[trk]) {
832         points->SetMarkerColor(kWhite); //PH This is the default color in TRD
833         points->SetMarkerSize(1); //PH Default size=1
834         points->SetPolyMarker(ntrkE[trk],coorE[trk],1); //PH Default style=1
835         delete [] coorE[trk];
836         coorE[trk] = 0;
837       }
838
839       // Set the TR photon points
840       if (ntrkT[trk]) {
841         points->SetTRpoints(ntrkT[trk],coorT[trk]);
842         delete [] coorT[trk];
843         coorT[trk] = 0;
844       }
845
846       fPoints->AddAt(points,trk);
847
848     }
849
850   }
851
852   delete [] coorE;
853   delete [] coorT;
854   delete [] ntrkE;
855   delete [] ntrkT;
856   delete [] limiE;
857   delete [] limiT;
858
859 }
860
861 //_____________________________________________________________________________
862 void AliTRD::MakeBranch(Option_t *option)
863 {
864   //
865   // Create Tree branches for the TRD digits.
866   //
867
868   Int_t  buffersize = 4000;
869   Char_t branchname[15];
870   sprintf(branchname,"%s",GetName());
871
872   const Char_t *cD = strstr(option,"D");
873
874   AliDetector::MakeBranch(option);
875
876   if (fDigits         && 
877       gAlice->TreeD() && 
878       cD) {
879     MakeBranchInTree(gAlice->TreeD(),branchname,&fDigits,buffersize,0);
880   }
881
882 }
883
884 //_____________________________________________________________________________
885 void AliTRD::ResetDigits()
886 {
887   //
888   // Reset number of digits and the digits array for this detector
889   //
890
891   fNdigits = 0;
892
893   if (fDigits) {
894     fDigits->Clear();
895   }
896
897 }
898
899 //_____________________________________________________________________________
900 void AliTRD::SetTreeAddress()
901 {
902   //
903   // Set the branch addresses for the trees.
904   //
905
906   if (fLoader->TreeH() && 
907       (fHits == 0x0)) {
908     fHits = new TClonesArray("AliTRDhit",405);
909   }
910   AliDetector::SetTreeAddress();
911
912 }
913
914 //_____________________________________________________________________________
915 Bool_t AliTRD::Raw2SDigits(AliRawReader *rawReader)
916 {
917   //
918   // Converts RAW data to SDigits
919   //
920
921   AliLoader *loader = fRunLoader->GetLoader("TRDLoader");
922   if (!loader) {
923     AliError("Can not get TRD loader from Run Loader");
924     return kFALSE;
925   }
926     
927   TTree *tree = 0;
928   tree = loader->TreeS();
929   if (!tree) {
930     loader->MakeTree("S");
931     tree = loader->TreeS();
932   }
933
934   AliTRDrawData       *rawdata        = new AliTRDrawData();
935   AliTRDdigitsManager *sdigitsManager = rawdata->Raw2Digits(rawReader);
936   if (sdigitsManager) {
937     sdigitsManager->SetSDigits(kTRUE);
938     sdigitsManager->MakeBranch(tree);
939     sdigitsManager->WriteDigits();
940     return kTRUE;
941   } 
942   else {
943     return kFALSE;
944   }
945
946 }
947
948 //_____________________________________________________________________________
949 AliLoader* AliTRD::MakeLoader(const char* topfoldername)
950 {
951  fLoader = new AliLoader(GetName(),topfoldername);
952
953  AliInfo("Adding Tracklets-loader");
954  AliDataLoader *dl = new AliDataLoader("TRD.Tracklets.root","tracklets", "tracklets");
955  fLoader->AddDataLoader(dl);
956
957  return fLoader;
958 }
959
960 //_____________________________________________________________________________
961 AliTRD &AliTRD::operator=(const AliTRD &trd)
962 {
963   //
964   // Assignment operator
965   //
966
967   if (this != &trd) {
968     ((AliTRD &) trd).Copy(*this);
969   }
970
971   return *this;
972
973 }