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