Several Changes:
[u/mrichter/AliRoot.git] / 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 <TRotMatrix.h>         // ROOT_TRotMatrix
90 #include <TTree.h>              // ROOT_TTree
91 #include <TBrowser.h>           // ROOT_TBrowser
92 #include <TVirtualMC.h>         // ROOT_TVirtualMC
93 #include <TVector2.h>           // ROOT_TVector2 
94 #include <TGeoManager.h>        // ROOT_TGeoManager
95
96 #include <AliRunDigitizer.h>    // ALIRUNDIGITIZER_H
97 #include <AliLoader.h>          // ALILOADER_H
98 #include <AliRun.h>             // ALIRUN_H
99 #include <AliMC.h>              // ALIMC_H
100 #include <AliMagF.h>            // ALIMAGF_H
101 // #include <AliLog.h>          // ALILOG_H
102 #include "AliFMDDebug.h" // Better debug macros
103 #include "AliFMD.h"             // ALIFMD_H
104 #include "AliFMDDigit.h"        // ALIFMDDIGIT_H
105 #include "AliFMDSDigit.h"       // ALIFMDSDIGIT_H
106 #include "AliFMDHit.h"          // ALIFMDHIT_H
107 #include "AliFMDGeometry.h"     // ALIFMDGEOMETRY_H
108 #include "AliFMDDetector.h"     // ALIFMDDETECTOR_H
109 #include "AliFMDRing.h"         // ALIFMDRING_H
110 #include "AliFMDDigitizer.h"    // ALIFMDDIGITIZER_H
111 #include "AliFMDHitDigitizer.h" // ALIFMDSDIGITIZER_H
112 // #define USE_SSDIGITIZER 
113 //#ifdef USE_SSDIGITIZER
114 //# include "AliFMDSSDigitizer.h"       // ALIFMDSDIGITIZER_H
115 //#endif
116 // #include "AliFMDGeometryBuilder.h"
117 #include "AliFMDRawWriter.h"    // ALIFMDRAWWRITER_H
118
119 //____________________________________________________________________
120 ClassImp(AliFMD)
121 #if 0
122   ; // This is to keep Emacs from indenting the next line 
123 #endif 
124
125 //____________________________________________________________________
126 AliFMD::AliFMD()
127   : AliDetector(),
128     fSDigits(0), 
129     fNsdigits(0),
130     fDetailed(kTRUE),
131     fUseOld(kFALSE),
132     fUseAssembly(kTRUE),
133     fBad(0) 
134 {
135   //
136   // Default constructor for class AliFMD
137   //
138   AliFMDDebug(10, ("\tDefault CTOR"));
139   fHits        = 0;
140   fDigits      = 0;
141   fIshunt      = 0;
142   fBad         = new TClonesArray("AliFMDHit");
143 }
144
145 //____________________________________________________________________
146 AliFMD::AliFMD(const char *name, const char *title)
147   : AliDetector (name, title),
148     fSDigits(0),
149     fNsdigits(0),
150     fDetailed(kTRUE),
151     fUseOld(kFALSE),
152     fUseAssembly(kFALSE),
153     fBad(0)
154 {
155   //
156   // Standard constructor for Forward Multiplicity Detector
157   //
158   AliFMDDebug(10, ("\tStandard CTOR"));
159   fBad         = new TClonesArray("AliFMDHit");
160   
161   // Initialise Hit array
162   HitsArray();
163   gAlice->GetMCApp()->AddHitList(fHits);
164
165   // (S)Digits for the detectors disk
166   DigitsArray();
167   SDigitsArray();
168   
169   // CHC: What is this?
170   fIshunt = 0;
171   //PH  SetMarkerColor(kRed);
172   //PH  SetLineColor(kYellow);
173 }
174
175 //____________________________________________________________________
176 AliFMD::~AliFMD ()
177 {
178   // Destructor for base class AliFMD
179   if (fHits) {
180     fHits->Delete();
181     delete fHits;
182     fHits = 0;
183   }
184   if (fDigits) {
185     fDigits->Delete();
186     delete fDigits;
187     fDigits = 0;
188   }
189   if (fSDigits) {
190     fSDigits->Delete();
191     delete fSDigits;
192     fSDigits = 0;
193   }
194   if (fBad) {
195     fBad->Delete();
196     delete fBad;
197     fBad = 0;
198   }
199 }
200
201
202 //====================================================================
203 //
204 // GEometry ANd Traking
205 //
206 //____________________________________________________________________
207 void 
208 AliFMD::CreateGeometry()
209 {
210   //
211   // Create the geometry of Forward Multiplicity Detector.  The actual
212   // construction of the geometry is delegated to the class
213   // AliFMDGeometryBuilder, invoked by the singleton manager
214   // AliFMDGeometry. 
215   //
216   AliFMDGeometry*  fmd = AliFMDGeometry::Instance();
217   fmd->SetDetailed(fDetailed);
218   fmd->UseAssembly(fUseAssembly);
219   fmd->Build();
220 }    
221
222 //____________________________________________________________________
223 void AliFMD::CreateMaterials() 
224 {
225   // Define the materials and tracking mediums needed by the FMD
226   // simulation.   These mediums are made by sending the messages
227   // AliMaterial, AliMixture, and AliMedium to the passed AliModule
228   // object module.   The defined mediums are 
229   // 
230   //    FMD Si$         Silicon (active medium in sensors)
231   //    FMD C$          Carbon fibre (support cone for FMD3 and vacuum pipe)
232   //    FMD Al$         Aluminium (honeycomb support plates)
233   //    FMD PCB$        Printed Circuit Board (FEE board with VA1_3)
234   //    FMD Chip$       Electronics chips (currently not used)
235   //    FMD Air$        Air (Air in the FMD)
236   //    FMD Plastic$    Plastic (Support legs for the hybrid cards)
237   //
238   // The geometry builder should really be the one that creates the
239   // materials, but the architecture of AliROOT makes that design
240   // akward.  What should happen, was that the AliFMDGeometryBuilder
241   // made the mediums, and that this class retrives pointers from the
242   // TGeoManager, and registers the mediums here.  Alas, it's not
243   // really that easy. 
244   //
245   AliFMDDebug(10, ("\tCreating materials"));
246   // Get pointer to geometry singleton object. 
247   AliFMDGeometry* geometry = AliFMDGeometry::Instance();
248   geometry->Init();
249 #if 0
250   if (gGeoManager && gGeoManager->GetMedium("FMD Si$")) {
251     // We need to figure out the some stuff about the geometry
252     fmd->ExtractGeomInfo();
253     return;
254   }
255 #endif  
256   Int_t    id;
257   Double_t a                = 0;
258   Double_t z                = 0;
259   Double_t density          = 0;
260   Double_t radiationLength  = 0;
261   Double_t absorbtionLength = 999;
262   Int_t    fieldType        = gAlice->Field()->Integ();     // Field type 
263   Double_t maxField         = gAlice->Field()->Max();     // Field max.
264   Double_t maxBending       = 0;     // Max Angle
265   Double_t maxStepSize      = 0.001; // Max step size 
266   Double_t maxEnergyLoss    = 1;     // Max Delta E
267   Double_t precision        = 0.001; // Precision
268   Double_t minStepSize      = 0.001; // Minimum step size 
269  
270   // Silicon 
271   a                = 28.0855;
272   z                = 14.;
273   density          = geometry->GetSiDensity();
274   radiationLength  = 9.36;
275   maxBending       = 1;
276   maxStepSize      = .001;
277   precision        = .001;
278   minStepSize      = .001;
279   id               = kSiId;
280   AliMaterial(id, "Si$", a, z, density, radiationLength, absorbtionLength);
281   AliMedium(kSiId, "Si$", id,1,fieldType,maxField,maxBending,
282             maxStepSize,maxEnergyLoss,precision,minStepSize);
283   
284
285   // Carbon 
286   a                = 12.011;
287   z                = 6.;
288   density          = 2.265;
289   radiationLength  = 18.8;
290   maxBending       = 10;
291   maxStepSize      = .01;
292   precision        = .003;
293   minStepSize      = .003;
294   id               = kCarbonId;
295   AliMaterial(id, "Carbon$", a, z, density, radiationLength, absorbtionLength);
296   AliMedium(kCarbonId, "Carbon$", id,0,fieldType,maxField,maxBending,
297                     maxStepSize,maxEnergyLoss,precision,minStepSize);
298
299   // Aluminum
300   a                = 26.981539;
301   z                = 13.;
302   density          = 2.7;
303   radiationLength  = 8.9;
304   id               = kAlId;
305   AliMaterial(id, "Aluminum$",a,z, density, radiationLength, absorbtionLength);
306   AliMedium(kAlId, "Aluminum$", id, 0, fieldType, maxField, maxBending,
307             maxStepSize, maxEnergyLoss, precision, minStepSize);
308   
309   
310   // Copper 
311   a                = 63.546;
312   z                = 29;
313   density          =  8.96;
314   radiationLength  =  1.43;
315   id               = kCopperId;
316   AliMaterial(id, "Copper$", 
317                       a, z, density, radiationLength, absorbtionLength);
318   AliMedium(kCopperId, "Copper$", id, 0, fieldType, maxField, maxBending,
319             maxStepSize, maxEnergyLoss, precision, minStepSize);
320   
321
322   // Silicon chip 
323   {
324     Float_t as[] = { 12.0107,      14.0067,      15.9994,
325                       1.00794,     28.0855,     107.8682 };
326     Float_t zs[] = {  6.,           7.,           8.,
327                       1.,          14.,          47. };
328     Float_t ws[] = {  0.039730642,  0.001396798,  0.01169634,
329                       0.004367771,  0.844665,     0.09814344903 };
330     density          = 2.36436;
331     maxBending       = 10;
332     maxStepSize      = .01;
333     precision        = .003;
334     minStepSize      = .003;
335     id               = kSiChipId;
336     AliMixture(id, "Si Chip$", as, zs, density, 6, ws);
337     AliMedium(kSiChipId, "Si Chip$",  id, 0, fieldType, maxField, maxBending, 
338               maxStepSize, maxEnergyLoss, precision, minStepSize);
339   }
340   
341   // Kaption
342   {
343     Float_t as[] = { 1.00794,  12.0107,  14.010,   15.9994};
344     Float_t zs[] = { 1.,        6.,       7.,       8.};
345     Float_t ws[] = { 0.026362,  0.69113,  0.07327,  0.209235};
346     density          = 1.42;
347     maxBending       = 1;
348     maxStepSize      = .001;
349     precision        = .001;
350     minStepSize      = .001;
351     id               = kKaptonId;
352     AliMixture(id, "Kaption$", as, zs, density, 4, ws);
353     AliMedium(kKaptonId, "Kaption$", id,0,fieldType,maxField,maxBending,
354               maxStepSize,maxEnergyLoss,precision,minStepSize);
355   }
356
357   // Air
358   {
359     Float_t as[] = { 12.0107, 14.0067,   15.9994,  39.948 };
360     Float_t zs[] = {  6.,      7.,       8.,       18. };
361     Float_t ws[] = { 0.000124, 0.755267, 0.231781, 0.012827 }; 
362     density      = .00120479;
363     maxBending   = 1;
364     maxStepSize  = .001;
365     precision    = .001;
366     minStepSize  = .001;
367     id           = kAirId;
368     AliMixture(id, "Air$", as, zs, density, 4, ws);
369     AliMedium(kAirId, "Air$", id,0,fieldType,maxField,maxBending,
370               maxStepSize,maxEnergyLoss,precision,minStepSize);
371   }
372   
373   // PCB
374   {
375     Float_t zs[] = { 14.,         20.,         13.,         12.,
376                       5.,         22.,         11.,         19.,
377                      26.,          9.,          8.,          6.,
378                       7.,          1.};
379     Float_t as[] = { 28.0855,     40.078,      26.981538,   24.305, 
380                      10.811,      47.867,      22.98977,    39.0983,
381                      55.845,      18.9984,     15.9994,     12.0107,
382                      14.0067,      1.00794};
383     Float_t ws[] = {  0.15144894,  0.08147477,  0.04128158,  0.00904554, 
384                       0.01397570,  0.00287685,  0.00445114,  0.00498089,
385                       0.00209828,  0.00420000,  0.36043788,  0.27529426,
386                       0.01415852,  0.03427566};
387     density      = 1.8;
388     maxBending   = 1;
389     maxStepSize  = .001;
390     precision    = .001;
391     minStepSize  = .001;
392     id           = kPcbId;
393     AliMixture(id, "PCB$", as, zs, density, 14, ws);
394     AliMedium(kPcbId, "PCB$", id,0,fieldType,maxField,maxBending,
395               maxStepSize,maxEnergyLoss,precision,minStepSize);
396   }
397   
398   // Stainless steel
399   {
400     Float_t as[] = { 55.847, 51.9961, 58.6934, 28.0855 };
401     Float_t zs[] = { 26.,    24.,     28.,     14.     };
402     Float_t ws[] = { .715,   .18,     .1,      .005    };
403     density      = 7.88;
404     id           = kSteelId;
405     AliMixture(id, "Steel$", as, zs, density, 4, ws);
406     AliMedium(kSteelId, "Steel$", id, 0, fieldType, maxField, maxBending, 
407               maxStepSize, maxEnergyLoss, precision, minStepSize);
408   }
409   // Plastic 
410   {
411     Float_t as[] = { 1.01, 12.01 };
412     Float_t zs[] = { 1.,   6.    };
413     Float_t ws[] = { 1.,   1.    };
414     density      = 1.03;
415     maxBending   = 10;
416     maxStepSize  = .01;
417     precision    = .003;
418     minStepSize  = .003;
419     id           = kPlasticId;
420     AliMixture(id, "Plastic$", as, zs, density, -2, ws);
421     AliMedium(kPlasticId, "Plastic$", id,0,fieldType,maxField,maxBending,
422               maxStepSize,maxEnergyLoss,precision,minStepSize);
423   }
424
425 }
426
427 //____________________________________________________________________
428 void  
429 AliFMD::Init()
430 {
431   // Initialize the detector 
432   // 
433   AliFMDDebug(1, ("Initialising FMD detector object"));
434   TVirtualMC*      mc     = TVirtualMC::GetMC();
435   AliFMDGeometry*  fmd    = AliFMDGeometry::Instance();
436   const TArrayI&   actGeo = fmd->ActiveIds();
437   TArrayI          actVmc(actGeo.fN);
438   for (Int_t i = 0; i < actGeo.fN; i++) {
439     TGeoVolume *sens = gGeoManager->GetVolume(actGeo[i]);
440     if (!sens) {
441       AliError(Form("No TGeo volume for sensitive volume ID=%d",actGeo[i]));
442       continue;
443     }   
444     actVmc[i] = mc->VolId(sens->GetName());
445     AliFMDDebug(1, ("Active vol id # %d: %d changed to %d", 
446                     i, actGeo[i], actVmc[i]));
447   }
448   fmd->SetActive(actVmc.fArray, actVmc.fN);
449   // fmd->InitTransformations();
450 }
451
452 //____________________________________________________________________
453 void
454 AliFMD::FinishEvent()
455 {
456   // Called at the end of the an event in simulations.  If the debug
457   // level is high enough, then the `bad' hits are printed.
458   // 
459   if (AliLog::GetDebugLevel("FMD", "AliFMD") < 10) return;
460   if (fBad && fBad->GetEntries() > 0) {
461     AliWarning((Form("EndEvent", "got %d 'bad' hits", fBad->GetEntries())));
462     TIter next(fBad);
463     AliFMDHit* hit;
464     while ((hit = static_cast<AliFMDHit*>(next()))) hit->Print("D");
465     fBad->Clear();
466   }
467 }
468
469
470
471 //====================================================================
472 //
473 // Hit and Digit managment 
474 //
475 //____________________________________________________________________
476 void 
477 AliFMD::MakeBranch(Option_t * option)
478 {
479   // Create Tree branches for the FMD.
480   //
481   // Options:
482   //
483   //    H          Make a branch of TClonesArray of AliFMDHit's
484   //    D          Make a branch of TClonesArray of AliFMDDigit's
485   //    S          Make a branch of TClonesArray of AliFMDSDigit's
486   // 
487   const Int_t kBufferSize = 16000;
488   TString branchname(GetName());
489   TString opt(option);
490   
491   if (opt.Contains("H", TString::kIgnoreCase)) {
492     HitsArray();
493     AliDetector::MakeBranch(option); 
494   }
495   if (opt.Contains("D", TString::kIgnoreCase)) { 
496     DigitsArray();
497     MakeBranchInTree(fLoader->TreeD(), branchname.Data(),
498                      &fDigits, kBufferSize, 0);
499   }
500   if (opt.Contains("S", TString::kIgnoreCase)) { 
501     SDigitsArray();
502     MakeBranchInTree(fLoader->TreeS(), branchname.Data(),
503                      &fSDigits, kBufferSize, 0);
504   }
505 }
506
507 //____________________________________________________________________
508 void 
509 AliFMD::SetTreeAddress()
510 {
511   // Set branch address for the Hits, Digits, and SDigits Tree.
512   if (fLoader->TreeH()) HitsArray();
513   AliDetector::SetTreeAddress();
514
515   TTree *treeD = fLoader->TreeD();
516   if (treeD) {
517     DigitsArray();
518     TBranch* branch = treeD->GetBranch ("FMD");
519     if (branch) branch->SetAddress(&fDigits);
520   }
521
522   TTree *treeS = fLoader->TreeS();
523   if (treeS) {
524     SDigitsArray();
525     TBranch* branch = treeS->GetBranch ("FMD");
526     if (branch) branch->SetAddress(&fSDigits);
527   }
528 }
529
530 //____________________________________________________________________
531 void 
532 AliFMD::SetHitsAddressBranch(TBranch *b)
533 {
534   // Set the TClonesArray to read hits into. 
535   b->SetAddress(&fHits);
536 }
537 //____________________________________________________________________
538 void 
539 AliFMD::SetSDigitsAddressBranch(TBranch *b)
540 {
541   // Set the TClonesArray to read hits into. 
542   b->SetAddress(&fSDigits);
543 }
544
545 //____________________________________________________________________
546 void 
547 AliFMD::AddHit(Int_t track, Int_t *vol, Float_t *hits) 
548 {
549   // Add a hit to the hits tree 
550   // 
551   // The information of the two arrays are decoded as 
552   // 
553   // Parameters
554   //    track                Track #
555   //    ivol[0]  [UShort_t ] Detector # 
556   //    ivol[1]  [Char_t   ] Ring ID 
557   //    ivol[2]  [UShort_t ] Sector #
558   //    ivol[3]  [UShort_t ] Strip # 
559   //    hits[0]  [Float_t  ] Track's X-coordinate at hit 
560   //    hits[1]  [Float_t  ] Track's Y-coordinate at hit
561   //    hits[3]  [Float_t  ] Track's Z-coordinate at hit
562   //    hits[4]  [Float_t  ] X-component of track's momentum             
563   //    hits[5]  [Float_t  ] Y-component of track's momentum             
564   //    hits[6]  [Float_t  ] Z-component of track's momentum            
565   //    hits[7]  [Float_t  ] Energy deposited by track                  
566   //    hits[8]  [Int_t    ] Track's particle Id # 
567   //    hits[9]  [Float_t  ] Time when the track hit
568   // 
569   // 
570   AddHitByFields(track, 
571                  UShort_t(vol[0]),  // Detector # 
572                  Char_t(vol[1]),    // Ring ID
573                  UShort_t(vol[2]),  // Sector # 
574                  UShort_t(vol[3]),  // Strip # 
575                  hits[0],           // X
576                  hits[1],           // Y
577                  hits[2],           // Z
578                  hits[3],           // Px
579                  hits[4],           // Py
580                  hits[5],           // Pz
581                  hits[6],           // Energy loss 
582                  Int_t(hits[7]),    // PDG 
583                  hits[8]);          // Time
584 }
585
586 //____________________________________________________________________
587 AliFMDHit*
588 AliFMD::AddHitByFields(Int_t    track, 
589                        UShort_t detector, 
590                        Char_t   ring, 
591                        UShort_t sector, 
592                        UShort_t strip, 
593                        Float_t  x, 
594                        Float_t  y, 
595                        Float_t  z,
596                        Float_t  px, 
597                        Float_t  py, 
598                        Float_t  pz,
599                        Float_t  edep,
600                        Int_t    pdg,
601                        Float_t  t, 
602                        Float_t  l, 
603                        Bool_t   stop)
604 {
605   // Add a hit to the list
606   //
607   // Parameters:
608   // 
609   //    track     Track #
610   //    detector  Detector # (1, 2, or 3)                      
611   //    ring      Ring ID ('I' or 'O')
612   //    sector    Sector # (For inner/outer rings: 0-19/0-39)
613   //    strip     Strip # (For inner/outer rings: 0-511/0-255)
614   //    x         Track's X-coordinate at hit
615   //    y         Track's Y-coordinate at hit
616   //    z         Track's Z-coordinate at hit
617   //    px        X-component of track's momentum 
618   //    py        Y-component of track's momentum
619   //    pz        Z-component of track's momentum
620   //    edep      Energy deposited by track
621   //    pdg       Track's particle Id #
622   //    t         Time when the track hit 
623   //    l         Track length through the material. 
624   //    stop      Whether track was stopped or disappeared
625   // 
626   TClonesArray& a = *(HitsArray());
627   // Search through the list of already registered hits, and see if we
628   // find a hit with the same parameters.  If we do, then don't create
629   // a new hit, but rather update the energy deposited in the hit.
630   // This is done, so that a FLUKA based simulation will get the
631   // number of hits right, not just the enerrgy deposition. 
632   AliFMDHit* hit = 0;
633   for (Int_t i = 0; i < fNhits; i++) {
634     if (!a.At(i)) continue;
635     hit = static_cast<AliFMDHit*>(a.At(i));
636     if (hit->Detector() == detector 
637         && hit->Ring() == ring
638         && hit->Sector() == sector 
639         && hit->Strip() == strip
640         && hit->Track() == track) {
641       AliFMDDebug(1, ("already had a hit in FMD%d%c[%2d,%3d] for track # %d,"
642                        " adding energy (%f) to that hit (%f) -> %f", 
643                        detector, ring, sector, strip, track, edep, hit->Edep(),
644                        hit->Edep() + edep));
645       hit->SetEdep(hit->Edep() + edep);
646       return hit;
647     }
648   }
649   // If hit wasn't already registered, do so know. 
650   hit = new (a[fNhits]) AliFMDHit(fIshunt, track, detector, ring, sector, 
651                                   strip, x, y, z, px, py, pz, edep, pdg, t, 
652                                   l, stop);
653   fNhits++;
654   return hit;
655 }
656
657 //____________________________________________________________________
658 void 
659 AliFMD::AddDigit(Int_t* digits, Int_t*)
660 {
661   // Add a digit to the Digit tree 
662   // 
663   // Paramters 
664   //
665   //    digits[0]  [UShort_t] Detector #
666   //    digits[1]  [Char_t]   Ring ID
667   //    digits[2]  [UShort_t] Sector #
668   //    digits[3]  [UShort_t] Strip #
669   //    digits[4]  [UShort_t] ADC Count 
670   //    digits[5]  [Short_t]  ADC Count, -1 if not used
671   //    digits[6]  [Short_t]  ADC Count, -1 if not used 
672   // 
673   AddDigitByFields(UShort_t(digits[0]),  // Detector #
674                    Char_t(digits[1]),    // Ring ID
675                    UShort_t(digits[2]),  // Sector #
676                    UShort_t(digits[3]),  // Strip #
677                    UShort_t(digits[4]),  // ADC Count1 
678                    Short_t(digits[5]),   // ADC Count2 
679                    Short_t(digits[6]),   // ADC Count3 
680                    Short_t(digits[7])); 
681 }
682
683 //____________________________________________________________________
684 void 
685 AliFMD::AddDigitByFields(UShort_t detector, 
686                          Char_t   ring, 
687                          UShort_t sector, 
688                          UShort_t strip, 
689                          UShort_t count1, 
690                          Short_t  count2,
691                          Short_t  count3, 
692                          Short_t  count4)
693 {
694   // add a real digit - as coming from data
695   // 
696   // Parameters 
697   //
698   //    detector  Detector # (1, 2, or 3)                      
699   //    ring      Ring ID ('I' or 'O')
700   //    sector    Sector # (For inner/outer rings: 0-19/0-39)
701   //    strip     Strip # (For inner/outer rings: 0-511/0-255)
702   //    count1    ADC count (a 10-bit word)
703   //    count2    ADC count (a 10-bit word), or -1 if not used
704   //    count3    ADC count (a 10-bit word), or -1 if not used
705   TClonesArray& a = *(DigitsArray());
706   
707   new (a[fNdigits++]) 
708     AliFMDDigit(detector, ring, sector, strip, count1, count2, count3, count4);
709   AliFMDDebug(15, ("Adding digit # %5d/%5d for FMD%d%c[%2d,%3d]=(%d,%d,%d,%d)",
710                    fNdigits-1, a.GetEntriesFast(),
711                    detector, ring, sector, strip, 
712                    count1, count2, count3, count4));
713   
714 }
715
716 //____________________________________________________________________
717 void 
718 AliFMD::AddSDigit(Int_t* digits)
719 {
720   // Add a digit to the SDigit tree 
721   // 
722   // Paramters 
723   //
724   //    digits[0]  [UShort_t] Detector #
725   //    digits[1]  [Char_t]   Ring ID
726   //    digits[2]  [UShort_t] Sector #
727   //    digits[3]  [UShort_t] Strip #
728   //    digits[4]  [Float_t]  Total energy deposited 
729   //    digits[5]  [UShort_t] ADC Count 
730   //    digits[6]  [Short_t]  ADC Count, -1 if not used
731   //    digits[7]  [Short_t]  ADC Count, -1 if not used 
732   // 
733   AddSDigitByFields(UShort_t(digits[0]),   // Detector #
734                     Char_t(digits[1]),     // Ring ID
735                     UShort_t(digits[2]),   // Sector #
736                     UShort_t(digits[3]),   // Strip #
737                     Float_t(digits[4]),    // Edep
738                     UShort_t(digits[5]),   // ADC Count1 
739                     Short_t(digits[6]),    // ADC Count2 
740                     Short_t(digits[7]),    // ADC Count3 
741                     Short_t(digits[8]),    // ADC Count4
742                     UShort_t(digits[9]),   // N particles
743                     UShort_t(digits[10])); // N primaries
744    
745 }
746
747 //____________________________________________________________________
748 void 
749 AliFMD::AddSDigitByFields(UShort_t detector, 
750                           Char_t   ring, 
751                           UShort_t sector, 
752                           UShort_t strip, 
753                           Float_t  edep,
754                           UShort_t count1, 
755                           Short_t  count2,
756                           Short_t  count3, 
757                           Short_t  count4, 
758                           UShort_t ntot, 
759                           UShort_t nprim)
760 {
761   // add a summable digit
762   // 
763   // Parameters 
764   //
765   //    detector  Detector # (1, 2, or 3)                      
766   //    ring      Ring ID ('I' or 'O')
767   //    sector    Sector # (For inner/outer rings: 0-19/0-39)
768   //    strip     Strip # (For inner/outer rings: 0-511/0-255)
769   //    edep      Total energy deposited
770   //    count1    ADC count (a 10-bit word)
771   //    count2    ADC count (a 10-bit word), or -1 if not used
772   //    count3    ADC count (a 10-bit word), or -1 if not used
773   //
774   TClonesArray& a = *(SDigitsArray());
775   // AliFMDDebug(0, ("Adding sdigit # %d", fNsdigits));
776   
777   new (a[fNsdigits++]) 
778     AliFMDSDigit(detector, ring, sector, strip, edep, 
779                  count1, count2, count3, count4, ntot, nprim);
780 }
781
782 //____________________________________________________________________
783 void 
784 AliFMD::ResetSDigits()
785 {
786   // Reset number of digits and the digits array for this detector. 
787   //
788   fNsdigits   = 0;
789   if (fSDigits) fSDigits->Clear();
790 }
791
792
793 //____________________________________________________________________
794 TClonesArray*
795 AliFMD::HitsArray() 
796 {
797   // Initialize hit array if not already, and return pointer to it. 
798   if (!fHits) { 
799     fHits = new TClonesArray("AliFMDHit", 1000);
800     fNhits = 0;
801   }
802   return fHits;
803 }
804
805 //____________________________________________________________________
806 TClonesArray*
807 AliFMD::DigitsArray() 
808 {
809   // Initialize digit array if not already, and return pointer to it. 
810   if (!fDigits) { 
811     fDigits = new TClonesArray("AliFMDDigit", 1000);
812     fNdigits = 0;
813   }
814   return fDigits;
815 }
816
817 //____________________________________________________________________
818 TClonesArray*
819 AliFMD::SDigitsArray() 
820 {
821   // Initialize digit array if not already, and return pointer to it. 
822   if (!fSDigits) { 
823     fSDigits = new TClonesArray("AliFMDSDigit", 1000);
824     fNsdigits = 0;
825   }
826   return fSDigits;
827 }
828
829 //====================================================================
830 //
831 // Digitization 
832 //
833 //____________________________________________________________________
834 void 
835 AliFMD::Hits2Digits() 
836 {
837   // Create AliFMDDigit's from AliFMDHit's.  This is done by making a
838   // AliFMDDigitizer, and executing that code.
839   // 
840   AliFMDHitDigitizer digitizer(this, AliFMDHitDigitizer::kDigits);
841   digitizer.Init();
842   digitizer.Exec("");
843 }
844
845 //____________________________________________________________________
846 void 
847 AliFMD::Hits2SDigits() 
848 {
849   // Create AliFMDSDigit's from AliFMDHit's.  This is done by creating
850   // an AliFMDSDigitizer object, and executing it. 
851   // 
852   AliFMDHitDigitizer digitizer(this, AliFMDHitDigitizer::kSDigits);
853   digitizer.Init();
854   digitizer.Exec("");
855 }
856
857   
858 //____________________________________________________________________
859 AliDigitizer* 
860 AliFMD::CreateDigitizer(AliRunDigitizer* manager) const
861 {
862   // Create a digitizer object 
863   
864   /* This is what we probably _should_ do */
865   AliFMDBaseDigitizer* digitizer = 0;
866   
867 #ifdef USE_SSDIGITIZER
868   digitizer = new AliFMDSSDigitizer(manager);
869 #else 
870   /* This is what we actually do, and will work */
871 #if 0
872   AliInfo("SDigit->Digit conversion not really supported, "
873           "doing Hit->Digit conversion instead");
874 #endif
875   digitizer = new AliFMDDigitizer(manager);
876 #endif
877   return digitizer;
878 }
879
880 //====================================================================
881 //
882 // Raw data simulation 
883 //
884 //__________________________________________________________________
885 void 
886 AliFMD::Digits2Raw() 
887 {
888   // Turn digits into raw data. 
889   // 
890   // This uses the class AliFMDRawWriter to do the job.   Please refer
891   // to that class for more information. 
892   AliFMDRawWriter writer(this);
893   writer.Exec();
894 }
895
896
897 //====================================================================
898 //
899 // Utility 
900 //
901 //__________________________________________________________________
902 void 
903 AliFMD::Browse(TBrowser* b) 
904 {
905   // Browse this object. 
906   //
907   AliFMDDebug(30, ("\tBrowsing the FMD"));
908   AliDetector::Browse(b);
909   b->Add(AliFMDGeometry::Instance());
910 }
911
912 //____________________________________________________________________  
913 void
914 AliFMD::AddAlignableVolumes() const
915 {
916   //
917   // Create entries for alignable volumes associating the symbolic volume
918   // name with the corresponding volume path. Needs to be syncronized with
919   // eventual changes in the geometry.
920   // 
921   // This code was made by Raffaele Grosso <rgrosso@mail.cern.ch>.  I
922   // (cholm) will probably want to change it.   For one, I think it
923   // should be the job of the geometry manager to deal with this. 
924   AliInfo("Add FMD alignable volumes");
925   AliFMDGeometry::Instance()->SetAlignableVolumes();
926 #if 0  
927   for(size_t f = 1; f <= 3; f++){ // Detector 1,2,3
928     for(size_t tb =  0; tb <2 ; tb++){ // Top/Bottom 
929       char     stb = tb == 0 ? 'T' : 'B';
930       unsigned min = tb == 0 ? 0   : 5;
931
932       TString halfVol(Form("/ALIC_1/F%dM%c_%d", f, stb, f));
933       TString halfSym(halfVol);
934       if(!gGeoManager->SetAlignableEntry(halfSym.Data(),halfVol.Data()))
935         AliFatal(Form("Alignable entry %s not created. "
936                       "Volume path %s not valid", 
937                       halfSym.Data(),halfVol.Data()));
938       for(size_t io = 0; io < 2; io++){ // inner, outer 
939         if (f==1 && io==1) continue; // Only one ring in FMD1 
940         if(tb == 1 && io==1) min=10;
941         char     sio = (io == 0 ? 'I' : 'O');
942         unsigned nio = (io == 0 ? 3   : 9);
943         unsigned max = (io == 0 ? 5   : 10) + min;
944         
945         for(size_t i = min; i < max; i++) { // Modules
946           TString modVol(Form("%s/F%c%cV_7%d/F%cSE_%d", halfVol.Data(), 
947                               sio, stb, nio, sio, i));
948           TString modSym(modVol);
949           if(!gGeoManager->SetAlignableEntry(modSym.Data(),modVol.Data()))
950             AliFatal(Form("Alignable entry %s not created. "
951                           "Volume path %s not valid", 
952                           modSym.Data(), modVol.Data()));
953         }
954       }
955     }
956   }
957 #endif
958 }
959 //___________________________________________________________________
960 //
961 // EOF
962 //