Removed some obsolete TOF classes
[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));
5f7ce82b 419 //Float_t zposition = TMath::Abs(zPos);
f117e3aa 420 if (zPos < fZPos)
b2b13ee2 421 {
422 // CPV
423 fDetNo = 1;
424 }
f117e3aa 425 else if (zPos > fZPos)
b2b13ee2 426 {
427 // PMD
428 fDetNo = 0;
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
f117e3aa 664 if (zPos < fZPos)
b2b13ee2 665 {
666 // CPV
667 fDetNo = 1;
668 }
f117e3aa 669 else if (zPos > fZPos)
b2b13ee2 670 {
671 // PMD
672 fDetNo = 0;
673 }
8599803d 674
24c5571f 675 Int_t smn = smnumber;
8599803d 676 Int_t ixx = xpad - 1;
677 Int_t iyy = ypad - 1;
678 if (fDetNo == 0)
b2b13ee2 679 {
36031625 680 fPRE[smn][ixx][iyy] += edep;
681 fPRECounter[smn][ixx][iyy]++;
8599803d 682
ebd83c56 683 AliPMDcell* cell = new AliPMDcell(mtrackno,smn,ixx,iyy,edep);
8599803d 684
ebd83c56 685 fCell.Add(cell);
b2b13ee2 686 }
8599803d 687 else if(fDetNo == 1)
b2b13ee2 688 {
8599803d 689 fCPV[smn][ixx][iyy] += edep;
4d59e381 690 fCPVCounter[smn][ixx][iyy]++;
691 AliPMDcell* cpvcell = new AliPMDcell(mtrackno,smn,ixx,iyy,edep);
692 fCPVCell.Add(cpvcell);
b2b13ee2 693 }
694 }
695 }
696 } // Track Loop ended
4d59e381 697 TrackAssignment2CPVCell();
b2b13ee2 698 TrackAssignment2Cell();
699 ResetCell();
700
50555ba1 701 Float_t gain1;
7088dbe9 702 Float_t adc;
b2b13ee2 703 Float_t deltaE = 0.;
704 Int_t detno = 0;
b2b13ee2 705 Int_t trno = 1;
920e13db 706 Int_t trpid = -99;
707
b2b13ee2 708 for (Int_t idet = 0; idet < 2; idet++)
5e2cc5bc 709 {
36031625 710 for (Int_t ism = 0; ism < fgkTotUM; ism++)
5e2cc5bc 711 {
36031625 712 for (Int_t jrow = 0; jrow < fgkRow; jrow++)
5e2cc5bc 713 {
36031625 714 for (Int_t kcol = 0; kcol < fgkCol; kcol++)
5e2cc5bc 715 {
b2b13ee2 716 if (idet == 0)
5e2cc5bc 717 {
718 deltaE = fPRE[ism][jrow][kcol];
36031625 719 trno = fPRETrackNo[ism][jrow][kcol];
920e13db 720 detno = 0;
5e2cc5bc 721 }
b2b13ee2 722 else if (idet == 1)
5e2cc5bc 723 {
724 deltaE = fCPV[ism][jrow][kcol];
8599803d 725 trno = fCPVTrackNo[ism][jrow][kcol];
920e13db 726 detno = 1;
5e2cc5bc 727 }
b2b13ee2 728 if (deltaE > 0.)
5e2cc5bc 729 {
7088dbe9 730 MeV2ADC(deltaE,adc);
5e2cc5bc 731
732 // To decalibrate the adc values
733 //
734 gain1 = Gain(idet,ism,jrow,kcol);
735 if (gain1 != 0.)
736 {
737 Int_t adcDecalib = (Int_t)(adc/gain1);
738 adc = (Float_t) adcDecalib;
739 }
740 else if(gain1 == 0.)
741 {
742 adc = 0.;
743 }
35535af7 744
745 // Pedestal Decalibration
746 Int_t pedmeanrms =
747 fCalibPed->GetPedMeanRms(idet,ism,jrow,kcol);
37edc588 748 Int_t pedrms1 = (Int_t) pedmeanrms%100;
35535af7 749 Float_t pedrms = (Float_t)pedrms1/10.;
750 Float_t pedmean =
751 (Float_t) (pedmeanrms - pedrms1)/1000.0;
35535af7 752 if (adc > 0.)
753 {
754 adc += (pedmean + 3.0*pedrms);
8c24e5f2 755 TParticle *mparticle
756 = gAlice->GetMCApp()->Particle(trno);
757 trpid = mparticle->GetPdgCode();
758
920e13db 759 AddDigit(trno,trpid,detno,ism,jrow,kcol,adc);
35535af7 760 }
5e2cc5bc 761 }
762 } // column loop
763 } // row loop
03965b13 764 treeD->Fill();
765 ResetDigit();
5e2cc5bc 766 } // supermodule loop
767 } // detector loop
03965b13 768
36031625 769 fPMDLoader->WriteDigits("OVERWRITE");
b2b13ee2 770 ResetCellADC();
920e13db 771
b2b13ee2 772}
7088dbe9 773//____________________________________________________________________________
b2b13ee2 774
775
776void AliPMDDigitizer::SDigits2Digits(Int_t ievt)
777{
36031625 778 // This reads the PMD sdigits tree and converts energy deposition
779 // in a cell to ADC and stores in the digits tree
780 //
ecee2a1a 781
b2b13ee2 782 fRunLoader->GetEvent(ievt);
783
ebd83c56 784 TTree* treeS = fPMDLoader->TreeS();
b2b13ee2 785 AliPMDsdigit *pmdsdigit;
ebd83c56 786 TBranch *branch = treeS->GetBranch("PMDSDigit");
ecee2a1a 787 if(!branch)
788 {
789 AliError("PMD Sdigit branch does not exist");
790 return;
791 }
ebd83c56 792 if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
b2b13ee2 793 branch->SetAddress(&fSDigits);
794
ebd83c56 795 TTree* treeD = fPMDLoader->TreeD();
796 if (treeD == 0x0)
b2b13ee2 797 {
36031625 798 fPMDLoader->MakeTree("D");
ebd83c56 799 treeD = fPMDLoader->TreeD();
b2b13ee2 800 }
801 Int_t bufsize = 16000;
ebd83c56 802 if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
5e2cc5bc 803 treeD->Branch("PMDDigit", &fDigits, bufsize);
b2b13ee2 804
920e13db 805 Int_t trno, trpid, det, smn;
7088dbe9 806 Int_t irow, icol;
b2b13ee2 807 Float_t edep, adc;
808
ebd83c56 809 Int_t nmodules = (Int_t) treeS->GetEntries();
ecee2a1a 810 AliDebug(1,Form("Number of modules = %d",nmodules));
b2b13ee2 811
812 for (Int_t imodule = 0; imodule < nmodules; imodule++)
813 {
5e2cc5bc 814 treeS->GetEntry(imodule);
b2b13ee2 815 Int_t nentries = fSDigits->GetLast();
ecee2a1a 816 AliDebug(2,Form("Number of entries per module = %d",nentries+1));
b2b13ee2 817 for (Int_t ient = 0; ient < nentries+1; ient++)
818 {
819 pmdsdigit = (AliPMDsdigit*)fSDigits->UncheckedAt(ient);
820 trno = pmdsdigit->GetTrackNumber();
920e13db 821 trpid = pmdsdigit->GetTrackPid();
b2b13ee2 822 det = pmdsdigit->GetDetector();
823 smn = pmdsdigit->GetSMNumber();
7088dbe9 824 irow = pmdsdigit->GetRow();
825 icol = pmdsdigit->GetColumn();
b2b13ee2 826 edep = pmdsdigit->GetCellEdep();
827
828 MeV2ADC(edep,adc);
5e2cc5bc 829
830 // To decalibrte the adc values
831 //
832 Float_t gain1 = Gain(det,smn,irow,icol);
833 if (gain1 != 0.)
834 {
04a35d8d 835 Int_t adcDecalib = (Int_t)(adc/gain1);
836 adc = (Float_t) adcDecalib;
5e2cc5bc 837 }
838 else if(gain1 == 0.)
839 {
840 adc = 0.;
841 }
35535af7 842 // Pedestal Decalibration
843 Int_t pedmeanrms = fCalibPed->GetPedMeanRms(det,smn,irow,icol);
37edc588 844 Int_t pedrms1 = (Int_t) pedmeanrms%100;
35535af7 845 Float_t pedrms = (Float_t)pedrms1/10.;
846 Float_t pedmean = (Float_t) (pedmeanrms - pedrms1)/1000.0;
35535af7 847 if(adc > 0.)
848 {
849 adc += (pedmean + 3.0*pedrms);
920e13db 850 AddDigit(trno,trpid,det,smn,irow,icol,adc);
35535af7 851 }
5e2cc5bc 852
b2b13ee2 853 }
ebd83c56 854 treeD->Fill();
b2b13ee2 855 ResetDigit();
856 }
36031625 857 fPMDLoader->WriteDigits("OVERWRITE");
ecee2a1a 858
b2b13ee2 859}
7088dbe9 860//____________________________________________________________________________
861void AliPMDDigitizer::Exec(Option_t *option)
862{
863 // Does the event merging and digitization
7088dbe9 864 const char *cdeb = strstr(option,"deb");
865 if(cdeb)
866 {
ecee2a1a 867 AliDebug(100," *** PMD Exec is called ***");
7088dbe9 868 }
ecee2a1a 869
7088dbe9 870 Int_t ninputs = fManager->GetNinputs();
ecee2a1a 871 AliDebug(1,Form("Number of files to be processed = %d",ninputs));
7088dbe9 872 ResetCellADC();
873
874 for (Int_t i = 0; i < ninputs; i++)
875 {
876 Int_t troffset = fManager->GetMask(i);
877 MergeSDigits(i, troffset);
878 }
b2b13ee2 879
7088dbe9 880 fRunLoader = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
881 fPMD = (AliPMD*)gAlice->GetDetector("PMD");
882 fPMDLoader = fRunLoader->GetLoader("PMDLoader");
883 if (fPMDLoader == 0x0)
884 {
ecee2a1a 885 AliError("Can not find PMD or PMDLoader");
7088dbe9 886 }
06f2c79d 887 fPMDLoader->LoadDigits("update");
ebd83c56 888 TTree* treeD = fPMDLoader->TreeD();
889 if (treeD == 0x0)
7088dbe9 890 {
891 fPMDLoader->MakeTree("D");
ebd83c56 892 treeD = fPMDLoader->TreeD();
7088dbe9 893 }
894 Int_t bufsize = 16000;
ebd83c56 895 if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
5e2cc5bc 896 treeD->Branch("PMDDigit", &fDigits, bufsize);
7088dbe9 897
898 Float_t adc;
899 Float_t deltaE = 0.;
920e13db 900 Int_t detno = 0;
901 Int_t trno = 1;
902 Int_t trpid = -99;
7088dbe9 903
904 for (Int_t idet = 0; idet < 2; idet++)
905 {
906 for (Int_t ism = 0; ism < fgkTotUM; ism++)
907 {
908 for (Int_t jrow = 0; jrow < fgkRow; jrow++)
909 {
910 for (Int_t kcol = 0; kcol < fgkCol; kcol++)
911 {
912 if (idet == 0)
913 {
914 deltaE = fPRE[ism][jrow][kcol];
915 trno = fPRETrackNo[ism][jrow][kcol];
920e13db 916 trpid = fPRETrackPid[ism][jrow][kcol];
917 detno = 0;
7088dbe9 918 }
919 else if (idet == 1)
920 {
921 deltaE = fCPV[ism][jrow][kcol];
922 trno = fCPVTrackNo[ism][jrow][kcol];
920e13db 923 trpid = fCPVTrackPid[ism][jrow][kcol];
924 detno = 1;
7088dbe9 925 }
926 if (deltaE > 0.)
927 {
928 MeV2ADC(deltaE,adc);
35535af7 929
5f566ada 930 //
35535af7 931 // Gain decalibration
5f566ada 932 //
933 Float_t gain1 = Gain(idet,ism,jrow,kcol);
35535af7 934
5f566ada 935 if (gain1 != 0.)
936 {
937 Int_t adcDecalib = (Int_t)(adc/gain1);
938 adc = (Float_t) adcDecalib;
939 }
940 else if(gain1 == 0.)
941 {
942 adc = 0.;
943 }
35535af7 944 // Pedestal Decalibration
945 Int_t pedmeanrms =
946 fCalibPed->GetPedMeanRms(idet,ism,jrow,kcol);
37edc588 947 Int_t pedrms1 = (Int_t) pedmeanrms%100;
35535af7 948 Float_t pedrms = (Float_t)pedrms1/10.;
949 Float_t pedmean =
950 (Float_t) (pedmeanrms - pedrms1)/1000.0;
35535af7 951 if (adc > 0.)
952 {
953 adc += (pedmean + 3.0*pedrms);
920e13db 954 AddDigit(trno,trpid,detno,ism,jrow,kcol,adc);
35535af7 955 }
5f566ada 956
7088dbe9 957 }
958 } // column loop
959 } // row loop
ebd83c56 960 treeD->Fill();
7088dbe9 961 ResetDigit();
962 } // supermodule loop
963 } // detector loop
5e2cc5bc 964 fPMDLoader->WriteDigits("OVERWRITE");
06f2c79d 965 fPMDLoader->UnloadDigits();
966 ResetCellADC();
7088dbe9 967}
968//____________________________________________________________________________
4d59e381 969void AliPMDDigitizer::TrackAssignment2CPVCell()
970{
971 // This block assigns the cell id when there are
972 // multiple tracks in a cell according to the
973 // energy deposition
974 // This method added by Ajay
975 Bool_t jsort = false;
976
977 Int_t i, j, k;
978
979 Int_t *status1;
980 Int_t *status2;
981 Int_t *trnarray;
982 Float_t *fracEdp;
983 Float_t *trEdp;
984
985 Int_t ****cpvTrack;
986 Float_t ****cpvEdep;
987
988 cpvTrack = new Int_t ***[fgkTotUM];
989 cpvEdep = new Float_t ***[fgkTotUM];
990 for (i=0; i<fgkTotUM; i++)
991 {
992 cpvTrack[i] = new Int_t **[fgkRow];
993 cpvEdep[i] = new Float_t **[fgkRow];
994 }
995
996 for (i = 0; i < fgkTotUM; i++)
997 {
998 for (j = 0; j < fgkRow; j++)
999 {
1000 cpvTrack[i][j] = new Int_t *[fgkCol];
1001 cpvEdep[i][j] = new Float_t *[fgkCol];
1002 }
1003 }
1004 for (i = 0; i < fgkTotUM; i++)
1005 {
1006 for (j = 0; j < fgkRow; j++)
1007 {
1008 for (k = 0; k < fgkCol; k++)
1009 {
1010 Int_t nn = fCPVCounter[i][j][k];
1011 if(nn > 0)
1012 {
1013 cpvTrack[i][j][k] = new Int_t[nn];
1014 cpvEdep[i][j][k] = new Float_t[nn];
1015 }
1016 else
1017 {
1018 nn = 1;
1019 cpvTrack[i][j][k] = new Int_t[nn];
1020 cpvEdep[i][j][k] = new Float_t[nn];
1021 }
1022 fCPVCounter[i][j][k] = 0;
1023 }
1024 }
1025 }
1026
1027
1028 Int_t nentries = fCPVCell.GetEntries();
1029
1030 Int_t mtrackno, ism, ixp, iyp;
1031 Float_t edep;
1032 for (i = 0; i < nentries; i++)
1033 {
1034 AliPMDcell* cpvcell = (AliPMDcell*)fCPVCell.UncheckedAt(i);
1035
1036 mtrackno = cpvcell->GetTrackNumber();
1037 ism = cpvcell->GetSMNumber();
1038 ixp = cpvcell->GetX();
1039 iyp = cpvcell->GetY();
1040 edep = cpvcell->GetEdep();
1041 Int_t nn = fCPVCounter[ism][ixp][iyp];
1042 cpvTrack[ism][ixp][iyp][nn] = (Int_t) mtrackno;
1043 cpvEdep[ism][ixp][iyp][nn] = edep;
1044 fCPVCounter[ism][ixp][iyp]++;
1045 }
1046
1047 Int_t iz, il;
1048 Int_t im, ix, iy;
1049 Int_t nn;
1050 for (im=0; im<fgkTotUM; im++)
1051 {
1052 for (ix=0; ix<fgkRow; ix++)
1053 {
1054 for (iy=0; iy<fgkCol; iy++)
1055 {
1056 nn = fCPVCounter[im][ix][iy];
1057 if (nn > 1)
1058 {
1059 // This block handles if a cell is fired
1060 // many times by many tracks
1061 status1 = new Int_t[nn];
1062 status2 = new Int_t[nn];
1063 trnarray = new Int_t[nn];
1064 for (iz = 0; iz < nn; iz++)
1065 {
1066 status1[iz] = cpvTrack[im][ix][iy][iz];
1067 }
1068 TMath::Sort(nn,status1,status2,jsort);
1069 Int_t trackOld = -99999;
1070 Int_t track, trCount = 0;
1071 for (iz = 0; iz < nn; iz++)
1072 {
1073 track = status1[status2[iz]];
1074 if (trackOld != track)
1075 {
1076 trnarray[trCount] = track;
1077 trCount++;
1078 }
1079 trackOld = track;
1080 }
1081 delete [] status1;
1082 delete [] status2;
1083 Float_t totEdp = 0.;
1084 trEdp = new Float_t[trCount];
1085 fracEdp = new Float_t[trCount];
1086 for (il = 0; il < trCount; il++)
1087 {
1088 trEdp[il] = 0.;
1089 track = trnarray[il];
1090 for (iz = 0; iz < nn; iz++)
1091 {
1092 if (track == cpvTrack[im][ix][iy][iz])
1093 {
1094 trEdp[il] += cpvEdep[im][ix][iy][iz];
1095 }
1096 }
1097 totEdp += trEdp[il];
1098 }
1099 Int_t ilOld = 0;
1100 Float_t fracOld = 0.;
1101
1102 for (il = 0; il < trCount; il++)
1103 {
1104 fracEdp[il] = trEdp[il]/totEdp;
1105 if (fracOld < fracEdp[il])
1106 {
1107 fracOld = fracEdp[il];
1108 ilOld = il;
1109 }
1110 }
1111 fCPVTrackNo[im][ix][iy] = trnarray[ilOld];
1112 delete [] fracEdp;
1113 delete [] trEdp;
1114 delete [] trnarray;
1115 }
1116 else if (nn == 1)
1117 {
1118 // This only handles if a cell is fired
1119 // by only one track
1120
1121 fCPVTrackNo[im][ix][iy] = cpvTrack[im][ix][iy][0];
1122
1123 }
1124 else if (nn ==0)
1125 {
1126 // This is if no cell is fired
1127 fCPVTrackNo[im][ix][iy] = -999;
1128 }
1129 } // end of iy
1130 } // end of ix
1131 } // end of im
1132
1133 // Delete all the pointers
1134
1135 for (i = 0; i < fgkTotUM; i++)
1136 {
1137 for (j = 0; j < fgkRow; j++)
1138 {
1139 for (k = 0; k < fgkCol; k++)
1140 {
1141 delete []cpvTrack[i][j][k];
1142 delete []cpvEdep[i][j][k];
1143 }
1144 }
1145 }
1146
1147 for (i = 0; i < fgkTotUM; i++)
1148 {
1149 for (j = 0; j < fgkRow; j++)
1150 {
1151 delete [] cpvTrack[i][j];
1152 delete [] cpvEdep[i][j];
1153 }
1154 }
1155
1156 for (i = 0; i < fgkTotUM; i++)
1157 {
1158 delete [] cpvTrack[i];
1159 delete [] cpvEdep[i];
1160 }
1161 delete [] cpvTrack;
1162 delete [] cpvEdep;
1163
1164 //
1165 // End of the cell id assignment
1166 //
1167}
1168//____________________________________________________________________________
7088dbe9 1169
1170void AliPMDDigitizer::MergeSDigits(Int_t filenumber, Int_t troffset)
1171{
1172 // merging sdigits
1173 fRunLoader = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(filenumber));
1174 fPMDLoader = fRunLoader->GetLoader("PMDLoader");
1175 fPMDLoader->LoadSDigits("read");
ebd83c56 1176 TTree* treeS = fPMDLoader->TreeS();
7088dbe9 1177 AliPMDsdigit *pmdsdigit;
ebd83c56 1178 TBranch *branch = treeS->GetBranch("PMDSDigit");
1179 if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
7088dbe9 1180 branch->SetAddress(&fSDigits);
1181
920e13db 1182 Int_t itrackno, itrackpid, idet, ism;
7088dbe9 1183 Int_t ixp, iyp;
1184 Float_t edep;
ebd83c56 1185 Int_t nmodules = (Int_t) treeS->GetEntries();
ecee2a1a 1186 AliDebug(1,Form("Number of Modules in the treeS = %d",nmodules));
1187 AliDebug(1,Form("Track Offset = %d",troffset));
7088dbe9 1188 for (Int_t imodule = 0; imodule < nmodules; imodule++)
1189 {
5e2cc5bc 1190 treeS->GetEntry(imodule);
7088dbe9 1191 Int_t nentries = fSDigits->GetLast();
ecee2a1a 1192 AliDebug(2,Form("Number of Entries per Module = %d",nentries));
7088dbe9 1193 for (Int_t ient = 0; ient < nentries+1; ient++)
1194 {
1195 pmdsdigit = (AliPMDsdigit*)fSDigits->UncheckedAt(ient);
1196 itrackno = pmdsdigit->GetTrackNumber();
920e13db 1197 itrackpid = pmdsdigit->GetTrackPid();
7088dbe9 1198 idet = pmdsdigit->GetDetector();
1199 ism = pmdsdigit->GetSMNumber();
1200 ixp = pmdsdigit->GetRow();
1201 iyp = pmdsdigit->GetColumn();
1202 edep = pmdsdigit->GetCellEdep();
7088dbe9 1203 if (idet == 0)
1204 {
1205 if (fPRE[ism][ixp][iyp] < edep)
1206 {
1207 fPRETrackNo[ism][ixp][iyp] = troffset + itrackno;
920e13db 1208 fPRETrackPid[ism][ixp][iyp] = itrackpid;
7088dbe9 1209 }
1210 fPRE[ism][ixp][iyp] += edep;
1211 }
1212 else if (idet == 1)
1213 {
1214 if (fCPV[ism][ixp][iyp] < edep)
1215 {
1216 fCPVTrackNo[ism][ixp][iyp] = troffset + itrackno;
920e13db 1217 fCPVTrackPid[ism][ixp][iyp] = itrackpid;
7088dbe9 1218 }
1219 fCPV[ism][ixp][iyp] += edep;
1220 }
1221 }
1222 }
1223
1224}
1225// ----------------------------------------------------------------------
b2b13ee2 1226void AliPMDDigitizer::TrackAssignment2Cell()
1227{
b2b13ee2 1228 //
75bb524e 1229 // This block assigns the cell id when there are
b2b13ee2 1230 // multiple tracks in a cell according to the
1231 // energy deposition
1232 //
75bb524e 1233 Bool_t jsort = false;
8599803d 1234
b2b13ee2 1235 Int_t i, j, k;
1236
4d59e381 1237 Int_t *status1;
1238 Int_t *status2;
1239 Int_t *trnarray;
36031625 1240 Float_t *fracEdp;
1241 Float_t *trEdp;
4d59e381 1242
36031625 1243 Int_t ****pmdTrack;
1244 Float_t ****pmdEdep;
b2b13ee2 1245
36031625 1246 pmdTrack = new Int_t ***[fgkTotUM];
1247 pmdEdep = new Float_t ***[fgkTotUM];
1248 for (i=0; i<fgkTotUM; i++)
b2b13ee2 1249 {
36031625 1250 pmdTrack[i] = new Int_t **[fgkRow];
1251 pmdEdep[i] = new Float_t **[fgkRow];
b2b13ee2 1252 }
1253
36031625 1254 for (i = 0; i < fgkTotUM; i++)
b2b13ee2 1255 {
36031625 1256 for (j = 0; j < fgkRow; j++)
b2b13ee2 1257 {
36031625 1258 pmdTrack[i][j] = new Int_t *[fgkCol];
1259 pmdEdep[i][j] = new Float_t *[fgkCol];
b2b13ee2 1260 }
1261 }
1262
36031625 1263 for (i = 0; i < fgkTotUM; i++)
b2b13ee2 1264 {
36031625 1265 for (j = 0; j < fgkRow; j++)
b2b13ee2 1266 {
36031625 1267 for (k = 0; k < fgkCol; k++)
b2b13ee2 1268 {
36031625 1269 Int_t nn = fPRECounter[i][j][k];
b2b13ee2 1270 if(nn > 0)
1271 {
36031625 1272 pmdTrack[i][j][k] = new Int_t[nn];
1273 pmdEdep[i][j][k] = new Float_t[nn];
b2b13ee2 1274 }
1275 else
1276 {
1277 nn = 1;
36031625 1278 pmdTrack[i][j][k] = new Int_t[nn];
1279 pmdEdep[i][j][k] = new Float_t[nn];
5e2cc5bc 1280 }
36031625 1281 fPRECounter[i][j][k] = 0;
b2b13ee2 1282 }
1283 }
1284 }
1285
1286
ebd83c56 1287 Int_t nentries = fCell.GetEntries();
b2b13ee2 1288
1289 Int_t mtrackno, ism, ixp, iyp;
1290 Float_t edep;
1291
1292 for (i = 0; i < nentries; i++)
1293 {
ebd83c56 1294 AliPMDcell* cell = (AliPMDcell*)fCell.UncheckedAt(i);
b2b13ee2 1295
ebd83c56 1296 mtrackno = cell->GetTrackNumber();
1297 ism = cell->GetSMNumber();
1298 ixp = cell->GetX();
1299 iyp = cell->GetY();
1300 edep = cell->GetEdep();
36031625 1301 Int_t nn = fPRECounter[ism][ixp][iyp];
36031625 1302 pmdTrack[ism][ixp][iyp][nn] = (Int_t) mtrackno;
1303 pmdEdep[ism][ixp][iyp][nn] = edep;
1304 fPRECounter[ism][ixp][iyp]++;
b2b13ee2 1305 }
1306
1307 Int_t iz, il;
1308 Int_t im, ix, iy;
1309 Int_t nn;
1310
36031625 1311 for (im=0; im<fgkTotUM; im++)
b2b13ee2 1312 {
36031625 1313 for (ix=0; ix<fgkRow; ix++)
b2b13ee2 1314 {
36031625 1315 for (iy=0; iy<fgkCol; iy++)
b2b13ee2 1316 {
36031625 1317 nn = fPRECounter[im][ix][iy];
b2b13ee2 1318 if (nn > 1)
1319 {
1320 // This block handles if a cell is fired
1321 // many times by many tracks
75bb524e 1322 status1 = new Int_t[nn];
1323 status2 = new Int_t[nn];
1324 trnarray = new Int_t[nn];
b2b13ee2 1325 for (iz = 0; iz < nn; iz++)
1326 {
36031625 1327 status1[iz] = pmdTrack[im][ix][iy][iz];
b2b13ee2 1328 }
75bb524e 1329 TMath::Sort(nn,status1,status2,jsort);
36031625 1330 Int_t trackOld = -99999;
1331 Int_t track, trCount = 0;
b2b13ee2 1332 for (iz = 0; iz < nn; iz++)
1333 {
75bb524e 1334 track = status1[status2[iz]];
36031625 1335 if (trackOld != track)
b2b13ee2 1336 {
36031625 1337 trnarray[trCount] = track;
1338 trCount++;
5e2cc5bc 1339 }
36031625 1340 trackOld = track;
b2b13ee2 1341 }
6eb80a7c 1342 delete [] status1;
1343 delete [] status2;
36031625 1344 Float_t totEdp = 0.;
1345 trEdp = new Float_t[trCount];
1346 fracEdp = new Float_t[trCount];
1347 for (il = 0; il < trCount; il++)
b2b13ee2 1348 {
36031625 1349 trEdp[il] = 0.;
75bb524e 1350 track = trnarray[il];
b2b13ee2 1351 for (iz = 0; iz < nn; iz++)
1352 {
36031625 1353 if (track == pmdTrack[im][ix][iy][iz])
b2b13ee2 1354 {
36031625 1355 trEdp[il] += pmdEdep[im][ix][iy][iz];
b2b13ee2 1356 }
1357 }
36031625 1358 totEdp += trEdp[il];
b2b13ee2 1359 }
36031625 1360 Int_t ilOld = 0;
1361 Float_t fracOld = 0.;
b2b13ee2 1362
36031625 1363 for (il = 0; il < trCount; il++)
b2b13ee2 1364 {
36031625 1365 fracEdp[il] = trEdp[il]/totEdp;
1366 if (fracOld < fracEdp[il])
b2b13ee2 1367 {
36031625 1368 fracOld = fracEdp[il];
1369 ilOld = il;
b2b13ee2 1370 }
1371 }
36031625 1372 fPRETrackNo[im][ix][iy] = trnarray[ilOld];
6eb80a7c 1373 delete [] fracEdp;
1374 delete [] trEdp;
1375 delete [] trnarray;
b2b13ee2 1376 }
1377 else if (nn == 1)
1378 {
1379 // This only handles if a cell is fired
1380 // by only one track
1381
36031625 1382 fPRETrackNo[im][ix][iy] = pmdTrack[im][ix][iy][0];
b2b13ee2 1383
1384 }
1385 else if (nn ==0)
1386 {
1387 // This is if no cell is fired
36031625 1388 fPRETrackNo[im][ix][iy] = -999;
b2b13ee2 1389 }
1390 } // end of iy
1391 } // end of ix
1392 } // end of im
1393
1394 // Delete all the pointers
1395
36031625 1396 for (i = 0; i < fgkTotUM; i++)
b2b13ee2 1397 {
36031625 1398 for (j = 0; j < fgkRow; j++)
b2b13ee2 1399 {
36031625 1400 for (k = 0; k < fgkCol; k++)
b2b13ee2 1401 {
36031625 1402 delete [] pmdTrack[i][j][k];
1403 delete [] pmdEdep[i][j][k];
b2b13ee2 1404 }
1405 }
1406 }
1407
36031625 1408 for (i = 0; i < fgkTotUM; i++)
b2b13ee2 1409 {
36031625 1410 for (j = 0; j < fgkRow; j++)
b2b13ee2 1411 {
36031625 1412 delete [] pmdTrack[i][j];
1413 delete [] pmdEdep[i][j];
b2b13ee2 1414 }
1415 }
1416
36031625 1417 for (i = 0; i < fgkTotUM; i++)
b2b13ee2 1418 {
36031625 1419 delete [] pmdTrack[i];
1420 delete [] pmdEdep[i];
b2b13ee2 1421 }
6eb80a7c 1422 delete [] pmdTrack;
1423 delete [] pmdEdep;
b2b13ee2 1424 //
1425 // End of the cell id assignment
1426 //
1427}
7088dbe9 1428//____________________________________________________________________________
a7d110b8 1429void AliPMDDigitizer::MeV2ADC(Float_t mev, Float_t & adc) const
b2b13ee2 1430{
36031625 1431 // This converts the simulated edep to ADC according to the
1432 // Test Beam Data
ffdc86b5 1433 //PS Test in May 2009, Voltage @ 1350 V
4d59e381 1434 // KeV - ADC conversion for 12bit ADC
7ba814f3 1435 // MPV data used for the fit and taken here
1436
ffdc86b5 1437 const Float_t kConstant = -0.1602;
1438 const Float_t kErConstant = 0.9914;
1439 const Float_t kSlope = 77.47;
1440 const Float_t kErSlope = 3.16;
4d59e381 1441
b8e69f1f 1442 Float_t cons = gRandom->Gaus(kConstant,kErConstant);
1443 Float_t slop = gRandom->Gaus(kSlope,kErSlope);
4d59e381 1444
1445 Float_t adc12bit = slop*mev*0.001 + cons;
1446
7ba814f3 1447 if (adc12bit < 0.) adc12bit = 0.;
04a35d8d 1448
4d59e381 1449 if(adc12bit < 1600.0)
b8e69f1f 1450 {
1451 adc = (Float_t) adc12bit;
1452 }
4d59e381 1453 else if (adc12bit >= 1600.0)
b8e69f1f 1454 {
4d59e381 1455 adc = 1600.0;
b8e69f1f 1456 }
b2b13ee2 1457}
7088dbe9 1458//____________________________________________________________________________
920e13db 1459void AliPMDDigitizer::AddSDigit(Int_t trnumber, Int_t trpid, Int_t det,
1460 Int_t smnumber, Int_t irow, Int_t icol,
1461 Float_t adc)
b2b13ee2 1462{
36031625 1463 // Add SDigit
1464 //
ebd83c56 1465 if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
b2b13ee2 1466 TClonesArray &lsdigits = *fSDigits;
920e13db 1467 new(lsdigits[fNsdigit++]) AliPMDsdigit(trnumber,trpid,det,smnumber,irow,icol,adc);
b2b13ee2 1468}
7088dbe9 1469//____________________________________________________________________________
b2b13ee2 1470
920e13db 1471void AliPMDDigitizer::AddDigit(Int_t trnumber, Int_t trpid, Int_t det,
1472 Int_t smnumber, Int_t irow, Int_t icol,
1473 Float_t adc)
b2b13ee2 1474{
36031625 1475 // Add Digit
1476 //
ebd83c56 1477 if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
b2b13ee2 1478 TClonesArray &ldigits = *fDigits;
920e13db 1479 new(ldigits[fNdigit++]) AliPMDdigit(trnumber,trpid, det,smnumber,irow,icol,adc);
b2b13ee2 1480}
7088dbe9 1481//____________________________________________________________________________
b2b13ee2 1482
b2b13ee2 1483void AliPMDDigitizer::SetZPosition(Float_t zpos)
1484{
1485 fZPos = zpos;
1486}
7088dbe9 1487//____________________________________________________________________________
b2b13ee2 1488Float_t AliPMDDigitizer::GetZPosition() const
1489{
1490 return fZPos;
1491}
7088dbe9 1492//____________________________________________________________________________
b2b13ee2 1493
1494void AliPMDDigitizer::ResetCell()
1495{
36031625 1496 // clears the cell array and also the counter
1497 // for each cell
1498 //
4d59e381 1499 fCPVCell.Delete();
ebd83c56 1500 fCell.Delete();
36031625 1501 for (Int_t i = 0; i < fgkTotUM; i++)
b2b13ee2 1502 {
36031625 1503 for (Int_t j = 0; j < fgkRow; j++)
b2b13ee2 1504 {
36031625 1505 for (Int_t k = 0; k < fgkCol; k++)
b2b13ee2 1506 {
4d59e381 1507 fCPVCounter[i][j][k] = 0;
5e2cc5bc 1508 fPRECounter[i][j][k] = 0;
b2b13ee2 1509 }
1510 }
1511 }
1512}
7088dbe9 1513//____________________________________________________________________________
b2b13ee2 1514void AliPMDDigitizer::ResetSDigit()
1515{
36031625 1516 // Clears SDigits
b2b13ee2 1517 fNsdigit = 0;
ebd83c56 1518 if (fSDigits) fSDigits->Delete();
b2b13ee2 1519}
7088dbe9 1520//____________________________________________________________________________
b2b13ee2 1521void AliPMDDigitizer::ResetDigit()
1522{
36031625 1523 // Clears Digits
b2b13ee2 1524 fNdigit = 0;
ebd83c56 1525 if (fDigits) fDigits->Delete();
b2b13ee2 1526}
7088dbe9 1527//____________________________________________________________________________
b2b13ee2 1528
1529void AliPMDDigitizer::ResetCellADC()
1530{
7088dbe9 1531 // Clears individual cells edep and track number
36031625 1532 for (Int_t i = 0; i < fgkTotUM; i++)
b2b13ee2 1533 {
36031625 1534 for (Int_t j = 0; j < fgkRow; j++)
b2b13ee2 1535 {
36031625 1536 for (Int_t k = 0; k < fgkCol; k++)
b2b13ee2 1537 {
920e13db 1538 fCPV[i][j][k] = 0.;
1539 fPRE[i][j][k] = 0.;
1540 fCPVTrackNo[i][j][k] = 0;
1541 fPRETrackNo[i][j][k] = 0;
1542 fCPVTrackPid[i][j][k] = -1;
1543 fPRETrackPid[i][j][k] = -1;
b2b13ee2 1544 }
1545 }
1546 }
1547}
7088dbe9 1548//____________________________________________________________________________
b2b13ee2 1549
1550void AliPMDDigitizer::UnLoad(Option_t *option)
1551{
36031625 1552 // Unloads all the root files
1553 //
b2b13ee2 1554 const char *cS = strstr(option,"S");
1555 const char *cD = strstr(option,"D");
1556
1557 fRunLoader->UnloadgAlice();
1558 fRunLoader->UnloadHeader();
1559 fRunLoader->UnloadKinematics();
1560
1561 if (cS)
1562 {
36031625 1563 fPMDLoader->UnloadHits();
b2b13ee2 1564 }
1565 if (cD)
1566 {
36031625 1567 fPMDLoader->UnloadHits();
1568 fPMDLoader->UnloadSDigits();
b2b13ee2 1569 }
1570}
09a06455 1571
1572//----------------------------------------------------------------------
1573Float_t AliPMDDigitizer::Gain(Int_t det, Int_t smn, Int_t row, Int_t col) const
1574{
1575 // returns of the gain of the cell
1576 // Added this method by ZA
1577
1578 //cout<<" I am here in gain "<<fCalibData<< "smn,row, col "<<smn
1579 //<<" "<<row<<" "<<col<<endl;
1580
35535af7 1581 if(!fCalibGain) {
09a06455 1582 AliError("No calibration data loaded from CDB!!!");
50555ba1 1583 return 1;
09a06455 1584 }
1585
1586 Float_t GainFact;
35535af7 1587 GainFact = fCalibGain->GetGainFact(det,smn,row,col);
09a06455 1588 return GainFact;
1589}
1590//----------------------------------------------------------------------
35535af7 1591AliPMDCalibData* AliPMDDigitizer::GetCalibGain() const
09a06455 1592{
1593 // The run number will be centralized in AliCDBManager,
1594 // you don't need to set it here!
1595 // Added this method by ZA
0dd3d6f9 1596 // Cleaned up by Alberto
35535af7 1597 AliCDBEntry *entry = AliCDBManager::Instance()->Get("PMD/Calib/Gain");
09a06455 1598
0dd3d6f9 1599 if(!entry) AliFatal("Calibration object retrieval failed!");
09a06455 1600
1601 AliPMDCalibData *calibdata=0;
1602 if (entry) calibdata = (AliPMDCalibData*) entry->GetObject();
1603
0dd3d6f9 1604 if (!calibdata) AliFatal("No calibration data from calibration database !");
09a06455 1605
1606 return calibdata;
1607}
35535af7 1608//----------------------------------------------------------------------
1609AliPMDPedestal* AliPMDDigitizer::GetCalibPed() const
1610{
1611 // The run number will be centralized in AliCDBManager,
1612 // you don't need to set it here!
1613
1614 AliCDBEntry *entry = AliCDBManager::Instance()->Get("PMD/Calib/Ped");
1615
1616 if(!entry) AliFatal("Pedestal object retrieval failed!");
1617
1618 AliPMDPedestal *pedestal=0;
1619 if (entry) pedestal = (AliPMDPedestal*) entry->GetObject();
1620
1621 if (!pedestal) AliFatal("No pedestal data from calibration database !");
1622
1623 return pedestal;
1624}