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 "AliMagF.h"
  44 #include "AliRun.h"
  45 #include "AliTRD.h"
  46 #include "AliTRDdigit.h"
  47 #include "AliTRDdigitizer.h"
  48 #include "AliTRDdigitsManager.h"
  49 #include "AliTRDgeometryFull.h"
  50 #include "AliTRDhit.h"
  51 #include "AliTRDpoints.h"
  52 #include "AliTRDrawData.h"
  53 #include "AliTrackReference.h"
  54 #include "AliMC.h"
  55
  56 ClassImp(AliTRD)
  57
  58 //_____________________________________________________________________________
  59 AliTRD::AliTRD()
  60 {
  61   //
  62   // Default constructor
  63   //
  64
  65   fIshunt        = 0;
  66   fGasMix        = 0;
  67   fHits          = 0;
  68   fDigits        = 0;
  69
  70   fGeometry      = 0;
  71
  72   fGasDensity    = 0;
  73   fFoilDensity   = 0;
  74
  75   fDrawTR        = 0;
  76   fDisplayType   = 0;
  77
  78 }
  79
  80 //_____________________________________________________________________________
  81 AliTRD::AliTRD(const char *name, const char *title)
  82        : AliDetector(name,title)
  83 {
  84   //
  85   // Standard constructor for the TRD
  86   //
  87
  88   // Check that FRAME is there otherwise we have no place where to
  89   // put TRD
  90   AliModule* frame = gAlice->GetModule("FRAME");
  91   if (!frame) {
  92     Error("Ctor","TRD needs FRAME to be present\n");
  93     exit(1);
  94   }
  95
  96   // Define the TRD geometry
  97   if ((frame->IsVersion() == 0) ||
  98       (frame->IsVersion() == 1)) {
  99     fGeometry = new AliTRDgeometryFull();
 100   }
 101   else {
 102     Error("Ctor","Could not find valid FRAME version\n");
 103     exit(1);
 104   }
 105
 106   // Save the geometry
 107   TDirectory* saveDir = gDirectory;
 108   gAlice->GetRunLoader()->CdGAFile();
 109   fGeometry->Write("TRDgeometry");
 110   saveDir->cd();
 111
 112   // Allocate the hit array
 113   fHits           = new TClonesArray("AliTRDhit"     ,405);
 114   gAlice->GetMCApp()->AddHitList(fHits);
 115
 116   // Allocate the digits array
 117   fDigits         = 0;
 118
 119   fIshunt        = 0;
 120   fGasMix        = 1;
 121
 122   fGasDensity    = 0;
 123   fFoilDensity   = 0;
 124
 125   fDrawTR        = 0;
 126   fDisplayType   = 0;
 127
 128   SetMarkerColor(kWhite);
 129
 130 }
 131
 132 //_____________________________________________________________________________
 133 AliTRD::AliTRD(const AliTRD &trd):AliDetector(trd)
 134 {
 135   //
 136   // Copy constructor
 137   //
 138
 139   ((AliTRD &) trd).Copy(*this);
 140
 141 }
 142
 143 //_____________________________________________________________________________
 144 AliTRD::~AliTRD()
 145 {
 146   //
 147   // TRD destructor
 148   //
 149
 150   fIshunt = 0;
 151
 152   if (fGeometry) {
 153     delete fGeometry;
 154     fGeometry  = 0;
 155   }
 156   if (fHits) {
 157     delete fHits;
 158     fHits      = 0;
 159   }
 160
 161 }
 162
 163 //_____________________________________________________________________________
 164 void AliTRD::Hits2Digits()
 165 {
 166   //
 167   // Create digits
 168   //
 169   AliTRDdigitizer digitizer("TRDdigitizer","TRD digitizer class");
 170   digitizer.SetDebug(GetDebug());
 171
 172   // Initialization
 173   digitizer.InitDetector();
 174
 175   if (!fLoader->TreeH()) fLoader->LoadHits("read");
 176   fLoader->LoadDigits("recreate");
 177   AliRunLoader* runLoader = fLoader->GetRunLoader();
 178
 179   for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
 180     digitizer.Open(runLoader->GetFileName().Data(), iEvent);
 181     digitizer.MakeDigits();
 182     digitizer.WriteDigits();
 183   }
 184
 185   fLoader->UnloadHits();
 186   fLoader->UnloadDigits();
 187
 188 }
 189
 190 //_____________________________________________________________________________
 191 void AliTRD::Hits2SDigits()
 192 {
 193   //
 194   // Create summable digits
 195   //
 196
 197   AliTRDdigitizer digitizer("TRDdigitizer","TRD digitizer class");
 198   // For the summable digits
 199   digitizer.SetSDigits(kTRUE);
 200   digitizer.SetDebug(GetDebug());
 201
 202   // Initialization
 203   digitizer.InitDetector();
 204
 205   if (!fLoader->TreeH()) fLoader->LoadHits("read");
 206   fLoader->LoadSDigits("recreate");
 207   AliRunLoader* runLoader = fLoader->GetRunLoader();
 208
 209   for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
 210     digitizer.Open(runLoader->GetFileName().Data(), iEvent);
 211     digitizer.MakeDigits();
 212     digitizer.WriteDigits();
 213   }
 214
 215   fLoader->UnloadHits();
 216   fLoader->UnloadSDigits();
 217
 218 }
 219
 220 //_____________________________________________________________________________
 221 AliDigitizer* AliTRD::CreateDigitizer(AliRunDigitizer* manager) const
 222 {
 223   //
 224   // Creates a new digitizer object
 225   //
 226
 227   return new AliTRDdigitizer(manager);
 228
 229 }
 230
 231 //_____________________________________________________________________________
 232 void AliTRD::SDigits2Digits()
 233 {
 234   //
 235   // Create final digits from summable digits
 236   //
 237
 238    // Create the TRD digitizer
 239   AliTRDdigitizer digitizer("TRDdigitizer","TRD digitizer class");
 240   digitizer.SetDebug(GetDebug());
 241
 242   // Set the parameter
 243   digitizer.SetEvent(gAlice->GetEvNumber());
 244
 245   // Initialization
 246   digitizer.InitDetector();
 247
 248   // Read the s-digits via digits manager
 249   AliTRDdigitsManager sdigitsManager;
 250
 251   sdigitsManager.SetDebug(GetDebug());
 252   sdigitsManager.SetSDigits(kTRUE);
 253   sdigitsManager.CreateArrays();
 254
 255   if (!fLoader->TreeS())
 256     if (fLoader->LoadSDigits("read"))
 257      {
 258        Error("SDigits2Digits","Error while reading SDigits for event %d",gAlice->GetEvNumber());
 259        return;
 260      }
 261   if (!fLoader->TreeS()) return;
 262
 263   sdigitsManager.ReadDigits(fLoader->TreeS());
 264
 265   // Add the s-digits to the input list
 266   digitizer.AddSDigitsManager(&sdigitsManager);
 267
 268   // Convert the s-digits to normal digits
 269   digitizer.SDigits2Digits();
 270
 271   // Store the digits
 272   if (!fLoader->TreeD()) fLoader->MakeTree("D");
 273   if (digitizer.MakeBranch(fLoader->TreeD())){
 274     digitizer.WriteDigits();
 275   }
 276
 277 }
 278
 279 //_____________________________________________________________________________
 280 void AliTRD::Digits2Raw()
 281 {
 282   //
 283   // convert digits of the current event to raw data
 284   //
 285
 286   fLoader->LoadDigits();
 287   TTree* digits = fLoader->TreeD();
 288   if (!digits) {
 289     Error("Digits2Raw", "no digits tree");
 290     return;
 291   }
 292
 293   AliTRDrawData rawWriter;
 294 //  rawWriter.SetDebug(2);
 295   if (!rawWriter.Digits2Raw(digits)) {
 296     Error("AliTRD::Digits2Raw","The raw writer could not load the digits tree");
 297   }
 298
 299   fLoader->UnloadDigits();
 300
 301 }
 302
 303 //_____________________________________________________________________________
 304 void AliTRD::AddHit(Int_t track, Int_t det, Float_t *hits, Int_t q
 305                   , Bool_t inDrift)
 306 {
 307   //
 308   // Add a hit for the TRD
 309   //
 310
 311   TClonesArray &lhits = *fHits;
 312   AliTRDhit *hit = new(lhits[fNhits++]) AliTRDhit(fIshunt,track,det,hits,q);
 313   if (inDrift) {
 314     hit->SetDrift();
 315   }
 316   else {
 317     hit->SetAmplification();
 318   }
 319   if (q < 0) {
 320     hit->SetTRphoton();
 321   }
 322
 323 }
 324
 325 //_____________________________________________________________________________
 326 void AliTRD::BuildGeometry()
 327 {
 328   //
 329   // Create the ROOT TNode geometry for the TRD
 330   //
 331
 332   TNode *node, *top;
 333   TPGON *pgon;
 334
 335   Float_t rmin, rmax;
 336   Float_t zmax1, zmax2;
 337
 338   Int_t   iPlan;
 339
 340   const Int_t kColorTRD = 46;
 341
 342   // Find the top node alice
 343   top = gAlice->GetGeometry()->GetNode("alice");
 344
 345   if      (fDisplayType == 0) {
 346
 347     pgon = new TPGON("S_TRD","TRD","void",0,360,AliTRDgeometry::Nsect(),4);
 348     rmin = AliTRDgeometry::Rmin();
 349     rmax = AliTRDgeometry::Rmax();
 350     pgon->DefineSection(0,-AliTRDgeometry::Zmax1(),rmax,rmax);
 351     pgon->DefineSection(1,-AliTRDgeometry::Zmax2(),rmin,rmax);
 352     pgon->DefineSection(2, AliTRDgeometry::Zmax2(),rmin,rmax);
 353     pgon->DefineSection(3, AliTRDgeometry::Zmax1(),rmax,rmax);
 354     top->cd();
 355     node = new TNode("TRD","TRD","S_TRD",0,0,0,"");
 356     node->SetLineColor(kColorTRD);
 357     fNodes->Add(node);
 358
 359   }
 360   else if (fDisplayType == 1) {
 361
 362     Char_t name[7];
 363
 364     Float_t slope = (AliTRDgeometry::Zmax1() - AliTRDgeometry::Zmax2())
 365                   / (AliTRDgeometry::Rmax()  - AliTRDgeometry::Rmin());
 366
 367     rmin  = AliTRDgeometry::Rmin() + AliTRDgeometry::CraHght();
 368     rmax  = rmin                   + AliTRDgeometry::CdrHght();
 369
 370     Float_t thickness = rmin - AliTRDgeometry::Rmin();
 371     zmax2 = AliTRDgeometry::Zmax2() + slope * thickness;
 372     zmax1 = zmax2 + slope * AliTRDgeometry::DrThick();
 373
 374     for (iPlan = 0; iPlan < AliTRDgeometry::Nplan(); iPlan++) {
 375
 376       sprintf(name,"S_TR1%d",iPlan);
 377       pgon  = new TPGON(name,"TRD","void",0,360,AliTRDgeometry::Nsect(),4);
 378       pgon->DefineSection(0,-zmax1,rmax,rmax);
 379       pgon->DefineSection(1,-zmax2,rmin,rmax);
 380       pgon->DefineSection(2, zmax2,rmin,rmax);
 381       pgon->DefineSection(3, zmax1,rmax,rmax);
 382       top->cd();
 383       node = new TNode("TRD","TRD",name,0,0,0,"");
 384       node->SetLineColor(kColorTRD);
 385       fNodes->Add(node);
 386
 387       Float_t height = AliTRDgeometry::Cheight() + AliTRDgeometry::Cspace();
 388       rmin  = rmin  + height;
 389       rmax  = rmax  + height;
 390       zmax1 = zmax1 + slope * height;
 391       zmax2 = zmax2 + slope * height;
 392
 393     }
 394
 395     thickness += AliTRDgeometry::DrThick();
 396     rmin  = AliTRDgeometry::Rmin() + thickness;
 397     rmax  = rmin + AliTRDgeometry::AmThick();
 398     zmax2 = AliTRDgeometry::Zmax2() + slope * thickness;
 399     zmax1 = zmax2 + slope * AliTRDgeometry::AmThick();
 400
 401     for (iPlan = 0; iPlan < AliTRDgeometry::Nplan(); iPlan++) {
 402
 403       sprintf(name,"S_TR2%d",iPlan);
 404       pgon  = new TPGON(name,"TRD","void",0,360,AliTRDgeometry::Nsect(),4);
 405       pgon->DefineSection(0,-zmax1,rmax,rmax);
 406       pgon->DefineSection(1,-zmax2,rmin,rmax);
 407       pgon->DefineSection(2, zmax2,rmin,rmax);
 408       pgon->DefineSection(3, zmax1,rmax,rmax);
 409       top->cd();
 410       node = new TNode("TRD","TRD",name,0,0,0,"");
 411       node->SetLineColor(kColorTRD);
 412       fNodes->Add(node);
 413
 414       Float_t height = AliTRDgeometry::Cheight() + AliTRDgeometry::Cspace();
 415       rmin  = rmin  + height;
 416       rmax  = rmax  + height;
 417       zmax1 = zmax1 + slope * height;
 418       zmax2 = zmax2 + slope * height;
 419
 420     }
 421
 422   }
 423
 424 }
 425
 426 //_____________________________________________________________________________
 427 void AliTRD::Copy(TObject &trd) const
 428 {
 429   //
 430   // Copy function
 431   //
 432
 433   ((AliTRD &) trd).fGasMix      = fGasMix;
 434   ((AliTRD &) trd).fGeometry    = fGeometry;
 435   ((AliTRD &) trd).fGasDensity  = fGasDensity;
 436   ((AliTRD &) trd).fFoilDensity = fFoilDensity;
 437   ((AliTRD &) trd).fDrawTR      = fDrawTR;
 438   ((AliTRD &) trd).fDisplayType = fDisplayType;
 439
 440   //AliDetector::Copy(trd);
 441
 442 }
 443
 444 //_____________________________________________________________________________
 445 void AliTRD::CreateGeometry()
 446 {
 447   //
 448   // Creates the volumes for the TRD chambers
 449   //
 450
 451   // Check that FRAME is there otherwise we have no place where to put the TRD
 452   AliModule* frame = gAlice->GetModule("FRAME");
 453   if (!frame) {
 454     printf(" The TRD needs the FRAME to be defined first\n");
 455     return;
 456   }
 457
 458   fGeometry->CreateGeometry(fIdtmed->GetArray() - 1299);
 459
 460 }
 461
 462 //_____________________________________________________________________________
 463 void AliTRD::CreateMaterials()
 464 {
 465   //
 466   // Create the materials for the TRD
 467   // Origin Y.Foka
 468   //
 469
 470   Int_t   isxfld = gAlice->Field()->Integ();
 471   Float_t sxmgmx = gAlice->Field()->Max();
 472
 473   // For polyethilene (CH2)
 474   Float_t ape[2] = { 12., 1. };
 475   Float_t zpe[2] = {  6., 1. };
 476   Float_t wpe[2] = {  1., 2. };
 477   Float_t dpe    = 0.95;
 478
 479   // For mylar (C5H4O2)
 480   Float_t amy[3] = { 12., 1., 16. };
 481   Float_t zmy[3] = {  6., 1.,  8. };
 482   Float_t wmy[3] = {  5., 4.,  2. };
 483   Float_t dmy    = 1.39;
 484
 485   // For CO2
 486   Float_t aco[2] = { 12., 16. };
 487   Float_t zco[2] = {  6.,  8. };
 488   Float_t wco[2] = {  1.,  2. };
 489   Float_t dco    = 0.001977;
 490
 491   // For water
 492   Float_t awa[2] = {  1., 16. };
 493   Float_t zwa[2] = {  1.,  8. };
 494   Float_t wwa[2] = {  2.,  1. };
 495   Float_t dwa    = 1.0;
 496
 497   // For isobutane (C4H10)
 498   Float_t ais[2] = { 12.,  1. };
 499   Float_t zis[2] = {  6.,  1. };
 500   Float_t wis[2] = {  4., 10. };
 501   Float_t dis    = 0.00267;
 502
 503   // For plexiglas (C5H8O2)
 504   Float_t apg[3] = { 12.011 ,  1.0    , 15.9994 };
 505   Float_t zpg[3] = {  6.0   ,  1.0    ,  8.0    };
 506   Float_t wpg[3] = {  5.0   ,  8.0    ,  2.0    };
 507   Float_t dpg    = 1.18;
 508
 509   // For epoxy (C18H19O3)
 510   Float_t aEpoxy[3] = { 15.9994 ,  1.00794 , 12.0107 };
 511   Float_t zEpoxy[3] = {  8.0    ,  1.0     ,  6.0    };
 512   Float_t wEpoxy[3] = {  3.0    , 19.0     , 18.0    };
 513   Float_t dEpoxy    = 1.8 ;
 514
 515   // For air
 516   Float_t aAir[4]   = { 12.0      , 14.0      , 16.0      , 36.0     };
 517   Float_t zAir[4]   = {  6.0      ,  7.0      ,  8.0      , 18.0     };
 518   Float_t wAir[4]   = {  0.000124 ,  0.755267 ,  0.231781 ,  0.012827};
 519   Float_t dAir      = 1.20479E-3;
 520
 521   // For G10
 522   Float_t aG10[4]   = {  1.0         , 12.011       , 15.9994      , 28.086       };
 523   Float_t zG10[4]   = {  1.0         ,  6.0         ,  8.0         , 14.0         };
 524   //PH  Float_t wG10[4]   = {  0.148648649 ,  0.104054054 ,  0.483499056 ,  0.241666667 };
 525   Float_t wG10[4]={0.15201,0.10641,0.49444,0.24714};
 526   Float_t dG10      = 1.7;
 527
 528   // For Xe/CO2-gas-mixture
 529   Float_t aXeCO2[3] = { 39.948, 12.0107, 15.9994 };
 530   Float_t zXeCO2[3] = { 18.0,    6.0,     8.0    };
 531   Float_t wXeCO2[3] = {  0.85,   0.0375,  0.1125 };
 532   // Xe-content of the Xe/CO2-mixture (85% / 15%)
 533   Float_t fxc    = .85;
 534   // Xe-content of the Xe/Isobutane-mixture (97% / 3%)
 535   Float_t fxi    = .97;
 536   Float_t dxe    = .005858;
 537
 538   // General tracking parameter
 539   Float_t tmaxfd = -10.;
 540   Float_t stemax = -1e10;
 541   Float_t deemax = -0.1;
 542   Float_t epsil  =  1e-4;
 543   Float_t stmin  = -0.001;
 544
 545   Float_t absl, radl, d, buf[1];
 546   Float_t agm[2], zgm[2], wgm[2];
 547   Float_t dgm1, dgm2;
 548   Int_t   nbuf;
 549
 550   //////////////////////////////////////////////////////////////////////////
 551   //     Define Materials
 552   //////////////////////////////////////////////////////////////////////////
 553
 554   AliMaterial( 1, "Al"   ,  26.98, 13.0, 2.7     ,     8.9 ,    37.2);
 555   //AliMaterial( 2, "Air"  ,  14.61,  7.3, 0.001205, 30420.0 , 67500.0);
 556   AliMaterial( 4, "Xe"   , 131.29, 54.0, dxe     ,  1447.59,     0.0);
 557   AliMaterial( 5, "Cu"   ,  63.54, 29.0, 8.96    ,     1.43,    14.8);
 558   AliMaterial( 6, "C"    ,  12.01,  6.0, 2.265   ,    18.8 ,    74.4);
 559   //AliMaterial(12, "G10"  ,  20.00, 10.0, 1.7     ,    19.4 ,   999.0);
 560   AliMaterial(15, "Sn"   , 118.71, 50.0, 7.31    ,     1.21,    14.8);
 561   AliMaterial(16, "Si"   ,  28.09, 14.0, 2.33    ,     9.36,    37.2);
 562   //AliMaterial(17, "Epoxy",  17.75,  8.9, 1.8     ,    21.82,   999.0);
 563
 564   // Mixtures
 565   AliMixture(2, "Air"         ,  aAir,  zAir,  dAir,  4,   wAir);
 566   AliMixture(3, "Polyethilene",   ape,   zpe,   dpe, -2,    wpe);
 567   AliMixture(7, "Mylar",          amy,   zmy,   dmy, -3,    wmy);
 568   AliMixture(8, "CO2",            aco,   zco,   dco, -2,    wco);
 569   AliMixture(9, "Isobutane",      ais,   zis,   dis, -2,    wis);
 570   AliMixture(12,"G10",           aG10,  zG10,  dG10,  4,   wG10);
 571   AliMixture(13,"Water",          awa,   zwa,   dwa, -2,    wwa);
 572   AliMixture(14,"Plexiglas",      apg,   zpg,   dpg, -3,    wpg);
 573   AliMixture(17,"Epoxy",       aEpoxy,zEpoxy,dEpoxy, -3, wEpoxy);
 574
 575   // Gas mixtures
 576   Char_t namate[21]="";
 577   // Xe/CO2-mixture
 578   // Get properties of Xe
 579   //gMC->Gfmate((*fIdmate)[4], namate, agm[0], zgm[0], d, radl, absl, buf, nbuf);
 580   // Get properties of CO2
 581   //gMC->Gfmate((*fIdmate)[8], namate, agm[1], zgm[1], d, radl, absl, buf, nbuf);
 582   // Create gas mixture
 583   wgm[0] = fxc;
 584   wgm[1] = 1. - fxc;
 585   dgm1   = wgm[0] * dxe + wgm[1] * dco;
 586   AliMixture(10, "Gas mixture 1", aXeCO2, zXeCO2, dgm1,  3, wXeCO2);
 587
 588   // Xe/Isobutane-mixture
 589   // Get properties of Xe
 590   gMC->Gfmate((*fIdmate)[4], namate, agm[0], zgm[0], d, radl, absl, buf, nbuf);
 591   // Get properties of Isobutane
 592   gMC->Gfmate((*fIdmate)[9], namate, agm[1], zgm[1], d, radl, absl, buf, nbuf);
 593   // Create gas mixture
 594   wgm[0] = fxi;
 595   wgm[1] = 1. - fxi;
 596   dgm2   = wgm[0] * dxe + wgm[1] * dis;
 597   AliMixture(11, "Gas mixture 2", agm, zgm, dgm2,  2, wgm);
 598
 599   //////////////////////////////////////////////////////////////////////////
 600   //     Tracking Media Parameters
 601   //////////////////////////////////////////////////////////////////////////
 602
 603   // Al Frame
 604   AliMedium(1, "Al Frame",   1, 0, isxfld, sxmgmx
 605                 , tmaxfd, stemax, deemax, epsil, stmin);
 606   // Air
 607   AliMedium(2, "Air",        2, 0, isxfld, sxmgmx
 608                 , tmaxfd, stemax, deemax, epsil, stmin);
 609   // Polyethilene
 610   AliMedium(3, "Radiator",   3, 0, isxfld, sxmgmx
 611                 , tmaxfd, stemax, deemax, epsil, stmin);
 612   // Xe
 613   AliMedium(4, "Xe",         4, 1, isxfld, sxmgmx
 614                 , tmaxfd, stemax, deemax, epsil, stmin);
 615   // Cu pads
 616   AliMedium(5, "Padplane",   5, 1, isxfld, sxmgmx
 617                 , tmaxfd, stemax, deemax, epsil, stmin);
 618   // Fee + cables
 619   AliMedium(6, "Readout",    1, 0, isxfld, sxmgmx
 620                 , tmaxfd, stemax, deemax, epsil, stmin);
 621   // C frame
 622   AliMedium(7, "C Frame",    6, 0, isxfld, sxmgmx
 623                 , tmaxfd, stemax, deemax, epsil, stmin);
 624   // Mylar foils
 625   AliMedium(8, "Mylar",      7, 0, isxfld, sxmgmx
 626                 , tmaxfd, stemax, deemax, epsil, stmin);
 627   if (fGasMix == 1) {
 628     // Gas-mixture (Xe/CO2)
 629     AliMedium(9, "Gas-mix",   10, 1, isxfld, sxmgmx
 630                   , tmaxfd, stemax, deemax, epsil, stmin);
 631   }
 632   else {
 633     // Gas-mixture (Xe/Isobutane)
 634     AliMedium(9, "Gas-mix",   11, 1, isxfld, sxmgmx
 635                   , tmaxfd, stemax, deemax, epsil, stmin);
 636   }
 637   // Nomex-honeycomb (use carbon for the time being)
 638   AliMedium(10, "Nomex",      6, 0, isxfld, sxmgmx
 639                 , tmaxfd, stemax, deemax, epsil, stmin);
 640   // Kapton foils (use Mylar for the time being)
 641   AliMedium(11, "Kapton",     7, 0, isxfld, sxmgmx
 642                 , tmaxfd, stemax, deemax, epsil, stmin);
 643   // Gas-filling of the radiator
 644   AliMedium(12, "CO2",        8, 0, isxfld, sxmgmx
 645                 , tmaxfd, stemax, deemax, epsil, stmin);
 646   // G10-plates
 647   AliMedium(13, "G10-plates",12, 0, isxfld, sxmgmx
 648                 , tmaxfd, stemax, deemax, epsil, stmin);
 649   // Cooling water
 650   AliMedium(14, "Water",     13, 0, isxfld, sxmgmx
 651                 , tmaxfd, stemax, deemax, epsil, stmin);
 652   // Rohacell (plexiglas) for the radiator
 653   AliMedium(15, "Rohacell",  14, 0, isxfld, sxmgmx
 654                 , tmaxfd, stemax, deemax, epsil, stmin);
 655   // Al layer in MCMs
 656   AliMedium(16, "MCM-Al"  ,   1, 0, isxfld, sxmgmx
 657                 , tmaxfd, stemax, deemax, epsil, stmin);
 658   // Sn layer in MCMs
 659   AliMedium(17, "MCM-Sn"  ,  15, 0, isxfld, sxmgmx
 660                 , tmaxfd, stemax, deemax, epsil, stmin);
 661   // Cu layer in MCMs
 662   AliMedium(18, "MCM-Cu"  ,   5, 0, isxfld, sxmgmx
 663                 , tmaxfd, stemax, deemax, epsil, stmin);
 664   // G10 layer in MCMs
 665   AliMedium(19, "MCM-G10" ,  12, 0, isxfld, sxmgmx
 666                 , tmaxfd, stemax, deemax, epsil, stmin);
 667   // Si in readout chips
 668   AliMedium(20, "Chip-Si" ,  16, 0, isxfld, sxmgmx
 669                 , tmaxfd, stemax, deemax, epsil, stmin);
 670   // Epoxy in readout chips
 671   AliMedium(21, "Chip-Ep" ,  17, 0, isxfld, sxmgmx
 672                 , tmaxfd, stemax, deemax, epsil, stmin);
 673   // PE in connectors
 674   AliMedium(22, "Conn-PE" ,   3, 0, isxfld, sxmgmx
 675                 , tmaxfd, stemax, deemax, epsil, stmin);
 676   // Cu in connectors
 677   AliMedium(23, "Chip-Cu" ,   5, 0, isxfld, sxmgmx
 678                 , tmaxfd, stemax, deemax, epsil, stmin);
 679   // Al of cooling pipes
 680   AliMedium(24, "Cooling" ,   1, 0, isxfld, sxmgmx
 681                 , tmaxfd, stemax, deemax, epsil, stmin);
 682   // Cu in services
 683   AliMedium(25, "Serv-Cu" ,   5, 0, isxfld, sxmgmx
 684                 , tmaxfd, stemax, deemax, epsil, stmin);
 685
 686   // Save the density values for the TRD absorbtion
 687   fFoilDensity = dmy;
 688   if (fGasMix == 1)
 689     fGasDensity = dgm1;
 690   else
 691     fGasDensity = dgm2;
 692
 693 }
 694
 695 //_____________________________________________________________________________
 696 void AliTRD::DrawModule() const
 697 {
 698   //
 699   // Draw a shaded view of the Transition Radiation Detector version 0
 700   //
 701
 702   // Set everything unseen
 703   gMC->Gsatt("*"   ,"SEEN",-1);
 704
 705   // Set ALIC mother transparent
 706   gMC->Gsatt("ALIC","SEEN", 0);
 707
 708   // Set the volumes visible
 709   if (fGeometry->IsVersion() == 0) {
 710     gMC->Gsatt("B071","SEEN", 0);
 711     gMC->Gsatt("B074","SEEN", 0);
 712     gMC->Gsatt("B075","SEEN", 0);
 713     gMC->Gsatt("B077","SEEN", 0);
 714     gMC->Gsatt("BTR1","SEEN", 0);
 715     gMC->Gsatt("BTR2","SEEN", 0);
 716     gMC->Gsatt("BTR3","SEEN", 0);
 717     gMC->Gsatt("UTR1","SEEN", 0);
 718     gMC->Gsatt("UTR2","SEEN", 0);
 719     gMC->Gsatt("UTR3","SEEN", 0);
 720   }
 721   else {
 722     gMC->Gsatt("B071","SEEN", 0);
 723     gMC->Gsatt("B074","SEEN", 0);
 724     gMC->Gsatt("B075","SEEN", 0);
 725     gMC->Gsatt("B077","SEEN", 0);
 726     gMC->Gsatt("BTR1","SEEN", 0);
 727     gMC->Gsatt("BTR2","SEEN", 0);
 728     gMC->Gsatt("BTR3","SEEN", 0);
 729     gMC->Gsatt("UTR1","SEEN", 0);
 730     if (fGeometry->GetPHOShole())
 731       gMC->Gsatt("UTR2","SEEN", 0);
 732     if (fGeometry->GetRICHhole())
 733       gMC->Gsatt("UTR3","SEEN", 0);
 734   }
 735 //   gMC->Gsatt("UCII","SEEN", 0);
 736 //   gMC->Gsatt("UCIM","SEEN", 0);
 737 //   gMC->Gsatt("UCIO","SEEN", 0);
 738 //   gMC->Gsatt("UL02","SEEN", 1);
 739 //   gMC->Gsatt("UL05","SEEN", 1);
 740 //   gMC->Gsatt("UL06","SEEN", 1);
 741
 742   gMC->Gdopt("hide", "on");
 743   gMC->Gdopt("shad", "on");
 744   gMC->Gsatt("*", "fill", 7);
 745   gMC->SetClipBox(".");
 746   gMC->SetClipBox("*", 0, 2000, -2000, 2000, -2000, 2000);
 747   gMC->DefaultRange();
 748   gMC->Gdraw("alic", 40, 30, 0, 12, 9.4, .021, .021);
 749   gMC->Gdhead(1111, "Transition Radiation Detector");
 750   gMC->Gdman(18, 4, "MAN");
 751
 752 }
 753
 754 //_____________________________________________________________________________
 755 Int_t AliTRD::DistancetoPrimitive(Int_t , Int_t )
 756 {
 757   //
 758   // Distance between the mouse and the TRD detector on the screen
 759   // Dummy routine
 760
 761   return 9999;
 762
 763 }
 764
 765 //_____________________________________________________________________________
 766 void AliTRD::Init()
 767 {
 768   //
 769   // Initialize the TRD detector after the geometry has been created
 770   //
 771
 772   Int_t i;
 773
 774   if (fDebug) {
 775     printf("\n%s: ",ClassName());
 776     for (i = 0; i < 35; i++) printf("*");
 777     printf(" TRD_INIT ");
 778     for (i = 0; i < 35; i++) printf("*");
 779     printf("\n");
 780   }
 781
 782   if (fGeometry->IsVersion() == 1) {
 783     printf("%s: Full geometry version initialized\n",ClassName());
 784     if (fGeometry->GetPHOShole())
 785       printf("%s: Leave space in front of PHOS free\n",ClassName());
 786     if (fGeometry->GetRICHhole())
 787       printf("%s: Leave space in front of RICH free\n",ClassName());
 788   }
 789   else {
 790     printf("%s: Not a valid geometry\n",ClassName());
 791   }
 792
 793   if (fGasMix == 1) {
 794     printf("%s: Gas Mixture: 85%% Xe + 15%% CO2\n",ClassName());
 795   }
 796   else {
 797     printf("%s: Gas Mixture: 97%% Xe + 3%% Isobutane\n",ClassName());
 798   }
 799
 800 }
 801
 802 //_____________________________________________________________________________
 803 void AliTRD::LoadPoints(Int_t )
 804 {
 805   //
 806   // Store x, y, z of all hits in memory.
 807   // Hit originating from TR photons are given a different color
 808   //
 809
 810   //if (!fDrawTR) {
 811   //  AliDetector::LoadPoints(track);
 812   //  return;
 813   //}
 814
 815   if (fHits == 0) return;
 816
 817   Int_t nhits = fHits->GetEntriesFast();
 818   if (nhits == 0) return;
 819
 820   Int_t tracks = gAlice->GetMCApp()->GetNtrack();
 821   if (fPoints == 0) fPoints = new TObjArray(tracks);
 822
 823   AliTRDhit *ahit;
 824
 825   Int_t    *ntrkE = new Int_t[tracks];
 826   Int_t    *ntrkT = new Int_t[tracks];
 827   Int_t    *limiE = new Int_t[tracks];
 828   Int_t    *limiT = new Int_t[tracks];
 829   Float_t **coorE = new Float_t*[tracks];
 830   Float_t **coorT = new Float_t*[tracks];
 831   for(Int_t i = 0; i < tracks; i++) {
 832     ntrkE[i] = 0;
 833     ntrkT[i] = 0;
 834     coorE[i] = 0;
 835     coorT[i] = 0;
 836     limiE[i] = 0;
 837     limiT[i] = 0;
 838   }
 839
 840   AliTRDpoints *points = 0;
 841   Float_t      *fp     = 0;
 842   Int_t         trk;
 843   Int_t         chunk  = nhits / 4 + 1;
 844
 845   // Loop over all the hits and store their position
 846   ahit = (AliTRDhit *) FirstHit(-1);
 847   while (ahit) {
 848
 849     // dEdx hits
 850     if (ahit->GetCharge() >= 0) {
 851
 852       trk = ahit->GetTrack();
 853       if (ntrkE[trk] == limiE[trk]) {
 854         // Initialise a new track
 855         fp = new Float_t[3*(limiE[trk]+chunk)];
 856         if (coorE[trk]) {
 857           memcpy(fp,coorE[trk],sizeof(Float_t)*3*limiE[trk]);
 858           delete [] coorE[trk];
 859         }
 860         limiE[trk] += chunk;
 861         coorE[trk]  = fp;
 862       }
 863       else {
 864         fp = coorE[trk];
 865       }
 866       fp[3*ntrkE[trk]  ] = ahit->X();
 867       fp[3*ntrkE[trk]+1] = ahit->Y();
 868       fp[3*ntrkE[trk]+2] = ahit->Z();
 869       ntrkE[trk]++;
 870
 871     }
 872     // TR photon hits
 873     else if ((ahit->GetCharge() < 0) && (fDrawTR)) {
 874
 875       trk = ahit->GetTrack();
 876       if (ntrkT[trk] == limiT[trk]) {
 877         // Initialise a new track
 878         fp = new Float_t[3*(limiT[trk]+chunk)];
 879         if (coorT[trk]) {
 880           memcpy(fp,coorT[trk],sizeof(Float_t)*3*limiT[trk]);
 881           delete [] coorT[trk];
 882         }
 883         limiT[trk] += chunk;
 884         coorT[trk]  = fp;
 885       }
 886       else {
 887         fp = coorT[trk];
 888       }
 889       fp[3*ntrkT[trk]  ] = ahit->X();
 890       fp[3*ntrkT[trk]+1] = ahit->Y();
 891       fp[3*ntrkT[trk]+2] = ahit->Z();
 892       ntrkT[trk]++;
 893
 894     }
 895
 896     ahit = (AliTRDhit *) NextHit();
 897
 898   }
 899
 900   for (trk = 0; trk < tracks; ++trk) {
 901
 902     if (ntrkE[trk] || ntrkT[trk]) {
 903
 904       points = new AliTRDpoints();
 905       points->SetDetector(this);
 906       points->SetParticle(trk);
 907
 908       // Set the dEdx points
 909       if (ntrkE[trk]) {
 910         points->SetMarkerColor(GetMarkerColor());
 911         points->SetMarkerSize(GetMarkerSize());
 912         points->SetPolyMarker(ntrkE[trk],coorE[trk],GetMarkerStyle());
 913         delete [] coorE[trk];
 914         coorE[trk] = 0;
 915       }
 916
 917       // Set the TR photon points
 918       if (ntrkT[trk]) {
 919         points->SetTRpoints(ntrkT[trk],coorT[trk]);
 920         delete [] coorT[trk];
 921         coorT[trk] = 0;
 922       }
 923
 924       fPoints->AddAt(points,trk);
 925
 926     }
 927
 928   }
 929
 930   delete [] coorE;
 931   delete [] coorT;
 932   delete [] ntrkE;
 933   delete [] ntrkT;
 934   delete [] limiE;
 935   delete [] limiT;
 936
 937 }
 938
 939 //_____________________________________________________________________________
 940 void AliTRD::MakeBranch(Option_t* option)
 941 {
 942   //
 943   // Create Tree branches for the TRD digits.
 944   //
 945
 946   Int_t  buffersize = 4000;
 947   Char_t branchname[15];
 948   sprintf(branchname,"%s",GetName());
 949
 950   const char *cD = strstr(option,"D");
 951
 952   AliDetector::MakeBranch(option);
 953
 954   if (fDigits && gAlice->TreeD() && cD) {
 955     MakeBranchInTree(gAlice->TreeD(),branchname,&fDigits,buffersize,0);
 956   }
 957
 958 }
 959
 960 //_____________________________________________________________________________
 961 void AliTRD::ResetDigits()
 962 {
 963   //
 964   // Reset number of digits and the digits array for this detector
 965   //
 966
 967   fNdigits = 0;
 968   if (fDigits) fDigits->Clear();
 969
 970 }
 971
 972 //_____________________________________________________________________________
 973 void AliTRD::SetTreeAddress()
 974 {
 975   //
 976   // Set the branch addresses for the trees.
 977   //
 978
 979   if ( fLoader->TreeH() && (fHits == 0x0)) {
 980     fHits = new TClonesArray("AliTRDhit",405);
 981   }
 982   AliDetector::SetTreeAddress();
 983
 984 }
 985
 986 //_____________________________________________________________________________
 987 void AliTRD::SetGasMix(Int_t imix)
 988 {
 989   //
 990   // Defines the gas mixture (imix=0:  Xe/Isobutane imix=1: Xe/CO2)
 991   //
 992
 993   if ((imix < 0) || (imix > 1)) {
 994     printf("Wrong input value: %d\n",imix);
 995     printf("Use standard setting\n");
 996     fGasMix = 1;
 997     return;
 998   }
 999
1000   fGasMix = imix;
1001
1002 }
1003
1004 //_____________________________________________________________________________
1005 void AliTRD::SetPHOShole()
1006 {
1007   //
1008   // Selects a geometry with a hole in front of the PHOS
1009   //
1010
1011   fGeometry->SetPHOShole();
1012
1013 }
1014
1015 //_____________________________________________________________________________
1016 void AliTRD::SetRICHhole()
1017 {
1018   //
1019   // Selects a geometry with a hole in front of the RICH
1020   //
1021
1022   fGeometry->SetRICHhole();
1023
1024 }
1025
1026 //_____________________________________________________________________________
1027 AliTRD &AliTRD::operator=(const AliTRD &trd)
1028 {
1029   //
1030   // Assignment operator
1031   //
1032
1033   if (this != &trd) ((AliTRD &) trd).Copy(*this);
1034   return *this;
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
1273
1274
1275
1276
1277
1278
1279
1280
1281
1282
1283
1284
1285
1286
1287
1288
1289
1290
1291
1292
1293
1294
1295
1296
1297
1298
1299
1300
1301
1302
1303
1304
1305
1306
1307
1308
1309
1310
1311
1312
1313
1314
1315
1316
1317
1318
1319
1320
1321
1322
1323
1324
1325
1326
1327
1328
1329
1330
1331
1332
1333
1334
1335
1336
1337
1338
1339
1340
1341
1342
1343
1344
1345
1346
1347
1348
1349
1350
1351
1352
1353
1354
1355
1356
1357
1358
1359
1360
1361
1362
1363
1364
1365
1366
1367
1368
1369
1370
1371