]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TRD/AliTRDdigitizer.cxx
Automatically find the library for the fpu exception
[u/mrichter/AliRoot.git] / TRD / AliTRDdigitizer.cxx
CommitLineData
f7336fa3 1/**************************************************************************
2 * Copyright(c) 1998-1999, 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/*
17$Log$
94de3818 18Revision 1.8 2000/06/09 11:10:07 cblume
19Compiler warnings and coding conventions, next round
20
dd9a6ee3 21Revision 1.7 2000/06/08 18:32:58 cblume
22Make code compliant to coding conventions
23
8230f242 24Revision 1.6 2000/06/07 16:27:32 cblume
25Try to remove compiler warnings on Sun and HP
26
9d0b222b 27Revision 1.5 2000/05/09 16:38:57 cblume
28Removed PadResponse(). Merge problem
29
f0a7bf65 30Revision 1.4 2000/05/08 15:53:45 cblume
31Resolved merge conflict
32
da581aea 33Revision 1.3 2000/04/28 14:49:27 cblume
34Only one declaration of iDict in MakeDigits()
35
36Revision 1.1.4.1 2000/05/08 14:42:04 cblume
37Introduced AliTRDdigitsManager
28329a48 38
1befd3b2 39Revision 1.1 2000/02/28 19:00:13 cblume
40Add new TRD classes
41
f7336fa3 42*/
43
44///////////////////////////////////////////////////////////////////////////////
45// //
46// Creates and handles digits from TRD hits //
47// //
48// The following effects are included: //
49// - Diffusion //
50// - ExB effects //
51// - Gas gain including fluctuations //
52// - Pad-response (simple Gaussian approximation) //
53// - Electronics noise //
54// - Electronics gain //
55// - Digitization //
56// - ADC threshold //
57// The corresponding parameter can be adjusted via the various //
58// Set-functions. If these parameters are not explicitly set, default //
59// values are used (see Init-function). //
60// To produce digits from a root-file with TRD-hits use the //
61// slowDigitsCreate.C macro. //
62// //
63///////////////////////////////////////////////////////////////////////////////
64
65#include <TMath.h>
66#include <TVector.h>
67#include <TRandom.h>
94de3818 68#include <TROOT.h>
69#include <TTree.h>
f7336fa3 70
71#include "AliTRD.h"
72#include "AliTRDdigitizer.h"
da581aea 73#include "AliTRDdataArrayI.h"
74#include "AliTRDdataArrayF.h"
75#include "AliTRDdigitsManager.h"
f7336fa3 76
77ClassImp(AliTRDdigitizer)
78
79//_____________________________________________________________________________
80AliTRDdigitizer::AliTRDdigitizer():TNamed()
81{
82 //
83 // AliTRDdigitizer default constructor
84 //
85
da581aea 86 fInputFile = NULL;
87 fDigits = NULL;
88 fTRD = NULL;
89 fGeo = NULL;
90 fPRF = NULL;
91
92 fEvent = 0;
93 fGasGain = 0.0;
94 fNoise = 0.0;
95 fChipGain = 0.0;
96 fADCoutRange = 0.0;
97 fADCinRange = 0.0;
98 fADCthreshold = 0;
99 fDiffusionOn = 0;
100 fDiffusionT = 0.0;
101 fDiffusionL = 0.0;
102 fElAttachOn = 0;
103 fElAttachProp = 0.0;
104 fExBOn = 0;
105 fLorentzAngle = 0.0;
f7336fa3 106
107}
108
109//_____________________________________________________________________________
110AliTRDdigitizer::AliTRDdigitizer(const Text_t *name, const Text_t *title)
111 :TNamed(name,title)
112{
113 //
114 // AliTRDdigitizer default constructor
115 //
116
da581aea 117 fInputFile = NULL;
118 fDigits = NULL;
119 fTRD = NULL;
120 fGeo = NULL;
f7336fa3 121
da581aea 122 fEvent = 0;
f7336fa3 123
124 Init();
125
126}
127
8230f242 128//_____________________________________________________________________________
dd9a6ee3 129AliTRDdigitizer::AliTRDdigitizer(const AliTRDdigitizer &d)
8230f242 130{
131 //
132 // AliTRDdigitizer copy constructor
133 //
134
dd9a6ee3 135 ((AliTRDdigitizer &) d).Copy(*this);
8230f242 136
137}
138
f7336fa3 139//_____________________________________________________________________________
140AliTRDdigitizer::~AliTRDdigitizer()
141{
8230f242 142 //
143 // AliTRDdigitizer destructor
144 //
f7336fa3 145
146 if (fInputFile) {
147 fInputFile->Close();
148 delete fInputFile;
149 }
150
da581aea 151 if (fDigits) {
152 delete fDigits;
f7336fa3 153 }
154
da581aea 155 if (fPRF) delete fPRF;
f7336fa3 156
157}
158
8230f242 159//_____________________________________________________________________________
dd9a6ee3 160AliTRDdigitizer &AliTRDdigitizer::operator=(const AliTRDdigitizer &d)
161{
162 //
163 // Assignment operator
164 //
165
166 if (this != &d) ((AliTRDdigitizer &) d).Copy(*this);
167 return *this;
168
169}
170
171//_____________________________________________________________________________
172void AliTRDdigitizer::Copy(TObject &d)
8230f242 173{
174 //
175 // Copy function
176 //
177
dd9a6ee3 178 ((AliTRDdigitizer &) d).fInputFile = NULL;
179 ((AliTRDdigitizer &) d).fDigits = NULL;
180 ((AliTRDdigitizer &) d).fTRD = NULL;
181 ((AliTRDdigitizer &) d).fGeo = NULL;
182
183 ((AliTRDdigitizer &) d).fEvent = 0;
184
185 ((AliTRDdigitizer &) d).fGasGain = fGasGain;
186 ((AliTRDdigitizer &) d).fNoise = fNoise;
187 ((AliTRDdigitizer &) d).fChipGain = fChipGain;
188 ((AliTRDdigitizer &) d).fADCoutRange = fADCoutRange;
189 ((AliTRDdigitizer &) d).fADCinRange = fADCinRange;
190 ((AliTRDdigitizer &) d).fADCthreshold = fADCthreshold;
191 ((AliTRDdigitizer &) d).fDiffusionOn = fDiffusionOn;
192 ((AliTRDdigitizer &) d).fDiffusionT = fDiffusionT;
193 ((AliTRDdigitizer &) d).fDiffusionL = fDiffusionL;
194 ((AliTRDdigitizer &) d).fElAttachOn = fElAttachOn;
195 ((AliTRDdigitizer &) d).fElAttachProp = fElAttachProp;
196 ((AliTRDdigitizer &) d).fExBOn = fExBOn;
197 ((AliTRDdigitizer &) d).fLorentzAngle = fLorentzAngle;
198 ((AliTRDdigitizer &) d).fLorentzFactor = fLorentzFactor;
199
200 fPRF->Copy(*((AliTRDdigitizer &) d).fPRF);
8230f242 201
202}
203
f7336fa3 204//_____________________________________________________________________________
205Int_t AliTRDdigitizer::Diffusion(Float_t driftlength, Float_t *xyz)
206{
207 //
208 // Applies the diffusion smearing to the position of a single electron
209 //
210
211 Float_t driftSqrt = TMath::Sqrt(driftlength);
212 Float_t sigmaT = driftSqrt * fDiffusionT;
213 Float_t sigmaL = driftSqrt * fDiffusionL;
214 xyz[0] = gRandom->Gaus(xyz[0], sigmaL * fLorentzFactor);
215 xyz[1] = gRandom->Gaus(xyz[1], sigmaT * fLorentzFactor);
216 xyz[2] = gRandom->Gaus(xyz[2], sigmaT);
217 return 1;
218
219}
220
221//_____________________________________________________________________________
222Int_t AliTRDdigitizer::ExB(Float_t driftlength, Float_t *xyz)
223{
224 //
225 // Applies E x B effects to the position of a single electron
226 //
227
228 xyz[0] = xyz[0];
229 xyz[1] = xyz[1] + fLorentzAngle * driftlength;
230 xyz[2] = xyz[2];
231
232 return 1;
233
234}
235
236//_____________________________________________________________________________
237void AliTRDdigitizer::Init()
238{
239 //
240 // Initializes the digitization procedure with standard values
241 //
242
243 // The default parameter for the digitization
244 fGasGain = 2.0E3;
245 fNoise = 3000.;
246 fChipGain = 10.;
247 fADCoutRange = 255.;
248 fADCinRange = 2000.;
249 fADCthreshold = 1;
250
251 // Transverse and longitudinal diffusion coefficients (Xe/Isobutane)
252 fDiffusionOn = 1;
253 fDiffusionT = 0.060;
254 fDiffusionL = 0.017;
255
256 // Propability for electron attachment
257 fElAttachOn = 0;
258 fElAttachProp = 0.0;
259
260 // E x B effects
261 fExBOn = 0;
262 // omega * tau. (tau ~ 12 * 10^-12, B = 0.2T)
263 fLorentzAngle = 17.6 * 12.0 * 0.2 * 0.01;
264
da581aea 265 // The pad response function
266 fPRF = new TF1("PRF","[0]*([1]+exp(-x*x/(2.0*[2])))",-2,2);
267 fPRF->SetParameter(0, 0.8872);
268 fPRF->SetParameter(1,-0.00573);
269 fPRF->SetParameter(2, 0.454 * 0.454);
270
f7336fa3 271}
272
273//_____________________________________________________________________________
274Bool_t AliTRDdigitizer::Open(const Char_t *name, Int_t nEvent)
275{
276 //
277 // Opens a ROOT-file with TRD-hits and reads in the hit-tree
278 //
279
280 // Connect the AliRoot file containing Geometry, Kine, and Hits
281 fInputFile = (TFile*) gROOT->GetListOfFiles()->FindObject(name);
282 if (!fInputFile) {
283 printf("AliTRDdigitizer::Open -- ");
284 printf("Open the ALIROOT-file %s.\n",name);
285 fInputFile = new TFile(name,"UPDATE");
286 }
287 else {
288 printf("AliTRDdigitizer::Open -- ");
289 printf("%s is already open.\n",name);
290 }
291
da581aea 292 gAlice = (AliRun*) fInputFile->Get("gAlice");
293 if (gAlice) {
294 printf("AliTRDdigitizer::Open -- ");
295 printf("AliRun object found on file.\n");
296 }
297 else {
298 printf("AliTRDdigitizer::Open -- ");
299 printf("Could not find AliRun object.\n");
300 return kFALSE;
301 }
f7336fa3 302
303 fEvent = nEvent;
304
305 // Import the Trees for the event nEvent in the file
306 Int_t nparticles = gAlice->GetEvent(fEvent);
307 if (nparticles <= 0) {
308 printf("AliTRDdigitizer::Open -- ");
309 printf("No entries in the trees for event %d.\n",fEvent);
310 return kFALSE;
311 }
312
dd9a6ee3 313 // Get the pointer to the detector class and check for version 1
314 fTRD = (AliTRD*) gAlice->GetDetector("TRD");
315 if (fTRD->IsVersion() != 1) {
316 printf("AliTRDdigitizer::Open -- ");
317 printf("TRD must be version 1 (slow simulator).\n");
318 exit(1);
319 }
320
321 // Get the geometry
322 fGeo = fTRD->GetGeometry();
323 printf("AliTRDdigitizer::Open -- ");
324 printf("Geometry version %d\n",fGeo->IsVersion());
325
f7336fa3 326 return kTRUE;
327
328}
329
f7336fa3 330//_____________________________________________________________________________
331Bool_t AliTRDdigitizer::MakeDigits()
332{
333 //
334 // Loops through the TRD-hits and creates the digits.
335 //
336
f7336fa3 337 ///////////////////////////////////////////////////////////////
338 // Parameter
339 ///////////////////////////////////////////////////////////////
340
341 // Converts number of electrons to fC
8230f242 342 const Float_t kEl2fC = 1.602E-19 * 1.0E15;
f7336fa3 343
344 ///////////////////////////////////////////////////////////////
345
346 Int_t iRow, iCol, iTime;
347 Int_t nBytes = 0;
28329a48 348 Int_t iDict;
f7336fa3 349
350 Int_t totalSizeDigits = 0;
351 Int_t totalSizeDict0 = 0;
352 Int_t totalSizeDict1 = 0;
353 Int_t totalSizeDict2 = 0;
354
8230f242 355 AliTRDdataArrayI *digits;
356 AliTRDdataArrayI *dictionary[kNDict];
da581aea 357
358 if (!fGeo) {
359 printf("AliTRDdigitizer::MakeDigits -- ");
360 printf("No geometry defined\n");
361 return kFALSE;
362 }
363
364 // Create a digits manager
365 fDigits = new AliTRDdigitsManager();
366
367 // Create detector arrays to keep the signal and track numbers
8230f242 368 AliTRDdataArrayF *signal = new AliTRDdataArrayF();
369 AliTRDdataArrayI *tracks[kNDict];
da581aea 370 for (iDict = 0; iDict < kNDict; iDict++) {
8230f242 371 tracks[iDict] = new AliTRDdataArrayI();
da581aea 372 }
f7336fa3 373
374 // Get the pointer to the hit tree
8230f242 375 TTree *hitTree = gAlice->TreeH();
da581aea 376
377 // Get the number of entries in the hit tree
378 // (Number of primary particles creating a hit somewhere)
8230f242 379 Int_t nTrack = (Int_t) hitTree->GetEntries();
da581aea 380
381 printf("AliTRDdigitizer::MakeDigits -- ");
382 printf("Start creating digits.\n");
f7336fa3 383
384 // The Lorentz factor
385 if (fExBOn) {
386 fLorentzFactor = 1.0 / (1.0 + fLorentzAngle*fLorentzAngle);
387 }
388 else {
389 fLorentzFactor = 1.0;
390 }
391
f7336fa3 392 Int_t chamBeg = 0;
393 Int_t chamEnd = kNcham;
9d0b222b 394 if (fTRD->GetSensChamber() >= 0) {
da581aea 395 chamBeg = fTRD->GetSensChamber();
396 chamEnd = chamBeg + 1;
f7336fa3 397 }
398 Int_t planBeg = 0;
399 Int_t planEnd = kNplan;
9d0b222b 400 if (fTRD->GetSensPlane() >= 0) {
da581aea 401 planBeg = fTRD->GetSensPlane();
f7336fa3 402 planEnd = planBeg + 1;
403 }
404 Int_t sectBeg = 0;
405 Int_t sectEnd = kNsect;
f7336fa3 406
8230f242 407 Int_t countHits = 0;
da581aea 408
f7336fa3 409 // Loop through all the chambers
410 for (Int_t iCham = chamBeg; iCham < chamEnd; iCham++) {
411 for (Int_t iPlan = planBeg; iPlan < planEnd; iPlan++) {
412 for (Int_t iSect = sectBeg; iSect < sectEnd; iSect++) {
413
9d0b222b 414 if (fTRD->GetSensSector() >= 0) {
415 Int_t sens1 = fTRD->GetSensSector();
416 Int_t sens2 = sens1 + fTRD->GetSensSectorRange();
417 sens2 -= ((Int_t) (sens2 / kNsect)) * kNsect;
dd9a6ee3 418 if (sens1 < sens2) {
9d0b222b 419 if ((iSect < sens1) || (iSect >= sens2)) continue;
dd9a6ee3 420 }
421 else {
9d0b222b 422 if ((iSect < sens1) && (iSect >= sens2)) continue;
dd9a6ee3 423 }
9d0b222b 424 }
425
f7336fa3 426 Int_t nDigits = 0;
427
f7336fa3 428 printf("AliTRDdigitizer::MakeDigits -- ");
429 printf("Digitizing chamber %d, plane %d, sector %d.\n"
430 ,iCham,iPlan,iSect);
431
da581aea 432 Int_t iDet = fGeo->GetDetector(iPlan,iCham,iSect);
433 Int_t nRowMax = fGeo->GetRowMax(iPlan,iCham,iSect);
434 Int_t nColMax = fGeo->GetColMax(iPlan);
435 Int_t nTimeMax = fGeo->GetTimeMax();
436 Float_t row0 = fGeo->GetRow0(iPlan,iCham,iSect);
437 Float_t col0 = fGeo->GetCol0(iPlan);
438 Float_t time0 = fGeo->GetTime0(iPlan);
439 Float_t rowPadSize = fGeo->GetRowPadSize();
440 Float_t colPadSize = fGeo->GetColPadSize();
441 Float_t timeBinSize = fGeo->GetTimeBinSize();
442
443 // Adjust the size of the detector arrays
8230f242 444 signal->Allocate(nRowMax,nColMax,nTimeMax);
da581aea 445 for (iDict = 0; iDict < kNDict; iDict++) {
8230f242 446 tracks[iDict]->Allocate(nRowMax,nColMax,nTimeMax);
da581aea 447 }
f7336fa3 448
449 // Loop through all entries in the tree
450 for (Int_t iTrack = 0; iTrack < nTrack; iTrack++) {
451
452 gAlice->ResetHits();
8230f242 453 nBytes += hitTree->GetEvent(iTrack);
f7336fa3 454
455 // Get the number of hits in the TRD created by this particle
da581aea 456 Int_t nHit = fTRD->Hits()->GetEntriesFast();
f7336fa3 457
458 // Loop through the TRD hits
459 for (Int_t iHit = 0; iHit < nHit; iHit++) {
460
8230f242 461 countHits++;
f7336fa3 462
8230f242 463 AliTRDhit *hit = (AliTRDhit *) fTRD->Hits()->UncheckedAt(iHit);
f7336fa3 464 Float_t pos[3];
94de3818 465 pos[0] = hit->X();
466 pos[1] = hit->Y();
467 pos[2] = hit->Z();
8230f242 468 Float_t q = hit->GetCharge();
94de3818 469 Int_t track = hit->Track();
8230f242 470 Int_t detector = hit->GetDetector();
da581aea 471 Int_t plane = fGeo->GetPlane(detector);
472 Int_t sector = fGeo->GetSector(detector);
473 Int_t chamber = fGeo->GetChamber(detector);
f7336fa3 474
475 if ((sector != iSect) ||
476 (plane != iPlan) ||
477 (chamber != iCham))
478 continue;
479
480 // Rotate the sectors on top of each other
481 Float_t rot[3];
da581aea 482 fGeo->Rotate(detector,pos,rot);
f7336fa3 483
484 // The hit position in pad coordinates (center pad)
485 // The pad row (z-direction)
486 Int_t rowH = (Int_t) ((rot[2] - row0) / rowPadSize);
487 // The pad column (rphi-direction)
488 Int_t colH = (Int_t) ((rot[1] - col0) / colPadSize);
489 // The time bucket
490 Int_t timeH = (Int_t) ((rot[0] - time0) / timeBinSize);
491
492 // Array to sum up the signal in a box surrounding the
493 // hit postition
8230f242 494 const Int_t kTimeBox = 7;
495 const Int_t kColBox = 9;
496 const Int_t kRowBox = 7;
497 Float_t signalSum[kRowBox][kColBox][kTimeBox];
498 for (iRow = 0; iRow < kRowBox; iRow++ ) {
499 for (iCol = 0; iCol < kColBox; iCol++ ) {
500 for (iTime = 0; iTime < kTimeBox; iTime++) {
f7336fa3 501 signalSum[iRow][iCol][iTime] = 0;
502 }
503 }
504 }
505
506 // Loop over all electrons of this hit
507 Int_t nEl = (Int_t) q;
508 for (Int_t iEl = 0; iEl < nEl; iEl++) {
509
510 // The driftlength
511 Float_t driftlength = rot[0] - time0;
512 if ((driftlength < 0) ||
513 (driftlength > kDrThick)) break;
514 Float_t driftlengthL = driftlength;
515 if (fExBOn) driftlengthL /= TMath::Sqrt(fLorentzFactor);
516 Float_t xyz[3];
517 xyz[0] = rot[0];
518 xyz[1] = rot[1];
519 xyz[2] = rot[2];
520
521 // Electron attachment
522 if (fElAttachOn) {
523 if (gRandom->Rndm() < (driftlengthL * fElAttachProp / 100.)) continue;
524 }
525
526 // Apply the diffusion smearing
527 if (fDiffusionOn) {
528 if (!(Diffusion(driftlengthL,xyz))) continue;
529 }
530
531 // Apply E x B effects
532 if (fExBOn) {
533 if (!(ExB(driftlength,xyz))) continue;
534 }
535
536 // The electron position and the distance to the hit position
537 // in pad units
538 // The pad row (z-direction)
539 Int_t rowE = (Int_t) ((xyz[2] - row0) / rowPadSize);
540 Int_t rowD = rowH - rowE;
541 // The pad column (rphi-direction)
542 Int_t colE = (Int_t) ((xyz[1] - col0) / colPadSize);
543 Int_t colD = colH - colE;
544 // The time bucket
545 Int_t timeE = (Int_t) ((xyz[0] - time0) / timeBinSize);
546 Int_t timeD = timeH - timeE;
547
548 // Apply the gas gain including fluctuations
549 Int_t signal = (Int_t) (-fGasGain * TMath::Log(gRandom->Rndm()));
550
551 // The distance of the electron to the center of the pad
552 // in units of pad width
553 Float_t dist = (xyz[1] - col0 - (colE + 0.5) * colPadSize)
554 / colPadSize;
555
556 // Sum up the signal in the different pixels
557 // and apply the pad response
8230f242 558 Int_t rowIdx = rowD + (Int_t) ( kRowBox / 2);
559 Int_t colIdx = colD + (Int_t) ( kColBox / 2);
560 Int_t timeIdx = timeD + (Int_t) (kTimeBox / 2);
f7336fa3 561
8230f242 562 if (( rowIdx < 0) || ( rowIdx > kRowBox)) {
f7336fa3 563 printf("AliTRDdigitizer::MakeDigits -- ");
8230f242 564 printf("Boundary error. rowIdx = %d (%d)\n", rowIdx, kRowBox);
f7336fa3 565 continue;
566 }
8230f242 567 if (( colIdx < 0) || ( colIdx > kColBox)) {
f7336fa3 568 printf("AliTRDdigitizer::MakeDigits -- ");
8230f242 569 printf("Boundary error. colIdx = %d (%d)\n", colIdx, kColBox);
f7336fa3 570 continue;
571 }
8230f242 572 if ((timeIdx < 0) || (timeIdx > kTimeBox)) {
f7336fa3 573 printf("AliTRDdigitizer::MakeDigits -- ");
8230f242 574 printf("Boundary error. timeIdx = %d (%d)\n",timeIdx,kTimeBox);
f7336fa3 575 continue;
576 }
da581aea 577 signalSum[rowIdx][colIdx-1][timeIdx] += fPRF->Eval(dist-1.0,0,0) * signal;
578 signalSum[rowIdx][colIdx ][timeIdx] += fPRF->Eval(dist ,0,0) * signal;
579 signalSum[rowIdx][colIdx+1][timeIdx] += fPRF->Eval(dist+1.0,0,0) * signal;
f7336fa3 580
581 }
582
583 // Add the padcluster to the detector matrix
8230f242 584 for (iRow = 0; iRow < kRowBox; iRow++ ) {
585 for (iCol = 0; iCol < kColBox; iCol++ ) {
586 for (iTime = 0; iTime < kTimeBox; iTime++) {
f7336fa3 587
8230f242 588 Int_t rowB = rowH + iRow - (Int_t) ( kRowBox / 2);
589 Int_t colB = colH + iCol - (Int_t) ( kColBox / 2);
590 Int_t timeB = timeH + iTime - (Int_t) (kTimeBox / 2);
f7336fa3 591 Float_t signalB = signalSum[iRow][iCol][iTime];
da581aea 592 if (( rowB < 0) || ( rowB >= nRowMax)) continue;
593 if (( colB < 0) || ( colB >= nColMax)) continue;
594 if ((timeB < 0) || (timeB >= nTimeMax)) continue;
f7336fa3 595 if (signalB > 0.0) {
da581aea 596
597 // Add the signal sum
8230f242 598 signalB += signal->GetData(rowB,colB,timeB);
599 signal->SetData(rowB,colB,timeB,signalB);
da581aea 600 // Store the track index in the dictionary
601 // Note: We store index+1 in order to allow the array to be compressed
602 for (iDict = 0; iDict < kNDict; iDict++) {
8230f242 603 Int_t oldTrack = tracks[iDict]->GetData(rowB,colB,timeB);
da581aea 604 if (oldTrack == track+1) break;
605 if (oldTrack == -1) break;
606 if (oldTrack == 0) {
8230f242 607 tracks[iDict]->SetData(rowB,colB,timeB,track+1);
da581aea 608 break;
609 }
610 }
611 if (iDict == kNDict) {
f7336fa3 612 printf("AliTRDdigitizer::MakeDigits -- ");
da581aea 613 printf("More than three tracks for one digit!\n");
614 }
615 }
f7336fa3 616
617 }
618 }
619 }
620
621 }
622
623 }
624
625 // Add a container for the digits of this detector
8230f242 626 digits = fDigits->GetDigits(iDet);
f7336fa3 627 // Allocate memory space for the digits buffer
8230f242 628 digits->Allocate(nRowMax,nColMax,nTimeMax);
f7336fa3 629
da581aea 630 // Do the same for the dictionary arrays
28329a48 631 for (iDict = 0; iDict < kNDict; iDict++) {
8230f242 632 dictionary[iDict] = fDigits->GetDictionary(iDet,iDict);
633 dictionary[iDict]->Allocate(nRowMax,nColMax,nTimeMax);
f7336fa3 634 }
635
da581aea 636 // Create the digits for this chamber
f7336fa3 637 for (iRow = 0; iRow < nRowMax; iRow++ ) {
638 for (iCol = 0; iCol < nColMax; iCol++ ) {
639 for (iTime = 0; iTime < nTimeMax; iTime++) {
640
8230f242 641 Float_t signalAmp = signal->GetData(iRow,iCol,iTime);
f7336fa3 642
643 // Add the noise
1befd3b2 644 signalAmp = TMath::Max((Float_t) gRandom->Gaus(signalAmp,fNoise)
645 ,(Float_t) 0.0);
f7336fa3 646 // Convert to fC
8230f242 647 signalAmp *= kEl2fC;
f7336fa3 648 // Convert to mV
649 signalAmp *= fChipGain;
650 // Convert to ADC counts
651 Int_t adc = (Int_t) (signalAmp * (fADCoutRange / fADCinRange));
652
da581aea 653 if (adc > fADCthreshold) {
f7336fa3 654
da581aea 655 nDigits++;
656
657 // Store the amplitude of the digit
8230f242 658 digits->SetData(iRow,iCol,iTime,adc);
da581aea 659
660 // Store the track index in the dictionary
661 // Note: We store index+1 in order to allow the array to be compressed
662 for (iDict = 0; iDict < kNDict; iDict++) {
8230f242 663 dictionary[iDict]->SetData(iRow,iCol,iTime
664 ,tracks[iDict]->GetData(iRow,iCol,iTime));
da581aea 665 }
f7336fa3 666
da581aea 667 }
f7336fa3 668
669 }
670 }
671 }
672
673 // Compress the arrays
8230f242 674 digits->Compress(1,0);
28329a48 675 for (iDict = 0; iDict < kNDict; iDict++) {
8230f242 676 dictionary[iDict]->Compress(1,0);
f7336fa3 677 }
678
8230f242 679 totalSizeDigits += digits->GetSize();
680 totalSizeDict0 += dictionary[0]->GetSize();
681 totalSizeDict1 += dictionary[1]->GetSize();
682 totalSizeDict2 += dictionary[2]->GetSize();
f7336fa3 683
684 printf("AliTRDdigitizer::MakeDigits -- ");
685 printf("Number of digits found: %d.\n",nDigits);
da581aea 686
687 // Reset the arrays
8230f242 688 signal->Reset();
da581aea 689 for (iDict = 0; iDict < kNDict; iDict++) {
8230f242 690 tracks[iDict]->Reset();
da581aea 691 }
f7336fa3 692
693 }
694 }
695 }
696
da581aea 697 printf("AliTRDdigitizer::MakeDigits -- ");
8230f242 698 printf("Total number of analyzed hits = %d\n",countHits);
da581aea 699
f7336fa3 700 printf("AliTRDdigitizer::MakeDigits -- ");
701 printf("Total digits data size = %d, %d, %d, %d\n",totalSizeDigits
702 ,totalSizeDict0
703 ,totalSizeDict1
704 ,totalSizeDict2);
705
706 return kTRUE;
707
708}
709
f7336fa3 710//_____________________________________________________________________________
711Bool_t AliTRDdigitizer::WriteDigits()
712{
713 //
714 // Writes out the TRD-digits and the dictionaries
715 //
716
717 // Create the branches
718 if (!(gAlice->TreeD()->GetBranch("TRDdigits"))) {
da581aea 719 if (!fDigits->MakeBranch()) return kFALSE;
f7336fa3 720 }
721
da581aea 722 // Store the digits and the dictionary in the tree
723 fDigits->WriteDigits();
f7336fa3 724
725 // Write the new tree into the input file (use overwrite option)
726 Char_t treeName[7];
727 sprintf(treeName,"TreeD%d",fEvent);
728 printf("AliTRDdigitizer::WriteDigits -- ");
729 printf("Write the digits tree %s for event %d.\n"
730 ,treeName,fEvent);
731 gAlice->TreeD()->Write(treeName,2);
732
733 return kTRUE;
734
735}