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