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 /**
  17  * @file    AliFMDDigitizer.cxx
  18  *
  19  * @author  Christian Holm Christensen <cholm@nbi.dk>
  20  * @date    Mon Mar 27 12:38:26 2006
  21  * @brief   FMD Digitizers implementation
  22  *
  23  * @ingroup FMD_sim
  24  */
  25 //////////////////////////////////////////////////////////////////////////////
  26 //
  27 //  This class contains the procedures simulation ADC  signal for the
  28 //  Forward Multiplicity detector  : SDigits->Digits
  29 //
  30 //  Digits consists of
  31 //   - Detector #
  32 //   - Ring ID
  33 //   - Sector #
  34 //   - Strip #
  35 //   - ADC count in this channel
  36 //
  37 //  Digits consists of
  38 //   - Detector #
  39 //   - Ring ID
  40 //   - Sector #
  41 //   - Strip #
  42 //   - Total energy deposited in the strip
  43 //   - ADC count in this channel
  44 //
  45 // As the Digits and SDigits have so much in common, the classes
  46 // AliFMDDigitizer and AliFMDSDigitizer are implemented via a base
  47 // class AliFMDBaseDigitizer.
  48 //
  49 //                 +---------------------+
  50 //                 | AliFMDBaseDigitizer |
  51 //                 +---------------------+
  52 //                           ^
  53 //                           |
  54 //                +----------+---------+
  55 //                |                    |
  56 //      +-----------------+     +------------------+
  57 //      | AliFMDDigitizer |     | AliFMDSDigitizer |
  58 //      +-----------------+     +------------------+
  59 //                |
  60 //     +-------------------+
  61 //     | AliFMDSSDigitizer |
  62 //     +-------------------+
  63 //
  64 // These classes has several paramters:
  65 //
  66 //     fPedestal
  67 //     fPedestalWidth
  68 //         (Only AliFMDDigitizer)
  69 //         Mean and width of the pedestal.  The pedestal is simulated
  70 //         by a Guassian, but derived classes my override MakePedestal
  71 //         to simulate it differently (or pick it up from a database).
  72 //
  73 //     fVA1MipRange
  74 //         The dymamic MIP range of the VA1_ALICE pre-amplifier chip
  75 //
  76 //     fAltroChannelSize
  77 //         The largest number plus one that can be stored in one
  78 //         channel in one time step in the ALTRO ADC chip.
  79 //
  80 //     fSampleRate
  81 //         How many times the ALTRO ADC chip samples the VA1_ALICE
  82 //         pre-amplifier signal.   The VA1_ALICE chip is read-out at
  83 //         10MHz, while it's possible to drive the ALTRO chip at
  84 //         25MHz.  That means, that the ALTRO chip can have time to
  85 //         sample each VA1_ALICE signal up to 2 times.  Although it's
  86 //         not certain this feature will be used in the production,
  87 //         we'd like have the option, and so it should be reflected in
  88 //         the code.
  89 //
  90 //
  91 // The shaping function of the VA1_ALICE is generally given by
  92 //
  93 //      f(x) = A(1 - exp(-Bx))
  94 //
  95 // where A is the total charge collected in the pre-amp., and B is a
  96 // paramter that depends on the shaping time of the VA1_ALICE circut.
  97 //
  98 // When simulating the shaping function of the VA1_ALICe
  99 // pre-amp. chip, we have to take into account, that the shaping
 100 // function depends on the previous value of read from the pre-amp.
 101 //
 102 // That results in the following algorithm:
 103 //
 104 //    last = 0;
 105 //    FOR charge IN pre-amp. charge train DO
 106 //      IF last < charge THEN
 107 //        f(t) = (charge - last) * (1 - exp(-B * t)) + last
 108 //      ELSE
 109 //        f(t) = (last - charge) * exp(-B * t) + charge)
 110 //      ENDIF
 111 //      FOR i IN # samples DO
 112 //        adc_i = f(i / (# samples))
 113 //      DONE
 114 //      last = charge
 115 //   DONE
 116 //
 117 // Here,
 118 //
 119 //   pre-amp. charge train
 120 //       is a series of 128 charges read from the VA1_ALICE chip
 121 //
 122 //   # samples
 123 //       is the number of times the ALTRO ADC samples each of the 128
 124 //       charges from the pre-amp.
 125 //
 126 // Where Q is the total charge collected by the VA1_ALICE
 127 // pre-amplifier.   Q is then given by
 128 //
 129 //           E S
 130 //      Q =  - -
 131 //           e R
 132 //
 133 // where E is the total energy deposited in a silicon strip, R is the
 134 // dynamic range of the VA1_ALICE pre-amp (fVA1MipRange), e is the
 135 // energy deposited by a single MIP, and S ALTRO channel size in each
 136 // time step (fAltroChannelSize).
 137 //
 138 // The energy deposited per MIP is given by
 139 //
 140 //      e = M * rho * w
 141 //
 142 // where M is the universal number 1.664, rho is the density of
 143 // silicon, and w is the depth of the silicon sensor.
 144 //
 145 // The final ADC count is given by
 146 //
 147 //      C' = C + P
 148 //
 149 // where P is the (randomized) pedestal (see MakePedestal)
 150 //
 151 // This class uses the class template AliFMDMap<Type> to make an
 152 // internal cache of the energy deposted of the hits.  The class
 153 // template is instantasized as
 154 //
 155 //  typedef AliFMDMap<std::pair<Float_t, UShort_t> > AliFMDEdepMap;
 156 //
 157 // The first member of the values is the summed energy deposition in a
 158 // given strip, while the second member of the values is the number of
 159 // hits in a given strip.  Using the second member, it's possible to
 160 // do some checks on just how many times a strip got hit, and what
 161 // kind of error we get in our reconstructed hits.  Note, that this
 162 // information is currently not written to the digits tree.  I think a
 163 // QA (Quality Assurance) digit tree is better suited for that task.
 164 // However, the information is there to be used in the future.
 165 //
 166 //
 167 // Latest changes by Christian Holm Christensen
 168 //
 169 //////////////////////////////////////////////////////////////////////////////
 170
 171 //      /1
 172 //      |           A(-1 + B + exp(-B))
 173 //      | f(x) dx = ------------------- = 1
 174 //      |                    B
 175 //      / 0
 176 //
 177 // and B is the a parameter defined by the shaping time (fShapingTime).
 178 //
 179 // Solving the above equation, for A gives
 180 //
 181 //                 B
 182 //      A = ----------------
 183 //          -1 + B + exp(-B)
 184 //
 185 // So, if we define the function g: [0,1] -> [0:1] by
 186 //
 187 //               / v
 188 //               |              Bu + exp(-Bu) - Bv - exp(-Bv)
 189 //      g(u,v) = | f(x) dx = -A -----------------------------
 190 //               |                            B
 191 //               / u
 192 //
 193 // we can evaluate the ALTRO sample of the VA1_ALICE pre-amp between
 194 // any two times (u, v), by
 195 //
 196 //
 197 //                                B         Bu + exp(-Bu) - Bv - exp(-Bv)
 198 //      C = Q g(u,v) = - Q ---------------- -----------------------------
 199 //                         -1 + B + exp(-B)              B
 200 //
 201 //               Bu + exp(-Bu) - Bv - exp(-Bv)
 202 //        = -  Q -----------------------------
 203 //                    -1 + B + exp(-B)
 204 //
 205
 206 #include <TTree.h>              // ROOT_TTree
 207 #include <TFile.h>
 208 #include "AliFMDDebug.h"        // Better debug macros
 209 #include "AliFMDDigitizer.h"    // ALIFMDSSDIGITIZER_H
 210 #include "AliFMD.h"             // ALIFMD_H
 211 #include "AliFMDSDigit.h"       // ALIFMDDIGIT_H
 212 #include "AliFMDDigit.h"        // ALIFMDDIGIT_H
 213 #include "AliFMDParameters.h"   // ALIFMDPARAMETERS_H
 214 #include <AliDigitizationInput.h>       // ALIRUNDIGITIZER_H
 215 #include <AliRun.h>             // ALIRUN_H
 216 #include <AliLoader.h>          // ALILOADER_H
 217 #include <AliRunLoader.h>       // ALIRUNLOADER_H
 218
 219 //====================================================================
 220 ClassImp(AliFMDDigitizer)
 221 #if 0
 222 ;
 223 #endif
 224
 225 //____________________________________________________________________
 226 Bool_t
 227 AliFMDDigitizer::Init()
 228 {
 229   //
 230   // Initialisation
 231   //
 232   if (!AliFMDBaseDigitizer::Init()) return kFALSE;
 233
 234 #if 0
 235   // Get the AliRun object
 236   AliRun* run = fRunLoader->GetAliRun();
 237   if (!run) {
 238     AliWarning("Loading gAlice");
 239     fRunLoader->LoadgAlice();
 240     if (!run) {
 241       AliError("Can not get Run from Run Loader");
 242       return kFALSE;
 243     }
 244   }
 245
 246   // Get the AliFMD object
 247   fFMD = static_cast<AliFMD*>(run->GetDetector("FMD"));
 248   if (!fFMD) {
 249     AliError("Can not get FMD from gAlice");
 250     return kFALSE;
 251   }
 252 #endif
 253   return kTRUE;
 254 }
 255
 256
 257 //____________________________________________________________________
 258 void
 259 AliFMDDigitizer::Digitize(Option_t*)
 260 {
 261   //
 262   // Execute this digitizer.
 263   // This member function will be called once per event by the passed
 264   // AliDigitizationInput* digInput object.
 265   //
 266   // Parameters:
 267   //    options Not used
 268   //
 269   if (!fDigInput) {
 270     AliError("No digitisation input defined");
 271     return;
 272   }
 273
 274   // Clear array of deposited energies
 275   fEdep.Reset();
 276
 277   AliRunLoader* runLoader = 0;
 278   if (!gAlice) {
 279     TString folderName(fDigInput->GetInputFolderName(0));
 280     runLoader = AliRunLoader::GetRunLoader(folderName.Data());
 281     if (!runLoader) {
 282       AliError(Form("Failed at getting run loader from %s",
 283                     folderName.Data()));
 284       return;
 285     }
 286     if (!runLoader->GetAliRun()) runLoader->LoadgAlice();
 287     runLoader->GetAliRun();
 288   }
 289   if (!gAlice) {
 290     AliError("Can not get Run from Run Loader");
 291     return;
 292   }
 293
 294   // Get the AliFMD object
 295   fFMD = static_cast<AliFMD*>(gAlice->GetDetector("FMD"));
 296   if (!fFMD) {
 297     AliError("Can not get FMD from gAlice");
 298     return;
 299   }
 300
 301
 302   // Loop over input files
 303   Int_t nFiles= fDigInput->GetNinputs();
 304   AliFMDDebug(1, (" Digitizing event number %d, got %d inputs",
 305                   fDigInput->GetOutputEventNr(), nFiles));
 306   for (Int_t inputFile = 0; inputFile < nFiles; inputFile++) {
 307     AliFMDDebug(5, ("Now reading input # %d", inputFile));
 308     // Get the current loader
 309     AliRunLoader* currentLoader =
 310       AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(inputFile));
 311     if (!currentLoader) {
 312       Error("Exec", "no run loader for input file # %d", inputFile);
 313       continue;
 314     }
 315
 316     // Cache contriutions
 317     AliFMDDebug(5, ("Now summing the contributions from input # %d",inputFile));
 318
 319     // Get the FMD loader
 320     AliLoader* inFMD = currentLoader->GetLoader("FMDLoader");
 321     // And load the summable digits
 322     inFMD->LoadSDigits("READ");
 323
 324     // Get the tree of summable digits
 325     TTree* sdigitsTree = inFMD->TreeS();
 326     if (!sdigitsTree)  {
 327       AliError("No sdigit tree from input");
 328       continue;
 329     }
 330     if (AliLog::GetDebugLevel("FMD","") >= 10) {
 331       TFile* file = sdigitsTree->GetCurrentFile();
 332       if (!file) {
 333         AliWarning("Input tree has no file!");
 334       }
 335       else {
 336         AliFMDDebug(10, ("Input tree file %s content:", file->GetName()));
 337         file->ls();
 338       }
 339       // AliFMDDebug(5, ("Input tree %s file structure:",
 340       //                 sdigitsTree->GetName()));
 341       // sdigitsTree->Print();
 342     }
 343
 344     // Get the FMD branch
 345     TBranch* sdigitsBranch = sdigitsTree->GetBranch("FMD");
 346     if (!sdigitsBranch) {
 347       AliError("Failed to get sdigit branch");
 348       return;
 349     }
 350
 351     // Set the branch addresses
 352     fFMD->SetTreeAddress();
 353
 354     // Sum contributions from the sdigits
 355     AliFMDDebug(3, ("Will now sum contributions from SDigits"));
 356     SumContributions(sdigitsBranch);
 357
 358     // Unload the sdigits
 359     inFMD->UnloadSDigits();
 360   }
 361
 362   TString       outFolder(fDigInput->GetOutputFolderName());
 363   AliRunLoader* out    = AliRunLoader::GetRunLoader(outFolder.Data());
 364   AliLoader*    outFMD = out->GetLoader("FMDLoader");
 365   if (!outFMD) {
 366     AliError("Cannot get the FMDLoader output folder");
 367     return;
 368   }
 369   TTree* outTree = MakeOutputTree(outFMD);
 370   if (!outTree) {
 371     AliError("Failed to get output tree");
 372     return;
 373   }
 374   // Set the branch address
 375   fFMD->SetTreeAddress();
 376
 377   // And digitize the cached data
 378   DigitizeHits();
 379
 380   // Fill the tree
 381   Int_t write = outTree->Fill();
 382   AliFMDDebug(5, ("Wrote %d bytes to digit tree", write));
 383
 384   // Store the digits
 385   StoreDigits(outFMD);
 386 }
 387
 388 //____________________________________________________________________
 389 void
 390 AliFMDDigitizer::SumContributions(TBranch* sdigitsBranch)
 391 {
 392   //
 393   // Sum contributions from SDigits
 394   //
 395   // Parameters:
 396   //    sdigitsBranch Branch of SDigit data
 397   //
 398   AliFMDDebug(3, ("Runnin our version of SumContributions"));
 399
 400   // Get a list of hits from the FMD manager
 401   TClonesArray *fmdSDigits = fFMD->SDigits();
 402
 403   // Get number of entries in the tree
 404   Int_t nevents  = Int_t(sdigitsBranch->GetEntries());
 405
 406   Int_t read = 0;
 407   // Loop over the events in the
 408   for (Int_t event = 0; event < nevents; event++)  {
 409     // Read in entry number `event'
 410     read += sdigitsBranch->GetEntry(event);
 411
 412     // Get the number of sdigits
 413     Int_t nsdigits = fmdSDigits->GetEntries ();
 414     AliFMDDebug(3, ("Got %5d SDigits", nsdigits));
 415     for (Int_t sdigit = 0; sdigit < nsdigits; sdigit++) {
 416       // Get the sdigit number `sdigit'
 417       AliFMDSDigit* fmdSDigit =
 418         static_cast<AliFMDSDigit*>(fmdSDigits->UncheckedAt(sdigit));
 419
 420       AliFMDDebug(5, ("Adding contribution of %d tracks",
 421                       fmdSDigit->GetNTrack()));
 422       AliFMDDebug(15, ("Contrib from FMD%d%c[%2d,%3d] (%s) from track %d",
 423                        fmdSDigit->Detector(),
 424                        fmdSDigit->Ring(),
 425                        fmdSDigit->Sector(),
 426                        fmdSDigit->Strip(),
 427                        fmdSDigit->GetName(),
 428                        fmdSDigit->GetTrack(0)));
 429
 430       // Extract parameters
 431       AddContribution(fmdSDigit->Detector(),
 432                       fmdSDigit->Ring(),
 433                       fmdSDigit->Sector(),
 434                       fmdSDigit->Strip(),
 435                       fmdSDigit->Edep(),
 436                       kTRUE,
 437                       fmdSDigit->GetNTrack(),
 438                       fmdSDigit->GetTracks());
 439     }  // sdigit loop
 440   } // event loop
 441
 442
 443   AliFMDDebug(3, ("Size of cache: %d bytes, read %d bytes",
 444                   int(sizeof(fEdep)), read));
 445 }
 446
 447 //____________________________________________________________________
 448 //
 449 // EOF
 450 //
 451
 452
 453
 454