FMD/AliFMD.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 /* $Id$ */
  16 /** @file    AliFMD.cxx
  17     @author  Christian Holm Christensen <cholm@nbi.dk>
  18     @date    Sun Mar 26 17:59:18 2006
  19     @brief   Implementation of AliFMD base class
  20 */
  21 //____________________________________________________________________
  22 //
  23 // Forward Multiplicity Detector based on Silicon wafers. This class
  24 // is the driver for especially simulation.
  25 //
  26 // The Forward Multiplicity Detector consists of 3 sub-detectors FMD1,
  27 // FMD2, and FMD3, each of which has 1 or 2 rings of silicon sensors.
  28 //
  29 // This is the base class for all FMD manager classes.
  30 //
  31 // The actual code is done by various separate classes.   Below is
  32 // diagram showing the relationship between the various FMD classes
  33 // that handles the simulation
  34 //
  35 //
  36 //       +----------+   +----------+
  37 //       | AliFMDv1 |   | AliFMDv0 |
  38 //       +----------+   +----------+
  39 //            |              |                    +-----------------+
  40 //       +----+--------------+                 +--| AliFMDDigitizer |
  41 //       |                                     |  +-----------------+
  42 //       |           +---------------------+   |
  43 //       |        +--| AliFMDBaseDigitizer |<--+
  44 //       V     1  |  +---------------------+   |
  45 //  +--------+<>--+                            |  +------------------+
  46 //  | AliFMD |                                 +--| AliFMDSDigitizer |
  47 //  +--------+<>--+                               +------------------+
  48 //             1  |  +---------------------+
  49 //                +--| AliFMDReconstructor |
  50 //                   +---------------------+
  51 //
  52 // *  AliFMD
  53 //    This defines the interface for the various parts of AliROOT that
  54 //    uses the FMD, like AliFMDSimulator, AliFMDDigitizer,
  55 //    AliFMDReconstructor, and so on.
  56 //
  57 // *  AliFMDv0
  58 //    This is a concrete implementation of the AliFMD interface.
  59 //    It is the responsibility of this class to create the FMD
  60 //    geometry.
  61 //
  62 // *  AliFMDv1
  63 //    This is a concrete implementation of the AliFMD interface.
  64 //    It is the responsibility of this class to create the FMD
  65 //    geometry, process hits in the FMD, and serve hits and digits to
  66 //    the various clients.
  67 //
  68 // *  AliFMDSimulator
  69 //    This is the base class for the FMD simulation tasks.   The
  70 //    simulator tasks are responsible to implment the geoemtry, and
  71 //    process hits.
  72 //
  73 // *  AliFMDReconstructor
  74 //    This is a concrete implementation of the AliReconstructor that
  75 //    reconstructs pseudo-inclusive-multiplicities from digits (raw or
  76 //    from simulation)
  77 //
  78 // Calibration and geometry parameters are managed by separate
  79 // singleton managers.  These are AliFMDGeometry and
  80 // AliFMDParameters.  Please refer to these classes for more
  81 // information on these.
  82 //
  83
  84 // These files are not in the same directory, so there's no reason to
  85 // ask the preprocessor to search in the current directory for these
  86 // files by including them with `#include "..."'
  87 #include <cmath>                // __CMATH__
  88 #include <TClonesArray.h>       // ROOT_TClonesArray
  89 #include <TGeometry.h>          // ROOT_TGeomtry
  90 #include <TNode.h>              // ROOT_TNode
  91 #include <TXTRU.h>              // ROOT_TXTRU
  92 #include <TRotMatrix.h>         // ROOT_TRotMatrix
  93 #include <TTUBE.h>              // ROOT_TTUBE
  94 #include <TTree.h>              // ROOT_TTree
  95 #include <TBrowser.h>           // ROOT_TBrowser
  96 #include <TVirtualMC.h>         // ROOT_TVirtualMC
  97 #include <TVector2.h>           // ROOT_TVector2
  98 #include <TGeoManager.h>        // ROOT_TGeoManager
  99
 100 #include <AliRunDigitizer.h>    // ALIRUNDIGITIZER_H
 101 #include <AliLoader.h>          // ALILOADER_H
 102 #include <AliRun.h>             // ALIRUN_H
 103 #include <AliMC.h>              // ALIMC_H
 104 #include <AliMagF.h>            // ALIMAGF_H
 105 // #include <AliLog.h>          // ALILOG_H
 106 #include "AliFMDDebug.h" // Better debug macros
 107 #include "AliFMD.h"             // ALIFMD_H
 108 #include "AliFMDDigit.h"        // ALIFMDDIGIT_H
 109 #include "AliFMDSDigit.h"       // ALIFMDSDIGIT_H
 110 #include "AliFMDHit.h"          // ALIFMDHIT_H
 111 #include "AliFMDGeometry.h"     // ALIFMDGEOMETRY_H
 112 #include "AliFMDDetector.h"     // ALIFMDDETECTOR_H
 113 #include "AliFMDRing.h"         // ALIFMDRING_H
 114 #include "AliFMDDigitizer.h"    // ALIFMDDIGITIZER_H
 115 #include "AliFMDHitDigitizer.h" // ALIFMDSDIGITIZER_H
 116 // #define USE_SSDIGITIZER
 117 //#ifdef USE_SSDIGITIZER
 118 //# include "AliFMDSSDigitizer.h"       // ALIFMDSDIGITIZER_H
 119 //#endif
 120 // #include "AliFMDGeometryBuilder.h"
 121 #include "AliFMDRawWriter.h"    // ALIFMDRAWWRITER_H
 122 #include "AliFMDPoints.h"       // ALIFMDPOINTS_H
 123
 124 //____________________________________________________________________
 125 ClassImp(AliFMD)
 126 #if 0
 127   ; // This is to keep Emacs from indenting the next line
 128 #endif
 129
 130 //____________________________________________________________________
 131 AliFMD::AliFMD()
 132   : AliDetector(),
 133     fSDigits(0),
 134     fNsdigits(0),
 135     fDetailed(kTRUE),
 136     fUseOld(kFALSE),
 137     fUseAssembly(kTRUE),
 138     fBad(0)
 139 {
 140   //
 141   // Default constructor for class AliFMD
 142   //
 143   AliFMDDebug(10, ("\tDefault CTOR"));
 144   fHits        = 0;
 145   fDigits      = 0;
 146   fIshunt      = 0;
 147   fBad         = new TClonesArray("AliFMDHit");
 148 }
 149
 150 //____________________________________________________________________
 151 AliFMD::AliFMD(const char *name, const char *title)
 152   : AliDetector (name, title),
 153     fSDigits(0),
 154     fNsdigits(0),
 155     fDetailed(kTRUE),
 156     fUseOld(kFALSE),
 157     fUseAssembly(kFALSE),
 158     fBad(0)
 159 {
 160   //
 161   // Standard constructor for Forward Multiplicity Detector
 162   //
 163   AliFMDDebug(10, ("\tStandard CTOR"));
 164   fBad         = new TClonesArray("AliFMDHit");
 165
 166   // Initialise Hit array
 167   HitsArray();
 168   gAlice->GetMCApp()->AddHitList(fHits);
 169
 170   // (S)Digits for the detectors disk
 171   DigitsArray();
 172   SDigitsArray();
 173
 174   // CHC: What is this?
 175   fIshunt = 0;
 176   //PH  SetMarkerColor(kRed);
 177   //PH  SetLineColor(kYellow);
 178 }
 179
 180 //____________________________________________________________________
 181 AliFMD::~AliFMD ()
 182 {
 183   // Destructor for base class AliFMD
 184   if (fHits) {
 185     fHits->Delete();
 186     delete fHits;
 187     fHits = 0;
 188   }
 189   if (fDigits) {
 190     fDigits->Delete();
 191     delete fDigits;
 192     fDigits = 0;
 193   }
 194   if (fSDigits) {
 195     fSDigits->Delete();
 196     delete fSDigits;
 197     fSDigits = 0;
 198   }
 199   if (fBad) {
 200     fBad->Delete();
 201     delete fBad;
 202     fBad = 0;
 203   }
 204 }
 205
 206
 207 //====================================================================
 208 //
 209 // GEometry ANd Traking
 210 //
 211 //____________________________________________________________________
 212 void
 213 AliFMD::CreateGeometry()
 214 {
 215   //
 216   // Create the geometry of Forward Multiplicity Detector.  The actual
 217   // construction of the geometry is delegated to the class
 218   // AliFMDGeometryBuilder, invoked by the singleton manager
 219   // AliFMDGeometry.
 220   //
 221   AliFMDGeometry*  fmd = AliFMDGeometry::Instance();
 222   fmd->SetDetailed(fDetailed);
 223   fmd->UseAssembly(fUseAssembly);
 224   fmd->Build();
 225 }
 226
 227 //____________________________________________________________________
 228 void AliFMD::CreateMaterials()
 229 {
 230   // Define the materials and tracking mediums needed by the FMD
 231   // simulation.   These mediums are made by sending the messages
 232   // AliMaterial, AliMixture, and AliMedium to the passed AliModule
 233   // object module.   The defined mediums are
 234   //
 235   //    FMD Si$         Silicon (active medium in sensors)
 236   //    FMD C$          Carbon fibre (support cone for FMD3 and vacuum pipe)
 237   //    FMD Al$         Aluminium (honeycomb support plates)
 238   //    FMD PCB$        Printed Circuit Board (FEE board with VA1_3)
 239   //    FMD Chip$       Electronics chips (currently not used)
 240   //    FMD Air$        Air (Air in the FMD)
 241   //    FMD Plastic$    Plastic (Support legs for the hybrid cards)
 242   //
 243   // The geometry builder should really be the one that creates the
 244   // materials, but the architecture of AliROOT makes that design
 245   // akward.  What should happen, was that the AliFMDGeometryBuilder
 246   // made the mediums, and that this class retrives pointers from the
 247   // TGeoManager, and registers the mediums here.  Alas, it's not
 248   // really that easy.
 249   //
 250   AliFMDDebug(10, ("\tCreating materials"));
 251   // Get pointer to geometry singleton object.
 252   AliFMDGeometry* geometry = AliFMDGeometry::Instance();
 253   geometry->Init();
 254 #if 0
 255   if (gGeoManager && gGeoManager->GetMedium("FMD Si$")) {
 256     // We need to figure out the some stuff about the geometry
 257     fmd->ExtractGeomInfo();
 258     return;
 259   }
 260 #endif
 261   Int_t    id;
 262   Double_t a                = 0;
 263   Double_t z                = 0;
 264   Double_t density          = 0;
 265   Double_t radiationLength  = 0;
 266   Double_t absorbtionLength = 999;
 267   Int_t    fieldType        = gAlice->Field()->Integ();     // Field type
 268   Double_t maxField         = gAlice->Field()->Max();     // Field max.
 269   Double_t maxBending       = 0;     // Max Angle
 270   Double_t maxStepSize      = 0.001; // Max step size
 271   Double_t maxEnergyLoss    = 1;     // Max Delta E
 272   Double_t precision        = 0.001; // Precision
 273   Double_t minStepSize      = 0.001; // Minimum step size
 274
 275   // Silicon
 276   a                = 28.0855;
 277   z                = 14.;
 278   density          = geometry->GetSiDensity();
 279   radiationLength  = 9.36;
 280   maxBending       = 1;
 281   maxStepSize      = .001;
 282   precision        = .001;
 283   minStepSize      = .001;
 284   id               = kSiId;
 285   AliMaterial(id, "Si$", a, z, density, radiationLength, absorbtionLength);
 286   AliMedium(kSiId, "Si$", id,1,fieldType,maxField,maxBending,
 287             maxStepSize,maxEnergyLoss,precision,minStepSize);
 288
 289
 290   // Carbon
 291   a                = 12.011;
 292   z                = 6.;
 293   density          = 2.265;
 294   radiationLength  = 18.8;
 295   maxBending       = 10;
 296   maxStepSize      = .01;
 297   precision        = .003;
 298   minStepSize      = .003;
 299   id               = kCarbonId;
 300   AliMaterial(id, "Carbon$", a, z, density, radiationLength, absorbtionLength);
 301   AliMedium(kCarbonId, "Carbon$", id,0,fieldType,maxField,maxBending,
 302                     maxStepSize,maxEnergyLoss,precision,minStepSize);
 303
 304   // Aluminum
 305   a                = 26.981539;
 306   z                = 13.;
 307   density          = 2.7;
 308   radiationLength  = 8.9;
 309   id               = kAlId;
 310   AliMaterial(id, "Aluminum$",a,z, density, radiationLength, absorbtionLength);
 311   AliMedium(kAlId, "Aluminum$", id, 0, fieldType, maxField, maxBending,
 312             maxStepSize, maxEnergyLoss, precision, minStepSize);
 313
 314
 315   // Copper
 316   a                = 63.546;
 317   z                = 29;
 318   density          =  8.96;
 319   radiationLength  =  1.43;
 320   id               = kCopperId;
 321   AliMaterial(id, "Copper$",
 322                       a, z, density, radiationLength, absorbtionLength);
 323   AliMedium(kCopperId, "Copper$", id, 0, fieldType, maxField, maxBending,
 324             maxStepSize, maxEnergyLoss, precision, minStepSize);
 325
 326
 327   // Silicon chip
 328   {
 329     Float_t as[] = { 12.0107,      14.0067,      15.9994,
 330                       1.00794,     28.0855,     107.8682 };
 331     Float_t zs[] = {  6.,           7.,           8.,
 332                       1.,          14.,          47. };
 333     Float_t ws[] = {  0.039730642,  0.001396798,  0.01169634,
 334                       0.004367771,  0.844665,     0.09814344903 };
 335     density          = 2.36436;
 336     maxBending       = 10;
 337     maxStepSize      = .01;
 338     precision        = .003;
 339     minStepSize      = .003;
 340     id               = kSiChipId;
 341     AliMixture(id, "Si Chip$", as, zs, density, 6, ws);
 342     AliMedium(kSiChipId, "Si Chip$",  id, 0, fieldType, maxField, maxBending,
 343               maxStepSize, maxEnergyLoss, precision, minStepSize);
 344   }
 345
 346   // Kaption
 347   {
 348     Float_t as[] = { 1.00794,  12.0107,  14.010,   15.9994};
 349     Float_t zs[] = { 1.,        6.,       7.,       8.};
 350     Float_t ws[] = { 0.026362,  0.69113,  0.07327,  0.209235};
 351     density          = 1.42;
 352     maxBending       = 1;
 353     maxStepSize      = .001;
 354     precision        = .001;
 355     minStepSize      = .001;
 356     id               = kKaptonId;
 357     AliMixture(id, "Kaption$", as, zs, density, 4, ws);
 358     AliMedium(kKaptonId, "Kaption$", id,0,fieldType,maxField,maxBending,
 359               maxStepSize,maxEnergyLoss,precision,minStepSize);
 360   }
 361
 362   // Air
 363   {
 364     Float_t as[] = { 12.0107, 14.0067,   15.9994,  39.948 };
 365     Float_t zs[] = {  6.,      7.,       8.,       18. };
 366     Float_t ws[] = { 0.000124, 0.755267, 0.231781, 0.012827 };
 367     density      = .00120479;
 368     maxBending   = 1;
 369     maxStepSize  = .001;
 370     precision    = .001;
 371     minStepSize  = .001;
 372     id           = kAirId;
 373     AliMixture(id, "Air$", as, zs, density, 4, ws);
 374     AliMedium(kAirId, "Air$", id,0,fieldType,maxField,maxBending,
 375               maxStepSize,maxEnergyLoss,precision,minStepSize);
 376   }
 377
 378   // PCB
 379   {
 380     Float_t zs[] = { 14.,         20.,         13.,         12.,
 381                       5.,         22.,         11.,         19.,
 382                      26.,          9.,          8.,          6.,
 383                       7.,          1.};
 384     Float_t as[] = { 28.0855,     40.078,      26.981538,   24.305,
 385                      10.811,      47.867,      22.98977,    39.0983,
 386                      55.845,      18.9984,     15.9994,     12.0107,
 387                      14.0067,      1.00794};
 388     Float_t ws[] = {  0.15144894,  0.08147477,  0.04128158,  0.00904554,
 389                       0.01397570,  0.00287685,  0.00445114,  0.00498089,
 390                       0.00209828,  0.00420000,  0.36043788,  0.27529426,
 391                       0.01415852,  0.03427566};
 392     density      = 1.8;
 393     maxBending   = 1;
 394     maxStepSize  = .001;
 395     precision    = .001;
 396     minStepSize  = .001;
 397     id           = kPcbId;
 398     AliMixture(id, "PCB$", as, zs, density, 14, ws);
 399     AliMedium(kPcbId, "PCB$", id,0,fieldType,maxField,maxBending,
 400               maxStepSize,maxEnergyLoss,precision,minStepSize);
 401   }
 402
 403   // Stainless steel
 404   {
 405     Float_t as[] = { 55.847, 51.9961, 58.6934, 28.0855 };
 406     Float_t zs[] = { 26.,    24.,     28.,     14.     };
 407     Float_t ws[] = { .715,   .18,     .1,      .005    };
 408     density      = 7.88;
 409     id           = kSteelId;
 410     AliMixture(id, "Steel$", as, zs, density, 4, ws);
 411     AliMedium(kSteelId, "Steel$", id, 0, fieldType, maxField, maxBending,
 412               maxStepSize, maxEnergyLoss, precision, minStepSize);
 413   }
 414   // Plastic
 415   {
 416     Float_t as[] = { 1.01, 12.01 };
 417     Float_t zs[] = { 1.,   6.    };
 418     Float_t ws[] = { 1.,   1.    };
 419     density      = 1.03;
 420     maxBending   = 10;
 421     maxStepSize  = .01;
 422     precision    = .003;
 423     minStepSize  = .003;
 424     id           = kPlasticId;
 425     AliMixture(id, "Plastic$", as, zs, density, -2, ws);
 426     AliMedium(kPlasticId, "Plastic$", id,0,fieldType,maxField,maxBending,
 427               maxStepSize,maxEnergyLoss,precision,minStepSize);
 428   }
 429
 430 }
 431
 432 //____________________________________________________________________
 433 void
 434 AliFMD::Init()
 435 {
 436   // Initialize the detector
 437   //
 438   AliFMDDebug(1, ("Initialising FMD detector object"));
 439   TVirtualMC*      mc     = TVirtualMC::GetMC();
 440   AliFMDGeometry*  fmd    = AliFMDGeometry::Instance();
 441   const TArrayI&   actGeo = fmd->ActiveIds();
 442   TArrayI          actVmc(actGeo.fN);
 443   for (Int_t i = 0; i < actGeo.fN; i++) {
 444     TGeoVolume *sens = gGeoManager->GetVolume(actGeo[i]);
 445     if (!sens) {
 446       AliError(Form("No TGeo volume for sensitive volume ID=%d",actGeo[i]));
 447       continue;
 448     }
 449     actVmc[i] = mc->VolId(sens->GetName());
 450     AliFMDDebug(1, ("Active vol id # %d: %d changed to %d",
 451                     i, actGeo[i], actVmc[i]));
 452   }
 453   fmd->SetActive(actVmc.fArray, actVmc.fN);
 454   // fmd->InitTransformations();
 455 }
 456
 457 //____________________________________________________________________
 458 void
 459 AliFMD::FinishEvent()
 460 {
 461   // Called at the end of the an event in simulations.  If the debug
 462   // level is high enough, then the `bad' hits are printed.
 463   //
 464   if (AliLog::GetDebugLevel("FMD", "AliFMD") < 10) return;
 465   if (fBad && fBad->GetEntries() > 0) {
 466     AliWarning((Form("EndEvent", "got %d 'bad' hits", fBad->GetEntries())));
 467     TIter next(fBad);
 468     AliFMDHit* hit;
 469     while ((hit = static_cast<AliFMDHit*>(next()))) hit->Print("D");
 470     fBad->Clear();
 471   }
 472 }
 473
 474
 475 //====================================================================
 476 //
 477 // Graphics and event display
 478 //
 479 //____________________________________________________________________
 480 void
 481 AliFMD::BuildGeometry()
 482 {
 483   //
 484   // Build simple ROOT TNode geometry for event display. With the new
 485   // geometry modeller, TGeoManager, this seems rather redundant.
 486   AliFMDDebug(10, ("\tCreating a simplified geometry"));
 487
 488   AliFMDGeometry* fmd = AliFMDGeometry::Instance();
 489
 490   static TXTRU*     innerShape = 0;
 491   static TXTRU*     outerShape = 0;
 492   static TObjArray* innerRot   = 0;
 493   static TObjArray* outerRot   = 0;
 494
 495   if (!innerShape || !outerShape) {
 496     // Make the shapes for the modules
 497     for (Int_t i = 0; i < 2; i++) {
 498       AliFMDRing* r = 0;
 499       switch (i) {
 500       case 0: r = fmd->GetRing('I'); break;
 501       case 1: r = fmd->GetRing('O'); break;
 502       }
 503       if (!r) {
 504         AliError(Form("no ring found for i=%d", i));
 505         return;
 506       }
 507       Double_t    siThick  = r->GetSiThickness();
 508       const Int_t knv      = r->GetNVerticies();
 509       Double_t    theta    = r->GetTheta();
 510       Int_t       nmod     = r->GetNModules();
 511
 512       TXTRU* shape = new TXTRU(r->GetName(), r->GetTitle(), "void", knv, 2);
 513       for (Int_t j = 0; j < knv; j++) {
 514         TVector2* vv = r->GetVertex(knv - 1 - j);
 515         shape->DefineVertex(j, vv->X(), vv->Y());
 516       }
 517       shape->DefineSection(0, -siThick / 2, 1, 0, 0);
 518       shape->DefineSection(1, +siThick / 2, 1, 0, 0);
 519       shape->SetLineColor(kYellow); //PH kYellow is the default line color in FMD
 520
 521       TObjArray* rots = new TObjArray(nmod);
 522       for (Int_t j = 0; j < nmod; j++) {
 523         Double_t th = (j + .5) * theta * 2;
 524         TString name(Form("FMD_ring_%c_rot_%02d", r->GetId(), j));
 525         TString title(Form("FMD Ring %c Rotation # %d", r->GetId(), j));
 526         TRotMatrix* rot = new TRotMatrix(name.Data(), title.Data(),
 527                                          90, th, 90, fmod(90+th,360), 0, 0);
 528         rots->AddAt(rot, j);
 529       }
 530
 531       switch (r->GetId()) {
 532       case 'i':
 533       case 'I': innerShape = shape; innerRot = rots; break;
 534       case 'o':
 535       case 'O': outerShape = shape; outerRot = rots; break;
 536       }
 537     }
 538   }
 539
 540   TNode* top = gAlice->GetGeometry()->GetNode("alice");
 541
 542   for (Int_t i = 1; i <= 3; i++) {
 543     AliFMDDetector* det = fmd->GetDetector(i);
 544     if (!det) {
 545       Warning("BuildGeometry", "FMD%d seems to be disabled", i);
 546       continue;
 547     }
 548     Double_t w  = 0;
 549     Double_t rh = det->GetRing('I')->GetHighR();
 550     Char_t   id = 'I';
 551     if (det->GetRing('O')) {
 552       w  = TMath::Abs(det->GetRingZ('O') - det->GetRingZ('I'));
 553       id = (TMath::Abs(det->GetRingZ('O'))
 554             > TMath::Abs(det->GetRingZ('I')) ? 'O' : 'I');
 555       rh = det->GetRing('O')->GetHighR();
 556     }
 557     w += (det->GetRing(id)->GetModuleSpacing() +
 558           det->GetRing(id)->GetSiThickness());
 559     TShape* shape = new TTUBE(det->GetName(), det->GetTitle(), "void",
 560                               det->GetRing('I')->GetLowR(), rh, w / 2);
 561     Double_t z = (det->GetRingZ('I') - w / 2);
 562     if (z > 0) z += det->GetRing(id)->GetModuleSpacing();
 563     top->cd();
 564     TNode* node = new TNode(det->GetName(), det->GetTitle(), shape,
 565                             0, 0, z, 0);
 566     fNodes->Add(node);
 567
 568     for (Int_t j = 0; j < 2; j++) {
 569       AliFMDRing* r      = 0;
 570       TShape*     rshape = 0;
 571       TObjArray*  rots   = 0;
 572       switch (j) {
 573       case 0:
 574         r = det->GetRing('I'); rshape = innerShape; rots = innerRot; break;
 575       case 1:
 576         r = det->GetRing('O'); rshape = outerShape; rots = outerRot; break;
 577       }
 578       if (!r) continue;
 579
 580       Double_t    siThick  = r->GetSiThickness();
 581       Int_t       nmod     = r->GetNModules();
 582       Double_t    modspace = r->GetModuleSpacing();
 583       Double_t    rz       = - (z - det->GetRingZ(r->GetId()));
 584
 585       for (Int_t k = 0; k < nmod; k++) {
 586         node->cd();
 587         Double_t    offz    = (k % 2 == 1 ? modspace : 0);
 588         TRotMatrix* rot     = static_cast<TRotMatrix*>(rots->At(k));
 589         TString name(Form("%s%c_module_%02d", det->GetName(), r->GetId(),k));
 590         TString title(Form("%s%c Module %d", det->GetName(), r->GetId(),k));
 591         TNode* mnod = new TNode(name.Data(), title.Data(), rshape,
 592                                 0, 0, rz - siThick / 2
 593                                 + TMath::Sign(offz,z), rot);
 594         mnod->SetLineColor(kYellow); //PH kYellow is the default line color in FMD
 595         fNodes->Add(mnod);
 596       } // for (Int_t k = 0 ; ...)
 597     } // for (Int_t j = 0 ; ...)
 598   } // for (Int_t i = 1 ; ...)
 599 }
 600
 601 //____________________________________________________________________
 602 void
 603 AliFMD::LoadPoints(Int_t /* track */)
 604 {
 605   // Store x, y, z of all hits in memory for display.
 606   //
 607   // Normally, the hits are drawn using TPolyMarker3D - however, that
 608   // is not very useful for the FMD.  Therefor, this member function
 609   // is overloaded to make TMarker3D, via the class AliFMDPoints.
 610   // AliFMDPoints is a local class.
 611   //
 612   if (!fHits) {
 613     AliError(Form("fHits == 0. Name is %s",GetName()));
 614     return;
 615   }
 616   Int_t nHits = fHits->GetEntriesFast();
 617   if (nHits == 0) {
 618     return;
 619   }
 620   Int_t tracks = gAlice->GetMCApp()->GetNtrack();
 621   if (fPoints == 0) fPoints = new TObjArray(2 * tracks);
 622
 623   // Get geometry
 624   AliFMDGeometry* geom = AliFMDGeometry::Instance();
 625   geom->Init();
 626   geom->InitTransformations();
 627
 628   // Now make markers for each hit
 629   // AliInfo(Form("Drawing %d hits (have %d points) for track %d",
 630   //              nHits, fPoints->GetEntriesFast(), track));
 631   for (Int_t ihit = 0; ihit < nHits; ihit++) {
 632     AliFMDHit* hit = static_cast<AliFMDHit*>(fHits->At(ihit));
 633     if (!hit) continue;
 634     Double_t edep    = hit->Edep();
 635     Double_t m       = hit->M();
 636     Double_t poverm  = (m == 0 ? 0 : hit->P());
 637     Double_t absQ    = TMath::Abs(hit->Q());
 638     Bool_t   bad     = kFALSE;
 639     // This `if' is to debug abnormal energy depositions.  We trigger on
 640     // p/m approx larger than or equal to a MIP, and a large edep - more
 641     // than 1 keV - a MIP is 100 eV.
 642     if (edep > absQ * absQ && poverm > 1) bad = kTRUE;
 643
 644     AliFMDPoints* p1 = new AliFMDPoints(hit, kRed); //PH kRed is the default marker color in FMD
 645     // AliPoints* p1 = new AliPoints();
 646     // p1->SetMarkerColor(GetMarkerColor());
 647     // p1->SetMarkerSize(GetMarkerSize());
 648     // p1->SetPoint(0, hit->X(), hit->Y(), hit->Z());
 649     p1->SetDetector(this);
 650     p1->SetParticle(hit->GetTrack());
 651     fPoints->AddAt(p1, hit->GetTrack());
 652     if (bad) {
 653       p1->SetMarkerColor(4);
 654       // p1->SetMarkerSize(2 * GetMarkerSize());
 655     }
 656
 657     Double_t x, y, z;
 658     geom->Detector2XYZ(hit->Detector(), hit->Ring(), hit->Sector(),
 659                        hit->Strip(), x, y, z);
 660     AliFMDPoints* p = new AliFMDPoints(hit, 3);
 661     // AliPoints* p = new AliPoints();
 662     // p->SetMarkerColor(3);
 663     // p->SetMarkerSize(GetMarkerSize());
 664     // p->SetPoint(0, x, y, z);
 665     p->SetDetector(this);
 666     p->SetParticle(hit->GetTrack());
 667     p->SetXYZ(x, y, z);
 668     p->SetMarkerColor(3);
 669     fPoints->AddAt(p, tracks+hit->GetTrack());
 670     if (bad) {
 671       p->SetMarkerColor(5);
 672       // p->SetMarkerSize(2 * GetMarkerSize());
 673     }
 674     // AliInfo(Form("Adding point at %d", tracks+hit->GetTrack()));
 675   }
 676 }
 677
 678 //____________________________________________________________________
 679 void
 680 AliFMD::DrawDetector()
 681 {
 682   // Draw a shaded view of the Forward multiplicity detector.  This
 683   // isn't really useful anymore.
 684   AliFMDDebug(10, ("\tDraw detector"));
 685 }
 686
 687 //____________________________________________________________________
 688 Int_t
 689 AliFMD::DistancetoPrimitive(Int_t, Int_t)
 690 {
 691   // Calculate the distance from the mouse to the FMD on the screen
 692   // Dummy routine.
 693   //
 694   return 9999;
 695 }
 696
 697 //====================================================================
 698 //
 699 // Hit and Digit managment
 700 //
 701 //____________________________________________________________________
 702 void
 703 AliFMD::MakeBranch(Option_t * option)
 704 {
 705   // Create Tree branches for the FMD.
 706   //
 707   // Options:
 708   //
 709   //    H          Make a branch of TClonesArray of AliFMDHit's
 710   //    D          Make a branch of TClonesArray of AliFMDDigit's
 711   //    S          Make a branch of TClonesArray of AliFMDSDigit's
 712   //
 713   const Int_t kBufferSize = 16000;
 714   TString branchname(GetName());
 715   TString opt(option);
 716
 717   if (opt.Contains("H", TString::kIgnoreCase)) {
 718     HitsArray();
 719     AliDetector::MakeBranch(option);
 720   }
 721   if (opt.Contains("D", TString::kIgnoreCase)) {
 722     DigitsArray();
 723     MakeBranchInTree(fLoader->TreeD(), branchname.Data(),
 724                      &fDigits, kBufferSize, 0);
 725   }
 726   if (opt.Contains("S", TString::kIgnoreCase)) {
 727     SDigitsArray();
 728     MakeBranchInTree(fLoader->TreeS(), branchname.Data(),
 729                      &fSDigits, kBufferSize, 0);
 730   }
 731 }
 732
 733 //____________________________________________________________________
 734 void
 735 AliFMD::SetTreeAddress()
 736 {
 737   // Set branch address for the Hits, Digits, and SDigits Tree.
 738   if (fLoader->TreeH()) HitsArray();
 739   AliDetector::SetTreeAddress();
 740
 741   TTree *treeD = fLoader->TreeD();
 742   if (treeD) {
 743     DigitsArray();
 744     TBranch* branch = treeD->GetBranch ("FMD");
 745     if (branch) branch->SetAddress(&fDigits);
 746   }
 747
 748   TTree *treeS = fLoader->TreeS();
 749   if (treeS) {
 750     SDigitsArray();
 751     TBranch* branch = treeS->GetBranch ("FMD");
 752     if (branch) branch->SetAddress(&fSDigits);
 753   }
 754 }
 755
 756 //____________________________________________________________________
 757 void
 758 AliFMD::SetHitsAddressBranch(TBranch *b)
 759 {
 760   // Set the TClonesArray to read hits into.
 761   b->SetAddress(&fHits);
 762 }
 763 //____________________________________________________________________
 764 void
 765 AliFMD::SetSDigitsAddressBranch(TBranch *b)
 766 {
 767   // Set the TClonesArray to read hits into.
 768   b->SetAddress(&fSDigits);
 769 }
 770
 771 //____________________________________________________________________
 772 void
 773 AliFMD::AddHit(Int_t track, Int_t *vol, Float_t *hits)
 774 {
 775   // Add a hit to the hits tree
 776   //
 777   // The information of the two arrays are decoded as
 778   //
 779   // Parameters
 780   //    track                Track #
 781   //    ivol[0]  [UShort_t ] Detector #
 782   //    ivol[1]  [Char_t   ] Ring ID
 783   //    ivol[2]  [UShort_t ] Sector #
 784   //    ivol[3]  [UShort_t ] Strip #
 785   //    hits[0]  [Float_t  ] Track's X-coordinate at hit
 786   //    hits[1]  [Float_t  ] Track's Y-coordinate at hit
 787   //    hits[3]  [Float_t  ] Track's Z-coordinate at hit
 788   //    hits[4]  [Float_t  ] X-component of track's momentum
 789   //    hits[5]  [Float_t  ] Y-component of track's momentum
 790   //    hits[6]  [Float_t  ] Z-component of track's momentum
 791   //    hits[7]  [Float_t  ] Energy deposited by track
 792   //    hits[8]  [Int_t    ] Track's particle Id #
 793   //    hits[9]  [Float_t  ] Time when the track hit
 794   //
 795   //
 796   AddHitByFields(track,
 797                  UShort_t(vol[0]),  // Detector #
 798                  Char_t(vol[1]),    // Ring ID
 799                  UShort_t(vol[2]),  // Sector #
 800                  UShort_t(vol[3]),  // Strip #
 801                  hits[0],           // X
 802                  hits[1],           // Y
 803                  hits[2],           // Z
 804                  hits[3],           // Px
 805                  hits[4],           // Py
 806                  hits[5],           // Pz
 807                  hits[6],           // Energy loss
 808                  Int_t(hits[7]),    // PDG
 809                  hits[8]);          // Time
 810 }
 811
 812 //____________________________________________________________________
 813 AliFMDHit*
 814 AliFMD::AddHitByFields(Int_t    track,
 815                        UShort_t detector,
 816                        Char_t   ring,
 817                        UShort_t sector,
 818                        UShort_t strip,
 819                        Float_t  x,
 820                        Float_t  y,
 821                        Float_t  z,
 822                        Float_t  px,
 823                        Float_t  py,
 824                        Float_t  pz,
 825                        Float_t  edep,
 826                        Int_t    pdg,
 827                        Float_t  t,
 828                        Float_t  l,
 829                        Bool_t   stop)
 830 {
 831   // Add a hit to the list
 832   //
 833   // Parameters:
 834   //
 835   //    track     Track #
 836   //    detector  Detector # (1, 2, or 3)
 837   //    ring      Ring ID ('I' or 'O')
 838   //    sector    Sector # (For inner/outer rings: 0-19/0-39)
 839   //    strip     Strip # (For inner/outer rings: 0-511/0-255)
 840   //    x         Track's X-coordinate at hit
 841   //    y         Track's Y-coordinate at hit
 842   //    z         Track's Z-coordinate at hit
 843   //    px        X-component of track's momentum
 844   //    py        Y-component of track's momentum
 845   //    pz        Z-component of track's momentum
 846   //    edep      Energy deposited by track
 847   //    pdg       Track's particle Id #
 848   //    t         Time when the track hit
 849   //    l         Track length through the material.
 850   //    stop      Whether track was stopped or disappeared
 851   //
 852   TClonesArray& a = *(HitsArray());
 853   // Search through the list of already registered hits, and see if we
 854   // find a hit with the same parameters.  If we do, then don't create
 855   // a new hit, but rather update the energy deposited in the hit.
 856   // This is done, so that a FLUKA based simulation will get the
 857   // number of hits right, not just the enerrgy deposition.
 858   AliFMDHit* hit = 0;
 859   for (Int_t i = 0; i < fNhits; i++) {
 860     if (!a.At(i)) continue;
 861     hit = static_cast<AliFMDHit*>(a.At(i));
 862     if (hit->Detector() == detector
 863         && hit->Ring() == ring
 864         && hit->Sector() == sector
 865         && hit->Strip() == strip
 866         && hit->Track() == track) {
 867       AliFMDDebug(1, ("already had a hit in FMD%d%c[%2d,%3d] for track # %d,"
 868                        " adding energy (%f) to that hit (%f) -> %f",
 869                        detector, ring, sector, strip, track, edep, hit->Edep(),
 870                        hit->Edep() + edep));
 871       hit->SetEdep(hit->Edep() + edep);
 872       return hit;
 873     }
 874   }
 875   // If hit wasn't already registered, do so know.
 876   hit = new (a[fNhits]) AliFMDHit(fIshunt, track, detector, ring, sector,
 877                                   strip, x, y, z, px, py, pz, edep, pdg, t,
 878                                   l, stop);
 879   fNhits++;
 880   return hit;
 881 }
 882
 883 //____________________________________________________________________
 884 void
 885 AliFMD::AddDigit(Int_t* digits, Int_t*)
 886 {
 887   // Add a digit to the Digit tree
 888   //
 889   // Paramters
 890   //
 891   //    digits[0]  [UShort_t] Detector #
 892   //    digits[1]  [Char_t]   Ring ID
 893   //    digits[2]  [UShort_t] Sector #
 894   //    digits[3]  [UShort_t] Strip #
 895   //    digits[4]  [UShort_t] ADC Count
 896   //    digits[5]  [Short_t]  ADC Count, -1 if not used
 897   //    digits[6]  [Short_t]  ADC Count, -1 if not used
 898   //
 899   AddDigitByFields(UShort_t(digits[0]),  // Detector #
 900                    Char_t(digits[1]),    // Ring ID
 901                    UShort_t(digits[2]),  // Sector #
 902                    UShort_t(digits[3]),  // Strip #
 903                    UShort_t(digits[4]),  // ADC Count1
 904                    Short_t(digits[5]),   // ADC Count2
 905                    Short_t(digits[6]),   // ADC Count3
 906                    Short_t(digits[7]));
 907 }
 908
 909 //____________________________________________________________________
 910 void
 911 AliFMD::AddDigitByFields(UShort_t detector,
 912                          Char_t   ring,
 913                          UShort_t sector,
 914                          UShort_t strip,
 915                          UShort_t count1,
 916                          Short_t  count2,
 917                          Short_t  count3,
 918                          Short_t  count4)
 919 {
 920   // add a real digit - as coming from data
 921   //
 922   // Parameters
 923   //
 924   //    detector  Detector # (1, 2, or 3)
 925   //    ring      Ring ID ('I' or 'O')
 926   //    sector    Sector # (For inner/outer rings: 0-19/0-39)
 927   //    strip     Strip # (For inner/outer rings: 0-511/0-255)
 928   //    count1    ADC count (a 10-bit word)
 929   //    count2    ADC count (a 10-bit word), or -1 if not used
 930   //    count3    ADC count (a 10-bit word), or -1 if not used
 931   TClonesArray& a = *(DigitsArray());
 932
 933   new (a[fNdigits++])
 934     AliFMDDigit(detector, ring, sector, strip, count1, count2, count3, count4);
 935   AliFMDDebug(15, ("Adding digit # %5d/%5d for FMD%d%c[%2d,%3d]=(%d,%d,%d,%d)",
 936                    fNdigits-1, a.GetEntriesFast(),
 937                    detector, ring, sector, strip,
 938                    count1, count2, count3, count4));
 939
 940 }
 941
 942 //____________________________________________________________________
 943 void
 944 AliFMD::AddSDigit(Int_t* digits)
 945 {
 946   // Add a digit to the SDigit tree
 947   //
 948   // Paramters
 949   //
 950   //    digits[0]  [UShort_t] Detector #
 951   //    digits[1]  [Char_t]   Ring ID
 952   //    digits[2]  [UShort_t] Sector #
 953   //    digits[3]  [UShort_t] Strip #
 954   //    digits[4]  [Float_t]  Total energy deposited
 955   //    digits[5]  [UShort_t] ADC Count
 956   //    digits[6]  [Short_t]  ADC Count, -1 if not used
 957   //    digits[7]  [Short_t]  ADC Count, -1 if not used
 958   //
 959   AddSDigitByFields(UShort_t(digits[0]),   // Detector #
 960                     Char_t(digits[1]),     // Ring ID
 961                     UShort_t(digits[2]),   // Sector #
 962                     UShort_t(digits[3]),   // Strip #
 963                     Float_t(digits[4]),    // Edep
 964                     UShort_t(digits[5]),   // ADC Count1
 965                     Short_t(digits[6]),    // ADC Count2
 966                     Short_t(digits[7]),    // ADC Count3
 967                     Short_t(digits[8]),    // ADC Count4
 968                     UShort_t(digits[9]),   // N particles
 969                     UShort_t(digits[10])); // N primaries
 970
 971 }
 972
 973 //____________________________________________________________________
 974 void
 975 AliFMD::AddSDigitByFields(UShort_t detector,
 976                           Char_t   ring,
 977                           UShort_t sector,
 978                           UShort_t strip,
 979                           Float_t  edep,
 980                           UShort_t count1,
 981                           Short_t  count2,
 982                           Short_t  count3,
 983                           Short_t  count4,
 984                           UShort_t ntot,
 985                           UShort_t nprim)
 986 {
 987   // add a summable digit
 988   //
 989   // Parameters
 990   //
 991   //    detector  Detector # (1, 2, or 3)
 992   //    ring      Ring ID ('I' or 'O')
 993   //    sector    Sector # (For inner/outer rings: 0-19/0-39)
 994   //    strip     Strip # (For inner/outer rings: 0-511/0-255)
 995   //    edep      Total energy deposited
 996   //    count1    ADC count (a 10-bit word)
 997   //    count2    ADC count (a 10-bit word), or -1 if not used
 998   //    count3    ADC count (a 10-bit word), or -1 if not used
 999   //
1000   TClonesArray& a = *(SDigitsArray());
1001   // AliFMDDebug(0, ("Adding sdigit # %d", fNsdigits));
1002
1003   new (a[fNsdigits++])
1004     AliFMDSDigit(detector, ring, sector, strip, edep,
1005                  count1, count2, count3, count4, ntot, nprim);
1006 }
1007
1008 //____________________________________________________________________
1009 void
1010 AliFMD::ResetSDigits()
1011 {
1012   // Reset number of digits and the digits array for this detector.
1013   //
1014   fNsdigits   = 0;
1015   if (fSDigits) fSDigits->Clear();
1016 }
1017
1018
1019 //____________________________________________________________________
1020 TClonesArray*
1021 AliFMD::HitsArray()
1022 {
1023   // Initialize hit array if not already, and return pointer to it.
1024   if (!fHits) {
1025     fHits = new TClonesArray("AliFMDHit", 1000);
1026     fNhits = 0;
1027   }
1028   return fHits;
1029 }
1030
1031 //____________________________________________________________________
1032 TClonesArray*
1033 AliFMD::DigitsArray()
1034 {
1035   // Initialize digit array if not already, and return pointer to it.
1036   if (!fDigits) {
1037     fDigits = new TClonesArray("AliFMDDigit", 1000);
1038     fNdigits = 0;
1039   }
1040   return fDigits;
1041 }
1042
1043 //____________________________________________________________________
1044 TClonesArray*
1045 AliFMD::SDigitsArray()
1046 {
1047   // Initialize digit array if not already, and return pointer to it.
1048   if (!fSDigits) {
1049     fSDigits = new TClonesArray("AliFMDSDigit", 1000);
1050     fNsdigits = 0;
1051   }
1052   return fSDigits;
1053 }
1054
1055 //====================================================================
1056 //
1057 // Digitization
1058 //
1059 //____________________________________________________________________
1060 void
1061 AliFMD::Hits2Digits()
1062 {
1063   // Create AliFMDDigit's from AliFMDHit's.  This is done by making a
1064   // AliFMDDigitizer, and executing that code.
1065   //
1066   AliFMDHitDigitizer digitizer(this, AliFMDHitDigitizer::kDigits);
1067   digitizer.Init();
1068   digitizer.Exec("");
1069 }
1070
1071 //____________________________________________________________________
1072 void
1073 AliFMD::Hits2SDigits()
1074 {
1075   // Create AliFMDSDigit's from AliFMDHit's.  This is done by creating
1076   // an AliFMDSDigitizer object, and executing it.
1077   //
1078   AliFMDHitDigitizer digitizer(this, AliFMDHitDigitizer::kSDigits);
1079   digitizer.Init();
1080   digitizer.Exec("");
1081 }
1082
1083
1084 //____________________________________________________________________
1085 AliDigitizer*
1086 AliFMD::CreateDigitizer(AliRunDigitizer* manager) const
1087 {
1088   // Create a digitizer object
1089
1090   /* This is what we probably _should_ do */
1091   AliFMDBaseDigitizer* digitizer = 0;
1092
1093 #ifdef USE_SSDIGITIZER
1094   digitizer = new AliFMDSSDigitizer(manager);
1095 #else
1096   /* This is what we actually do, and will work */
1097 #if 0
1098   AliInfo("SDigit->Digit conversion not really supported, "
1099           "doing Hit->Digit conversion instead");
1100 #endif
1101   digitizer = new AliFMDDigitizer(manager);
1102 #endif
1103   return digitizer;
1104 }
1105
1106 //====================================================================
1107 //
1108 // Raw data simulation
1109 //
1110 //__________________________________________________________________
1111 void
1112 AliFMD::Digits2Raw()
1113 {
1114   // Turn digits into raw data.
1115   //
1116   // This uses the class AliFMDRawWriter to do the job.   Please refer
1117   // to that class for more information.
1118   AliFMDRawWriter writer(this);
1119   writer.Exec();
1120 }
1121
1122
1123 //====================================================================
1124 //
1125 // Utility
1126 //
1127 //__________________________________________________________________
1128 void
1129 AliFMD::Browse(TBrowser* b)
1130 {
1131   // Browse this object.
1132   //
1133   AliFMDDebug(30, ("\tBrowsing the FMD"));
1134   AliDetector::Browse(b);
1135   b->Add(AliFMDGeometry::Instance());
1136 }
1137
1138 //____________________________________________________________________
1139 void
1140 AliFMD::AddAlignableVolumes() const
1141 {
1142   //
1143   // Create entries for alignable volumes associating the symbolic volume
1144   // name with the corresponding volume path. Needs to be syncronized with
1145   // eventual changes in the geometry.
1146   //
1147   // This code was made by Raffaele Grosso <rgrosso@mail.cern.ch>.  I
1148   // (cholm) will probably want to change it.   For one, I think it
1149   // should be the job of the geometry manager to deal with this.
1150   AliInfo("Add FMD alignable volumes");
1151   AliFMDGeometry::Instance()->SetAlignableVolumes();
1152 #if 0
1153   for(size_t f = 1; f <= 3; f++){ // Detector 1,2,3
1154     for(size_t tb =  0; tb <2 ; tb++){ // Top/Bottom
1155       char     stb = tb == 0 ? 'T' : 'B';
1156       unsigned min = tb == 0 ? 0   : 5;
1157
1158       TString halfVol(Form("/ALIC_1/F%dM%c_%d", f, stb, f));
1159       TString halfSym(halfVol);
1160       if(!gGeoManager->SetAlignableEntry(halfSym.Data(),halfVol.Data()))
1161         AliFatal(Form("Alignable entry %s not created. "
1162                       "Volume path %s not valid",
1163                       halfSym.Data(),halfVol.Data()));
1164       for(size_t io = 0; io < 2; io++){ // inner, outer
1165         if (f==1 && io==1) continue; // Only one ring in FMD1
1166         if(tb == 1 && io==1) min=10;
1167         char     sio = (io == 0 ? 'I' : 'O');
1168         unsigned nio = (io == 0 ? 3   : 9);
1169         unsigned max = (io == 0 ? 5   : 10) + min;
1170
1171         for(size_t i = min; i < max; i++) { // Modules
1172           TString modVol(Form("%s/F%c%cV_7%d/F%cSE_%d", halfVol.Data(),
1173                               sio, stb, nio, sio, i));
1174           TString modSym(modVol);
1175           if(!gGeoManager->SetAlignableEntry(modSym.Data(),modVol.Data()))
1176             AliFatal(Form("Alignable entry %s not created. "
1177                           "Volume path %s not valid",
1178                           modSym.Data(), modVol.Data()));
1179         }
1180       }
1181     }
1182   }
1183 #endif
1184 }
1185 //___________________________________________________________________
1186 //
1187 // EOF
1188 //