Take the absolute time of the hits into account in the digitization. Needed for pileu...
[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
50 #include "AliTRD.h"
51 #include "AliTRDdigit.h"
52 #include "AliTRDdigitizer.h"
53 #include "AliTRDdigitsManager.h"
54 #include "AliTRDgeometry.h"
55 #include "AliTRDhit.h"
56 #include "AliTRDpoints.h"
57 #include "AliTRDrawData.h"
58 #include "AliTRDSimParam.h"
59 #include "AliTRDRecParam.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   ,fDrawTR(0)
72   ,fDisplayType(0)
73 {
74   //
75   // Default constructor
76   //
77  
78 }
79  
80 //_____________________________________________________________________________
81 AliTRD::AliTRD(const char *name, const char *title)
82   :AliDetector(name,title)
83   ,fGeometry(0)
84   ,fGasDensity(0)
85   ,fFoilDensity(0)
86   ,fDrawTR(0)
87   ,fDisplayType(0)
88 {
89   //
90   // Standard constructor for the TRD
91   //
92
93   // Check that FRAME is there otherwise we have no place where to put TRD
94   AliModule *frame = gAlice->GetModule("FRAME");
95   if (!frame) {
96     AliError("TRD needs FRAME to be present\n");
97     exit(1);
98   } 
99
100   // Define the TRD geometry
101   if ((frame->IsVersion() == 0) ||
102       (frame->IsVersion() == 1)) {
103     fGeometry = new AliTRDgeometry();
104   }
105   else {
106     AliError("Could not find valid FRAME version\n");
107     exit(1);
108   }
109
110   // Save the geometry
111   TDirectory *saveDir = gDirectory;
112   gAlice->GetRunLoader()->CdGAFile();
113   fGeometry->Write("TRDgeometry");
114   saveDir->cd();
115
116   // Allocate the hit array
117   fHits = new TClonesArray("AliTRDhit",405);
118   gAlice->GetMCApp()->AddHitList(fHits);
119
120   //PH SetMarkerColor(kWhite);   
121
122 }
123
124 //_____________________________________________________________________________
125 AliTRD::~AliTRD()
126 {
127   //
128   // TRD destructor
129   //
130
131   if (fGeometry) {
132     delete fGeometry;
133     fGeometry = 0;
134   }
135
136   if (fHits) {
137     delete fHits;
138     fHits     = 0;
139   }
140
141 }
142
143 //_____________________________________________________________________________
144 void AliTRD::Hits2Digits()
145 {
146   //
147   // Create digits
148   //
149
150   AliTRDdigitizer digitizer("TRDdigitizer","TRD digitizer class");
151   AliLog::SetClassDebugLevel("TRDdigitizer",AliDebugLevel());
152   
153   // Initialization
154   digitizer.InitDetector();
155     
156   if (!fLoader->TreeH()) {
157     fLoader->LoadHits("read");
158   }
159   fLoader->LoadDigits("recreate");
160
161   AliRunLoader *runLoader = fLoader->GetRunLoader(); 
162
163   for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
164     runLoader->GetEvent(iEvent);
165     digitizer.Open(runLoader,iEvent);
166     digitizer.MakeDigits();
167     digitizer.WriteDigits();
168   }
169
170   fLoader->UnloadHits();
171   fLoader->UnloadDigits();
172
173 }
174
175 //_____________________________________________________________________________
176 void AliTRD::Hits2SDigits()
177 {
178   //
179   // Create summable digits
180   //
181
182   AliTRDdigitizer digitizer("TRDdigitizer","TRD digitizer class");
183   // For the summable digits
184   digitizer.SetSDigits(kTRUE);
185   AliLog::SetClassDebugLevel("TRDdigitizer",AliDebugLevel());
186
187   // Initialization
188   digitizer.InitDetector();
189     
190   if (!fLoader->TreeH()) {
191     fLoader->LoadHits("read");
192   }
193   fLoader->LoadSDigits("recreate");
194
195   AliRunLoader *runLoader = fLoader->GetRunLoader(); 
196
197   for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
198     runLoader->GetEvent(iEvent);
199     digitizer.Open(runLoader,iEvent);
200     digitizer.MakeDigits();
201     digitizer.WriteDigits();
202   }
203
204   fLoader->UnloadHits();
205   fLoader->UnloadSDigits();
206   
207 }
208
209 //_____________________________________________________________________________
210 AliDigitizer *AliTRD::CreateDigitizer(AliRunDigitizer *manager) const
211 {
212   //
213   // Creates a new digitizer object
214   //
215
216   return new AliTRDdigitizer(manager);
217
218 }
219
220 //_____________________________________________________________________________
221 void AliTRD::SDigits2Digits()
222 {
223   //
224   // Create final digits from summable digits
225   //
226
227   // Create the TRD digitizer
228   AliTRDdigitizer digitizer("TRDdigitizer","TRD digitizer class");  
229   AliLog::SetClassDebugLevel("TRDdigitizer",AliDebugLevel());
230
231   // Set the parameter
232   digitizer.SetEvent(gAlice->GetEvNumber());
233
234   // Initialization
235   digitizer.InitDetector();
236
237   // Read the s-digits via digits manager
238   AliTRDdigitsManager sdigitsManager;
239  
240   AliLog::SetClassDebugLevel("TRDdigitisManager",AliDebugLevel());
241   sdigitsManager.SetSDigits(kTRUE);
242   sdigitsManager.CreateArrays();
243   
244   if (!fLoader->TreeS()) { 
245     if (fLoader->LoadSDigits("read")) {
246       return;
247     }
248   }
249   if (!fLoader->TreeS()) {
250     AliError(Form("Error while reading SDigits for event %d",gAlice->GetEvNumber()));
251     return;
252   }
253   
254   sdigitsManager.ReadDigits(fLoader->TreeS());
255
256   // Add the s-digits to the input list 
257   digitizer.AddSDigitsManager(&sdigitsManager);
258
259   // Convert the s-digits to normal digits
260   digitizer.SDigits2Digits();
261
262   // Store the digits
263   if (!fLoader->TreeD()) {
264     fLoader->MakeTree("D");
265   }
266   if (digitizer.MakeBranch(fLoader->TreeD())){
267     digitizer.WriteDigits();
268   }
269
270 }
271
272 //_____________________________________________________________________________
273 void AliTRD::Digits2Raw() 
274 {
275   //
276   // Convert digits of the current event to raw data
277   //
278
279   fLoader->LoadDigits();
280   TTree *digits = fLoader->TreeD();
281   if (!digits) {
282     AliError("No digits tree");
283     return;
284   }
285
286   AliTRDrawData rawWriter;
287   if (!rawWriter.Digits2Raw(digits)) {
288     AliError("The raw writer could not load the digits tree");
289   }
290
291   fLoader->UnloadDigits();
292
293 }
294
295 //_____________________________________________________________________________
296 void AliTRD::AddHit(Int_t track, Int_t det, Float_t *hits, Int_t q
297                   , Float_t time, Bool_t inDrift)
298 {
299   //
300   // Add a hit for the TRD
301   // 
302
303   TClonesArray &lhits = *fHits;
304   AliTRDhit *hit = new(lhits[fNhits++]) AliTRDhit(fIshunt
305                                                  ,track
306                                                  ,det
307                                                  ,hits
308                                                  ,q
309                                                  ,time);
310
311   printf("AliTRD::AddHit --- time = %f\n",time);
312
313   if (inDrift) {
314     hit->SetDrift();
315   }
316   else {
317     hit->SetAmplification();
318   }
319
320   if (q < 0) {
321     hit->SetTRphoton();
322   }
323
324 }
325
326 //_____________________________________________________________________________
327 void AliTRD::BuildGeometry()
328 {
329   //
330   // Create the ROOT TNode geometry for the TRD
331   //
332
333   TNode *node;
334   TNode *top;
335   TPGON *pgon;
336
337   // The dimensions of the TRD super module
338   const Float_t kRmin  = 291.0;
339   const Float_t kRmax  = 369.0;
340   const Float_t kZmax1 = 378.35;
341   const Float_t kZmax2 = 302.0;
342
343   Float_t rmin;
344   Float_t rmax;
345   Float_t zmax1;
346   Float_t zmax2;
347
348   Int_t   iPlan;
349  
350   const Int_t kColorTRD = 46;
351   
352   // Find the top node alice
353   top = gAlice->GetGeometry()->GetNode("alice");
354   
355   if      (fDisplayType == 0) {
356
357     pgon = new TPGON("S_TRD","TRD","void",0,360,AliTRDgeometry::Nsect(),4);
358     rmin = kRmin;
359     rmax = kRmax;
360     pgon->DefineSection(0,-kZmax1,rmax,rmax);
361     pgon->DefineSection(1,-kZmax2,rmin,rmax);
362     pgon->DefineSection(2, kZmax2,rmin,rmax);
363     pgon->DefineSection(3, kZmax1,rmax,rmax);
364     top->cd();
365     node = new TNode("TRD","TRD","S_TRD",0,0,0,"");
366     node->SetLineColor(kColorTRD);
367     fNodes->Add(node);
368
369   }
370   else if (fDisplayType == 1) {
371
372     Char_t name[7];
373
374     Float_t slope = (kZmax1 - kZmax2) / (kRmax  - kRmin);
375
376     rmin  = kRmin + AliTRDgeometry::CraHght();
377     rmax  = rmin  + AliTRDgeometry::CdrHght();
378
379     Float_t thickness = rmin - kRmin;
380     zmax2 = kZmax2 + slope * thickness;
381     zmax1 = zmax2 + slope * AliTRDgeometry::DrThick();
382
383     for (iPlan = 0; iPlan < AliTRDgeometry::Nplan(); iPlan++) {
384
385       sprintf(name,"S_TR1%d",iPlan);
386       pgon  = new TPGON(name,"TRD","void",0,360,AliTRDgeometry::Nsect(),4);
387       pgon->DefineSection(0,-zmax1,rmax,rmax);
388       pgon->DefineSection(1,-zmax2,rmin,rmax);
389       pgon->DefineSection(2, zmax2,rmin,rmax);
390       pgon->DefineSection(3, zmax1,rmax,rmax);
391       top->cd();
392       node = new TNode("TRD","TRD",name,0,0,0,"");
393       node->SetLineColor(kColorTRD);
394       fNodes->Add(node);
395
396       Float_t height = AliTRDgeometry::Cheight() + AliTRDgeometry::Cspace(); 
397       rmin  = rmin  + height;
398       rmax  = rmax  + height;
399       zmax1 = zmax1 + slope * height;
400       zmax2 = zmax2 + slope * height;
401
402     }
403
404     thickness += AliTRDgeometry::DrThick();
405     rmin       = kRmin  + thickness;
406     rmax       = rmin   + AliTRDgeometry::AmThick();
407     zmax2      = kZmax2 + slope * thickness;
408     zmax1      = zmax2  + slope * AliTRDgeometry::AmThick();
409
410     for (iPlan = 0; iPlan < AliTRDgeometry::Nplan(); iPlan++) {
411
412       sprintf(name,"S_TR2%d",iPlan);
413       pgon  = new TPGON(name,"TRD","void",0,360,AliTRDgeometry::Nsect(),4);
414       pgon->DefineSection(0,-zmax1,rmax,rmax);
415       pgon->DefineSection(1,-zmax2,rmin,rmax);
416       pgon->DefineSection(2, zmax2,rmin,rmax);
417       pgon->DefineSection(3, zmax1,rmax,rmax);
418       top->cd();
419       node = new TNode("TRD","TRD",name,0,0,0,"");
420       node->SetLineColor(kColorTRD);
421       fNodes->Add(node);
422
423       Float_t height = AliTRDgeometry::Cheight() + AliTRDgeometry::Cspace(); 
424       rmin  = rmin  + height;
425       rmax  = rmax  + height;
426       zmax1 = zmax1 + slope * height;
427       zmax2 = zmax2 + slope * height;
428
429     }
430
431   }
432
433 }
434  
435 //_____________________________________________________________________________
436 void AliTRD::CreateGeometry()
437 {
438   //
439   // Creates the volumes for the TRD chambers
440   //
441
442   // Check that FRAME is there otherwise we have no place where to put the TRD
443   AliModule *frame = gAlice->GetModule("FRAME");
444   if (!frame) {
445     AliFatal("The TRD needs the FRAME to be defined first");
446   }
447
448   fGeometry->CreateGeometry(fIdtmed->GetArray() - 1299);
449
450 }
451  
452 //_____________________________________________________________________________
453 void AliTRD::CreateMaterials()
454 {
455   //
456   // Create the materials for the TRD
457   //
458
459   Int_t   isxfld = gAlice->Field()->Integ();
460   Float_t sxmgmx = gAlice->Field()->Max();
461   
462   // For polyethilene (CH2) 
463   Float_t ape[2]    = { 12.011 ,  1.0079 };
464   Float_t zpe[2]    = {  6.0   ,  1.0    };
465   Float_t wpe[2]    = {  1.0   ,  2.0    };
466   Float_t dpe       = 0.95;
467
468   // For CO2 
469   Float_t aco[2]    = { 12.011 , 15.9994 };
470   Float_t zco[2]    = {  6.0   ,  8.0    };
471   Float_t wco[2]    = {  1.0   ,  2.0    };
472   Float_t dco       = 0.00186;
473
474   // For water
475   Float_t awa[2]    = {  1.0079, 15.9994 };
476   Float_t zwa[2]    = {  1.0   ,  8.0    };
477   Float_t wwa[2]    = {  2.0   ,  1.0    };
478   Float_t dwa       = 1.0;
479
480   // For isobutane (C4H10)
481   Float_t ais[2]    = { 12.011 ,  1.0079 };
482   Float_t zis[2]    = {  6.0   ,  1.0    };
483   Float_t wis[2]    = {  4.0   , 10.0    };
484   Float_t dis       = 0.00267;
485
486   // For plexiglas (C5H8O2)
487   Float_t apg[3]    = { 12.011 ,  1.0079, 15.9994 };
488   Float_t zpg[3]    = {  6.0   ,  1.0   ,  8.0    };
489   Float_t wpg[3]    = {  5.0   ,  8.0   ,  2.0    };
490   Float_t dpg       = 1.18; 
491   
492   // For epoxy (C18H19O3)
493   Float_t aEpoxy[3] = { 15.9994,  1.0079, 12.011  }; 
494   Float_t zEpoxy[3] = {  8.0   ,  1.0   ,  6.0    }; 
495   Float_t wEpoxy[3] = {  3.0   , 19.0   , 18.0    }; 
496   Float_t dEpoxy    = 1.8 ; 
497
498   // For Araldite, low density epoxy (C18H19O3)
499   Float_t aAral[3]  = { 15.9994,  1.0079, 12.011  }; 
500   Float_t zAral[3]  = {  8.0   ,  1.0   ,  6.0    }; 
501   Float_t wAral[3]  = {  3.0   , 19.0   , 18.0    }; 
502   Float_t dAral     = 1.05; 
503
504   // For air  
505   Float_t aAir[4]   = { 12.011   , 14.0     , 15.9994  , 36.0      };
506   Float_t zAir[4]   = {  6.0     ,  7.0     ,  8.0     , 18.0      };
507   Float_t wAir[4]   = {  0.000124,  0.755267,  0.231781,  0.012827 };
508   Float_t dAir      = 1.20479e-03;
509
510   // For G10
511   Float_t aG10[4]   = {  1.0079  , 12.011   , 15.9994  , 28.086    };
512   Float_t zG10[4]   = {  1.0     ,  6.0     ,  8.0     , 14.0      };
513   Float_t wG10[4]   = {  0.15201 ,  0.10641 ,  0.49444 ,  0.24714  };
514   Float_t dG10      = 1.7;
515
516   // For Xe/CO2-gas-mixture 
517   Float_t aXeCO2[3] = { 131.29   ,  12.0107 ,  15.9994  };
518   Float_t zXeCO2[3] = {  54.0    ,   6.0    ,   8.0     };
519   Float_t wXeCO2[3] = {   8.5    ,   1.5    ,   3.0     }; 
520   // Xe-content of the Xe/CO2-mixture (85% / 15%) 
521   Float_t fxc       = 0.85;
522   Float_t dxe       = 0.00549;
523   Float_t dgm       = fxc * dxe + (1.0 - fxc) * dco;
524   
525   // General tracking parameter
526   Float_t tmaxfd    = -10.0;
527   Float_t stemax    = -1.0e10;
528   Float_t deemax    = -0.1;
529   Float_t epsil     =  1.0e-4;
530   Float_t stmin     = -0.001;
531   
532   //////////////////////////////////////////////////////////////////////////
533   //     Define Materials 
534   //////////////////////////////////////////////////////////////////////////
535
536   AliMaterial( 1, "Al"   ,  26.98, 13.0, 2.7     ,     8.9 ,    37.2);
537   AliMaterial( 4, "Xe"   , 131.29, 54.0, dxe     ,  1546.16,     0.0);
538   AliMaterial( 5, "Cu"   ,  63.54, 29.0, 8.96    ,     1.43,    14.8);
539   AliMaterial( 6, "C"    ,  12.01,  6.0, 2.265   ,    18.8 ,    74.4);
540   AliMaterial(15, "Sn"   , 118.71, 50.0, 7.31    ,     1.21,    14.8);
541   AliMaterial(16, "Si"   ,  28.09, 14.0, 2.33    ,     9.36,    37.2);
542   AliMaterial(18, "Fe"   ,  55.85, 26.0, 7.87    ,     1.76,    14.8);
543
544   // Mixtures 
545   AliMixture(2, "Air"         , aAir,   zAir,   dAir,    4, wAir  );
546   AliMixture(3, "Polyethilene", ape,    zpe,    dpe,    -2, wpe   );
547   AliMixture(8, "CO2",          aco,    zco,    dco,    -2, wco   );
548   AliMixture(9, "Isobutane",    ais,    zis,    dis,    -2, wis   );
549   AliMixture(10,"Gas mixture",  aXeCO2, zXeCO2, dgm,    -3, wXeCO2);
550   AliMixture(12,"G10",          aG10,   zG10,   dG10,    4, wG10  );
551   AliMixture(13,"Water",        awa,    zwa,    dwa,    -2, wwa   );
552   AliMixture(14,"Plexiglas",    apg,    zpg,    dpg,    -3, wpg   );
553   AliMixture(17,"Epoxy",        aEpoxy, zEpoxy, dEpoxy, -3, wEpoxy);
554   AliMixture(19,"Araldite",     aAral,  zAral,  dAral,  -3, wAral );
555
556   //////////////////////////////////////////////////////////////////////////
557   //     Tracking Media Parameters 
558   //////////////////////////////////////////////////////////////////////////
559
560   // Al Frame 
561   AliMedium( 1,"Al Frame"   , 1,0,isxfld,sxmgmx
562               ,tmaxfd,stemax,deemax,epsil,stmin);
563   // Air 
564   AliMedium( 2,"Air"        , 2,0,isxfld,sxmgmx
565               ,tmaxfd,stemax,deemax,epsil,stmin);
566   // Wires
567   AliMedium( 3,"Wires"      , 5,0,isxfld,sxmgmx
568               ,tmaxfd,stemax,deemax,epsil,stmin);
569   // All other ROB materials (caps, etc.)
570   AliMedium( 4,"ROB Other"  , 5,0,isxfld,sxmgmx
571               ,tmaxfd,stemax,deemax,epsil,stmin);
572   // Cu pads 
573   AliMedium( 5,"Padplane"   , 5,1,isxfld,sxmgmx
574               ,tmaxfd,stemax,deemax,epsil,stmin);
575   // Fee + cables 
576   AliMedium( 6,"Readout"    , 5,0,isxfld,sxmgmx
577               ,tmaxfd,stemax,deemax,epsil,stmin);
578   // C frame 
579   AliMedium( 7,"C Frame"    , 6,0,isxfld,sxmgmx
580               ,tmaxfd,stemax,deemax,epsil,stmin);
581   // INOX of cooling bus bars
582   AliMedium( 8,"Cooling bus",18,0,isxfld,sxmgmx
583               ,tmaxfd,stemax,deemax,epsil,stmin);
584   // Gas-mixture (Xe/CO2) 
585   AliMedium( 9,"Gas-mix"    ,10,1,isxfld,sxmgmx
586               ,tmaxfd,stemax,deemax,epsil,stmin);
587   // Nomex-honeycomb
588   AliMedium(10,"Nomex"      ,12,0,isxfld,sxmgmx
589               ,tmaxfd,stemax,deemax,epsil,stmin);
590   // Araldite glue
591   AliMedium(11,"Glue"       ,19,0,isxfld,sxmgmx
592               ,tmaxfd,stemax,deemax,epsil,stmin);
593   // G10-plates
594   AliMedium(13,"G10-plates" ,12,0,isxfld,sxmgmx
595               ,tmaxfd,stemax,deemax,epsil,stmin);
596   // Cooling water
597   AliMedium(14,"Water"      ,13,0,isxfld,sxmgmx
598               ,tmaxfd,stemax,deemax,epsil,stmin);
599   // Rohacell (plexiglas) for the radiator
600   AliMedium(15,"Rohacell"   ,14,0,isxfld,sxmgmx
601               ,tmaxfd,stemax,deemax,epsil,stmin);
602   // Al layer in MCMs
603   AliMedium(16,"MCM-Al"     , 1,0,isxfld,sxmgmx
604               ,tmaxfd,stemax,deemax,epsil,stmin);
605   // Sn layer in MCMs
606   AliMedium(17,"MCM-Sn"     ,15,0,isxfld,sxmgmx
607               ,tmaxfd,stemax,deemax,epsil,stmin);
608   // Cu layer in MCMs
609   AliMedium(18,"MCM-Cu"     , 5,0,isxfld,sxmgmx
610               ,tmaxfd,stemax,deemax,epsil,stmin);
611   // G10 layer in MCMs
612   AliMedium(19,"MCM-G10"    ,12,0,isxfld,sxmgmx
613               ,tmaxfd,stemax,deemax,epsil,stmin);
614   // Si in readout chips
615   AliMedium(20,"Chip-Si"    ,16,0,isxfld,sxmgmx
616               ,tmaxfd,stemax,deemax,epsil,stmin);
617   // Epoxy in readout chips
618   AliMedium(21,"Chip-Ep"    ,17,0,isxfld,sxmgmx
619               ,tmaxfd,stemax,deemax,epsil,stmin);
620   // PE in connectors
621   AliMedium(22,"Conn-PE"    , 3,0,isxfld,sxmgmx
622               ,tmaxfd,stemax,deemax,epsil,stmin);
623   // Cu in connectors
624   AliMedium(23,"Chip-Cu"    , 5,0,isxfld,sxmgmx
625               ,tmaxfd,stemax,deemax,epsil,stmin);
626   // Al of cooling pipes
627   AliMedium(24,"Cooling"    , 1,0,isxfld,sxmgmx
628               ,tmaxfd,stemax,deemax,epsil,stmin);
629   // Cu in services
630   AliMedium(25,"Serv-Cu"    , 5,0,isxfld,sxmgmx
631               ,tmaxfd,stemax,deemax,epsil,stmin);
632
633   // Save the density values for the TRD absorbtion
634   Float_t dmy  = 1.39;
635   fFoilDensity = dmy;
636   fGasDensity  = dgm;
637
638 }
639
640 //_____________________________________________________________________________
641 void AliTRD::DrawModule() const
642 {
643   //
644   // Draw a shaded view of the Transition Radiation Detector version 0
645   //
646
647   // Set everything unseen
648   gMC->Gsatt("*"   ,"SEEN",-1);
649   
650   // Set ALIC mother transparent
651   gMC->Gsatt("ALIC","SEEN", 0);
652   
653   // Set the volumes visible
654   if (fGeometry->IsVersion() == 0) {
655     gMC->Gsatt("B071","SEEN", 0);
656     gMC->Gsatt("B074","SEEN", 0);
657     gMC->Gsatt("B075","SEEN", 0);
658     gMC->Gsatt("B077","SEEN", 0);
659     gMC->Gsatt("BTR1","SEEN", 0);
660     gMC->Gsatt("BTR2","SEEN", 0);
661     gMC->Gsatt("BTR3","SEEN", 0);
662     gMC->Gsatt("UTR1","SEEN", 0);
663     gMC->Gsatt("UTR2","SEEN", 0);
664     gMC->Gsatt("UTR3","SEEN", 0);
665   }
666   else {
667     gMC->Gsatt("B071","SEEN", 0);
668     gMC->Gsatt("B074","SEEN", 0);
669     gMC->Gsatt("B075","SEEN", 0);
670     gMC->Gsatt("B077","SEEN", 0);
671     gMC->Gsatt("BTR1","SEEN", 0);
672     gMC->Gsatt("BTR2","SEEN", 0);
673     gMC->Gsatt("BTR3","SEEN", 0);
674     gMC->Gsatt("UTR1","SEEN", 0);
675   }
676   
677   gMC->Gdopt("hide", "on");
678   gMC->Gdopt("shad", "on");
679   gMC->Gsatt("*", "fill", 7);
680   gMC->SetClipBox(".");
681   gMC->SetClipBox("*", 0, 2000, -2000, 2000, -2000, 2000);
682   gMC->DefaultRange();
683   gMC->Gdraw("alic", 40, 30, 0, 12, 9.4, .021, .021);
684   gMC->Gdhead(1111, "Transition Radiation Detector");
685   gMC->Gdman(18, 4, "MAN");
686
687 }
688
689 //_____________________________________________________________________________
690 Int_t AliTRD::DistancetoPrimitive(Int_t , Int_t )
691 {
692   //
693   // Distance between the mouse and the TRD detector on the screen
694   // Dummy routine
695   //
696   
697   return 9999;
698
699 }
700  
701 //_____________________________________________________________________________
702 void AliTRD::Init()
703 {
704   //
705   // Initialize the TRD detector after the geometry has been created
706   //
707
708   AliDebug(1,"++++++++++++++++++++++++++++++++++++++++++++++");
709
710   if (fGeometry->IsVersion() != 1) {
711     AliError("Not a valid geometry");
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   
717 }
718
719 //_____________________________________________________________________________
720 void AliTRD::LoadPoints(Int_t )
721 {
722   //
723   // Store x, y, z of all hits in memory.
724   // Hit originating from TR photons are given a different color
725   //
726
727   if (fHits == 0) {
728     return;
729   }
730
731   Int_t nhits  = fHits->GetEntriesFast();
732   if (nhits == 0) {
733     return;
734   }
735
736   Int_t tracks = gAlice->GetMCApp()->GetNtrack();
737   if (fPoints == 0) {
738     fPoints = new TObjArray(tracks);
739   }
740
741   AliTRDhit *ahit;
742   
743   Int_t    *ntrkE = new Int_t[tracks];
744   Int_t    *ntrkT = new Int_t[tracks];
745   Int_t    *limiE = new Int_t[tracks];
746   Int_t    *limiT = new Int_t[tracks];
747   Float_t **coorE = new Float_t*[tracks];
748   Float_t **coorT = new Float_t*[tracks];
749   for(Int_t i = 0; i < tracks; i++) {
750     ntrkE[i] = 0;
751     ntrkT[i] = 0;
752     coorE[i] = 0;
753     coorT[i] = 0;
754     limiE[i] = 0;
755     limiT[i] = 0;
756   }
757   
758   AliTRDpoints *points = 0;
759   Float_t      *fp     = 0;
760   Int_t         trk;
761   Int_t         chunk  = nhits / 4 + 1;
762
763   // Loop over all the hits and store their position
764   ahit = (AliTRDhit *) FirstHit(-1);
765   while (ahit) {
766
767     // dEdx hits
768     if (ahit->GetCharge() >= 0) {
769
770       trk = ahit->GetTrack();
771       if (ntrkE[trk] == limiE[trk]) {
772         // Initialise a new track
773         fp = new Float_t[3*(limiE[trk]+chunk)];
774         if (coorE[trk]) {
775           memcpy(fp,coorE[trk],sizeof(Float_t)*3*limiE[trk]);
776           delete [] coorE[trk];
777         }
778         limiE[trk] += chunk;
779         coorE[trk]  = fp;
780       } 
781       else {
782         fp = coorE[trk];
783       }
784       fp[3*ntrkE[trk]  ] = ahit->X();
785       fp[3*ntrkE[trk]+1] = ahit->Y();
786       fp[3*ntrkE[trk]+2] = ahit->Z();
787       ntrkE[trk]++;
788
789     }
790     // TR photon hits
791     else if ((ahit->GetCharge() < 0) && 
792              (fDrawTR)) {
793
794       trk = ahit->GetTrack();
795       if (ntrkT[trk] == limiT[trk]) {
796         // Initialise a new track
797         fp = new Float_t[3*(limiT[trk]+chunk)];
798         if (coorT[trk]) {
799           memcpy(fp,coorT[trk],sizeof(Float_t)*3*limiT[trk]);
800           delete [] coorT[trk];
801         }
802         limiT[trk] += chunk;
803         coorT[trk]  = fp;
804       } 
805       else {
806         fp = coorT[trk];
807       }
808       fp[3*ntrkT[trk]  ] = ahit->X();
809       fp[3*ntrkT[trk]+1] = ahit->Y();
810       fp[3*ntrkT[trk]+2] = ahit->Z();
811       ntrkT[trk]++;
812
813     }
814
815     ahit = (AliTRDhit *) NextHit();
816
817   }
818
819   for (trk = 0; trk < tracks; ++trk) {
820
821     if (ntrkE[trk] || ntrkT[trk]) {
822
823       points = new AliTRDpoints();
824       points->SetDetector(this);
825       points->SetParticle(trk);
826
827       // Set the dEdx points
828       if (ntrkE[trk]) {
829         points->SetMarkerColor(kWhite); //PH This is the default color in TRD
830         points->SetMarkerSize(1); //PH Default size=1
831         points->SetPolyMarker(ntrkE[trk],coorE[trk],1); //PH Default style=1
832         delete [] coorE[trk];
833         coorE[trk] = 0;
834       }
835
836       // Set the TR photon points
837       if (ntrkT[trk]) {
838         points->SetTRpoints(ntrkT[trk],coorT[trk]);
839         delete [] coorT[trk];
840         coorT[trk] = 0;
841       }
842
843       fPoints->AddAt(points,trk);
844
845     }
846
847   }
848
849   delete [] coorE;
850   delete [] coorT;
851   delete [] ntrkE;
852   delete [] ntrkT;
853   delete [] limiE;
854   delete [] limiT;
855
856 }
857
858 //_____________________________________________________________________________
859 void AliTRD::MakeBranch(Option_t *option)
860 {
861   //
862   // Create Tree branches for the TRD digits.
863   //
864
865   Int_t  buffersize = 4000;
866   Char_t branchname[15];
867   sprintf(branchname,"%s",GetName());
868
869   const Char_t *cD = strstr(option,"D");
870
871   AliDetector::MakeBranch(option);
872
873   if (fDigits         && 
874       gAlice->TreeD() && 
875       cD) {
876     MakeBranchInTree(gAlice->TreeD(),branchname,&fDigits,buffersize,0);
877   }
878
879 }
880
881 //_____________________________________________________________________________
882 void AliTRD::ResetDigits()
883 {
884   //
885   // Reset number of digits and the digits array for this detector
886   //
887
888   fNdigits = 0;
889
890   if (fDigits) {
891     fDigits->Clear();
892   }
893
894 }
895
896 //_____________________________________________________________________________
897 void AliTRD::SetTreeAddress()
898 {
899   //
900   // Set the branch addresses for the trees.
901   //
902
903   if (fLoader->TreeH() && 
904       (fHits == 0x0)) {
905     fHits = new TClonesArray("AliTRDhit",405);
906   }
907   AliDetector::SetTreeAddress();
908
909 }
910
911 //_____________________________________________________________________________
912 AliTRD &AliTRD::operator=(const AliTRD &trd)
913 {
914   //
915   // Assignment operator
916   //
917
918   if (this != &trd) {
919     ((AliTRD &) trd).Copy(*this);
920   }
921
922   return *this;
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
1001
1002
1003
1004
1005
1006
1007
1008
1009
1010
1011
1012
1013
1014
1015
1016
1017
1018
1019
1020
1021
1022
1023
1024
1025
1026
1027
1028
1029
1030
1031
1032
1033
1034
1035
1036
1037
1038
1039
1040
1041
1042
1043
1044
1045
1046
1047
1048
1049
1050
1051
1052
1053
1054
1055
1056
1057
1058
1059
1060
1061
1062
1063
1064
1065
1066
1067
1068
1069
1070
1071
1072
1073
1074
1075
1076
1077
1078
1079
1080
1081
1082
1083
1084
1085
1086
1087
1088
1089
1090
1091
1092
1093
1094
1095
1096
1097
1098
1099
1100
1101
1102
1103
1104
1105
1106
1107
1108
1109
1110
1111
1112
1113
1114
1115
1116
1117
1118
1119
1120
1121
1122
1123
1124
1125
1126
1127
1128
1129
1130
1131
1132
1133
1134
1135
1136
1137
1138
1139
1140
1141
1142
1143
1144
1145
1146
1147
1148
1149
1150
1151
1152
1153
1154
1155
1156
1157
1158
1159
1160
1161
1162
1163
1164
1165
1166
1167
1168
1169
1170
1171
1172
1173
1174
1175
1176
1177
1178
1179
1180
1181
1182
1183
1184
1185
1186
1187
1188
1189
1190
1191
1192
1193
1194
1195
1196
1197
1198
1199
1200
1201
1202
1203
1204
1205
1206
1207
1208
1209
1210
1211
1212
1213
1214
1215
1216
1217
1218
1219
1220
1221
1222
1223
1224
1225
1226
1227
1228
1229
1230
1231
1232
1233
1234
1235
1236
1237
1238
1239
1240
1241
1242
1243
1244
1245
1246
1247
1248
1249
1250
1251
1252
1253
1254
1255
1256
1257
1258
1259