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