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