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