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