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