AliRoot Core  edcc906 (edcc906)
AliFMDDigitizer.cxx
Go to the documentation of this file.
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$ */
25 //
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 //
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
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");
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
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 =
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(); // RS: this sets address to 1st streem always
353 
354  // Sum contributions from the sdigits
355  AliFMDDebug(3, ("Will now sum contributions from SDigits"));
356  SumContributions(sdigitsBranch,fDigInput->GetMask(inputFile));
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, int offset)
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  sdigitsBranch->SetAddress(&fmdSDigits); //RS Impose reading from supplied branch
404 
405  // Get number of entries in the tree
406  Int_t nevents = Int_t(sdigitsBranch->GetEntries());
407 
408  Int_t read = 0;
409  // Loop over the events in the
410  for (Int_t event = 0; event < nevents; event++) {
411  // Read in entry number `event'
412  read += sdigitsBranch->GetEntry(event);
413 
414  // Get the number of sdigits
415  Int_t nsdigits = fmdSDigits->GetEntries ();
416  AliFMDDebug(3, ("Got %5d SDigits", nsdigits));
417  for (Int_t sdigit = 0; sdigit < nsdigits; sdigit++) {
418  // Get the sdigit number `sdigit'
419  AliFMDSDigit* fmdSDigit =
420  static_cast<AliFMDSDigit*>(fmdSDigits->UncheckedAt(sdigit));
421 
422  AliFMDDebug(5, ("Adding contribution of %d tracks",
423  fmdSDigit->GetNTrack()));
424  AliFMDDebug(15, ("Contrib from FMD%d%c[%2d,%3d] (%s) from track %d",
425  fmdSDigit->Detector(),
426  fmdSDigit->Ring(),
427  fmdSDigit->Sector(),
428  fmdSDigit->Strip(),
429  fmdSDigit->GetName(),
430  fmdSDigit->GetTrack(0)));
431 
432  // Extract parameters
433  AddContribution(fmdSDigit->Detector(),
434  fmdSDigit->Ring(),
435  fmdSDigit->Sector(),
436  fmdSDigit->Strip(),
437  fmdSDigit->Edep(),
438  kTRUE,
439  fmdSDigit->GetNTrack(),
440  fmdSDigit->GetTracks(),
441  offset);
442  } // sdigit loop
443  } // event loop
444 
445 
446  AliFMDDebug(3, ("Size of cache: %d bytes, read %d bytes",
447  int(sizeof(fEdep)), read));
448 }
449 
450 //____________________________________________________________________
451 //
452 // EOF
453 //
454 
455 
456 
457 
Int_t LoadgAlice()
virtual void StoreDigits(const AliLoader *loader)
Definition: AliRun.h:27
AliDigitizationInput * fDigInput
Definition: AliDigitizer.h:43
const char * GetOutputFolderName()
Digits for the FMD.
Int_t GetMask(Int_t i) const
UShort_t Sector() const
Char_t Ring() const
Manager of FMD parameters.
virtual Int_t * GetTracks()
Definition: AliDigit.h:20
virtual Bool_t Init()
AliLoader * GetLoader(const char *detname) const
Int_t GetOutputEventNr() const
Forward Multiplicity Detector based on Silicon wafers. This class is the driver for especially simula...
Definition: AliFMD.h:306
TTree * TreeS() const
Definition: AliLoader.h:85
virtual TTree * MakeOutputTree(AliLoader *loader)
const TString & GetInputFolderName(Int_t i) const
#define AliWarning(message)
Definition: AliLog.h:541
FMD Digitizers declaration.
Concrete digitizer to make digits from hits. See also AliFMDBaseDigitizer documentation.
Int_t LoadSDigits(Option_t *opt="")
Definition: AliLoader.h:101
virtual void Digitize(Option_t *option="")
Float_t Edep() const
Definition: AliFMDSDigit.h:91
virtual void AddContribution(UShort_t detector, Char_t ring, UShort_t sector, UShort_t strip, Float_t edep, Bool_t isPrimary, Int_t nTrackno, Int_t *tracknos, Int_t offset=0)
virtual Int_t GetTrack(Int_t i) const
Definition: AliDigit.h:21
virtual void DigitizeHits() const
virtual TClonesArray * SDigits()
Definition: AliFMD.h:427
Declaration of AliFMD detector driver.
virtual void SetTreeAddress()
Definition: AliFMD.cxx:589
UShort_t Strip() const
#define AliFMDDebug(N, A)
Definition: AliFMDDebug.h:39
AliRun * gAlice
Definition: AliRun.cxx:62
const char * GetName() const
Digits for the FMD.
virtual void Reset()
void SumContributions(TBranch *sdigitsBranch, int offset=0)
class for summable digits
Definition: AliFMDSDigit.h:27
void UnloadSDigits() const
Definition: AliLoader.h:130
UShort_t Detector() const
static Int_t GetDebugLevel(const char *module, const char *className)
Definition: AliLog.cxx:843
AliDetector * GetDetector(const char *name) const
Definition: AliRun.cxx:200
#define AliError(message)
Definition: AliLog.h:591
static AliRunLoader * GetRunLoader(const char *eventfoldername)
AliRun * GetAliRun() const
AliFMDEdepMap fEdep
Run loader.
UShort_t GetNTrack() const