Line data Source code
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 : // Source File : PMDDigitizer.cxx, Version 00 //
18 : // //
19 : // Date : September 20 2002 //
20 : // //
21 : //-----------------------------------------------------//
22 :
23 : #include <Riostream.h>
24 : #include <TBRIK.h>
25 : #include <TNode.h>
26 : #include <TTree.h>
27 : #include <TGeometry.h>
28 : #include <TObjArray.h>
29 : #include <TClonesArray.h>
30 : #include <TFile.h>
31 : #include <TNtuple.h>
32 : #include <TParticle.h>
33 : #include <TRandom.h>
34 :
35 : #include "AliLog.h"
36 : #include "AliRun.h"
37 : #include "AliHit.h"
38 : #include "AliDetector.h"
39 : #include "AliRunLoader.h"
40 : #include "AliLoader.h"
41 : #include "AliConfig.h"
42 : #include "AliMagF.h"
43 : #include "AliDigitizationInput.h"
44 : #include "AliDigitizer.h"
45 : #include "AliHeader.h"
46 : #include "AliCDBManager.h"
47 : #include "AliCDBStorage.h"
48 : #include "AliCDBEntry.h"
49 : #include "AliMC.h"
50 :
51 : #include "AliPMD.h"
52 : #include "AliPMDhit.h"
53 : #include "AliPMDcell.h"
54 : #include "AliPMDsdigit.h"
55 : #include "AliPMDdigit.h"
56 : #include "AliPMDCalibData.h"
57 : #include "AliPMDPedestal.h"
58 : #include "AliPMDDigitizer.h"
59 :
60 :
61 12 : ClassImp(AliPMDDigitizer)
62 :
63 1 : AliPMDDigitizer::AliPMDDigitizer() :
64 1 : fRunLoader(0),
65 1 : fPMDHit(0),
66 1 : fPMD(0),
67 1 : fPMDLoader(0),
68 2 : fCalibGain(GetCalibGain()),
69 2 : fCalibPed(GetCalibPed()),
70 1 : fSDigits(0),
71 1 : fDigits(0),
72 1 : fCPVCell(0),
73 1 : fCell(0),
74 1 : fNsdigit(0),
75 1 : fNdigit(0),
76 1 : fDetNo(0),
77 1 : fZPos(361.5) // in units of cm, default position of PMD
78 5 : {
79 : // Default Constructor
80 : //
81 50 : for (Int_t i = 0; i < fgkTotUM; i++)
82 : {
83 2352 : for (Int_t j = 0; j < fgkRow; j++)
84 : {
85 223488 : for (Int_t k = 0; k < fgkCol; k++)
86 : {
87 110592 : fCPV[i][j][k] = 0.;
88 110592 : fPRE[i][j][k] = 0.;
89 110592 : fCPVCounter[i][j][k] = 0;
90 110592 : fPRECounter[i][j][k] = 0;
91 110592 : fCPVTrackNo[i][j][k] = -1;
92 110592 : fPRETrackNo[i][j][k] = -1;
93 110592 : fCPVTrackPid[i][j][k] = -1;
94 110592 : fPRETrackPid[i][j][k] = -1;
95 : }
96 : }
97 : }
98 :
99 2 : }
100 : //____________________________________________________________________________
101 : AliPMDDigitizer::AliPMDDigitizer(const AliPMDDigitizer& digitizer):
102 0 : AliDigitizer(digitizer),
103 0 : fRunLoader(0),
104 0 : fPMDHit(0),
105 0 : fPMD(0),
106 0 : fPMDLoader(0),
107 0 : fCalibGain(GetCalibGain()),
108 0 : fCalibPed(GetCalibPed()),
109 0 : fSDigits(0),
110 0 : fDigits(0),
111 0 : fCPVCell(0),
112 0 : fCell(0),
113 0 : fNsdigit(0),
114 0 : fNdigit(0),
115 0 : fDetNo(0),
116 0 : fZPos(361.5) // in units of cm, default position of PMD
117 0 : {
118 : // copy constructor
119 0 : AliError("Copy constructor not allowed ");
120 :
121 0 : }
122 : //____________________________________________________________________________
123 : AliPMDDigitizer & AliPMDDigitizer::operator=(const AliPMDDigitizer& /*digitizer*/)
124 : {
125 : // Assignment operator
126 0 : AliError("Assignement operator not allowed ");
127 :
128 0 : return *this;
129 : }
130 : //____________________________________________________________________________
131 : AliPMDDigitizer::AliPMDDigitizer(AliDigitizationInput* digInput):
132 1 : AliDigitizer(digInput),
133 1 : fRunLoader(0),
134 1 : fPMDHit(0),
135 1 : fPMD(0),
136 1 : fPMDLoader(0),
137 2 : fCalibGain(GetCalibGain()),
138 2 : fCalibPed(GetCalibPed()),
139 3 : fSDigits(new TClonesArray("AliPMDsdigit", 1000)),
140 3 : fDigits(new TClonesArray("AliPMDdigit", 1000)),
141 1 : fCPVCell(0),
142 1 : fCell(0),
143 1 : fNsdigit(0),
144 1 : fNdigit(0),
145 1 : fDetNo(0),
146 1 : fZPos(361.5)// in units of cm, This is the default position of PMD
147 5 : {
148 : // ctor which should be used
149 :
150 :
151 50 : for (Int_t i = 0; i < fgkTotUM; i++)
152 : {
153 2352 : for (Int_t j = 0; j < fgkRow; j++)
154 : {
155 223488 : for (Int_t k = 0; k < fgkCol; k++)
156 : {
157 110592 : fCPV[i][j][k] = 0.;
158 110592 : fPRE[i][j][k] = 0.;
159 110592 : fCPVCounter[i][j][k] = 0;
160 110592 : fPRECounter[i][j][k] = 0;
161 110592 : fCPVTrackNo[i][j][k] = -1;
162 110592 : fPRETrackNo[i][j][k] = -1;
163 110592 : fCPVTrackPid[i][j][k] = -1;
164 110592 : fPRETrackPid[i][j][k] = -1;
165 : }
166 : }
167 : }
168 2 : }
169 :
170 : //____________________________________________________________________________
171 : AliPMDDigitizer::~AliPMDDigitizer()
172 12 : {
173 : // Default Destructor
174 : //
175 2 : if (fSDigits) {
176 2 : fSDigits->Delete();
177 4 : delete fSDigits;
178 2 : fSDigits=0;
179 2 : }
180 2 : if (fDigits) {
181 1 : fDigits->Delete();
182 2 : delete fDigits;
183 1 : fDigits=0;
184 1 : }
185 2 : fCPVCell.Delete();
186 2 : fCell.Delete();
187 6 : }
188 : //
189 : // Member functions
190 : //
191 : //____________________________________________________________________________
192 : void AliPMDDigitizer::OpengAliceFile(const char *file, Option_t *option)
193 : {
194 : // Loads galice.root file and corresponding header, kinematics
195 : // hits and sdigits or digits depending on the option
196 : //
197 :
198 2 : TString evfoldname = AliConfig::GetDefaultEventFolderName();
199 3 : fRunLoader = AliRunLoader::GetRunLoader(evfoldname);
200 1 : if (!fRunLoader)
201 0 : fRunLoader = AliRunLoader::Open(file,AliConfig::GetDefaultEventFolderName(), "UPDATE");
202 :
203 1 : if (!fRunLoader)
204 : {
205 0 : AliError(Form("Can not open session for file %s.",file));
206 : }
207 :
208 1 : const char *cHS = strstr(option,"HS");
209 1 : const char *cHD = strstr(option,"HD");
210 1 : const char *cSD = strstr(option,"SD");
211 :
212 1 : if(cHS || cHD)
213 : {
214 2 : if (!fRunLoader->GetAliRun()) fRunLoader->LoadgAlice();
215 2 : if (!fRunLoader->TreeE()) fRunLoader->LoadHeader();
216 2 : if (!fRunLoader->TreeK()) fRunLoader->LoadKinematics();
217 :
218 2 : gAlice = fRunLoader->GetAliRun();
219 :
220 1 : if (gAlice)
221 : {
222 5 : AliDebug(1,"Alirun object found");
223 : }
224 : else
225 : {
226 0 : AliError("Could not found Alirun object");
227 : }
228 :
229 2 : fPMD = (AliPMD*)gAlice->GetDetector("PMD");
230 1 : }
231 :
232 2 : fPMDLoader = fRunLoader->GetLoader("PMDLoader");
233 1 : if (fPMDLoader == 0x0)
234 : {
235 0 : AliError("Can not find PMDLoader");
236 : }
237 :
238 :
239 1 : if (cHS)
240 : {
241 1 : fPMDLoader->LoadHits("READ");
242 1 : fPMDLoader->LoadSDigits("recreate");
243 : }
244 0 : else if (cHD)
245 : {
246 0 : fPMDLoader->LoadHits("READ");
247 0 : fPMDLoader->LoadDigits("recreate");
248 : }
249 0 : else if (cSD)
250 : {
251 0 : fPMDLoader->LoadSDigits("READ");
252 0 : fPMDLoader->LoadDigits("recreate");
253 : }
254 1 : }
255 : //____________________________________________________________________________
256 : void AliPMDDigitizer::Hits2SDigits(Int_t ievt)
257 : {
258 : // This reads the PMD Hits tree and assigns the right track number
259 : // to a cell and stores in the summable digits tree
260 : //
261 :
262 : const Int_t kPi0 = 111;
263 : const Int_t kGamma = 22;
264 : Int_t npmd = 0;
265 : Int_t trackno = 0;
266 : Int_t smnumber = 0;
267 : Int_t trackpid = 0;
268 : Int_t mtrackno = 0;
269 : Int_t mtrackpid = 0;
270 :
271 : Float_t xPos = 0., yPos = 0., zPos = 0.;
272 : Int_t xpad = -1, ypad = -1;
273 : Float_t edep = 0.;
274 : // Float_t vx = -999.0, vy = -999.0, vz = -999.0; //coverity (8443) fix satya (1/9/2014)
275 :
276 :
277 10 : if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
278 4 : ResetSDigit();
279 :
280 12 : AliDebug(1,Form("Event Number = %d",ievt));
281 4 : Int_t nparticles = fRunLoader->GetHeader()->GetNtrack();
282 12 : AliDebug(1,Form("Number of Particles = %d",nparticles));
283 :
284 : //
285 :
286 4 : fRunLoader->GetEvent(ievt);
287 : // ------------------------------------------------------- //
288 : // Pointer to specific detector hits.
289 : // Get pointers to Alice detectors and Hits containers
290 :
291 4 : TTree* treeH = fPMDLoader->TreeH();
292 :
293 4 : Int_t ntracks = (Int_t) treeH->GetEntries();
294 12 : AliDebug(1,Form("Number of Tracks in the TreeH = %d", ntracks));
295 4 : TTree* treeS = fPMDLoader->TreeS();
296 4 : if (treeS == 0x0)
297 : {
298 4 : fPMDLoader->MakeTree("S");
299 4 : treeS = fPMDLoader->TreeS();
300 4 : }
301 : Int_t bufsize = 16000;
302 4 : treeS->Branch("PMDSDigit", &fSDigits, bufsize);
303 :
304 : TClonesArray* hits = 0;
305 8 : if (fPMD) hits = fPMD->Hits();
306 :
307 : // Start loop on tracks in the hits containers
308 :
309 232 : for (Int_t track=0; track<ntracks;track++)
310 : {
311 112 : gAlice->GetMCApp()->ResetHits();
312 112 : treeH->GetEvent(track);
313 112 : if (fPMD)
314 : {
315 112 : npmd = hits->GetEntriesFast();
316 1168 : for (Int_t ipmd = 0; ipmd < npmd; ipmd++)
317 : {
318 472 : fPMDHit = (AliPMDhit*) hits->UncheckedAt(ipmd);
319 472 : trackno = fPMDHit->GetTrack();
320 : // get kinematics of the particles
321 :
322 472 : TParticle* mparticle = gAlice->GetMCApp()->Particle(trackno);
323 472 : trackpid = mparticle->GetPdgCode();
324 472 : Int_t ks = mparticle->GetStatusCode();
325 : Int_t imo;
326 : Int_t tracknoOld=0, trackpidOld=0;
327 : // Int_t statusOld = 0; //coverity (8443) fix satya (1/9/2014)
328 :
329 472 : if (mparticle->GetFirstMother() == -1)
330 : {
331 : tracknoOld = trackno;
332 : trackpidOld = trackpid;
333 : //statusOld = -1;
334 19 : }
335 : // Int_t igstatus = 0;
336 :
337 : Int_t trnotemp = trackno;
338 944 : if(ks==1||(imo = mparticle->GetFirstMother())<0 ){
339 : // vx = mparticle->Vx(); //coverity (8443) fix satya (1/9/2014)
340 : // vy = mparticle->Vy();
341 : // vz = mparticle->Vz();
342 :
343 : // if(trackpid==kGamma||trackpid==11||trackpid==-11||
344 : // trackpid==kPi0)igstatus=1;
345 : }
346 :
347 :
348 4224 : while(((imo = mparticle->GetFirstMother()) >= 0)&&
349 1876 : (ks = mparticle->GetStatusCode() <1) )
350 : {
351 1876 : mparticle = gAlice->GetMCApp()->Particle(imo);
352 1876 : trackpid = mparticle->GetPdgCode();
353 1876 : ks = mparticle->GetStatusCode();
354 : // vx = mparticle->Vx();//coverity (8443) fix satya (1/9/2014)
355 : // vy = mparticle->Vy();
356 : // vz = mparticle->Vz();
357 :
358 :
359 : trnotemp = trackno;
360 1876 : if(trackpid == 111)
361 : {
362 : trackno = trnotemp;
363 42 : }
364 1876 : if(trackpid != 111)
365 : {
366 : trackno=imo;
367 1834 : }
368 : // end of modification on 25th Nov 2009
369 : }
370 :
371 : // if(trackpid==kGamma||trackpid==11||trackpid==-11||
372 : // trackpid==kPi0)igstatus=1;
373 : mtrackpid=trackpid;
374 : mtrackno=trackno;
375 : trackpid=trackpidOld;
376 : trackno=tracknoOld;
377 :
378 : //-----------------end of modification----------------
379 472 : Float_t ptime = fPMDHit->GetTime()*1e6; // time in microsec
380 955 : if (ptime < 0. || ptime > 1.2) continue;
381 :
382 461 : xPos = fPMDHit->X();
383 461 : yPos = fPMDHit->Y();
384 461 : zPos = fPMDHit->Z();
385 :
386 461 : edep = fPMDHit->GetEnergy();
387 461 : Int_t vol1 = fPMDHit->GetVolume(1); // Column
388 461 : Int_t vol2 = fPMDHit->GetVolume(2); // Row
389 461 : Int_t vol7 = fPMDHit->GetVolume(4); // Serial Module No
390 :
391 :
392 : // -----------------------------------------//
393 : // In new geometry after adding electronics //
394 : // For Super Module 1 & 2 //
395 : // nrow = 48, ncol = 96 //
396 : // For Super Module 3 & 4 //
397 : // nrow = 96, ncol = 48 //
398 : // -----------------------------------------//
399 :
400 461 : if (vol7 < 24)
401 : {
402 : smnumber = vol7;
403 382 : }
404 : else
405 : {
406 79 : smnumber = vol7 - 24;
407 : }
408 461 : Int_t vol8 = smnumber/6 + 1; // fake supermodule
409 :
410 461 : if (vol8 == 1 || vol8 == 2)
411 : {
412 : xpad = vol2;
413 : ypad = vol1;
414 400 : }
415 61 : else if (vol8 == 3 || vol8 == 4)
416 : {
417 : xpad = vol1;
418 : ypad = vol2;
419 61 : }
420 :
421 1383 : AliDebug(2,Form("Zposition = %f Edeposition = %f",zPos,edep));
422 :
423 461 : if (vol7 < 24)
424 : {
425 : // PRE
426 843 : fDetNo = 0;
427 382 : }
428 : else
429 : {
430 : // CPV
431 79 : fDetNo = 1;
432 : }
433 :
434 : Int_t smn = smnumber;
435 461 : Int_t ixx = xpad - 1;
436 461 : Int_t iyy = ypad - 1;
437 461 : if (fDetNo == 0)
438 : {
439 382 : fPRE[smn][ixx][iyy] += edep;
440 382 : fPRECounter[smn][ixx][iyy]++;
441 :
442 382 : AliPMDcell* cell = new AliPMDcell(mtrackno,smn,ixx,iyy,edep);
443 382 : fCell.Add(cell);
444 382 : }
445 79 : else if(fDetNo == 1)
446 : {
447 79 : fCPV[smn][ixx][iyy] += edep;
448 79 : fCPVCounter[smn][ixx][iyy]++;
449 79 : AliPMDcell* cpvcell = new AliPMDcell(mtrackno,smn,ixx,iyy,edep);
450 79 : fCPVCell.Add(cpvcell);
451 79 : }
452 461 : }
453 112 : }
454 : } // Track Loop ended
455 4 : TrackAssignment2CPVCell();
456 4 : TrackAssignment2Cell();
457 4 : ResetCell();
458 :
459 : Float_t deltaE = 0.;
460 : Int_t detno = 0;
461 : Int_t trno = -1;
462 : Int_t trpid = -99;
463 :
464 24 : for (Int_t idet = 0; idet < 2; idet++)
465 : {
466 400 : for (Int_t ism = 0; ism < fgkTotUM; ism++)
467 : {
468 18816 : for (Int_t jrow = 0; jrow < fgkRow; jrow++)
469 : {
470 1787904 : for (Int_t kcol = 0; kcol < fgkCol; kcol++)
471 : {
472 884736 : if (idet == 0)
473 : {
474 442368 : deltaE = fPRE[ism][jrow][kcol];
475 442368 : trno = fPRETrackNo[ism][jrow][kcol];
476 : detno = 0;
477 442368 : }
478 442368 : else if (idet == 1)
479 : {
480 442368 : deltaE = fCPV[ism][jrow][kcol];
481 442368 : trno = fCPVTrackNo[ism][jrow][kcol];
482 : detno = 1;
483 442368 : }
484 884736 : if (deltaE > 0.)
485 : {
486 : // Natasha
487 221 : TParticle *mparticle = gAlice->GetMCApp()->Particle(trno);
488 221 : trpid = mparticle->GetPdgCode();
489 221 : AddSDigit(trno,trpid,detno,ism,jrow,kcol,deltaE);
490 221 : }
491 : }
492 : }
493 192 : treeS->Fill();
494 192 : ResetSDigit();
495 : }
496 : }
497 4 : fPMDLoader->WriteSDigits("OVERWRITE");
498 4 : ResetCellADC();
499 4 : }
500 : //____________________________________________________________________________
501 :
502 : void AliPMDDigitizer::Hits2Digits(Int_t ievt)
503 : {
504 : // This reads the PMD Hits tree and assigns the right track number
505 : // to a cell and stores in the digits tree
506 : //
507 : const Int_t kPi0 = 111;
508 : const Int_t kGamma = 22;
509 : Int_t npmd = 0;
510 : Int_t trackno = 0;
511 : Int_t smnumber = 0;
512 : Int_t trackpid = 0;
513 : Int_t mtrackno = 0;
514 : Int_t mtrackpid = 0;
515 :
516 : Float_t xPos = 0., yPos = 0., zPos = 0.;
517 : Int_t xpad = -1, ypad = -1;
518 : Float_t edep = 0.;
519 : // Float_t vx = -999.0, vy = -999.0, vz = -999.0; //coverity (8443) fix satya (1/9/2014)
520 :
521 0 : if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
522 0 : ResetDigit();
523 :
524 0 : AliDebug(1,Form("Event Number = %d",ievt));
525 0 : Int_t nparticles = fRunLoader->GetHeader()->GetNtrack();
526 0 : AliDebug(1,Form("Number of Particles = %d", nparticles));
527 :
528 0 : fRunLoader->GetEvent(ievt);
529 : // ------------------------------------------------------- //
530 : // Pointer to specific detector hits.
531 : // Get pointers to Alice detectors and Hits containers
532 :
533 0 : fPMD = (AliPMD*)gAlice->GetDetector("PMD");
534 0 : fPMDLoader = fRunLoader->GetLoader("PMDLoader");
535 :
536 0 : if (fPMDLoader == 0x0)
537 : {
538 0 : AliError("Can not find PMD or PMDLoader");
539 0 : }
540 0 : TTree* treeH = fPMDLoader->TreeH();
541 0 : Int_t ntracks = (Int_t) treeH->GetEntries();
542 0 : AliDebug(1,Form("Number of Tracks in the TreeH = %d", ntracks));
543 0 : fPMDLoader->LoadDigits("recreate");
544 0 : TTree* treeD = fPMDLoader->TreeD();
545 0 : if (treeD == 0x0)
546 : {
547 0 : fPMDLoader->MakeTree("D");
548 0 : treeD = fPMDLoader->TreeD();
549 0 : }
550 : Int_t bufsize = 16000;
551 0 : treeD->Branch("PMDDigit", &fDigits, bufsize);
552 :
553 : TClonesArray* hits = 0;
554 0 : if (fPMD) hits = fPMD->Hits();
555 :
556 : // Start loop on tracks in the hits containers
557 :
558 0 : for (Int_t track=0; track<ntracks;track++)
559 : {
560 0 : gAlice->GetMCApp()->ResetHits();
561 0 : treeH->GetEvent(track);
562 :
563 0 : if (fPMD)
564 : {
565 0 : npmd = hits->GetEntriesFast();
566 0 : for (Int_t ipmd = 0; ipmd < npmd; ipmd++)
567 : {
568 0 : fPMDHit = (AliPMDhit*) hits->UncheckedAt(ipmd);
569 0 : trackno = fPMDHit->GetTrack();
570 :
571 : // get kinematics of the particles
572 :
573 0 : TParticle* mparticle = gAlice->GetMCApp()->Particle(trackno);
574 0 : trackpid = mparticle->GetPdgCode();
575 0 : Int_t ks = mparticle->GetStatusCode();
576 : Int_t imo;
577 : Int_t tracknoOld=0, trackpidOld=0;
578 : //Int_t statusOld = 0; //coverity (8443) fix satya (1/9/2014)
579 0 : if (mparticle->GetFirstMother() == -1)
580 : {
581 : tracknoOld = trackno;
582 : trackpidOld = trackpid;
583 : // statusOld = -1;
584 0 : }
585 :
586 : // Int_t igstatus = 0;
587 :
588 : Int_t trnotemp = trackno;
589 0 : if(ks==1||(imo = mparticle->GetFirstMother())<0 ){
590 : // vx = mparticle->Vx(); //coverity (8443) fix satya (1/9/2014)
591 : // vy = mparticle->Vy();
592 : // vz = mparticle->Vz();
593 :
594 : // if(trackpid==kGamma||trackpid==11||trackpid==-11||trackpid==kPi0)
595 : // igstatus=1;
596 : }
597 :
598 :
599 0 : while(((imo = mparticle->GetFirstMother()) >= 0)&&
600 0 : (ks = mparticle->GetStatusCode() <1) )
601 : {
602 0 : mparticle = gAlice->GetMCApp()->Particle(imo);
603 0 : trackpid = mparticle->GetPdgCode();
604 0 : ks = mparticle->GetStatusCode();
605 : // vx = mparticle->Vx(); //coverity (8443) fix satya (1/9/2014)
606 : // vy = mparticle->Vy();
607 : // vz = mparticle->Vz();
608 :
609 :
610 : trnotemp = trackno;
611 0 : if(trackpid == 111)
612 : {
613 : trackno = trnotemp;
614 0 : }
615 0 : if(trackpid != 111)
616 : {
617 : trackno=imo;
618 0 : }
619 : }
620 :
621 : // if(trackpid==kGamma||trackpid==11||trackpid==-11||trackpid==kPi0)
622 : // igstatus=1; //coverity (8443) fix satya (1/9/2014)
623 : mtrackpid=trackpid;
624 : mtrackno=trackno;
625 : trackpid=trackpidOld;
626 : trackno=tracknoOld;
627 :
628 0 : Float_t ptime = fPMDHit->GetTime()*1e6;
629 0 : if (ptime < 0. || ptime > 1.2) continue;
630 :
631 0 : xPos = fPMDHit->X();
632 0 : yPos = fPMDHit->Y();
633 0 : zPos = fPMDHit->Z();
634 0 : edep = fPMDHit->GetEnergy();
635 0 : Int_t vol1 = fPMDHit->GetVolume(1); // Column
636 0 : Int_t vol2 = fPMDHit->GetVolume(2); // Row
637 0 : Int_t vol7 = fPMDHit->GetVolume(4); // Serial Module No
638 :
639 : // -----------------------------------------//
640 : // In new geometry after adding electronics //
641 : // For Super Module 1 & 2 //
642 : // nrow = 48, ncol = 96 //
643 : // For Super Module 3 & 4 //
644 : // nrow = 96, ncol = 48 //
645 : // -----------------------------------------//
646 :
647 0 : if (vol7 < 24)
648 : {
649 : smnumber = vol7;
650 0 : }
651 : else
652 : {
653 0 : smnumber = vol7 - 24;
654 : }
655 0 : Int_t vol8 = smnumber/6 + 1; // fake supermodule
656 :
657 0 : if (vol8 == 1 || vol8 == 2)
658 : {
659 : xpad = vol2;
660 : ypad = vol1;
661 0 : }
662 0 : else if (vol8 == 3 || vol8 == 4)
663 : {
664 : xpad = vol1;
665 : ypad = vol2;
666 0 : }
667 :
668 0 : AliDebug(2,Form("ZPosition = %f Edeposition = %f",zPos,edep));
669 :
670 0 : if (vol7 < 24)
671 : {
672 : // PRE
673 0 : fDetNo = 0;
674 0 : }
675 : else
676 : {
677 0 : fDetNo = 1;
678 : }
679 :
680 : Int_t smn = smnumber;
681 0 : Int_t ixx = xpad - 1;
682 0 : Int_t iyy = ypad - 1;
683 0 : if (fDetNo == 0)
684 : {
685 0 : fPRE[smn][ixx][iyy] += edep;
686 0 : fPRECounter[smn][ixx][iyy]++;
687 :
688 0 : AliPMDcell* cell = new AliPMDcell(mtrackno,smn,ixx,iyy,edep);
689 :
690 0 : fCell.Add(cell);
691 0 : }
692 0 : else if(fDetNo == 1)
693 : {
694 0 : fCPV[smn][ixx][iyy] += edep;
695 0 : fCPVCounter[smn][ixx][iyy]++;
696 0 : AliPMDcell* cpvcell = new AliPMDcell(mtrackno,smn,ixx,iyy,edep);
697 0 : fCPVCell.Add(cpvcell);
698 0 : }
699 0 : }
700 0 : }
701 : } // Track Loop ended
702 0 : TrackAssignment2CPVCell();
703 0 : TrackAssignment2Cell();
704 0 : ResetCell();
705 :
706 : Float_t gain1 = 1.;
707 0 : Float_t adc = 0. ;
708 : Float_t deltaE = 0.;
709 : Int_t detno = 0;
710 : Int_t trno = 1;
711 : Int_t trpid = -99;
712 :
713 0 : for (Int_t idet = 0; idet < 2; idet++)
714 : {
715 0 : for (Int_t ism = 0; ism < fgkTotUM; ism++)
716 : {
717 0 : for (Int_t jrow = 0; jrow < fgkRow; jrow++)
718 : {
719 0 : for (Int_t kcol = 0; kcol < fgkCol; kcol++)
720 : {
721 0 : if (idet == 0)
722 : {
723 0 : deltaE = fPRE[ism][jrow][kcol];
724 0 : trno = fPRETrackNo[ism][jrow][kcol];
725 : detno = 0;
726 0 : }
727 0 : else if (idet == 1)
728 : {
729 0 : deltaE = fCPV[ism][jrow][kcol];
730 0 : trno = fCPVTrackNo[ism][jrow][kcol];
731 : detno = 1;
732 0 : }
733 0 : if (deltaE > 0.)
734 : {
735 0 : MeV2ADC(deltaE,adc);
736 :
737 : // To decalibrate the adc values
738 : //
739 0 : gain1 = Gain(idet,ism,jrow,kcol);
740 0 : if (gain1 != 0.)
741 : {
742 0 : Int_t adcDecalib = (Int_t)(adc/gain1);
743 0 : adc = (Float_t) adcDecalib;
744 0 : }
745 0 : else if(gain1 == 0.)
746 : {
747 0 : adc = 0.;
748 0 : }
749 :
750 : // Pedestal Decalibration
751 : Int_t pedmeanrms =
752 0 : fCalibPed->GetPedMeanRms(idet,ism,jrow,kcol);
753 0 : Int_t pedrms1 = (Int_t) pedmeanrms%100;
754 0 : Float_t pedrms = (Float_t)pedrms1/10.;
755 : Float_t pedmean =
756 0 : (Float_t) (pedmeanrms - pedrms1)/1000.0;
757 0 : if (adc > 0.)
758 : {
759 0 : adc += (pedmean + 3.0*pedrms);
760 : TParticle *mparticle
761 0 : = gAlice->GetMCApp()->Particle(trno);
762 0 : trpid = mparticle->GetPdgCode();
763 :
764 0 : AddDigit(trno,trpid,detno,ism,jrow,kcol,adc);
765 0 : }
766 0 : }
767 : } // column loop
768 : } // row loop
769 0 : treeD->Fill();
770 0 : ResetDigit();
771 : } // supermodule loop
772 : } // detector loop
773 :
774 0 : fPMDLoader->WriteDigits("OVERWRITE");
775 0 : ResetCellADC();
776 :
777 0 : }
778 : //____________________________________________________________________________
779 :
780 :
781 : void AliPMDDigitizer::SDigits2Digits(Int_t ievt)
782 : {
783 : // This reads the PMD sdigits tree and converts energy deposition
784 : // in a cell to ADC and stores in the digits tree
785 : //
786 :
787 0 : fRunLoader->GetEvent(ievt);
788 :
789 0 : TTree* treeS = fPMDLoader->TreeS();
790 : AliPMDsdigit *pmdsdigit;
791 0 : TBranch *branch = treeS->GetBranch("PMDSDigit");
792 0 : if(!branch)
793 : {
794 0 : AliError("PMD Sdigit branch does not exist");
795 0 : return;
796 : }
797 0 : if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
798 0 : branch->SetAddress(&fSDigits);
799 :
800 0 : TTree* treeD = fPMDLoader->TreeD();
801 0 : if (treeD == 0x0)
802 : {
803 0 : fPMDLoader->MakeTree("D");
804 0 : treeD = fPMDLoader->TreeD();
805 0 : }
806 : Int_t bufsize = 16000;
807 0 : if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
808 0 : treeD->Branch("PMDDigit", &fDigits, bufsize);
809 :
810 : Int_t trno = 1, trpid = 0, det = 0, smn = 0;
811 : Int_t irow = 0, icol = 0;
812 0 : Float_t edep = 0., adc = 0.;
813 :
814 0 : Int_t nmodules = (Int_t) treeS->GetEntries();
815 0 : AliDebug(1,Form("Number of modules = %d",nmodules));
816 :
817 0 : for (Int_t imodule = 0; imodule < nmodules; imodule++)
818 : {
819 0 : treeS->GetEntry(imodule);
820 0 : Int_t nentries = fSDigits->GetLast();
821 0 : AliDebug(2,Form("Number of entries per module = %d",nentries+1));
822 0 : for (Int_t ient = 0; ient < nentries+1; ient++)
823 : {
824 0 : pmdsdigit = (AliPMDsdigit*)fSDigits->UncheckedAt(ient);
825 0 : trno = pmdsdigit->GetTrackNumber();
826 0 : trpid = pmdsdigit->GetTrackPid();
827 0 : det = pmdsdigit->GetDetector();
828 0 : smn = pmdsdigit->GetSMNumber();
829 0 : irow = pmdsdigit->GetRow();
830 0 : icol = pmdsdigit->GetColumn();
831 0 : edep = pmdsdigit->GetCellEdep();
832 :
833 0 : MeV2ADC(edep,adc);
834 :
835 : // To decalibrte the adc values
836 : //
837 0 : Float_t gain1 = Gain(det,smn,irow,icol);
838 0 : if (gain1 != 0.)
839 : {
840 0 : Int_t adcDecalib = (Int_t)(adc/gain1);
841 0 : adc = (Float_t) adcDecalib;
842 0 : }
843 0 : else if(gain1 == 0.)
844 : {
845 0 : adc = 0.;
846 0 : }
847 : // Pedestal Decalibration
848 0 : Int_t pedmeanrms = fCalibPed->GetPedMeanRms(det,smn,irow,icol);
849 0 : Int_t pedrms1 = (Int_t) pedmeanrms%100;
850 0 : Float_t pedrms = (Float_t)pedrms1/10.;
851 0 : Float_t pedmean = (Float_t) (pedmeanrms - pedrms1)/1000.0;
852 0 : if(adc > 0.)
853 : {
854 0 : adc += (pedmean + 3.0*pedrms);
855 0 : AddDigit(trno,trpid,det,smn,irow,icol,adc);
856 0 : }
857 :
858 : }
859 0 : treeD->Fill();
860 0 : ResetDigit();
861 : }
862 0 : fPMDLoader->WriteDigits("OVERWRITE");
863 :
864 0 : }
865 : //____________________________________________________________________________
866 : void AliPMDDigitizer::Digitize(Option_t *option)
867 : {
868 : // Does the event merging and digitization
869 8 : const char *cdeb = strstr(option,"deb");
870 4 : if(cdeb)
871 : {
872 0 : AliDebug(100," *** PMD Exec is called ***");
873 : }
874 :
875 4 : Int_t ninputs = fDigInput->GetNinputs();
876 12 : AliDebug(1,Form("Number of files to be processed = %d",ninputs));
877 4 : ResetCellADC();
878 :
879 16 : for (Int_t i = 0; i < ninputs; i++)
880 : {
881 4 : Int_t troffset = fDigInput->GetMask(i);
882 4 : MergeSDigits(i, troffset);
883 : }
884 :
885 4 : fRunLoader = AliRunLoader::GetRunLoader(fDigInput->GetOutputFolderName());
886 4 : fPMD = (AliPMD*)gAlice->GetDetector("PMD");
887 4 : fPMDLoader = fRunLoader->GetLoader("PMDLoader");
888 4 : if (fPMDLoader == 0x0)
889 : {
890 0 : AliError("Can not find PMD or PMDLoader");
891 0 : }
892 4 : fPMDLoader->LoadDigits("update");
893 4 : TTree* treeD = fPMDLoader->TreeD();
894 4 : if (treeD == 0x0)
895 : {
896 4 : fPMDLoader->MakeTree("D");
897 4 : treeD = fPMDLoader->TreeD();
898 4 : }
899 : Int_t bufsize = 16000;
900 4 : if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
901 4 : treeD->Branch("PMDDigit", &fDigits, bufsize);
902 :
903 4 : Float_t adc = 0.;
904 : Float_t deltaE = 0.;
905 : Int_t detno = 0;
906 : Int_t trno = 1;
907 : Int_t trpid = -99;
908 :
909 24 : for (Int_t idet = 0; idet < 2; idet++)
910 : {
911 400 : for (Int_t ism = 0; ism < fgkTotUM; ism++)
912 : {
913 18816 : for (Int_t jrow = 0; jrow < fgkRow; jrow++)
914 : {
915 1787904 : for (Int_t kcol = 0; kcol < fgkCol; kcol++)
916 : {
917 884736 : if (idet == 0)
918 : {
919 442368 : deltaE = fPRE[ism][jrow][kcol];
920 442368 : trno = fPRETrackNo[ism][jrow][kcol];
921 442368 : trpid = fPRETrackPid[ism][jrow][kcol];
922 : detno = 0;
923 442368 : }
924 442368 : else if (idet == 1)
925 : {
926 442368 : deltaE = fCPV[ism][jrow][kcol];
927 442368 : trno = fCPVTrackNo[ism][jrow][kcol];
928 442368 : trpid = fCPVTrackPid[ism][jrow][kcol];
929 : detno = 1;
930 442368 : }
931 884736 : if (deltaE > 0.)
932 : {
933 221 : MeV2ADC(deltaE,adc);
934 :
935 : //
936 : // Gain decalibration
937 : //
938 221 : Float_t gain1 = Gain(idet,ism,jrow,kcol);
939 :
940 221 : if (gain1 != 0.)
941 : {
942 180 : Int_t adcDecalib = (Int_t)(adc/gain1);
943 180 : adc = (Float_t) adcDecalib;
944 180 : }
945 41 : else if(gain1 == 0.)
946 : {
947 41 : adc = 0.;
948 41 : }
949 : // Pedestal Decalibration
950 : Int_t pedmeanrms =
951 221 : fCalibPed->GetPedMeanRms(idet,ism,jrow,kcol);
952 221 : Int_t pedrms1 = (Int_t) pedmeanrms%100;
953 221 : Float_t pedrms = (Float_t)pedrms1/10.;
954 : Float_t pedmean =
955 221 : (Float_t) (pedmeanrms - pedrms1)/1000.0;
956 221 : if (adc > 0.)
957 : {
958 180 : adc += (pedmean + 3.0*pedrms);
959 180 : AddDigit(trno,trpid,detno,ism,jrow,kcol,adc);
960 180 : }
961 :
962 221 : }
963 : } // column loop
964 : } // row loop
965 192 : treeD->Fill();
966 192 : ResetDigit();
967 : } // supermodule loop
968 : } // detector loop
969 4 : fPMDLoader->WriteDigits("OVERWRITE");
970 4 : fPMDLoader->UnloadDigits();
971 4 : ResetCellADC();
972 4 : }
973 : //____________________________________________________________________________
974 : void AliPMDDigitizer::TrackAssignment2CPVCell()
975 : {
976 : // This block assigns the cell id when there are
977 : // multiple tracks in a cell according to the
978 : // energy deposition
979 : // This method added by Ajay
980 : Bool_t jsort = false;
981 :
982 : Int_t i = 0, j = 0, k = 0;
983 :
984 : Int_t *status1;
985 : Int_t *status2;
986 : Int_t *trnarray;
987 : Float_t *fracEdp;
988 : Float_t *trEdp;
989 :
990 : Int_t ****cpvTrack;
991 : Float_t ****cpvEdep;
992 :
993 8 : cpvTrack = new Int_t ***[fgkTotUM];
994 4 : cpvEdep = new Float_t ***[fgkTotUM];
995 200 : for (i=0; i<fgkTotUM; i++)
996 : {
997 96 : cpvTrack[i] = new Int_t **[fgkRow];
998 96 : cpvEdep[i] = new Float_t **[fgkRow];
999 : }
1000 :
1001 200 : for (i = 0; i < fgkTotUM; i++)
1002 : {
1003 9408 : for (j = 0; j < fgkRow; j++)
1004 : {
1005 4608 : cpvTrack[i][j] = new Int_t *[fgkCol];
1006 4608 : cpvEdep[i][j] = new Float_t *[fgkCol];
1007 : }
1008 : }
1009 200 : for (i = 0; i < fgkTotUM; i++)
1010 : {
1011 9408 : for (j = 0; j < fgkRow; j++)
1012 : {
1013 893952 : for (k = 0; k < fgkCol; k++)
1014 : {
1015 442368 : Int_t nn = fCPVCounter[i][j][k];
1016 442368 : if(nn > 0)
1017 : {
1018 72 : cpvTrack[i][j][k] = new Int_t[nn];
1019 72 : cpvEdep[i][j][k] = new Float_t[nn];
1020 72 : }
1021 : else
1022 : {
1023 : nn = 1;
1024 442296 : cpvTrack[i][j][k] = new Int_t[nn];
1025 442296 : cpvEdep[i][j][k] = new Float_t[nn];
1026 : }
1027 442368 : fCPVCounter[i][j][k] = 0;
1028 : }
1029 : }
1030 : }
1031 :
1032 :
1033 4 : Int_t nentries = fCPVCell.GetEntries();
1034 :
1035 : Int_t mtrackno = 0, ism = 0, ixp = 0, iyp = 0;
1036 : Float_t edep = 0.;
1037 166 : for (i = 0; i < nentries; i++)
1038 : {
1039 79 : AliPMDcell* cpvcell = (AliPMDcell*)fCPVCell.UncheckedAt(i);
1040 :
1041 79 : mtrackno = cpvcell->GetTrackNumber();
1042 79 : ism = cpvcell->GetSMNumber();
1043 79 : ixp = cpvcell->GetX();
1044 79 : iyp = cpvcell->GetY();
1045 79 : edep = cpvcell->GetEdep();
1046 79 : Int_t nn = fCPVCounter[ism][ixp][iyp];
1047 79 : cpvTrack[ism][ixp][iyp][nn] = (Int_t) mtrackno;
1048 79 : cpvEdep[ism][ixp][iyp][nn] = edep;
1049 79 : fCPVCounter[ism][ixp][iyp]++;
1050 : }
1051 :
1052 : Int_t iz = 0, il = 0;
1053 : Int_t im = 0, ix = 0, iy = 0;
1054 : Int_t nn = 0;
1055 200 : for (im=0; im<fgkTotUM; im++)
1056 : {
1057 9408 : for (ix=0; ix<fgkRow; ix++)
1058 : {
1059 893952 : for (iy=0; iy<fgkCol; iy++)
1060 : {
1061 442368 : nn = fCPVCounter[im][ix][iy];
1062 442368 : if (nn > 1)
1063 : {
1064 : // This block handles if a cell is fired
1065 : // many times by many tracks
1066 7 : status1 = new Int_t[nn];
1067 7 : status2 = new Int_t[2*nn];
1068 7 : trnarray = new Int_t[nn];
1069 42 : for (iz = 0; iz < nn; iz++)
1070 : {
1071 14 : status1[iz] = cpvTrack[im][ix][iy][iz];
1072 : }
1073 7 : TMath::Sort(nn,status1,status2,jsort);
1074 : Int_t trackOld = -99999;
1075 : Int_t track, trCount = 0;
1076 42 : for (iz = 0; iz < nn; iz++)
1077 : {
1078 14 : track = status1[status2[iz]];
1079 14 : if (trackOld != track)
1080 : {
1081 7 : trnarray[trCount] = track;
1082 7 : trCount++;
1083 7 : }
1084 : trackOld = track;
1085 : }
1086 14 : delete [] status1;
1087 14 : delete [] status2;
1088 : Float_t totEdp = 0.;
1089 7 : trEdp = new Float_t[trCount];
1090 7 : fracEdp = new Float_t[trCount];
1091 28 : for (il = 0; il < trCount; il++)
1092 : {
1093 7 : trEdp[il] = 0.;
1094 7 : track = trnarray[il];
1095 42 : for (iz = 0; iz < nn; iz++)
1096 : {
1097 14 : if (track == cpvTrack[im][ix][iy][iz])
1098 : {
1099 14 : trEdp[il] += cpvEdep[im][ix][iy][iz];
1100 14 : }
1101 : }
1102 7 : totEdp += trEdp[il];
1103 : }
1104 : Int_t ilOld = 0;
1105 : Float_t fracOld = 0.;
1106 :
1107 28 : for (il = 0; il < trCount; il++)
1108 : {
1109 7 : fracEdp[il] = trEdp[il]/totEdp;
1110 7 : if (fracOld < fracEdp[il])
1111 : {
1112 : fracOld = fracEdp[il];
1113 : ilOld = il;
1114 7 : }
1115 : }
1116 7 : fCPVTrackNo[im][ix][iy] = trnarray[ilOld];
1117 14 : delete [] fracEdp;
1118 14 : delete [] trEdp;
1119 14 : delete [] trnarray;
1120 7 : }
1121 442361 : else if (nn == 1)
1122 : {
1123 : // This only handles if a cell is fired
1124 : // by only one track
1125 :
1126 65 : fCPVTrackNo[im][ix][iy] = cpvTrack[im][ix][iy][0];
1127 :
1128 65 : }
1129 442296 : else if (nn ==0)
1130 : {
1131 : // This is if no cell is fired
1132 442296 : fCPVTrackNo[im][ix][iy] = -999;
1133 442296 : }
1134 : } // end of iy
1135 : } // end of ix
1136 : } // end of im
1137 :
1138 : // Delete all the pointers
1139 :
1140 200 : for (i = 0; i < fgkTotUM; i++)
1141 : {
1142 9408 : for (j = 0; j < fgkRow; j++)
1143 : {
1144 893952 : for (k = 0; k < fgkCol; k++)
1145 : {
1146 884736 : delete []cpvTrack[i][j][k];
1147 884736 : delete []cpvEdep[i][j][k];
1148 : }
1149 : }
1150 : }
1151 :
1152 200 : for (i = 0; i < fgkTotUM; i++)
1153 : {
1154 9408 : for (j = 0; j < fgkRow; j++)
1155 : {
1156 9216 : delete [] cpvTrack[i][j];
1157 9216 : delete [] cpvEdep[i][j];
1158 : }
1159 : }
1160 :
1161 200 : for (i = 0; i < fgkTotUM; i++)
1162 : {
1163 192 : delete [] cpvTrack[i];
1164 192 : delete [] cpvEdep[i];
1165 : }
1166 8 : delete [] cpvTrack;
1167 8 : delete [] cpvEdep;
1168 :
1169 : //
1170 : // End of the cell id assignment
1171 : //
1172 4 : }
1173 : //____________________________________________________________________________
1174 :
1175 : void AliPMDDigitizer::MergeSDigits(Int_t filenumber, Int_t troffset)
1176 : {
1177 : // merging sdigits
1178 8 : fRunLoader = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(filenumber));
1179 4 : fPMDLoader = fRunLoader->GetLoader("PMDLoader");
1180 4 : fPMDLoader->LoadSDigits("read");
1181 4 : TTree* treeS = fPMDLoader->TreeS();
1182 : AliPMDsdigit *pmdsdigit;
1183 4 : TBranch *branch = treeS->GetBranch("PMDSDigit");
1184 4 : if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
1185 4 : branch->SetAddress(&fSDigits);
1186 :
1187 : Int_t itrackno = 1, itrackpid = 0, idet = 0, ism = 0;
1188 : Int_t ixp = 0, iyp = 0;
1189 : Float_t edep = 0.;
1190 4 : Int_t nmodules = (Int_t) treeS->GetEntries();
1191 12 : AliDebug(1,Form("Number of Modules in the treeS = %d",nmodules));
1192 12 : AliDebug(1,Form("Track Offset = %d",troffset));
1193 392 : for (Int_t imodule = 0; imodule < nmodules; imodule++)
1194 : {
1195 192 : treeS->GetEntry(imodule);
1196 192 : Int_t nentries = fSDigits->GetLast();
1197 576 : AliDebug(2,Form("Number of Entries per Module = %d",nentries));
1198 826 : for (Int_t ient = 0; ient < nentries+1; ient++)
1199 : {
1200 221 : pmdsdigit = (AliPMDsdigit*)fSDigits->UncheckedAt(ient);
1201 221 : itrackno = pmdsdigit->GetTrackNumber();
1202 221 : itrackpid = pmdsdigit->GetTrackPid();
1203 221 : idet = pmdsdigit->GetDetector();
1204 221 : ism = pmdsdigit->GetSMNumber();
1205 221 : ixp = pmdsdigit->GetRow();
1206 221 : iyp = pmdsdigit->GetColumn();
1207 221 : edep = pmdsdigit->GetCellEdep();
1208 221 : if (idet == 0)
1209 : {
1210 149 : if (fPRE[ism][ixp][iyp] < edep)
1211 : {
1212 149 : fPRETrackNo[ism][ixp][iyp] = troffset + itrackno;
1213 149 : fPRETrackPid[ism][ixp][iyp] = itrackpid;
1214 149 : }
1215 149 : fPRE[ism][ixp][iyp] += edep;
1216 149 : }
1217 72 : else if (idet == 1)
1218 : {
1219 72 : if (fCPV[ism][ixp][iyp] < edep)
1220 : {
1221 72 : fCPVTrackNo[ism][ixp][iyp] = troffset + itrackno;
1222 72 : fCPVTrackPid[ism][ixp][iyp] = itrackpid;
1223 72 : }
1224 72 : fCPV[ism][ixp][iyp] += edep;
1225 72 : }
1226 : }
1227 : }
1228 :
1229 4 : }
1230 : // ----------------------------------------------------------------------
1231 : void AliPMDDigitizer::TrackAssignment2Cell()
1232 : {
1233 : //
1234 : // This block assigns the cell id when there are
1235 : // multiple tracks in a cell according to the
1236 : // energy deposition
1237 : //
1238 : Bool_t jsort = false;
1239 :
1240 : Int_t i = 0, j = 0, k = 0;
1241 :
1242 : Int_t *status1;
1243 : Int_t *status2;
1244 : Int_t *trnarray;
1245 : Float_t *fracEdp;
1246 : Float_t *trEdp;
1247 :
1248 : Int_t ****pmdTrack;
1249 : Float_t ****pmdEdep;
1250 :
1251 8 : pmdTrack = new Int_t ***[fgkTotUM];
1252 4 : pmdEdep = new Float_t ***[fgkTotUM];
1253 200 : for (i=0; i<fgkTotUM; i++)
1254 : {
1255 96 : pmdTrack[i] = new Int_t **[fgkRow];
1256 96 : pmdEdep[i] = new Float_t **[fgkRow];
1257 : }
1258 :
1259 200 : for (i = 0; i < fgkTotUM; i++)
1260 : {
1261 9408 : for (j = 0; j < fgkRow; j++)
1262 : {
1263 4608 : pmdTrack[i][j] = new Int_t *[fgkCol];
1264 4608 : pmdEdep[i][j] = new Float_t *[fgkCol];
1265 : }
1266 : }
1267 :
1268 200 : for (i = 0; i < fgkTotUM; i++)
1269 : {
1270 9408 : for (j = 0; j < fgkRow; j++)
1271 : {
1272 893952 : for (k = 0; k < fgkCol; k++)
1273 : {
1274 442368 : Int_t nn = fPRECounter[i][j][k];
1275 442368 : if(nn > 0)
1276 : {
1277 149 : pmdTrack[i][j][k] = new Int_t[nn];
1278 149 : pmdEdep[i][j][k] = new Float_t[nn];
1279 149 : }
1280 : else
1281 : {
1282 : nn = 1;
1283 442219 : pmdTrack[i][j][k] = new Int_t[nn];
1284 442219 : pmdEdep[i][j][k] = new Float_t[nn];
1285 : }
1286 442368 : fPRECounter[i][j][k] = 0;
1287 : }
1288 : }
1289 : }
1290 :
1291 :
1292 4 : Int_t nentries = fCell.GetEntries();
1293 :
1294 : Int_t mtrackno, ism, ixp, iyp;
1295 : Float_t edep;
1296 :
1297 772 : for (i = 0; i < nentries; i++)
1298 : {
1299 382 : AliPMDcell* cell = (AliPMDcell*)fCell.UncheckedAt(i);
1300 :
1301 382 : mtrackno = cell->GetTrackNumber();
1302 382 : ism = cell->GetSMNumber();
1303 382 : ixp = cell->GetX();
1304 382 : iyp = cell->GetY();
1305 382 : edep = cell->GetEdep();
1306 382 : Int_t nn = fPRECounter[ism][ixp][iyp];
1307 382 : pmdTrack[ism][ixp][iyp][nn] = (Int_t) mtrackno;
1308 382 : pmdEdep[ism][ixp][iyp][nn] = edep;
1309 382 : fPRECounter[ism][ixp][iyp]++;
1310 : }
1311 :
1312 : Int_t iz = 0, il = 0;
1313 : Int_t im = 0, ix = 0, iy = 0;
1314 : Int_t nn = 0;
1315 :
1316 200 : for (im=0; im<fgkTotUM; im++)
1317 : {
1318 9408 : for (ix=0; ix<fgkRow; ix++)
1319 : {
1320 893952 : for (iy=0; iy<fgkCol; iy++)
1321 : {
1322 442368 : nn = fPRECounter[im][ix][iy];
1323 442368 : if (nn > 1)
1324 : {
1325 : // This block handles if a cell is fired
1326 : // many times by many tracks
1327 42 : status1 = new Int_t[nn];
1328 42 : status2 = new Int_t[2*nn];
1329 42 : trnarray = new Int_t[nn];
1330 634 : for (iz = 0; iz < nn; iz++)
1331 : {
1332 275 : status1[iz] = pmdTrack[im][ix][iy][iz];
1333 : }
1334 42 : TMath::Sort(nn,status1,status2,jsort);
1335 : Int_t trackOld = -99999;
1336 : Int_t track, trCount = 0;
1337 634 : for (iz = 0; iz < nn; iz++)
1338 : {
1339 275 : track = status1[status2[iz]];
1340 275 : if (trackOld != track)
1341 : {
1342 43 : trnarray[trCount] = track;
1343 43 : trCount++;
1344 43 : }
1345 : trackOld = track;
1346 : }
1347 84 : delete [] status1;
1348 84 : delete [] status2;
1349 : Float_t totEdp = 0.;
1350 42 : trEdp = new Float_t[trCount];
1351 42 : fracEdp = new Float_t[trCount];
1352 170 : for (il = 0; il < trCount; il++)
1353 : {
1354 43 : trEdp[il] = 0.;
1355 43 : track = trnarray[il];
1356 642 : for (iz = 0; iz < nn; iz++)
1357 : {
1358 278 : if (track == pmdTrack[im][ix][iy][iz])
1359 : {
1360 275 : trEdp[il] += pmdEdep[im][ix][iy][iz];
1361 275 : }
1362 : }
1363 43 : totEdp += trEdp[il];
1364 : }
1365 : Int_t ilOld = 0;
1366 : Float_t fracOld = 0.;
1367 :
1368 170 : for (il = 0; il < trCount; il++)
1369 : {
1370 43 : fracEdp[il] = trEdp[il]/totEdp;
1371 43 : if (fracOld < fracEdp[il])
1372 : {
1373 : fracOld = fracEdp[il];
1374 : ilOld = il;
1375 42 : }
1376 : }
1377 42 : fPRETrackNo[im][ix][iy] = trnarray[ilOld];
1378 84 : delete [] fracEdp;
1379 84 : delete [] trEdp;
1380 84 : delete [] trnarray;
1381 42 : }
1382 442326 : else if (nn == 1)
1383 : {
1384 : // This only handles if a cell is fired
1385 : // by only one track
1386 :
1387 107 : fPRETrackNo[im][ix][iy] = pmdTrack[im][ix][iy][0];
1388 :
1389 107 : }
1390 442219 : else if (nn ==0)
1391 : {
1392 : // This is if no cell is fired
1393 442219 : fPRETrackNo[im][ix][iy] = -999;
1394 442219 : }
1395 : } // end of iy
1396 : } // end of ix
1397 : } // end of im
1398 :
1399 : // Delete all the pointers
1400 :
1401 200 : for (i = 0; i < fgkTotUM; i++)
1402 : {
1403 9408 : for (j = 0; j < fgkRow; j++)
1404 : {
1405 893952 : for (k = 0; k < fgkCol; k++)
1406 : {
1407 884736 : delete [] pmdTrack[i][j][k];
1408 884736 : delete [] pmdEdep[i][j][k];
1409 : }
1410 : }
1411 : }
1412 :
1413 200 : for (i = 0; i < fgkTotUM; i++)
1414 : {
1415 9408 : for (j = 0; j < fgkRow; j++)
1416 : {
1417 9216 : delete [] pmdTrack[i][j];
1418 9216 : delete [] pmdEdep[i][j];
1419 : }
1420 : }
1421 :
1422 200 : for (i = 0; i < fgkTotUM; i++)
1423 : {
1424 192 : delete [] pmdTrack[i];
1425 192 : delete [] pmdEdep[i];
1426 : }
1427 8 : delete [] pmdTrack;
1428 8 : delete [] pmdEdep;
1429 : //
1430 : // End of the cell id assignment
1431 : //
1432 4 : }
1433 : //____________________________________________________________________________
1434 : void AliPMDDigitizer::MeV2ADC(Float_t mev, Float_t & adc) const
1435 : {
1436 : // This converts the simulated edep to ADC according to the
1437 : // Test Beam Data
1438 : // PS Test in June 2010, Voltage @ 1300 V
1439 : // KeV - ADC conversion for 12bit ADC
1440 : // MPV data used for the fit and taken here
1441 :
1442 : // constants are from Test Beam 2010
1443 :
1444 : const Float_t kConstant = 0.612796;
1445 : const Float_t kSlope = 130.158;
1446 :
1447 442 : Float_t adc12bit = kSlope*mev*0.001 + kConstant;
1448 221 : if (adc12bit < 0.) adc12bit = 0.;
1449 :
1450 : //Introducing Readout Resolution for ALICE-PMD
1451 :
1452 221 : Float_t sigrr = 0.605016 - 0.000273*adc12bit + 6.54e-8*adc12bit*adc12bit;
1453 221 : Float_t adcwithrr = gRandom->Gaus(adc12bit,sigrr);
1454 :
1455 221 : if(adcwithrr < 0.)
1456 : {
1457 0 : adc = 0.;
1458 0 : }
1459 442 : else if(adcwithrr >= 0. && adcwithrr < 1600.0)
1460 : {
1461 212 : adc = adcwithrr;
1462 212 : }
1463 9 : else if (adcwithrr >= 1600.0)
1464 : {
1465 9 : adc = 1600.0;
1466 9 : }
1467 :
1468 221 : }
1469 : //____________________________________________________________________________
1470 : void AliPMDDigitizer::AddSDigit(Int_t trnumber, Int_t trpid, Int_t det,
1471 : Int_t smnumber, Int_t irow, Int_t icol,
1472 : Float_t adc)
1473 : {
1474 : // Add SDigit
1475 : //
1476 442 : if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
1477 221 : TClonesArray &lsdigits = *fSDigits;
1478 221 : new(lsdigits[fNsdigit++]) AliPMDsdigit(trnumber,trpid,det,smnumber,irow,icol,adc);
1479 221 : }
1480 : //____________________________________________________________________________
1481 :
1482 : void AliPMDDigitizer::AddDigit(Int_t trnumber, Int_t trpid, Int_t det,
1483 : Int_t smnumber, Int_t irow, Int_t icol,
1484 : Float_t adc)
1485 : {
1486 : // Add Digit
1487 : //
1488 360 : if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
1489 180 : TClonesArray &ldigits = *fDigits;
1490 180 : new(ldigits[fNdigit++]) AliPMDdigit(trnumber,trpid, det,smnumber,irow,icol,adc);
1491 180 : }
1492 : //____________________________________________________________________________
1493 :
1494 : void AliPMDDigitizer::SetZPosition(Float_t zpos)
1495 : {
1496 2 : fZPos = zpos;
1497 1 : }
1498 : //____________________________________________________________________________
1499 : Float_t AliPMDDigitizer::GetZPosition() const
1500 : {
1501 0 : return fZPos;
1502 : }
1503 : //____________________________________________________________________________
1504 :
1505 : void AliPMDDigitizer::ResetCell()
1506 : {
1507 : // clears the cell array and also the counter
1508 : // for each cell
1509 : //
1510 8 : fCPVCell.Delete();
1511 4 : fCell.Delete();
1512 200 : for (Int_t i = 0; i < fgkTotUM; i++)
1513 : {
1514 9408 : for (Int_t j = 0; j < fgkRow; j++)
1515 : {
1516 893952 : for (Int_t k = 0; k < fgkCol; k++)
1517 : {
1518 442368 : fCPVCounter[i][j][k] = 0;
1519 442368 : fPRECounter[i][j][k] = 0;
1520 : }
1521 : }
1522 : }
1523 4 : }
1524 : //____________________________________________________________________________
1525 : void AliPMDDigitizer::ResetSDigit()
1526 : {
1527 : // Clears SDigits
1528 392 : fNsdigit = 0;
1529 392 : if (fSDigits) fSDigits->Delete();
1530 196 : }
1531 : //____________________________________________________________________________
1532 : void AliPMDDigitizer::ResetDigit()
1533 : {
1534 : // Clears Digits
1535 384 : fNdigit = 0;
1536 384 : if (fDigits) fDigits->Delete();
1537 192 : }
1538 : //____________________________________________________________________________
1539 :
1540 : void AliPMDDigitizer::ResetCellADC()
1541 : {
1542 : // Clears individual cells edep and track number
1543 612 : for (Int_t i = 0; i < fgkTotUM; i++)
1544 : {
1545 28224 : for (Int_t j = 0; j < fgkRow; j++)
1546 : {
1547 2681856 : for (Int_t k = 0; k < fgkCol; k++)
1548 : {
1549 1327104 : fCPV[i][j][k] = 0.;
1550 1327104 : fPRE[i][j][k] = 0.;
1551 1327104 : fCPVTrackNo[i][j][k] = 0;
1552 1327104 : fPRETrackNo[i][j][k] = 0;
1553 1327104 : fCPVTrackPid[i][j][k] = -1;
1554 1327104 : fPRETrackPid[i][j][k] = -1;
1555 : }
1556 : }
1557 : }
1558 12 : }
1559 : //____________________________________________________________________________
1560 :
1561 : void AliPMDDigitizer::UnLoad(Option_t *option)
1562 : {
1563 : // Unloads all the root files
1564 : //
1565 0 : const char *cS = strstr(option,"S");
1566 0 : const char *cD = strstr(option,"D");
1567 :
1568 0 : fRunLoader->UnloadgAlice();
1569 0 : fRunLoader->UnloadHeader();
1570 0 : fRunLoader->UnloadKinematics();
1571 :
1572 0 : if (cS)
1573 : {
1574 0 : fPMDLoader->UnloadHits();
1575 0 : }
1576 0 : if (cD)
1577 : {
1578 0 : fPMDLoader->UnloadHits();
1579 0 : fPMDLoader->UnloadSDigits();
1580 0 : }
1581 0 : }
1582 :
1583 : //----------------------------------------------------------------------
1584 : Float_t AliPMDDigitizer::Gain(Int_t det, Int_t smn, Int_t row, Int_t col) const
1585 : {
1586 : // returns of the gain of the cell
1587 : // Added this method by ZA
1588 :
1589 : //cout<<" I am here in gain "<<fCalibData<< "smn,row, col "<<smn
1590 : //<<" "<<row<<" "<<col<<endl;
1591 :
1592 442 : if(!fCalibGain) {
1593 0 : AliError("No calibration data loaded from CDB!!!");
1594 0 : return 1;
1595 : }
1596 :
1597 : Float_t GainFact;
1598 221 : GainFact = fCalibGain->GetGainFact(det,smn,row,col);
1599 : return GainFact;
1600 221 : }
1601 : //----------------------------------------------------------------------
1602 : AliPMDCalibData* AliPMDDigitizer::GetCalibGain() const
1603 : {
1604 : // The run number will be centralized in AliCDBManager,
1605 : // you don't need to set it here!
1606 : // Added this method by ZA
1607 : // Cleaned up by Alberto
1608 6 : AliCDBEntry *entry = AliCDBManager::Instance()->Get("PMD/Calib/Gain");
1609 :
1610 2 : if(!entry) AliFatal("Calibration object retrieval failed!");
1611 :
1612 : AliPMDCalibData *calibdata=0;
1613 4 : if (entry) calibdata = (AliPMDCalibData*) entry->GetObject();
1614 :
1615 2 : if (!calibdata) AliFatal("No calibration data from calibration database !");
1616 :
1617 2 : return calibdata;
1618 0 : }
1619 : //----------------------------------------------------------------------
1620 : AliPMDPedestal* AliPMDDigitizer::GetCalibPed() const
1621 : {
1622 : // The run number will be centralized in AliCDBManager,
1623 : // you don't need to set it here!
1624 :
1625 6 : AliCDBEntry *entry = AliCDBManager::Instance()->Get("PMD/Calib/Ped");
1626 :
1627 2 : if(!entry) AliFatal("Pedestal object retrieval failed!");
1628 :
1629 : AliPMDPedestal *pedestal=0;
1630 4 : if (entry) pedestal = (AliPMDPedestal*) entry->GetObject();
1631 :
1632 2 : if (!pedestal) AliFatal("No pedestal data from calibration database !");
1633 :
1634 2 : return pedestal;
1635 0 : }
|