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