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