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