Replace AliTRDCalibra
[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                   , 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,track,det,hits,q);
305
306   if (inDrift) {
307     hit->SetDrift();
308   }
309   else {
310     hit->SetAmplification();
311   }
312
313   if (q < 0) {
314     hit->SetTRphoton();
315   }
316
317 }
318
319 //_____________________________________________________________________________
320 void AliTRD::BuildGeometry()
321 {
322   //
323   // Create the ROOT TNode geometry for the TRD
324   //
325
326   TNode *node;
327   TNode *top;
328   TPGON *pgon;
329
330   // The dimensions of the TRD super module
331   const Float_t kRmin  = 291.0;
332   const Float_t kRmax  = 369.0;
333   const Float_t kZmax1 = 378.35;
334   const Float_t kZmax2 = 302.0;
335
336   Float_t rmin;
337   Float_t rmax;
338   Float_t zmax1;
339   Float_t zmax2;
340
341   Int_t   iPlan;
342  
343   const Int_t kColorTRD = 46;
344   
345   // Find the top node alice
346   top = gAlice->GetGeometry()->GetNode("alice");
347   
348   if      (fDisplayType == 0) {
349
350     pgon = new TPGON("S_TRD","TRD","void",0,360,AliTRDgeometry::Nsect(),4);
351     rmin = kRmin;
352     rmax = kRmax;
353     pgon->DefineSection(0,-kZmax1,rmax,rmax);
354     pgon->DefineSection(1,-kZmax2,rmin,rmax);
355     pgon->DefineSection(2, kZmax2,rmin,rmax);
356     pgon->DefineSection(3, kZmax1,rmax,rmax);
357     top->cd();
358     node = new TNode("TRD","TRD","S_TRD",0,0,0,"");
359     node->SetLineColor(kColorTRD);
360     fNodes->Add(node);
361
362   }
363   else if (fDisplayType == 1) {
364
365     Char_t name[7];
366
367     Float_t slope = (kZmax1 - kZmax2) / (kRmax  - kRmin);
368
369     rmin  = kRmin + AliTRDgeometry::CraHght();
370     rmax  = rmin  + AliTRDgeometry::CdrHght();
371
372     Float_t thickness = rmin - kRmin;
373     zmax2 = kZmax2 + slope * thickness;
374     zmax1 = zmax2 + slope * AliTRDgeometry::DrThick();
375
376     for (iPlan = 0; iPlan < AliTRDgeometry::Nplan(); iPlan++) {
377
378       sprintf(name,"S_TR1%d",iPlan);
379       pgon  = new TPGON(name,"TRD","void",0,360,AliTRDgeometry::Nsect(),4);
380       pgon->DefineSection(0,-zmax1,rmax,rmax);
381       pgon->DefineSection(1,-zmax2,rmin,rmax);
382       pgon->DefineSection(2, zmax2,rmin,rmax);
383       pgon->DefineSection(3, zmax1,rmax,rmax);
384       top->cd();
385       node = new TNode("TRD","TRD",name,0,0,0,"");
386       node->SetLineColor(kColorTRD);
387       fNodes->Add(node);
388
389       Float_t height = AliTRDgeometry::Cheight() + AliTRDgeometry::Cspace(); 
390       rmin  = rmin  + height;
391       rmax  = rmax  + height;
392       zmax1 = zmax1 + slope * height;
393       zmax2 = zmax2 + slope * height;
394
395     }
396
397     thickness += AliTRDgeometry::DrThick();
398     rmin       = kRmin  + thickness;
399     rmax       = rmin   + AliTRDgeometry::AmThick();
400     zmax2      = kZmax2 + slope * thickness;
401     zmax1      = zmax2  + slope * AliTRDgeometry::AmThick();
402
403     for (iPlan = 0; iPlan < AliTRDgeometry::Nplan(); iPlan++) {
404
405       sprintf(name,"S_TR2%d",iPlan);
406       pgon  = new TPGON(name,"TRD","void",0,360,AliTRDgeometry::Nsect(),4);
407       pgon->DefineSection(0,-zmax1,rmax,rmax);
408       pgon->DefineSection(1,-zmax2,rmin,rmax);
409       pgon->DefineSection(2, zmax2,rmin,rmax);
410       pgon->DefineSection(3, zmax1,rmax,rmax);
411       top->cd();
412       node = new TNode("TRD","TRD",name,0,0,0,"");
413       node->SetLineColor(kColorTRD);
414       fNodes->Add(node);
415
416       Float_t height = AliTRDgeometry::Cheight() + AliTRDgeometry::Cspace(); 
417       rmin  = rmin  + height;
418       rmax  = rmax  + height;
419       zmax1 = zmax1 + slope * height;
420       zmax2 = zmax2 + slope * height;
421
422     }
423
424   }
425
426 }
427  
428 //_____________________________________________________________________________
429 void AliTRD::CreateGeometry()
430 {
431   //
432   // Creates the volumes for the TRD chambers
433   //
434
435   // Check that FRAME is there otherwise we have no place where to put the TRD
436   AliModule *frame = gAlice->GetModule("FRAME");
437   if (!frame) {
438     AliFatal("The TRD needs the FRAME to be defined first");
439   }
440
441   fGeometry->CreateGeometry(fIdtmed->GetArray() - 1299);
442
443 }
444  
445 //_____________________________________________________________________________
446 void AliTRD::CreateMaterials()
447 {
448   //
449   // Create the materials for the TRD
450   //
451
452   Int_t   isxfld = gAlice->Field()->Integ();
453   Float_t sxmgmx = gAlice->Field()->Max();
454   
455   // For polyethilene (CH2) 
456   Float_t ape[2]    = { 12.011 ,  1.0079 };
457   Float_t zpe[2]    = {  6.0   ,  1.0    };
458   Float_t wpe[2]    = {  1.0   ,  2.0    };
459   Float_t dpe       = 0.95;
460
461   // For CO2 
462   Float_t aco[2]    = { 12.011 , 15.9994 };
463   Float_t zco[2]    = {  6.0   ,  8.0    };
464   Float_t wco[2]    = {  1.0   ,  2.0    };
465   Float_t dco       = 0.00186;
466
467   // For water
468   Float_t awa[2]    = {  1.0079, 15.9994 };
469   Float_t zwa[2]    = {  1.0   ,  8.0    };
470   Float_t wwa[2]    = {  2.0   ,  1.0    };
471   Float_t dwa       = 1.0;
472
473   // For isobutane (C4H10)
474   Float_t ais[2]    = { 12.011 ,  1.0079 };
475   Float_t zis[2]    = {  6.0   ,  1.0    };
476   Float_t wis[2]    = {  4.0   , 10.0    };
477   Float_t dis       = 0.00267;
478
479   // For plexiglas (C5H8O2)
480   Float_t apg[3]    = { 12.011 ,  1.0079, 15.9994 };
481   Float_t zpg[3]    = {  6.0   ,  1.0   ,  8.0    };
482   Float_t wpg[3]    = {  5.0   ,  8.0   ,  2.0    };
483   Float_t dpg       = 1.18; 
484   
485   // For epoxy (C18H19O3)
486   Float_t aEpoxy[3] = { 15.9994,  1.0079, 12.011  }; 
487   Float_t zEpoxy[3] = {  8.0   ,  1.0   ,  6.0    }; 
488   Float_t wEpoxy[3] = {  3.0   , 19.0   , 18.0    }; 
489   Float_t dEpoxy    = 1.8 ; 
490
491   // For Araldite, low density epoxy (C18H19O3)
492   Float_t aAral[3]  = { 15.9994,  1.0079, 12.011  }; 
493   Float_t zAral[3]  = {  8.0   ,  1.0   ,  6.0    }; 
494   Float_t wAral[3]  = {  3.0   , 19.0   , 18.0    }; 
495   Float_t dAral     = 1.05; 
496
497   // For air  
498   Float_t aAir[4]   = { 12.011   , 14.0     , 15.9994  , 36.0      };
499   Float_t zAir[4]   = {  6.0     ,  7.0     ,  8.0     , 18.0      };
500   Float_t wAir[4]   = {  0.000124,  0.755267,  0.231781,  0.012827 };
501   Float_t dAir      = 1.20479e-03;
502
503   // For G10
504   Float_t aG10[4]   = {  1.0079  , 12.011   , 15.9994  , 28.086    };
505   Float_t zG10[4]   = {  1.0     ,  6.0     ,  8.0     , 14.0      };
506   Float_t wG10[4]   = {  0.15201 ,  0.10641 ,  0.49444 ,  0.24714  };
507   Float_t dG10      = 1.7;
508
509   // For Xe/CO2-gas-mixture 
510   Float_t aXeCO2[3] = { 131.29   ,  12.0107 ,  15.9994  };
511   Float_t zXeCO2[3] = {  54.0    ,   6.0    ,   8.0     };
512   Float_t wXeCO2[3] = {   8.5    ,   1.5    ,   3.0     }; 
513   // Xe-content of the Xe/CO2-mixture (85% / 15%) 
514   Float_t fxc       = 0.85;
515   Float_t dxe       = 0.00549;
516   Float_t dgm       = fxc * dxe + (1.0 - fxc) * dco;
517   
518   // General tracking parameter
519   Float_t tmaxfd    = -10.0;
520   Float_t stemax    = -1.0e10;
521   Float_t deemax    = -0.1;
522   Float_t epsil     =  1.0e-4;
523   Float_t stmin     = -0.001;
524   
525   //////////////////////////////////////////////////////////////////////////
526   //     Define Materials 
527   //////////////////////////////////////////////////////////////////////////
528
529   AliMaterial( 1, "Al"   ,  26.98, 13.0, 2.7     ,     8.9 ,    37.2);
530   AliMaterial( 4, "Xe"   , 131.29, 54.0, dxe     ,  1546.16,     0.0);
531   AliMaterial( 5, "Cu"   ,  63.54, 29.0, 8.96    ,     1.43,    14.8);
532   AliMaterial( 6, "C"    ,  12.01,  6.0, 2.265   ,    18.8 ,    74.4);
533   AliMaterial(15, "Sn"   , 118.71, 50.0, 7.31    ,     1.21,    14.8);
534   AliMaterial(16, "Si"   ,  28.09, 14.0, 2.33    ,     9.36,    37.2);
535   AliMaterial(18, "Fe"   ,  55.85, 26.0, 7.87    ,     1.76,    14.8);
536
537   // Mixtures 
538   AliMixture(2, "Air"         , aAir,   zAir,   dAir,    4, wAir  );
539   AliMixture(3, "Polyethilene", ape,    zpe,    dpe,    -2, wpe   );
540   AliMixture(8, "CO2",          aco,    zco,    dco,    -2, wco   );
541   AliMixture(9, "Isobutane",    ais,    zis,    dis,    -2, wis   );
542   AliMixture(10,"Gas mixture",  aXeCO2, zXeCO2, dgm,    -3, wXeCO2);
543   AliMixture(12,"G10",          aG10,   zG10,   dG10,    4, wG10  );
544   AliMixture(13,"Water",        awa,    zwa,    dwa,    -2, wwa   );
545   AliMixture(14,"Plexiglas",    apg,    zpg,    dpg,    -3, wpg   );
546   AliMixture(17,"Epoxy",        aEpoxy, zEpoxy, dEpoxy, -3, wEpoxy);
547   AliMixture(19,"Araldite",     aAral,  zAral,  dAral,  -3, wAral );
548
549   //////////////////////////////////////////////////////////////////////////
550   //     Tracking Media Parameters 
551   //////////////////////////////////////////////////////////////////////////
552
553   // Al Frame 
554   AliMedium( 1,"Al Frame"   , 1,0,isxfld,sxmgmx
555               ,tmaxfd,stemax,deemax,epsil,stmin);
556   // Air 
557   AliMedium( 2,"Air"        , 2,0,isxfld,sxmgmx
558               ,tmaxfd,stemax,deemax,epsil,stmin);
559   // Wires
560   AliMedium( 3,"Wires"      , 5,0,isxfld,sxmgmx
561               ,tmaxfd,stemax,deemax,epsil,stmin);
562   // All other ROB materials (caps, etc.)
563   AliMedium( 4,"ROB Other"  , 5,0,isxfld,sxmgmx
564               ,tmaxfd,stemax,deemax,epsil,stmin);
565   // Cu pads 
566   AliMedium( 5,"Padplane"   , 5,1,isxfld,sxmgmx
567               ,tmaxfd,stemax,deemax,epsil,stmin);
568   // Fee + cables 
569   AliMedium( 6,"Readout"    , 5,0,isxfld,sxmgmx
570               ,tmaxfd,stemax,deemax,epsil,stmin);
571   // C frame 
572   AliMedium( 7,"C Frame"    , 6,0,isxfld,sxmgmx
573               ,tmaxfd,stemax,deemax,epsil,stmin);
574   // INOX of cooling bus bars
575   AliMedium( 8,"Cooling bus",18,0,isxfld,sxmgmx
576               ,tmaxfd,stemax,deemax,epsil,stmin);
577   // Gas-mixture (Xe/CO2) 
578   AliMedium( 9,"Gas-mix"    ,10,1,isxfld,sxmgmx
579               ,tmaxfd,stemax,deemax,epsil,stmin);
580   // Nomex-honeycomb
581   AliMedium(10,"Nomex"      ,12,0,isxfld,sxmgmx
582               ,tmaxfd,stemax,deemax,epsil,stmin);
583   // Araldite glue
584   AliMedium(11,"Glue"       ,19,0,isxfld,sxmgmx
585               ,tmaxfd,stemax,deemax,epsil,stmin);
586   // G10-plates
587   AliMedium(13,"G10-plates" ,12,0,isxfld,sxmgmx
588               ,tmaxfd,stemax,deemax,epsil,stmin);
589   // Cooling water
590   AliMedium(14,"Water"      ,13,0,isxfld,sxmgmx
591               ,tmaxfd,stemax,deemax,epsil,stmin);
592   // Rohacell (plexiglas) for the radiator
593   AliMedium(15,"Rohacell"   ,14,0,isxfld,sxmgmx
594               ,tmaxfd,stemax,deemax,epsil,stmin);
595   // Al layer in MCMs
596   AliMedium(16,"MCM-Al"     , 1,0,isxfld,sxmgmx
597               ,tmaxfd,stemax,deemax,epsil,stmin);
598   // Sn layer in MCMs
599   AliMedium(17,"MCM-Sn"     ,15,0,isxfld,sxmgmx
600               ,tmaxfd,stemax,deemax,epsil,stmin);
601   // Cu layer in MCMs
602   AliMedium(18,"MCM-Cu"     , 5,0,isxfld,sxmgmx
603               ,tmaxfd,stemax,deemax,epsil,stmin);
604   // G10 layer in MCMs
605   AliMedium(19,"MCM-G10"    ,12,0,isxfld,sxmgmx
606               ,tmaxfd,stemax,deemax,epsil,stmin);
607   // Si in readout chips
608   AliMedium(20,"Chip-Si"    ,16,0,isxfld,sxmgmx
609               ,tmaxfd,stemax,deemax,epsil,stmin);
610   // Epoxy in readout chips
611   AliMedium(21,"Chip-Ep"    ,17,0,isxfld,sxmgmx
612               ,tmaxfd,stemax,deemax,epsil,stmin);
613   // PE in connectors
614   AliMedium(22,"Conn-PE"    , 3,0,isxfld,sxmgmx
615               ,tmaxfd,stemax,deemax,epsil,stmin);
616   // Cu in connectors
617   AliMedium(23,"Chip-Cu"    , 5,0,isxfld,sxmgmx
618               ,tmaxfd,stemax,deemax,epsil,stmin);
619   // Al of cooling pipes
620   AliMedium(24,"Cooling"    , 1,0,isxfld,sxmgmx
621               ,tmaxfd,stemax,deemax,epsil,stmin);
622   // Cu in services
623   AliMedium(25,"Serv-Cu"    , 5,0,isxfld,sxmgmx
624               ,tmaxfd,stemax,deemax,epsil,stmin);
625
626   // Save the density values for the TRD absorbtion
627   Float_t dmy  = 1.39;
628   fFoilDensity = dmy;
629   fGasDensity  = dgm;
630
631 }
632
633 //_____________________________________________________________________________
634 void AliTRD::DrawModule() const
635 {
636   //
637   // Draw a shaded view of the Transition Radiation Detector version 0
638   //
639
640   // Set everything unseen
641   gMC->Gsatt("*"   ,"SEEN",-1);
642   
643   // Set ALIC mother transparent
644   gMC->Gsatt("ALIC","SEEN", 0);
645   
646   // Set the volumes visible
647   if (fGeometry->IsVersion() == 0) {
648     gMC->Gsatt("B071","SEEN", 0);
649     gMC->Gsatt("B074","SEEN", 0);
650     gMC->Gsatt("B075","SEEN", 0);
651     gMC->Gsatt("B077","SEEN", 0);
652     gMC->Gsatt("BTR1","SEEN", 0);
653     gMC->Gsatt("BTR2","SEEN", 0);
654     gMC->Gsatt("BTR3","SEEN", 0);
655     gMC->Gsatt("UTR1","SEEN", 0);
656     gMC->Gsatt("UTR2","SEEN", 0);
657     gMC->Gsatt("UTR3","SEEN", 0);
658   }
659   else {
660     gMC->Gsatt("B071","SEEN", 0);
661     gMC->Gsatt("B074","SEEN", 0);
662     gMC->Gsatt("B075","SEEN", 0);
663     gMC->Gsatt("B077","SEEN", 0);
664     gMC->Gsatt("BTR1","SEEN", 0);
665     gMC->Gsatt("BTR2","SEEN", 0);
666     gMC->Gsatt("BTR3","SEEN", 0);
667     gMC->Gsatt("UTR1","SEEN", 0);
668   }
669   
670   gMC->Gdopt("hide", "on");
671   gMC->Gdopt("shad", "on");
672   gMC->Gsatt("*", "fill", 7);
673   gMC->SetClipBox(".");
674   gMC->SetClipBox("*", 0, 2000, -2000, 2000, -2000, 2000);
675   gMC->DefaultRange();
676   gMC->Gdraw("alic", 40, 30, 0, 12, 9.4, .021, .021);
677   gMC->Gdhead(1111, "Transition Radiation Detector");
678   gMC->Gdman(18, 4, "MAN");
679
680 }
681
682 //_____________________________________________________________________________
683 Int_t AliTRD::DistancetoPrimitive(Int_t , Int_t )
684 {
685   //
686   // Distance between the mouse and the TRD detector on the screen
687   // Dummy routine
688   //
689   
690   return 9999;
691
692 }
693  
694 //_____________________________________________________________________________
695 void AliTRD::Init()
696 {
697   //
698   // Initialize the TRD detector after the geometry has been created
699   //
700
701   AliDebug(1,"++++++++++++++++++++++++++++++++++++++++++++++");
702
703   if (fGeometry->IsVersion() != 1) {
704     AliError("Not a valid geometry");
705   }
706   // Special tracking options for charged particles for XeCO2
707   gMC->Gstpar((* fIdtmed)[9],"DRAY",1.0);
708   gMC->Gstpar((* fIdtmed)[9],"STRA",1.0); 
709   
710 }
711
712 //_____________________________________________________________________________
713 void AliTRD::LoadPoints(Int_t )
714 {
715   //
716   // Store x, y, z of all hits in memory.
717   // Hit originating from TR photons are given a different color
718   //
719
720   if (fHits == 0) {
721     return;
722   }
723
724   Int_t nhits  = fHits->GetEntriesFast();
725   if (nhits == 0) {
726     return;
727   }
728
729   Int_t tracks = gAlice->GetMCApp()->GetNtrack();
730   if (fPoints == 0) {
731     fPoints = new TObjArray(tracks);
732   }
733
734   AliTRDhit *ahit;
735   
736   Int_t    *ntrkE = new Int_t[tracks];
737   Int_t    *ntrkT = new Int_t[tracks];
738   Int_t    *limiE = new Int_t[tracks];
739   Int_t    *limiT = new Int_t[tracks];
740   Float_t **coorE = new Float_t*[tracks];
741   Float_t **coorT = new Float_t*[tracks];
742   for(Int_t i = 0; i < tracks; i++) {
743     ntrkE[i] = 0;
744     ntrkT[i] = 0;
745     coorE[i] = 0;
746     coorT[i] = 0;
747     limiE[i] = 0;
748     limiT[i] = 0;
749   }
750   
751   AliTRDpoints *points = 0;
752   Float_t      *fp     = 0;
753   Int_t         trk;
754   Int_t         chunk  = nhits / 4 + 1;
755
756   // Loop over all the hits and store their position
757   ahit = (AliTRDhit *) FirstHit(-1);
758   while (ahit) {
759
760     // dEdx hits
761     if (ahit->GetCharge() >= 0) {
762
763       trk = ahit->GetTrack();
764       if (ntrkE[trk] == limiE[trk]) {
765         // Initialise a new track
766         fp = new Float_t[3*(limiE[trk]+chunk)];
767         if (coorE[trk]) {
768           memcpy(fp,coorE[trk],sizeof(Float_t)*3*limiE[trk]);
769           delete [] coorE[trk];
770         }
771         limiE[trk] += chunk;
772         coorE[trk]  = fp;
773       } 
774       else {
775         fp = coorE[trk];
776       }
777       fp[3*ntrkE[trk]  ] = ahit->X();
778       fp[3*ntrkE[trk]+1] = ahit->Y();
779       fp[3*ntrkE[trk]+2] = ahit->Z();
780       ntrkE[trk]++;
781
782     }
783     // TR photon hits
784     else if ((ahit->GetCharge() < 0) && 
785              (fDrawTR)) {
786
787       trk = ahit->GetTrack();
788       if (ntrkT[trk] == limiT[trk]) {
789         // Initialise a new track
790         fp = new Float_t[3*(limiT[trk]+chunk)];
791         if (coorT[trk]) {
792           memcpy(fp,coorT[trk],sizeof(Float_t)*3*limiT[trk]);
793           delete [] coorT[trk];
794         }
795         limiT[trk] += chunk;
796         coorT[trk]  = fp;
797       } 
798       else {
799         fp = coorT[trk];
800       }
801       fp[3*ntrkT[trk]  ] = ahit->X();
802       fp[3*ntrkT[trk]+1] = ahit->Y();
803       fp[3*ntrkT[trk]+2] = ahit->Z();
804       ntrkT[trk]++;
805
806     }
807
808     ahit = (AliTRDhit *) NextHit();
809
810   }
811
812   for (trk = 0; trk < tracks; ++trk) {
813
814     if (ntrkE[trk] || ntrkT[trk]) {
815
816       points = new AliTRDpoints();
817       points->SetDetector(this);
818       points->SetParticle(trk);
819
820       // Set the dEdx points
821       if (ntrkE[trk]) {
822         points->SetMarkerColor(kWhite); //PH This is the default color in TRD
823         points->SetMarkerSize(1); //PH Default size=1
824         points->SetPolyMarker(ntrkE[trk],coorE[trk],1); //PH Default style=1
825         delete [] coorE[trk];
826         coorE[trk] = 0;
827       }
828
829       // Set the TR photon points
830       if (ntrkT[trk]) {
831         points->SetTRpoints(ntrkT[trk],coorT[trk]);
832         delete [] coorT[trk];
833         coorT[trk] = 0;
834       }
835
836       fPoints->AddAt(points,trk);
837
838     }
839
840   }
841
842   delete [] coorE;
843   delete [] coorT;
844   delete [] ntrkE;
845   delete [] ntrkT;
846   delete [] limiE;
847   delete [] limiT;
848
849 }
850
851 //_____________________________________________________________________________
852 void AliTRD::MakeBranch(Option_t *option)
853 {
854   //
855   // Create Tree branches for the TRD digits.
856   //
857
858   Int_t  buffersize = 4000;
859   Char_t branchname[15];
860   sprintf(branchname,"%s",GetName());
861
862   const Char_t *cD = strstr(option,"D");
863
864   AliDetector::MakeBranch(option);
865
866   if (fDigits         && 
867       gAlice->TreeD() && 
868       cD) {
869     MakeBranchInTree(gAlice->TreeD(),branchname,&fDigits,buffersize,0);
870   }
871
872 }
873
874 //_____________________________________________________________________________
875 void AliTRD::ResetDigits()
876 {
877   //
878   // Reset number of digits and the digits array for this detector
879   //
880
881   fNdigits = 0;
882
883   if (fDigits) {
884     fDigits->Clear();
885   }
886
887 }
888
889 //_____________________________________________________________________________
890 void AliTRD::SetTreeAddress()
891 {
892   //
893   // Set the branch addresses for the trees.
894   //
895
896   if (fLoader->TreeH() && 
897       (fHits == 0x0)) {
898     fHits = new TClonesArray("AliTRDhit",405);
899   }
900   AliDetector::SetTreeAddress();
901
902 }
903
904 //_____________________________________________________________________________
905 AliTRD &AliTRD::operator=(const AliTRD &trd)
906 {
907   //
908   // Assignment operator
909   //
910
911   if (this != &trd) {
912     ((AliTRD &) trd).Copy(*this);
913   }
914
915   return *this;
916
917
918
919
920
921
922
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