Fixes, and extra debug
[u/mrichter/AliRoot.git] / FMD / AliFMDDigitizer.cxx
1 /**************************************************************************
2  * Copyright(c) 2004, 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
16 /* $Id$ */
17
18 //////////////////////////////////////////////////////////////////////////////
19 //
20 //  This class contains the procedures simulation ADC  signal for the
21 //  Forward Multiplicity detector  : Hits->Digits and Hits->SDigits
22 // 
23 //  Digits consists of
24 //   - Detector #
25 //   - Ring ID                                             
26 //   - Sector #     
27 //   - Strip #
28 //   - ADC count in this channel                                  
29 //
30 //  Digits consists of
31 //   - Detector #
32 //   - Ring ID                                             
33 //   - Sector #     
34 //   - Strip #
35 //   - Total energy deposited in the strip
36 //   - ADC count in this channel                                  
37 //
38 // As the Digits and SDigits have so much in common, the classes
39 // AliFMDDigitizer and AliFMDSDigitizer are implemented via a base
40 // class AliFMDBaseDigitizer.
41 //
42 //                 +---------------------+
43 //                 | AliFMDBaseDigitizer |
44 //                 +---------------------+
45 //                           ^
46 //                           |
47 //                +----------+---------+
48 //                |                    |
49 //      +-----------------+     +------------------+
50 //      | AliFMDDigitizer |     | AliFMDSDigitizer |
51 //      +-----------------+     +------------------+
52 //
53 // These classes has several paramters: 
54 //
55 //     fPedestal
56 //     fPedestalWidth
57 //         (Only AliFMDDigitizer)
58 //         Mean and width of the pedestal.  The pedestal is simulated
59 //         by a Guassian, but derived classes my override MakePedestal
60 //         to simulate it differently (or pick it up from a database).
61 //
62 //     fVA1MipRange
63 //         The dymamic MIP range of the VA1_ALICE pre-amplifier chip 
64 //
65 //     fAltroChannelSize
66 //         The largest number plus one that can be stored in one
67 //         channel in one time step in the ALTRO ADC chip. 
68 //
69 //     fSampleRate
70 //         How many times the ALTRO ADC chip samples the VA1_ALICE
71 //         pre-amplifier signal.   The VA1_ALICE chip is read-out at
72 //         10MHz, while it's possible to drive the ALTRO chip at
73 //         25MHz.  That means, that the ALTRO chip can have time to
74 //         sample each VA1_ALICE signal up to 2 times.  Although it's
75 //         not certain this feature will be used in the production,
76 //         we'd like have the option, and so it should be reflected in
77 //         the code.
78 //
79 //
80 // The shaping function of the VA1_ALICE is generally given by 
81 //
82 //      f(x) = A(1 - exp(-Bx))
83 //
84 // where A is the total charge collected in the pre-amp., and B is a
85 // paramter that depends on the shaping time of the VA1_ALICE circut.
86 // 
87 // When simulating the shaping function of the VA1_ALICe
88 // pre-amp. chip, we have to take into account, that the shaping
89 // function depends on the previous value of read from the pre-amp. 
90 //
91 // That results in the following algorithm:
92 //
93 //    last = 0;
94 //    FOR charge IN pre-amp. charge train DO 
95 //      IF last < charge THEN 
96 //        f(t) = (charge - last) * (1 - exp(-B * t)) + last
97 //      ELSE
98 //        f(t) = (last - charge) * exp(-B * t) + charge)
99 //      ENDIF
100 //      FOR i IN # samples DO 
101 //        adc_i = f(i / (# samples))
102 //      DONE
103 //      last = charge
104 //   DONE
105 //
106 // Here, 
107 //
108 //   pre-amp. charge train 
109 //       is a series of 128 charges read from the VA1_ALICE chip
110 //
111 //   # samples
112 //       is the number of times the ALTRO ADC samples each of the 128
113 //       charges from the pre-amp. 
114 //
115 // Where Q is the total charge collected by the VA1_ALICE
116 // pre-amplifier.   Q is then given by 
117 //
118 //           E S 
119 //      Q =  - -
120 //           e R
121 //
122 // where E is the total energy deposited in a silicon strip, R is the
123 // dynamic range of the VA1_ALICE pre-amp (fVA1MipRange), e is the
124 // energy deposited by a single MIP, and S ALTRO channel size in each
125 // time step (fAltroChannelSize).  
126 //
127 // The energy deposited per MIP is given by 
128 //
129 //      e = M * rho * w 
130 //
131 // where M is the universal number 1.664, rho is the density of
132 // silicon, and w is the depth of the silicon sensor. 
133 //
134 // The final ADC count is given by 
135 //
136 //      C' = C + P
137 //
138 // where P is the (randomized) pedestal (see MakePedestal)
139 //
140 // This class uses the class template AliFMDMap<Type> to make an
141 // internal cache of the energy deposted of the hits.  The class
142 // template is instantasized as 
143 //
144 //  typedef AliFMDMap<std::pair<Float_t, UShort_t> > AliFMDEdepMap;
145 //
146 // The first member of the values is the summed energy deposition in a
147 // given strip, while the second member of the values is the number of
148 // hits in a given strip.  Using the second member, it's possible to
149 // do some checks on just how many times a strip got hit, and what
150 // kind of error we get in our reconstructed hits.  Note, that this
151 // information is currently not written to the digits tree.  I think a
152 // QA (Quality Assurance) digit tree is better suited for that task.
153 // However, the information is there to be used in the future. 
154 //
155 //
156 // Latest changes by Christian Holm Christensen
157 //
158 //////////////////////////////////////////////////////////////////////////////
159
160 //      /1
161 //      |           A(-1 + B + exp(-B))
162 //      | f(x) dx = ------------------- = 1
163 //      |                    B
164 //      / 0
165 //
166 // and B is the a parameter defined by the shaping time (fShapingTime).  
167 //
168 // Solving the above equation, for A gives
169 //
170 //                 B
171 //      A = ----------------
172 //          -1 + B + exp(-B)
173 //
174 // So, if we define the function g: [0,1] -> [0:1] by 
175 //
176 //               / v
177 //               |              Bu + exp(-Bu) - Bv - exp(-Bv) 
178 //      g(u,v) = | f(x) dx = -A -----------------------------
179 //               |                            B
180 //               / u
181 //
182 // we can evaluate the ALTRO sample of the VA1_ALICE pre-amp between
183 // any two times (u, v), by 
184 //       
185 //
186 //                                B         Bu + exp(-Bu) - Bv - exp(-Bv)
187 //      C = Q g(u,v) = - Q ---------------- -----------------------------
188 //                         -1 + B + exp(-B)              B                  
189 //
190 //               Bu + exp(-Bu) - Bv - exp(-Bv) 
191 //        = -  Q -----------------------------
192 //                    -1 + B + exp(-B)
193 //
194
195 #include <TTree.h>              // ROOT_TTree
196 #include <TRandom.h>            // ROOT_TRandom
197 #include <AliLog.h>             // ALILOG_H
198 #include "AliFMDDigitizer.h"    // ALIFMDDIGITIZER_H
199 #include "AliFMD.h"             // ALIFMD_H
200 #include "AliFMDGeometry.h"     // ALIFMDGEOMETRY_H
201 #include "AliFMDDetector.h"     // ALIFMDDETECTOR_H
202 #include "AliFMDRing.h"         // ALIFMDRING_H
203 #include "AliFMDHit.h"          // ALIFMDHIT_H
204 #include "AliFMDDigit.h"        // ALIFMDDIGIT_H
205 #include <AliRunDigitizer.h>    // ALIRUNDIGITIZER_H
206 #include <AliRun.h>             // ALIRUN_H
207 #include <AliLoader.h>          // ALILOADER_H
208 #include <AliRunLoader.h>       // ALIRUNLOADER_H
209     
210 //====================================================================
211 ClassImp(AliFMDBaseDigitizer)
212 #if 0
213   ; // This is here to keep Emacs for indenting the next line
214 #endif
215
216 //____________________________________________________________________
217 AliFMDBaseDigitizer::AliFMDBaseDigitizer()  
218   : fRunLoader(0)
219 {
220   // Default ctor - don't use it
221 }
222
223 //____________________________________________________________________
224 AliFMDBaseDigitizer::AliFMDBaseDigitizer(AliRunDigitizer* manager) 
225   : AliDigitizer(manager, "AliFMDBaseDigitizer", "FMD Digitizer base class"), 
226     fRunLoader(0),
227     fEdep(AliFMDMap::kMaxDetectors, 
228           AliFMDMap::kMaxRings, 
229           AliFMDMap::kMaxSectors, 
230           AliFMDMap::kMaxStrips)
231 {
232   // Normal CTOR
233   AliDebug(1," processed");
234   SetVA1MipRange();
235   SetAltroChannelSize();
236   SetSampleRate();
237   SetShapingTime();
238 }
239
240 //____________________________________________________________________
241 AliFMDBaseDigitizer::AliFMDBaseDigitizer(const Char_t* name, 
242                                          const Char_t* title) 
243   : AliDigitizer(name, title),
244     fRunLoader(0),
245     fEdep(AliFMDMap::kMaxDetectors, 
246           AliFMDMap::kMaxRings, 
247           AliFMDMap::kMaxSectors, 
248           AliFMDMap::kMaxStrips)
249 {
250   // Normal CTOR
251   AliDebug(1," processed");
252   SetVA1MipRange();
253   SetAltroChannelSize();
254   SetSampleRate();
255   SetShapingTime();
256 }
257
258 //____________________________________________________________________
259 AliFMDBaseDigitizer::~AliFMDBaseDigitizer()
260 {
261   // Destructor
262 }
263
264 //____________________________________________________________________
265 Bool_t 
266 AliFMDBaseDigitizer::Init()
267 {
268   // Initialization
269   return kTRUE;
270 }
271  
272
273 //____________________________________________________________________
274 void
275 AliFMDBaseDigitizer::SumContributions(AliFMD* fmd) 
276 {
277   // Sum energy deposited contributions from each hit in a cache
278   // (fEdep).  
279   if (!fRunLoader) 
280     Fatal("SumContributions", "no run loader");
281   
282   // Clear array of deposited energies 
283   fEdep.Reset();
284   
285   // Get the FMD loader 
286   AliLoader* inFMD = fRunLoader->GetLoader("FMDLoader");
287   // And load the hits 
288   inFMD->LoadHits("READ");
289   
290   // Get the tree of hits 
291   TTree* hitsTree = inFMD->TreeH();
292   if (!hitsTree)  {
293     // Try again 
294     inFMD->LoadHits("READ");
295     hitsTree = inFMD->TreeH();
296   }
297   
298   // Get the FMD branch 
299   TBranch* hitsBranch = hitsTree->GetBranch("FMD");
300   if (hitsBranch) fmd->SetHitsAddressBranch(hitsBranch);
301   else            AliFatal("Branch FMD hit not found");
302   
303   // Get a list of hits from the FMD manager 
304   TClonesArray *fmdHits = fmd->Hits();
305   
306   // Get number of entries in the tree 
307   Int_t ntracks  = Int_t(hitsTree->GetEntries());
308   
309   // Loop over the tracks in the 
310   for (Int_t track = 0; track < ntracks; track++)  {
311     // Read in entry number `track' 
312     hitsBranch->GetEntry(track);
313     
314     // Get the number of hits 
315     Int_t nhits = fmdHits->GetEntries ();
316     for (Int_t hit = 0; hit < nhits; hit++) {
317       // Get the hit number `hit'
318       AliFMDHit* fmdHit = 
319         static_cast<AliFMDHit*>(fmdHits->UncheckedAt(hit));
320       
321       // Extract parameters 
322       UShort_t detector = fmdHit->Detector();
323       Char_t   ring     = fmdHit->Ring();
324       UShort_t sector   = fmdHit->Sector();
325       UShort_t strip    = fmdHit->Strip();
326       Float_t  edep     = fmdHit->Edep();
327       if (fEdep(detector, ring, sector, strip).fEdep != 0)
328         AliDebug(1, Form("Double hit in %d%c(%d,%d)", 
329                          detector, ring, sector, strip));
330       
331       fEdep(detector, ring, sector, strip).fEdep  += edep;
332       fEdep(detector, ring, sector, strip).fN     += 1;
333       // Add this to the energy deposited for this strip
334     }  // hit loop
335   } // track loop
336 }
337
338 //____________________________________________________________________
339 void
340 AliFMDBaseDigitizer::DigitizeHits(AliFMD* fmd) const
341 {
342   // For the stored energy contributions in the cache (fEdep), convert
343   // the energy signal to ADC counts, and store the created digit in
344   // the digits array (AliFMD::fDigits)
345   //
346   AliFMDGeometry* geometry = AliFMDGeometry::Instance();
347   
348   TArrayI counts(3);
349   for (UShort_t detector=1; detector <= 3; detector++) {
350     // Get pointer to subdetector 
351     AliFMDDetector* det = geometry->GetDetector(detector);
352     if (!det) continue;
353     for (UShort_t ringi = 0; ringi <= 1; ringi++) {
354       // Get pointer to Ring
355       AliFMDRing* r = det->GetRing((ringi == 0 ? 'I' : 'O'));
356       if (!r) continue;
357       
358       // Get number of sectors 
359       UShort_t nSectors = UShort_t(360. / r->GetTheta());
360       // Loop over the number of sectors 
361       for (UShort_t sector = 0; sector < nSectors; sector++) {
362         // Get number of strips 
363         UShort_t nStrips = r->GetNStrips();
364         // Loop over the stips 
365         Float_t last = 0;
366         for (UShort_t strip = 0; strip < nStrips; strip++) {
367           // Reset the counter array to the invalid value -1 
368           counts.Reset(-1);
369           // Reset the last `ADC' value when we've get to the end of a
370           // VA1_ALICE channel. 
371           if (strip % 128 == 0) last = 0;
372           
373           Float_t edep = fEdep(detector, r->GetId(), sector, strip).fEdep;
374           ConvertToCount(edep, last, r->GetSiThickness(), 
375                          geometry->GetSiDensity(), counts);
376           last = edep;
377           AddDigit(fmd, detector, r->GetId(), sector, strip, 
378                    edep, UShort_t(counts[0]), 
379                    Short_t(counts[1]), Short_t(counts[2]));
380 #if 0
381           // This checks if the digit created will give the `right'
382           // number of particles when reconstructed, using a naiive
383           // approach.  It's here only as a quality check - nothing
384           // else. 
385           CheckDigit(fEdep(detector, r->GetId(), sector, strip).fEdep,
386                      fEdep(detector, r->GetId(), sector, strip).fN,
387                      counts);
388 #endif
389         } // Strip
390       } // Sector 
391     } // Ring 
392   } // Detector 
393 }
394
395 //____________________________________________________________________
396 void
397 AliFMDBaseDigitizer::ConvertToCount(Float_t   edep, 
398                                     Float_t   last,
399                                     Float_t   siThickness, 
400                                     Float_t   siDensity, 
401                                     TArrayI&  counts) const
402 {
403   // Convert the total energy deposited to a (set of) ADC count(s). 
404   // 
405   // This is done by 
406   // 
407   //               Energy_Deposited      ALTRO_Channel_Size
408   //    ADC = -------------------------- ------------------- + pedestal
409   //          Energy_Deposition_Of_1_MIP VA1_ALICE_MIP_Range
410   //
411   //               Energy_Deposited             fAltroChannelSize
412   //        = --------------------------------- ----------------- + pedestal 
413   //          1.664 * Si_Thickness * Si_Density   fVA1MipRange   
414   //          
415   // 
416   //        = Energy_Deposited * ConversionFactor + pedestal
417   // 
418   // However, this is modified by the response function of the
419   // VA1_ALICE pre-amp. chip in case we are doing oversampling of the
420   // VA1_ALICE output. 
421   // 
422   // In that case, we get N=fSampleRate values of the ADC, and the
423   // `EnergyDeposited' is a function of which sample where are
424   // calculating the ADC for 
425   // 
426   //     ADC_i = f(EnergyDeposited, i/N, Last) * ConversionFactor + pedestal 
427   // 
428   // where Last is the Energy deposited in the previous strip. 
429   // 
430   // Here, f is the shaping function of the VA1_ALICE.   This is given
431   // by 
432   //                       
433   //                    |   (E - l) * (1 - exp(-B * t) + l   if E > l
434   //       f(E, t, l) = <
435   //                    |   (l - E) * exp(-B * t) + E        otherwise
436   //                       
437   // 
438   //                  = E + (l - E) * ext(-B * t)
439   // 
440   const Float_t mipEnergy = 1.664 * siThickness * siDensity;
441   const Float_t convf     = (1 / mipEnergy * Float_t(fAltroChannelSize) 
442                              / fVA1MipRange);
443   UShort_t ped            = MakePedestal();
444
445   // In case we don't oversample, just return the end value. 
446   if (fSampleRate == 1) {
447     counts[0] = UShort_t(TMath::Min(edep * convf + ped, 
448                                     Float_t(fAltroChannelSize)));
449     return;
450   }
451
452   // Create a pedestal 
453   Int_t   n = fSampleRate;
454   Float_t B = fShapingTime;
455   for (Ssiz_t i = 0; i < n;  i++) {
456     Float_t t = Float_t(i) / n;
457     Float_t s = edep + (last - edep) * TMath::Exp(-B * t);
458     counts[i] = UShort_t(TMath::Min(s * convf + ped, 
459                                     Float_t(fAltroChannelSize))); 
460   }
461 }
462
463
464 //====================================================================
465 ClassImp(AliFMDDigitizer)
466
467 //____________________________________________________________________
468 AliFMDDigitizer::AliFMDDigitizer()  
469   : AliFMDBaseDigitizer()
470 {
471   // Default ctor - don't use it
472 }
473
474 //____________________________________________________________________
475 AliFMDDigitizer::AliFMDDigitizer(AliRunDigitizer* manager) 
476   : AliFMDBaseDigitizer(manager)
477 {
478   // Normal CTOR
479   AliDebug(1," processed");
480   SetPedestal();
481 }
482
483 //____________________________________________________________________
484 void
485 AliFMDDigitizer::Exec(Option_t*) 
486 {
487   // Get the output manager 
488   TString outFolder(fManager->GetOutputFolderName());
489   AliRunLoader* out = 
490     AliRunLoader::GetRunLoader(outFolder.Data());
491   // Get the FMD output manager 
492   AliLoader* outFMD = out->GetLoader("FMDLoader");
493
494   // Get the input loader 
495   TString inFolder(fManager->GetInputFolderName(0));
496   fRunLoader = 
497     AliRunLoader::GetRunLoader(inFolder.Data());
498   if (!fRunLoader) {
499     AliError("Can not find Run Loader for input stream 0");
500     return;
501   }
502   // Get the AliRun object 
503   if (!fRunLoader->GetAliRun()) fRunLoader->LoadgAlice();
504
505   // Get the AliFMD object 
506   AliFMD* fmd = static_cast<AliFMD*>(fRunLoader->GetAliRun()->GetDetector("FMD"));
507   if (!fmd) {
508     AliError("Can not get FMD from gAlice");
509     return;
510   }
511
512   Int_t nFiles= fManager->GetNinputs();
513   for (Int_t inputFile = 0; inputFile < nFiles; inputFile++) {
514     AliDebug(1,Form(" Digitizing event number %d",
515                     fManager->GetOutputEventNr()));
516     // Get the current loader 
517     fRunLoader = 
518       AliRunLoader::GetRunLoader(fManager->GetInputFolderName(inputFile));
519     if (!fRunLoader) Fatal("Exec", "no run loader");
520     // Cache contriutions 
521     SumContributions(fmd);
522   }
523   // Digitize the event 
524   DigitizeHits(fmd);
525
526   // Load digits from the tree 
527   outFMD->LoadDigits("update");
528   // Get the tree of digits 
529   TTree* digitTree = outFMD->TreeD();
530   if (!digitTree) {
531     outFMD->MakeTree("D");
532     digitTree = outFMD->TreeD();
533   }
534   digitTree->Reset();
535   // Make a branch in the tree 
536   TClonesArray* digits = fmd->Digits();
537   fmd->MakeBranchInTree(digitTree, fmd->GetName(), &(digits), 4000, 0);
538   // TBranch* digitBranch = digitTree->GetBranch(fmd->GetName());
539   // Fill the tree 
540   digitTree->Fill();
541   
542   // Write the digits to disk 
543   outFMD->WriteDigits("OVERWRITE");
544   outFMD->UnloadHits();
545   outFMD->UnloadDigits();
546
547   // Reset the digits in the AliFMD object 
548   fmd->ResetDigits();
549 }
550
551
552 //____________________________________________________________________
553 UShort_t
554 AliFMDDigitizer::MakePedestal() const 
555 {
556   return UShort_t(TMath::Max(gRandom->Gaus(fPedestal, fPedestalWidth), 0.));
557 }
558
559 //____________________________________________________________________
560 void
561 AliFMDDigitizer::AddDigit(AliFMD*  fmd,
562                           UShort_t detector, 
563                           Char_t   ring,
564                           UShort_t sector, 
565                           UShort_t strip, 
566                           Float_t  /* edep */, 
567                           UShort_t count1, 
568                           Short_t  count2, 
569                           Short_t  count3) const
570 {
571   fmd->AddDigitByFields(detector, ring, sector, strip, count1, count2, count3);
572 }
573
574 //____________________________________________________________________
575 void
576 AliFMDDigitizer::CheckDigit(Float_t         /* edep */, 
577                             UShort_t        nhits,
578                             const TArrayI&  counts) 
579 {
580   Int_t integral = counts[0];
581   if (counts[1] >= 0) integral += counts[1];
582   if (counts[2] >= 0) integral += counts[2];
583   integral -= Int_t(fPedestal + 2 * fPedestalWidth);
584   if (integral < 0) integral = 0;
585   
586   Float_t convf = Float_t(fVA1MipRange) / fAltroChannelSize;
587   Float_t mips  = integral * convf;
588   if (mips > Float_t(nhits) + .5 || mips < Float_t(nhits) - .5) 
589     Warning("CheckDigit", "Digit -> %4.2f MIPS != %d +/- .5 hits", 
590             mips, nhits);
591 }
592
593 //====================================================================
594 ClassImp(AliFMDSDigitizer)
595
596 //____________________________________________________________________
597 AliFMDSDigitizer::AliFMDSDigitizer()  
598 {
599   // Default ctor - don't use it
600 }
601
602 //____________________________________________________________________
603 AliFMDSDigitizer::AliFMDSDigitizer(const Char_t* headerFile, 
604                                    const Char_t* /* sdigfile */)
605   : AliFMDBaseDigitizer("FMDSDigitizer", "FMD SDigitizer")
606 {
607   // Normal CTOR
608   AliDebug(1," processed");
609
610   fRunLoader = AliRunLoader::GetRunLoader(); // Open(headerFile);
611   if (!fRunLoader) 
612     Fatal("AliFMDSDigitizer", "cannot open session, header file '%s'",
613           headerFile);
614   AliLoader* loader = fRunLoader->GetLoader("FMDLoader");
615   if (!loader) 
616     Fatal("AliFMDSDigitizer", "cannot find FMD loader in specified event");
617
618   // Add task to tasks folder 
619   loader->PostSDigitizer(this);
620
621 }
622
623 //____________________________________________________________________
624 AliFMDSDigitizer::~AliFMDSDigitizer() 
625 {
626   AliLoader* loader = fRunLoader->GetLoader("FMDLoader");
627   loader->CleanSDigitizer();
628 }
629
630 //____________________________________________________________________
631 void
632 AliFMDSDigitizer::Exec(Option_t*) 
633 {
634   // Get the output manager 
635   if (!fRunLoader) {
636     Error("Exec", "Run loader is not set");
637     return;
638   }
639   if (!fRunLoader->GetAliRun()) fRunLoader->LoadgAlice();
640   if (!fRunLoader->TreeE())     fRunLoader->LoadHeader();
641   
642   AliLoader* fmdLoader = fRunLoader->GetLoader("FMDLoader");
643   if (!fmdLoader) Fatal("Exec", "no FMD loader");
644   
645   // Get the AliFMD object 
646   AliFMD* fmd = 
647     static_cast<AliFMD*>(fRunLoader->GetAliRun()->GetDetector("FMD"));
648   if (!fmd) {
649     AliError("Can not get FMD from gAlice");
650     return;
651   }
652
653   Int_t nEvents = Int_t(fRunLoader->TreeE()->GetEntries());
654   for (Int_t event = 0; event < nEvents; event++) {
655     AliDebug(1,Form(" Digitizing event number %d", event));
656     // Get the current loader 
657     fRunLoader->GetEvent(event);
658
659     if (!fmdLoader->TreeS()) fmdLoader->MakeTree("S");
660     // Make a branch
661     fmd->MakeBranch("S");
662     
663     // Cache contriutions 
664     SumContributions(fmd);
665
666     // Digitize the event 
667     DigitizeHits(fmd);
668
669     fmdLoader->TreeS()->Reset();
670     fmdLoader->TreeS()->Fill();
671     fmdLoader->WriteSDigits("OVERWRITE");
672   }
673 }
674
675 //____________________________________________________________________
676 void
677 AliFMDSDigitizer::AddDigit(AliFMD*  fmd,
678                            UShort_t detector, 
679                            Char_t   ring,
680                            UShort_t sector, 
681                            UShort_t strip, 
682                            Float_t  edep, 
683                            UShort_t count1, 
684                            Short_t  count2, 
685                            Short_t  count3) const
686 {
687   fmd->AddSDigitByFields(detector, ring, sector, strip, edep, 
688                          count1, count2, count3); 
689 }
690
691
692
693 //____________________________________________________________________
694 //
695 // EOF
696 // 
697
698
699
700