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