Corrected raw->sdigits method. Using default altro raw-data reader class. Avoid alloc...
[u/mrichter/AliRoot.git] / TPC / AliTPC.cxx
CommitLineData
4c039060 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
88cb7938 16/* $Id$ */
4c039060 17
fe4da5cc 18///////////////////////////////////////////////////////////////////////////////
19// //
20// Time Projection Chamber //
21// This class contains the basic functions for the Time Projection Chamber //
22// detector. Functions specific to one particular geometry are //
23// contained in the derived classes //
24// //
25//Begin_Html
26/*
1439f98e 27<img src="picts/AliTPCClass.gif">
fe4da5cc 28*/
29//End_Html
30// //
8c555625 31// //
fe4da5cc 32///////////////////////////////////////////////////////////////////////////////
33
73042f01 34//
35
88cb7938 36#include <Riostream.h>
37#include <stdlib.h>
38
f7a1cc68 39#include <TF2.h>
88cb7938 40#include <TFile.h>
f7a1cc68 41#include <TGeoGlobalMagField.h>
88cb7938 42#include <TInterpreter.h>
fe4da5cc 43#include <TMath.h>
e8d02863 44#include <TMatrixF.h>
fe4da5cc 45#include <TObjectTable.h>
88cb7938 46#include <TParticle.h>
afc42102 47#include <TROOT.h>
88cb7938 48#include <TRandom.h>
f7a1cc68 49#include <TStopwatch.h>
50#include <TString.h>
afc42102 51#include <TSystem.h>
88cb7938 52#include <TTree.h>
f7a1cc68 53#include <TVector.h>
88cb7938 54#include <TVirtualMC.h>
f546a4b4 55#include <TParameter.h>
cc80f89e 56
cc80f89e 57#include "AliDigits.h"
88cb7938 58#include "AliMagF.h"
88cb7938 59#include "AliRun.h"
60#include "AliRunLoader.h"
cc80f89e 61#include "AliSimDigits.h"
88cb7938 62#include "AliTPC.h"
63#include "AliTPC.h"
88cb7938 64#include "AliTPCDigitsArray.h"
65#include "AliTPCLoader.h"
66#include "AliTPCPRF2D.h"
67#include "AliTPCParamSR.h"
68#include "AliTPCRF1D.h"
792bb11c 69#include "AliTPCTrackHitsV2.h"
88cb7938 70#include "AliTrackReference.h"
5d12ce38 71#include "AliMC.h"
5bf6eb16 72#include "AliStack.h"
85a5290f 73#include "AliTPCDigitizer.h"
0421c3d1 74#include "AliTPCBuffer.h"
75#include "AliTPCDDLRawData.h"
8c2b3fd7 76#include "AliLog.h"
79fa3dc8 77#include "AliTPCcalibDB.h"
78#include "AliTPCCalPad.h"
79#include "AliTPCCalROC.h"
80#include "AliTPCExB.h"
bc807150 81#include "AliRawReader.h"
dbfbf54e 82#include "AliTPCRawStreamV3.h"
23effa8e 83#include "TTreeStream.h"
39c8eb58 84
fe4da5cc 85ClassImp(AliTPC)
fe4da5cc 86//_____________________________________________________________________________
179c6296 87 AliTPC::AliTPC():AliDetector(),
88 fDefaults(0),
89 fSens(0),
90 fNsectors(0),
91 fDigitsArray(0),
92 fTPCParam(0),
93 fTrackHits(0),
94 fHitType(0),
95 fDigitsSwitch(0),
96 fSide(0),
b97617a6 97 fPrimaryIonisation(0),
179c6296 98 fNoiseDepth(0),
99 fNoiseTable(0),
100 fCurrentNoise(0),
8261c673 101 fActiveSectors(0),
23effa8e 102 fGainFactor(1.),
f546a4b4 103 fDebugStreamer(0),
104 fLHCclockPhaseSw(0)
179c6296 105
fe4da5cc 106{
107 //
108 // Default constructor
109 //
110 fIshunt = 0;
179c6296 111
be5ffbfe 112 // fTrackHitsOld = 0;
8c2b3fd7 113#if ROOT_VERSION_CODE >= ROOT_VERSION(4,0,1)
114 fHitType = 4; // ROOT containers
115#else
116 fHitType = 2; //default CONTAINERS - based on ROOT structure
117#endif
fe4da5cc 118}
119
120//_____________________________________________________________________________
121AliTPC::AliTPC(const char *name, const char *title)
179c6296 122 : AliDetector(name,title),
123 fDefaults(0),
124 fSens(0),
125 fNsectors(0),
126 fDigitsArray(0),
127 fTPCParam(0),
128 fTrackHits(0),
129 fHitType(0),
130 fDigitsSwitch(0),
131 fSide(0),
b97617a6 132 fPrimaryIonisation(0),
179c6296 133 fNoiseDepth(0),
134 fNoiseTable(0),
135 fCurrentNoise(0),
8261c673 136 fActiveSectors(0),
f546a4b4 137 fGainFactor(1.),
138 fDebugStreamer(0),
139 fLHCclockPhaseSw(0)
b97617a6 140
fe4da5cc 141{
142 //
143 // Standard constructor
144 //
145
146 //
147 // Initialise arrays of hits and digits
148 fHits = new TClonesArray("AliTPChit", 176);
5d12ce38 149 gAlice->GetMCApp()->AddHitList(fHits);
fe4da5cc 150 //
792bb11c 151 fTrackHits = new AliTPCTrackHitsV2;
39c8eb58 152 fTrackHits->SetHitPrecision(0.002);
153 fTrackHits->SetStepPrecision(0.003);
792bb11c 154 fTrackHits->SetMaxDistance(100);
155
be5ffbfe 156 //fTrackHitsOld = new AliTPCTrackHits; //MI - 13.09.2000
157 //fTrackHitsOld->SetHitPrecision(0.002);
158 //fTrackHitsOld->SetStepPrecision(0.003);
159 //fTrackHitsOld->SetMaxDistance(100);
792bb11c 160
39c8eb58 161
8c2b3fd7 162#if ROOT_VERSION_CODE >= ROOT_VERSION(4,0,1)
163 fHitType = 4; // ROOT containers
164#else
da33556f 165 fHitType = 2;
8c2b3fd7 166#endif
179c6296 167
168
cc80f89e 169
fe4da5cc 170 //
171 fIshunt = 0;
172 //
173 // Initialise color attributes
e939a978 174 //PH SetMarkerColor(kYellow);
73042f01 175
176 //
177 // Set TPC parameters
178 //
179
afc42102 180
181 if (!strcmp(title,"Default")) {
db2fdcfb 182 //fTPCParam = new AliTPCParamSR;
183 fTPCParam = AliTPCcalibDB::Instance()->GetParameters();
73042f01 184 } else {
8c2b3fd7 185 AliWarning("In Config.C you must set non-default parameters.");
73042f01 186 fTPCParam=0;
187 }
23effa8e 188
fe4da5cc 189}
0f2ff24c 190void AliTPC::CreateDebugStremer(){
191 //
192 // Create Debug streamer to check simulation
193 //
194 fDebugStreamer = new TTreeSRedirector("TPCSimdebug.root");
195}
fe4da5cc 196//_____________________________________________________________________________
197AliTPC::~AliTPC()
198{
199 //
200 // TPC destructor
201 //
73042f01 202
fe4da5cc 203 fIshunt = 0;
204 delete fHits;
205 delete fDigits;
db2fdcfb 206 //delete fTPCParam;
39c8eb58 207 delete fTrackHits; //MI 15.09.2000
be5ffbfe 208 // delete fTrackHitsOld; //MI 10.12.2001
9ff4af81 209
210 fDigitsArray = 0x0;
211 delete [] fNoiseTable;
212 delete [] fActiveSectors;
23effa8e 213 if (fDebugStreamer) delete fDebugStreamer;
fe4da5cc 214}
215
fe4da5cc 216//_____________________________________________________________________________
217void AliTPC::AddHit(Int_t track, Int_t *vol, Float_t *hits)
218{
219 //
220 // Add a hit to the list
221 //
39c8eb58 222 if (fHitType&1){
223 TClonesArray &lhits = *fHits;
224 new(lhits[fNhits++]) AliTPChit(fIshunt,track,vol,hits);
225 }
792bb11c 226 if (fHitType>1)
8c2b3fd7 227 AddHit2(track,vol,hits);
fe4da5cc 228}
88cb7938 229
fe4da5cc 230//_____________________________________________________________________________
231void AliTPC::CreateMaterials()
232{
8c555625 233 //-----------------------------------------------
37831078 234 // Create Materials for for TPC simulations
8c555625 235 //-----------------------------------------------
236
237 //-----------------------------------------------------------------
238 // Origin: Marek Kowalski IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
239 //-----------------------------------------------------------------
1283eee5 240
f7a1cc68 241 Int_t iSXFLD=((AliMagF*)TGeoGlobalMagField::Instance()->GetField())->Integ();
242 Float_t sXMGMX=((AliMagF*)TGeoGlobalMagField::Instance()->GetField())->Max();
1283eee5 243
244 Float_t amat[5]; // atomic numbers
245 Float_t zmat[5]; // z
246 Float_t wmat[5]; // proportions
247
248 Float_t density;
37831078 249
1283eee5 250
1283eee5 251
c1637e41 252 //***************** Gases *************************
1283eee5 253
37831078 254
1283eee5 255 //--------------------------------------------------------------
c1637e41 256 // gases - air and CO2
1283eee5 257 //--------------------------------------------------------------
258
37831078 259 // CO2
1283eee5 260
261 amat[0]=12.011;
262 amat[1]=15.9994;
263
264 zmat[0]=6.;
265 zmat[1]=8.;
266
c1637e41 267 wmat[0]=0.2729;
268 wmat[1]=0.7271;
1283eee5 269
2556fc5b 270 density=0.001754609;
1283eee5 271
1283eee5 272
c1637e41 273 AliMixture(10,"CO2",amat,zmat,density,2,wmat);
274 //
275 // Air
276 //
277 amat[0]=15.9994;
278 amat[1]=14.007;
279 //
280 zmat[0]=8.;
281 zmat[1]=7.;
282 //
283 wmat[0]=0.233;
284 wmat[1]=0.767;
285 //
286 density=0.001205;
1283eee5 287
c1637e41 288 AliMixture(11,"Air",amat,zmat,density,2,wmat);
289
1283eee5 290 //----------------------------------------------------------------
c1637e41 291 // drift gases
1283eee5 292 //----------------------------------------------------------------
37831078 293
4f1c6f0b 294 //
37831078 295 // Drift gases 1 - nonsensitive, 2 - sensitive
4f1c6f0b 296 // Ne-CO2-N2 (90-10-5) (volume) values at 20deg and 1 atm.
297 // rho(Ne) = 0.839 g/cm^3, rho(CO2) = 1.842 g/cm^3, rho(N2) = 1.165 g/cm^3
298 // for the calculation - everything is normalized to 1
e4e763b9 299
300 amat[0]= 20.18;
301 amat[1]=12.011;
302 amat[2]=15.9994;
c1637e41 303 amat[3]=14.007;
e4e763b9 304
305 zmat[0]= 10.;
306 zmat[1]=6.;
307 zmat[2]=8.;
c1637e41 308 zmat[3]=7.;
e4e763b9 309
4f1c6f0b 310 wmat[0]=0.756992632;
311 wmat[1]=0.056235789;
312 wmat[2]=0.128469474;
313 wmat[3]=0.058395789;
c1637e41 314
2556fc5b 315 density=0.000904929;
e4e763b9 316
c1637e41 317 AliMixture(12,"Ne-CO2-N-1",amat,zmat,density,4,wmat);
318 AliMixture(13,"Ne-CO2-N-2",amat,zmat,density,4,wmat);
0f5dc26f 319 AliMixture(30,"Ne-CO2-N-3",amat,zmat,density,4,wmat);
1283eee5 320 //----------------------------------------------------------------------
321 // solid materials
322 //----------------------------------------------------------------------
323
1283eee5 324
37831078 325 // Kevlar C14H22O2N2
1283eee5 326
37831078 327 amat[0] = 12.011;
328 amat[1] = 1.;
329 amat[2] = 15.999;
330 amat[3] = 14.006;
1283eee5 331
37831078 332 zmat[0] = 6.;
333 zmat[1] = 1.;
334 zmat[2] = 8.;
335 zmat[3] = 7.;
336
337 wmat[0] = 14.;
338 wmat[1] = 22.;
339 wmat[2] = 2.;
340 wmat[3] = 2.;
341
342 density = 1.45;
343
c1637e41 344 AliMixture(14,"Kevlar",amat,zmat,density,-4,wmat);
37831078 345
346 // NOMEX
1283eee5 347
37831078 348 amat[0] = 12.011;
349 amat[1] = 1.;
350 amat[2] = 15.999;
351 amat[3] = 14.006;
352
353 zmat[0] = 6.;
354 zmat[1] = 1.;
355 zmat[2] = 8.;
356 zmat[3] = 7.;
357
358 wmat[0] = 14.;
359 wmat[1] = 22.;
360 wmat[2] = 2.;
361 wmat[3] = 2.;
362
c1637e41 363 density = 0.029;
364
365 AliMixture(15,"NOMEX",amat,zmat,density,-4,wmat);
37831078 366
367 // Makrolon C16H18O3
368
369 amat[0] = 12.011;
370 amat[1] = 1.;
371 amat[2] = 15.999;
1283eee5 372
37831078 373 zmat[0] = 6.;
374 zmat[1] = 1.;
375 zmat[2] = 8.;
376
377 wmat[0] = 16.;
378 wmat[1] = 18.;
379 wmat[2] = 3.;
380
381 density = 1.2;
382
c1637e41 383 AliMixture(16,"Makrolon",amat,zmat,density,-3,wmat);
384
385 // Tedlar C2H3F
386
387 amat[0] = 12.011;
388 amat[1] = 1.;
389 amat[2] = 18.998;
390
391 zmat[0] = 6.;
392 zmat[1] = 1.;
393 zmat[2] = 9.;
394
395 wmat[0] = 2.;
396 wmat[1] = 3.;
397 wmat[2] = 1.;
398
399 density = 1.71;
400
401 AliMixture(17, "Tedlar",amat,zmat,density,-3,wmat);
37831078 402
1283eee5 403 // Mylar C5H4O2
404
405 amat[0]=12.011;
406 amat[1]=1.;
407 amat[2]=15.9994;
408
409 zmat[0]=6.;
410 zmat[1]=1.;
411 zmat[2]=8.;
412
413 wmat[0]=5.;
414 wmat[1]=4.;
415 wmat[2]=2.;
416
417 density = 1.39;
fe4da5cc 418
c1637e41 419 AliMixture(18, "Mylar",amat,zmat,density,-3,wmat);
420 // material for "prepregs"
421 // Epoxy - C14 H20 O3
422 // Quartz SiO2
423 // Carbon C
424 // prepreg1 60% C-fiber, 40% epoxy (vol)
425 amat[0]=12.011;
426 amat[1]=1.;
427 amat[2]=15.994;
1283eee5 428
c1637e41 429 zmat[0]=6.;
430 zmat[1]=1.;
431 zmat[2]=8.;
1283eee5 432
c1637e41 433 wmat[0]=0.923;
434 wmat[1]=0.023;
435 wmat[2]=0.054;
1283eee5 436
c1637e41 437 density=1.859;
1283eee5 438
c1637e41 439 AliMixture(19, "Prepreg1",amat,zmat,density,3,wmat);
1283eee5 440
c1637e41 441 //prepreg2 60% glass-fiber, 40% epoxy
1283eee5 442
c1637e41 443 amat[0]=12.01;
444 amat[1]=1.;
445 amat[2]=15.994;
446 amat[3]=28.086;
447
448 zmat[0]=6.;
449 zmat[1]=1.;
450 zmat[2]=8.;
451 zmat[3]=14.;
452
453 wmat[0]=0.194;
454 wmat[1]=0.023;
455 wmat[2]=0.443;
efbf94c5 456 wmat[3]=0.34;
c1637e41 457
458 density=1.82;
459
460 AliMixture(20, "Prepreg2",amat,zmat,density,4,wmat);
461
462 //prepreg3 50% glass-fiber, 50% epoxy
463
464 amat[0]=12.01;
465 amat[1]=1.;
466 amat[2]=15.994;
467 amat[3]=28.086;
468
469 zmat[0]=6.;
470 zmat[1]=1.;
471 zmat[2]=8.;
472 zmat[3]=14.;
1283eee5 473
efbf94c5 474 wmat[0]=0.257;
c1637e41 475 wmat[1]=0.03;
efbf94c5 476 wmat[2]=0.412;
c1637e41 477 wmat[3]=0.3;
478
0f788f87 479 density=1.725;
c1637e41 480
481 AliMixture(21, "Prepreg3",amat,zmat,density,4,wmat);
482
483 // G10 60% SiO2 40% epoxy
484
485 amat[0]=12.01;
486 amat[1]=1.;
487 amat[2]=15.994;
488 amat[3]=28.086;
489
490 zmat[0]=6.;
491 zmat[1]=1.;
492 zmat[2]=8.;
493 zmat[3]=14.;
494
495 wmat[0]=0.194;
496 wmat[1]=0.023;
497 wmat[2]=0.443;
498 wmat[3]=0.340;
499
500 density=1.7;
501
502 AliMixture(22, "G10",amat,zmat,density,4,wmat);
503
37831078 504 // Al
1283eee5 505
37831078 506 amat[0] = 26.98;
507 zmat[0] = 13.;
1283eee5 508
37831078 509 density = 2.7;
1283eee5 510
c1637e41 511 AliMaterial(23,"Al",amat[0],zmat[0],density,999.,999.);
1283eee5 512
c1637e41 513 // Si (for electronics
1283eee5 514
37831078 515 amat[0] = 28.086;
9a3667bd 516 zmat[0] = 14.;
1283eee5 517
37831078 518 density = 2.33;
1283eee5 519
c1637e41 520 AliMaterial(24,"Si",amat[0],zmat[0],density,999.,999.);
1283eee5 521
37831078 522 // Cu
1283eee5 523
37831078 524 amat[0] = 63.546;
525 zmat[0] = 29.;
1283eee5 526
37831078 527 density = 8.96;
1283eee5 528
c1637e41 529 AliMaterial(25,"Cu",amat[0],zmat[0],density,999.,999.);
1283eee5 530
0f788f87 531 // brass
532
533 amat[0] = 63.546;
534 zmat[0] = 29.;
535 //
536 amat[1]= 65.409;
537 zmat[1]= 30.;
538 //
539 wmat[0]= 0.6;
540 wmat[1]= 0.4;
541
542 //
543 density = 8.23;
544
545
546 //
547 AliMixture(33,"Brass",amat,zmat,density,2,wmat);
548
c1637e41 549 // Epoxy - C14 H20 O3
550
37831078 551 amat[0]=12.011;
552 amat[1]=1.;
553 amat[2]=15.9994;
1283eee5 554
37831078 555 zmat[0]=6.;
556 zmat[1]=1.;
557 zmat[2]=8.;
1283eee5 558
c1637e41 559 wmat[0]=14.;
560 wmat[1]=20.;
561 wmat[2]=3.;
1283eee5 562
c1637e41 563 density=1.25;
1283eee5 564
c1637e41 565 AliMixture(26,"Epoxy",amat,zmat,density,-3,wmat);
efbf94c5 566 //
567 // epoxy film - 90% epoxy, 10% glass fiber
568 //
569 amat[0]=12.01;
570 amat[1]=1.;
571 amat[2]=15.994;
572 amat[3]=28.086;
573
574 zmat[0]=6.;
575 zmat[1]=1.;
576 zmat[2]=8.;
577 zmat[3]=14.;
578
579 wmat[0]=0.596;
580 wmat[1]=0.071;
581 wmat[2]=0.257;
582 wmat[3]=0.076;
583
584
585 density=1.345;
586
587 AliMixture(34, "Epoxy-film",amat,zmat,density,4,wmat);
1283eee5 588
c1637e41 589 // Plexiglas C5H8O2
1283eee5 590
5a28a08f 591 amat[0]=12.011;
592 amat[1]=1.;
593 amat[2]=15.9994;
594
595 zmat[0]=6.;
596 zmat[1]=1.;
597 zmat[2]=8.;
598
c1637e41 599 wmat[0]=5.;
600 wmat[1]=8.;
601 wmat[2]=2.;
5a28a08f 602
c1637e41 603 density=1.18;
5a28a08f 604
c1637e41 605 AliMixture(27,"Plexiglas",amat,zmat,density,-3,wmat);
37831078 606
ba1d05f9 607 // Carbon
608
609 amat[0]=12.011;
610 zmat[0]=6.;
611 density= 2.265;
612
c1637e41 613 AliMaterial(28,"C",amat[0],zmat[0],density,999.,999.);
ba1d05f9 614
c1637e41 615 // Fe (steel for the inner heat screen)
e4e763b9 616
c1637e41 617 amat[0]=55.845;
e4e763b9 618
c1637e41 619 zmat[0]=26.;
e4e763b9 620
c1637e41 621 density=7.87;
ba1d05f9 622
c1637e41 623 AliMaterial(29,"Fe",amat[0],zmat[0],density,999.,999.);
f4441b15 624 //
625 // Peek - (C6H4-O-OC6H4-O-C6H4-CO)n
626 amat[0]=12.011;
627 amat[1]=1.;
628 amat[2]=15.9994;
629
630 zmat[0]=6.;
631 zmat[1]=1.;
632 zmat[2]=8.;
633
634 wmat[0]=19.;
635 wmat[1]=12.;
636 wmat[2]=3.;
637 //
638 density=1.3;
639 //
640 AliMixture(30,"Peek",amat,zmat,density,-3,wmat);
641 //
642 // Ceramics - Al2O3
643 //
644 amat[0] = 26.98;
645 amat[1]= 15.9994;
646
647 zmat[0] = 13.;
648 zmat[1]=8.;
649
650 wmat[0]=2.;
651 wmat[1]=3.;
652
653 density = 3.97;
654
655 AliMixture(31,"Alumina",amat,zmat,density,-2,wmat);
656
657 //
658 // liquids
659 //
660
661 // water
662
663 amat[0]=1.;
664 amat[1]=15.9994;
665
666 zmat[0]=1.;
667 zmat[1]=8.;
668
669 wmat[0]=2.;
670 wmat[1]=1.;
671
672 density=1.;
673
674 AliMixture(32,"Water",amat,zmat,density,-2,wmat);
efbf94c5 675
ba1d05f9 676
37831078 677 //----------------------------------------------------------
678 // tracking media for gases
679 //----------------------------------------------------------
680
c1637e41 681 AliMedium(0, "Air", 11, 0, iSXFLD, sXMGMX, 10., 999., .1, .01, .1);
682 AliMedium(1, "Ne-CO2-N-1", 12, 0, iSXFLD, sXMGMX, 10., 999.,.1,.001, .001);
683 AliMedium(2, "Ne-CO2-N-2", 13, 1, iSXFLD, sXMGMX, 10., 999.,.1,.001, .001);
37831078 684 AliMedium(3,"CO2",10,0, iSXFLD, sXMGMX, 10., 999.,.1, .001, .001);
0f5dc26f 685 AliMedium(20, "Ne-CO2-N-3", 30, 1, iSXFLD, sXMGMX, 10., 999.,.1,.001, .001);
37831078 686 //-----------------------------------------------------------
687 // tracking media for solids
688 //-----------------------------------------------------------
689
c1637e41 690 AliMedium(4,"Al",23,0, iSXFLD, sXMGMX, 10., 999., .1, .0005, .001);
691 AliMedium(5,"Kevlar",14,0, iSXFLD, sXMGMX, 10., 999., .1, .0005, .001);
692 AliMedium(6,"Nomex",15,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
693 AliMedium(7,"Makrolon",16,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
694 AliMedium(8,"Mylar",18,0, iSXFLD, sXMGMX, 10., 999., .1, .0005, .001);
695 AliMedium(9,"Tedlar",17,0, iSXFLD, sXMGMX, 10., 999., .1, .0005, .001);
696 //
697 AliMedium(10,"Prepreg1",19,0, iSXFLD, sXMGMX, 10., 999., .1, .0005, .001);
698 AliMedium(11,"Prepreg2",20,0, iSXFLD, sXMGMX, 10., 999., .1, .0005, .001);
699 AliMedium(12,"Prepreg3",21,0, iSXFLD, sXMGMX, 10., 999., .1, .0005, .001);
700 AliMedium(13,"Epoxy",26,0, iSXFLD, sXMGMX, 10., 999., .1, .0005, .001);
701
702 AliMedium(14,"Cu",25,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
703 AliMedium(15,"Si",24,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
704 AliMedium(16,"G10",22,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
705 AliMedium(17,"Plexiglas",27,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
706 AliMedium(18,"Steel",29,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
f4441b15 707 AliMedium(19,"Peek",30,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
708 AliMedium(21,"Alumina",31,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
0f788f87 709 AliMedium(22,"Water",32,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
efbf94c5 710 AliMedium(23,"Brass",33,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
711 AliMedium(24,"Epoxyfm",34,0, iSXFLD, sXMGMX, 10., 999., .1, .0005, .001);
fe4da5cc 712}
713
407ff276 714void AliTPC::GenerNoise(Int_t tablesize)
715{
716 //
717 //Generate table with noise
718 //
719 if (fTPCParam==0) {
720 // error message
721 fNoiseDepth=0;
722 return;
723 }
724 if (fNoiseTable) delete[] fNoiseTable;
725 fNoiseTable = new Float_t[tablesize];
726 fNoiseDepth = tablesize;
727 fCurrentNoise =0; //!index of the noise in the noise table
728
729 Float_t norm = fTPCParam->GetNoise()*fTPCParam->GetNoiseNormFac();
730 for (Int_t i=0;i<tablesize;i++) fNoiseTable[i]= gRandom->Gaus(0,norm);
731}
732
733Float_t AliTPC::GetNoise()
734{
735 // get noise from table
736 // if ((fCurrentNoise%10)==0)
737 // fCurrentNoise= gRandom->Rndm()*fNoiseDepth;
738 if (fCurrentNoise>=fNoiseDepth) fCurrentNoise=0;
739 return fNoiseTable[fCurrentNoise++];
740 //gRandom->Gaus(0, fTPCParam->GetNoise()*fTPCParam->GetNoiseNormFac());
741}
742
743
9bdd974b 744Bool_t AliTPC::IsSectorActive(Int_t sec) const
792bb11c 745{
746 //
747 // check if the sector is active
748 if (!fActiveSectors) return kTRUE;
749 else return fActiveSectors[sec];
750}
751
752void AliTPC::SetActiveSectors(Int_t * sectors, Int_t n)
753{
754 // activate interesting sectors
88cb7938 755 SetTreeAddress();//just for security
7ed5ec98 756 if (!fActiveSectors) fActiveSectors = new Bool_t[fTPCParam->GetNSector()];
792bb11c 757 for (Int_t i=0;i<fTPCParam->GetNSector();i++) fActiveSectors[i]=kFALSE;
758 for (Int_t i=0;i<n;i++)
759 if ((sectors[i]>=0) && sectors[i]<fTPCParam->GetNSector()) fActiveSectors[sectors[i]]=kTRUE;
760
761}
762
12d8d654 763void AliTPC::SetActiveSectors(Int_t flag)
792bb11c 764{
765 //
766 // activate sectors which were hitted by tracks
767 //loop over tracks
88cb7938 768 SetTreeAddress();//just for security
792bb11c 769 if (fHitType==0) return; // if Clones hit - not short volume ID information
7ed5ec98 770 if (!fActiveSectors) fActiveSectors = new Bool_t[fTPCParam->GetNSector()];
12d8d654 771 if (flag) {
772 for (Int_t i=0;i<fTPCParam->GetNSector();i++) fActiveSectors[i]=kTRUE;
773 return;
774 }
792bb11c 775 for (Int_t i=0;i<fTPCParam->GetNSector();i++) fActiveSectors[i]=kFALSE;
776 TBranch * branch=0;
cb5b8b21 777 if (fLoader->TreeH() == 0x0)
88cb7938 778 {
8c2b3fd7 779 AliFatal("Can not find TreeH in folder");
88cb7938 780 return;
781 }
cb5b8b21 782 if (fHitType>1) branch = fLoader->TreeH()->GetBranch("TPC2");
783 else branch = fLoader->TreeH()->GetBranch("TPC");
784 Stat_t ntracks = fLoader->TreeH()->GetEntries();
792bb11c 785 // loop over all hits
7d85e147 786 AliDebug(1,Form("Got %d tracks", (Int_t) ntracks));
88cb7938 787
8c2b3fd7 788 for(Int_t track=0;track<ntracks;track++) {
792bb11c 789 ResetHits();
790 //
791 if (fTrackHits && fHitType&4) {
cb5b8b21 792 TBranch * br1 = fLoader->TreeH()->GetBranch("fVolumes");
793 TBranch * br2 = fLoader->TreeH()->GetBranch("fNVolumes");
792bb11c 794 br1->GetEvent(track);
795 br2->GetEvent(track);
796 Int_t *volumes = fTrackHits->GetVolumes();
7ed5ec98 797 for (Int_t j=0;j<fTrackHits->GetNVolumes(); j++) {
e8631c7d 798 if (volumes[j]>-1 && volumes[j]<fTPCParam->GetNSector()) {
7ed5ec98 799 fActiveSectors[volumes[j]]=kTRUE;
800 }
801 else {
802 AliError(Form("Volume %d -> sector number %d is outside (0..%d)",
803 j,
804 volumes[j],
805 fTPCParam->GetNSector()));
806 }
807 }
792bb11c 808 }
809
810 //
be5ffbfe 811// if (fTrackHitsOld && fHitType&2) {
cb5b8b21 812// TBranch * br = fLoader->TreeH()->GetBranch("fTrackHitsInfo");
be5ffbfe 813// br->GetEvent(track);
814// AliObjectArray * ar = fTrackHitsOld->fTrackHitsInfo;
815// for (UInt_t j=0;j<ar->GetSize();j++){
816// fActiveSectors[((AliTrackHitsInfo*)ar->At(j))->fVolumeID] =kTRUE;
817// }
818// }
792bb11c 819 }
792bb11c 820}
821
822
823
fe4da5cc 824
0421c3d1 825//_____________________________________________________________________________
826void AliTPC::Digits2Raw()
827{
828// convert digits of the current event to raw data
829
830 static const Int_t kThreshold = 0;
0421c3d1 831
832 fLoader->LoadDigits();
833 TTree* digits = fLoader->TreeD();
834 if (!digits) {
8c2b3fd7 835 AliError("No digits tree");
0421c3d1 836 return;
837 }
838
839 AliSimDigits digarr;
840 AliSimDigits* digrow = &digarr;
841 digits->GetBranch("Segment")->SetAddress(&digrow);
842
843 const char* fileName = "AliTPCDDL.dat";
844 AliTPCBuffer* buffer = new AliTPCBuffer(fileName);
845 //Verbose level
846 // 0: Silent
847 // 1: cout messages
848 // 2: txt files with digits
849 //BE CAREFUL, verbose level 2 MUST be used only for debugging and
850 //it is highly suggested to use this mode only for debugging digits files
851 //reasonably small, because otherwise the size of the txt files can reach
852 //quickly several MB wasting time and disk space.
853 buffer->SetVerbose(0);
854
855 Int_t nEntries = Int_t(digits->GetEntries());
856 Int_t previousSector = -1;
857 Int_t subSector = 0;
858 for (Int_t i = 0; i < nEntries; i++) {
859 digits->GetEntry(i);
860 Int_t sector, row;
861 fTPCParam->AdjustSectorRow(digarr.GetID(), sector, row);
862 if(previousSector != sector) {
863 subSector = 0;
864 previousSector = sector;
865 }
866
867 if (sector < 36) { //inner sector [0;35]
868 if (row != 30) {
869 //the whole row is written into the output file
870 buffer->WriteRowBinary(kThreshold, digrow, 0, 0, 0,
871 sector, subSector, row);
872 } else {
873 //only the pads in the range [37;48] are written into the output file
874 buffer->WriteRowBinary(kThreshold, digrow, 37, 48, 1,
875 sector, subSector, row);
876 subSector = 1;
877 //only the pads outside the range [37;48] are written into the output file
878 buffer->WriteRowBinary(kThreshold, digrow, 37, 48, 2,
879 sector, subSector, row);
880 }//end else
881
882 } else { //outer sector [36;71]
883 if (row == 54) subSector = 2;
884 if ((row != 27) && (row != 76)) {
885 buffer->WriteRowBinary(kThreshold, digrow, 0, 0, 0,
886 sector, subSector, row);
887 } else if (row == 27) {
8c2b3fd7 888 //only the pads outside the range [43;46] are written into the output file
889 buffer->WriteRowBinary(kThreshold, digrow, 43, 46, 2,
0421c3d1 890 sector, subSector, row);
8c2b3fd7 891 subSector = 1;
892 //only the pads in the range [43;46] are written into the output file
893 buffer->WriteRowBinary(kThreshold, digrow, 43, 46, 1,
0421c3d1 894 sector, subSector, row);
895 } else if (row == 76) {
8c2b3fd7 896 //only the pads outside the range [33;88] are written into the output file
897 buffer->WriteRowBinary(kThreshold, digrow, 33, 88, 2,
898 sector, subSector, row);
899 subSector = 3;
900 //only the pads in the range [33;88] are written into the output file
901 buffer->WriteRowBinary(kThreshold, digrow, 33, 88, 1,
0421c3d1 902 sector, subSector, row);
903 }
904 }//end else
905 }//end for
906
907 delete buffer;
908 fLoader->UnloadDigits();
909
910 AliTPCDDLRawData rawWriter;
911 rawWriter.SetVerbose(0);
912
913 rawWriter.RawData(fileName);
914 gSystem->Unlink(fileName);
915
0421c3d1 916}
917
918
bc807150 919//_____________________________________________________________________________
920Bool_t AliTPC::Raw2SDigits(AliRawReader* rawReader){
921 // Converts the TPC raw data into summable digits
922 // The method is used for merging simulated and
923 // real data events
924 if (fLoader->TreeS() == 0x0 ) {
925 fLoader->MakeTree("S");
926 }
927
928 if(fDefaults == 0) SetDefaults(); // check if the parameters are set
929
930 //setup TPCDigitsArray
931 if(GetDigitsArray()) delete GetDigitsArray();
932
933 AliTPCDigitsArray *arr = new AliTPCDigitsArray;
934 arr->SetClass("AliSimDigits");
935 arr->Setup(fTPCParam);
936 arr->MakeTree(fLoader->TreeS());
937
938 SetDigitsArray(arr);
939
940 // set zero suppression to "0"
941 fTPCParam->SetZeroSup(0);
942
943 // Loop over sectors
d899e254 944 const Int_t kmaxTime = fTPCParam->GetMaxTBin();
bc807150 945 const Int_t kNIS = fTPCParam->GetNInnerSector();
946 const Int_t kNOS = fTPCParam->GetNOuterSector();
947 const Int_t kNS = kNIS + kNOS;
948
dbfbf54e 949 // Setup storage
950 AliTPCROC * roc = AliTPCROC::Instance();
951 Int_t nRowsMax = roc->GetNRows(roc->GetNSector()-1);
952 Int_t nPadsMax = roc->GetNPads(roc->GetNSector()-1,nRowsMax-1);
953 Short_t** allBins = new Short_t*[nRowsMax];
954 for (Int_t iRow = 0; iRow < nRowsMax; iRow++) {
955 Int_t maxBin = kmaxTime*nPadsMax;
956 allBins[iRow] = new Short_t[maxBin];
957 memset(allBins[iRow],0,sizeof(Short_t)*maxBin);
958 }
959
bc807150 960 for(Int_t iSector = 0; iSector < kNS; iSector++) {
961
962 Int_t nRows = fTPCParam->GetNRow(iSector);
963 Int_t nDDLs = 0, indexDDL = 0;
964 if (iSector < kNIS) {
965 nDDLs = 2;
966 indexDDL = iSector * 2;
967 }
968 else {
969 nDDLs = 4;
970 indexDDL = (iSector-kNIS) * 4 + kNIS * 2;
971 }
972
dbfbf54e 973 // Load the raw data for corresponding DDLs
bc807150 974 rawReader->Reset();
dbfbf54e 975
976 AliTPCAltroMapping** mapping =AliTPCcalibDB::Instance()->GetMapping();
977 AliTPCRawStreamV3 input(rawReader,(AliAltroMapping**)mapping);
bc807150 978 rawReader->Select("TPC",indexDDL,indexDDL+nDDLs-1);
979
dbfbf54e 980 // Clean storage
981 for (Int_t iRow = 0; iRow < nRowsMax; iRow++) {
982 Int_t maxBin = kmaxTime*nPadsMax;
bc807150 983 memset(allBins[iRow],0,sizeof(Short_t)*maxBin);
984 }
985
986 // Begin loop over altro data
dbfbf54e 987 while (input.NextDDL()) {
bc807150 988
989 if (input.GetSector() != iSector)
990 AliFatal(Form("Sector index mismatch ! Expected (%d), but got (%d) !",iSector,input.GetSector()));
991
dbfbf54e 992 //loop over pads
993 while ( input.NextChannel() ) {
994
995 Int_t iRow = input.GetRow();
996 if (iRow < 0 || iRow >= nRows)
997 AliFatal(Form("Pad-row index (%d) outside the range (%d -> %d) !",
998 iRow, 0, nRows -1));
999 Int_t iPad = input.GetPad();
1000
1001 Int_t maxPad = fTPCParam->GetNPads(iSector,iRow);
1002
1003 if (iPad < 0 || iPad >= maxPad)
1004 AliFatal(Form("Pad index (%d) outside the range (%d -> %d) !",
1005 iPad, 0, maxPad -1));
1006
1007 //loop over bunches
1008 while ( input.NextBunch() ){
1009 Int_t startTbin = (Int_t)input.GetStartTimeBin();
1010 Int_t bunchlength = (Int_t)input.GetBunchLength();
1011 const UShort_t *sig = input.GetSignals();
1012 for (Int_t iTime = 0; iTime<bunchlength; iTime++){
1013 Int_t iTimeBin=startTbin-iTime;
1014 if ( iTimeBin < 0 || iTimeBin >= kmaxTime) {
1015 continue;
1016 //AliFatal(Form("Timebin index (%d) outside the range (%d -> %d) !",
1017 // iTimeBin, 0, kmaxTime -1));
1018 }
1019
1020 Int_t maxBin = kmaxTime*maxPad;
1021 if (((iPad*kmaxTime+iTimeBin) >= maxBin) ||
1022 ((iPad*kmaxTime+iTimeBin) < 0))
1023 AliFatal(Form("Index outside the allowed range"
1024 " Sector=%d Row=%d Pad=%d Timebin=%d"
1025 " (Max.index=%d)",iSector,iRow,iPad,iTimeBin,maxBin));
1026 allBins[iRow][iPad*kmaxTime+iTimeBin] = sig[iTime];
1027 }
1028 }
1029 } // End loop over altro data
1030 }
bc807150 1031
bc807150 1032 // Now fill the digits array
1033 if (fDigitsArray->GetTree()==0) {
1034 AliFatal("Tree not set in fDigitsArray");
1035 }
1036
1037 for (Int_t iRow = 0; iRow < nRows; iRow++) {
1038 AliDigits * dig = fDigitsArray->CreateRow(iSector,iRow);
bc807150 1039 Int_t maxPad = fTPCParam->GetNPads(iSector,iRow);
1040 for(Int_t iPad = 0; iPad < maxPad; iPad++) {
d899e254 1041 for(Int_t iTimeBin = 0; iTimeBin < kmaxTime; iTimeBin++) {
1042 Short_t q = allBins[iRow][iPad*kmaxTime + iTimeBin];
bc807150 1043 if (q <= 0) continue;
1044 q *= 16;
1045 dig->SetDigitFast((Short_t)q,iTimeBin,iPad);
1046 }
1047 }
1048 fDigitsArray->StoreRow(iSector,iRow);
1049 Int_t ndig = dig->GetDigitSize();
1050
1051 AliDebug(10,
1052 Form("*** Sector, row, compressed digits %d %d %d ***\n",
1053 iSector,iRow,ndig));
1054
1055 fDigitsArray->ClearRow(iSector,iRow);
1056
1057 } // end of the sector digitization
bc807150 1058 }
1059
1060 fLoader->WriteSDigits("OVERWRITE");
1061
1062 if(GetDigitsArray()) delete GetDigitsArray();
1063 SetDigitsArray(0x0);
1064
dbfbf54e 1065 // cleanup storage
1066 for (Int_t iRow = 0; iRow < nRowsMax; iRow++)
1067 delete [] allBins[iRow];
1068 delete [] allBins;
1069
bc807150 1070 return kTRUE;
1071}
0a61bf9d 1072
85a5290f 1073//______________________________________________________________________
9bdd974b 1074AliDigitizer* AliTPC::CreateDigitizer(AliRunDigitizer* manager) const
85a5290f 1075{
1076 return new AliTPCDigitizer(manager);
1077}
0a61bf9d 1078//__
176aff27 1079void AliTPC::SDigits2Digits2(Int_t /*eventnumber*/)
0a61bf9d 1080{
1081 //create digits from summable digits
407ff276 1082 GenerNoise(500000); //create teble with noise
0a61bf9d 1083
1084 //conect tree with sSDigits
88cb7938 1085 TTree *t = fLoader->TreeS();
1086
8c2b3fd7 1087 if (t == 0x0) {
1088 fLoader->LoadSDigits("READ");
1089 t = fLoader->TreeS();
1090 if (t == 0x0) {
1091 AliError("Can not get input TreeS");
1092 return;
1093 }
1094 }
88cb7938 1095
1096 if (fLoader->TreeD() == 0x0) fLoader->MakeTree("D");
1097
0a61bf9d 1098 AliSimDigits digarr, *dummy=&digarr;
88cb7938 1099 TBranch* sdb = t->GetBranch("Segment");
8c2b3fd7 1100 if (sdb == 0x0) {
1101 AliError("Can not find branch with segments in TreeS.");
1102 return;
1103 }
88cb7938 1104
1105 sdb->SetAddress(&dummy);
1106
0a61bf9d 1107 Stat_t nentries = t->GetEntries();
1108
5f16d237 1109 // set zero suppression
1110
1111 fTPCParam->SetZeroSup(2);
1112
1113 // get zero suppression
1114
1115 Int_t zerosup = fTPCParam->GetZeroSup();
1116
1117 //make tree with digits
1118
0a61bf9d 1119 AliTPCDigitsArray *arr = new AliTPCDigitsArray;
1120 arr->SetClass("AliSimDigits");
1121 arr->Setup(fTPCParam);
88cb7938 1122 arr->MakeTree(fLoader->TreeD());
0a61bf9d 1123
88cb7938 1124 AliTPCParam * par = fTPCParam;
5f16d237 1125
0a61bf9d 1126 //Loop over segments of the TPC
5f16d237 1127
0a61bf9d 1128 for (Int_t n=0; n<nentries; n++) {
1129 t->GetEvent(n);
1130 Int_t sec, row;
1131 if (!par->AdjustSectorRow(digarr.GetID(),sec,row)) {
8c2b3fd7 1132 AliWarning(Form("Invalid segment ID ! %d",digarr.GetID()));
0a61bf9d 1133 continue;
1134 }
8c2b3fd7 1135 if (!IsSectorActive(sec)) continue;
1136
0a61bf9d 1137 AliSimDigits * digrow =(AliSimDigits*) arr->CreateRow(sec,row);
1138 Int_t nrows = digrow->GetNRows();
1139 Int_t ncols = digrow->GetNCols();
1140
1141 digrow->ExpandBuffer();
1142 digarr.ExpandBuffer();
1143 digrow->ExpandTrackBuffer();
1144 digarr.ExpandTrackBuffer();
1145
5f16d237 1146
407ff276 1147 Short_t * pamp0 = digarr.GetDigits();
1148 Int_t * ptracks0 = digarr.GetTracks();
1149 Short_t * pamp1 = digrow->GetDigits();
1150 Int_t * ptracks1 = digrow->GetTracks();
1151 Int_t nelems =nrows*ncols;
e93a083a 1152 Int_t saturation = fTPCParam->GetADCSat() - 1;
407ff276 1153 //use internal structure of the AliDigits - for speed reason
1154 //if you cahnge implementation
1155 //of the Alidigits - it must be rewriten -
1156 for (Int_t i= 0; i<nelems; i++){
407ff276 1157 Float_t q = TMath::Nint(Float_t(*pamp0)/16.+GetNoise());
1158 if (q>zerosup){
1159 if (q>saturation) q=saturation;
1160 *pamp1=(Short_t)q;
8c2b3fd7 1161
407ff276 1162 ptracks1[0]=ptracks0[0];
1163 ptracks1[nelems]=ptracks0[nelems];
1164 ptracks1[2*nelems]=ptracks0[2*nelems];
1165 }
1166 pamp0++;
1167 pamp1++;
1168 ptracks0++;
1169 ptracks1++;
1170 }
5f16d237 1171
5f16d237 1172 arr->StoreRow(sec,row);
1173 arr->ClearRow(sec,row);
5f16d237 1174 }
0a61bf9d 1175
0a61bf9d 1176
5f16d237 1177 //write results
88cb7938 1178 fLoader->WriteDigits("OVERWRITE");
5f16d237 1179
88cb7938 1180 delete arr;
afc42102 1181}
1182//__________________________________________________________________
1183void AliTPC::SetDefaults(){
9bdd974b 1184 //
1185 // setting the defaults
1186 //
afc42102 1187
afc42102 1188 // Set response functions
1189
88cb7938 1190 //
024a7e64 1191 AliRunLoader* rl = (AliRunLoader*)fLoader->GetEventFolder()->FindObject(AliRunLoader::GetRunLoaderName());
88cb7938 1192 rl->CdGAFile();
2ab0c725 1193 AliTPCParamSR *param=(AliTPCParamSR*)gDirectory->Get("75x40_100x60");
db2fdcfb 1194 // if(param){
1195// AliInfo("You are using 2 pad-length geom hits with 3 pad-lenght geom digits...");
1196// delete param;
1197// param = new AliTPCParamSR();
1198// }
1199// else {
1200// param=(AliTPCParamSR*)gDirectory->Get("75x40_100x60_150x60");
1201// }
1202 param = (AliTPCParamSR*)AliTPCcalibDB::Instance()->GetParameters();
1203 if (!param->IsGeoRead()){
1204 //
1205 // read transformation matrices for gGeoManager
1206 //
1207 param->ReadGeoMatrices();
1208 }
7a09f434 1209 if(!param){
8c2b3fd7 1210 AliFatal("No TPC parameters found");
7a09f434 1211 }
1212
1213
2ab0c725 1214 AliTPCPRF2D * prfinner = new AliTPCPRF2D;
f03e3423 1215 AliTPCPRF2D * prfouter1 = new AliTPCPRF2D;
1216 AliTPCPRF2D * prfouter2 = new AliTPCPRF2D;
bfb1cfbd 1217
afc42102 1218
bfb1cfbd 1219 //AliTPCRF1D * rf = new AliTPCRF1D(kTRUE);
1220 //rf->SetGauss(param->GetZSigma(),param->GetZWidth(),1.);
1221 //rf->SetOffset(3*param->GetZSigma());
1222 //rf->Update();
1223 //
1224 // Use gamma 4
1225 //
1226 char strgamma4[1000];
04c5c48d 1227 sprintf(strgamma4,"AliTPCRF1D::Gamma4((x-0.135+%f)*%f,55,160)",3*param->GetZSigma(), 1000000000*param->GetTSample()/param->GetZWidth());
bfb1cfbd 1228
1229 TF1 * fgamma4 = new TF1("fgamma4",strgamma4, -1,1);
1230 AliTPCRF1D * rf = new AliTPCRF1D(kTRUE,1000);
1231 rf->SetParam(fgamma4,param->GetZWidth(), 1,0.2);
1232 rf->SetOffset(3*param->GetZSigma());
1233 rf->Update();
1234
afc42102 1235 TDirectory *savedir=gDirectory;
2ab0c725 1236 TFile *f=TFile::Open("$ALICE_ROOT/TPC/AliTPCprf2d.root");
8c2b3fd7 1237 if (!f->IsOpen())
1238 AliFatal("Can't open $ALICE_ROOT/TPC/AliTPCprf2d.root !");
2a336e15 1239
1240 TString s;
2ab0c725 1241 prfinner->Read("prf_07504_Gati_056068_d02");
2a336e15 1242 //PH Set different names
1243 s=prfinner->GetGRF()->GetName();
1244 s+="in";
1245 prfinner->GetGRF()->SetName(s.Data());
1246
f03e3423 1247 prfouter1->Read("prf_10006_Gati_047051_d03");
2a336e15 1248 s=prfouter1->GetGRF()->GetName();
1249 s+="out1";
1250 prfouter1->GetGRF()->SetName(s.Data());
1251
f03e3423 1252 prfouter2->Read("prf_15006_Gati_047051_d03");
2a336e15 1253 s=prfouter2->GetGRF()->GetName();
1254 s+="out2";
1255 prfouter2->GetGRF()->SetName(s.Data());
1256
2ab0c725 1257 f->Close();
afc42102 1258 savedir->cd();
1259
2ab0c725 1260 param->SetInnerPRF(prfinner);
f03e3423 1261 param->SetOuter1PRF(prfouter1);
1262 param->SetOuter2PRF(prfouter2);
2ab0c725 1263 param->SetTimeRF(rf);
1264
afc42102 1265 // set fTPCParam
1266
2ab0c725 1267 SetParam(param);
1268
afc42102 1269
1270 fDefaults = 1;
1271
1272}
1273//__________________________________________________________________
85a5290f 1274void AliTPC::Hits2Digits()
1275{
9bdd974b 1276 //
1277 // creates digits from hits
1278 //
506e60a6 1279 if (!fTPCParam->IsGeoRead()){
1280 //
1281 // read transformation matrices for gGeoManager
1282 //
1283 fTPCParam->ReadGeoMatrices();
1284 }
9bdd974b 1285
85a5290f 1286 fLoader->LoadHits("read");
1287 fLoader->LoadDigits("recreate");
1288 AliRunLoader* runLoader = fLoader->GetRunLoader();
1289
1290 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
7ed5ec98 1291 //PH runLoader->GetEvent(iEvent);
85a5290f 1292 Hits2Digits(iEvent);
1293 }
1294
1295 fLoader->UnloadHits();
1296 fLoader->UnloadDigits();
1297}
1298//__________________________________________________________________
afc42102 1299void AliTPC::Hits2Digits(Int_t eventnumber)
1300{
1301 //----------------------------------------------------
1302 // Loop over all sectors for a single event
1303 //----------------------------------------------------
024a7e64 1304 AliRunLoader* rl = (AliRunLoader*)fLoader->GetEventFolder()->FindObject(AliRunLoader::GetRunLoaderName());
88cb7938 1305 rl->GetEvent(eventnumber);
6db60b71 1306 SetActiveSectors();
8c2b3fd7 1307 if (fLoader->TreeH() == 0x0) {
1308 if(fLoader->LoadHits()) {
1309 AliError("Can not load hits.");
1310 }
1311 }
88cb7938 1312 SetTreeAddress();
1313
8c2b3fd7 1314 if (fLoader->TreeD() == 0x0 ) {
1315 fLoader->MakeTree("D");
1316 if (fLoader->TreeD() == 0x0 ) {
1317 AliError("Can not get TreeD");
1318 return;
1319 }
1320 }
afc42102 1321
1322 if(fDefaults == 0) SetDefaults(); // check if the parameters are set
407ff276 1323 GenerNoise(500000); //create teble with noise
2ab0c725 1324
1325 //setup TPCDigitsArray
afc42102 1326
1327 if(GetDigitsArray()) delete GetDigitsArray();
8c2b3fd7 1328
2ab0c725 1329 AliTPCDigitsArray *arr = new AliTPCDigitsArray;
1330 arr->SetClass("AliSimDigits");
afc42102 1331 arr->Setup(fTPCParam);
88cb7938 1332
1333 arr->MakeTree(fLoader->TreeD());
2ab0c725 1334 SetDigitsArray(arr);
1335
afc42102 1336 fDigitsSwitch=0; // standard digits
f546a4b4 1337 // here LHC clock phase
1338 Float_t lhcph = 0.;
1339 switch (fLHCclockPhaseSw){
1340 case 0:
1341 // no phase
1342 lhcph=0.;
1343 break;
1344 case 1:
1345 // random phase
1346 lhcph = (Int_t)(gRandom->Rndm()/0.25);
1347 break;
1348 case 2:
1349 lhcph=0.;
1350 // not implemented yet
1351 break;
1352 }
1353 // adding phase to the TreeD user info
1354 fLoader->TreeD()->GetUserInfo()->Add(new TParameter<float>("lhcphase0",lhcph));
1355 //
8c2b3fd7 1356 for(Int_t isec=0;isec<fTPCParam->GetNSector();isec++)
1357 if (IsSectorActive(isec)) {
7d85e147 1358 AliDebug(1,Form("Hits2Digits: Sector %d is active.",isec));
8c2b3fd7 1359 Hits2DigitsSector(isec);
1360 }
1361 else {
7d85e147 1362 AliDebug(1,Form("Hits2Digits: Sector %d is NOT active.",isec));
8c2b3fd7 1363 }
1364
88cb7938 1365 fLoader->WriteDigits("OVERWRITE");
f2a509af 1366
1367//this line prevents the crash in the similar one
1368//on the beginning of this method
1369//destructor attempts to reset the tree, which is deleted by the loader
1370//need to be redesign
8c2b3fd7 1371 if(GetDigitsArray()) delete GetDigitsArray();
1372 SetDigitsArray(0x0);
f2a509af 1373
8c555625 1374}
1375
f8cf550c 1376//__________________________________________________________________
0a61bf9d 1377void AliTPC::Hits2SDigits2(Int_t eventnumber)
f8cf550c 1378{
1379
1380 //-----------------------------------------------------------
1381 // summable digits - 16 bit "ADC", no noise, no saturation
1382 //-----------------------------------------------------------
1383
8c2b3fd7 1384 //----------------------------------------------------
1385 // Loop over all sectors for a single event
1386 //----------------------------------------------------
88cb7938 1387
1388 AliRunLoader* rl = fLoader->GetRunLoader();
1389
1390 rl->GetEvent(eventnumber);
8c2b3fd7 1391 if (fLoader->TreeH() == 0x0) {
1392 if(fLoader->LoadHits()) {
1393 AliError("Can not load hits.");
1394 return;
1395 }
1396 }
88cb7938 1397 SetTreeAddress();
f8cf550c 1398
afc42102 1399
8c2b3fd7 1400 if (fLoader->TreeS() == 0x0 ) {
1401 fLoader->MakeTree("S");
1402 }
88cb7938 1403
afc42102 1404 if(fDefaults == 0) SetDefaults();
88cb7938 1405
407ff276 1406 GenerNoise(500000); //create table with noise
afc42102 1407 //setup TPCDigitsArray
1408
1409 if(GetDigitsArray()) delete GetDigitsArray();
1410
88cb7938 1411
afc42102 1412 AliTPCDigitsArray *arr = new AliTPCDigitsArray;
1413 arr->SetClass("AliSimDigits");
1414 arr->Setup(fTPCParam);
88cb7938 1415 arr->MakeTree(fLoader->TreeS());
1416
afc42102 1417 SetDigitsArray(arr);
1418
afc42102 1419 fDigitsSwitch=1; // summable digits
5f16d237 1420
1421 // set zero suppression to "0"
f546a4b4 1422 // here LHC clock phase
1423 Float_t lhcph = 0.;
1424 switch (fLHCclockPhaseSw){
1425 case 0:
1426 // no phase
1427 lhcph=0.;
1428 break;
1429 case 1:
1430 // random phase
1431 lhcph = (Int_t)(gRandom->Rndm()/0.25);
1432 break;
1433 case 2:
1434 lhcph=0.;
1435 // not implemented yet
1436 break;
1437 }
1438 // adding phase to the TreeS user info
1439
1440 fLoader->TreeS()->GetUserInfo()->Add(new TParameter<float>("lhcphase0",lhcph));
5f16d237 1441
1442 fTPCParam->SetZeroSup(0);
f8cf550c 1443
8c2b3fd7 1444 for(Int_t isec=0;isec<fTPCParam->GetNSector();isec++)
1445 if (IsSectorActive(isec)) {
1446 Hits2DigitsSector(isec);
1447 }
88cb7938 1448
8c2b3fd7 1449 fLoader->WriteSDigits("OVERWRITE");
88cb7938 1450
1451//this line prevents the crash in the similar one
1452//on the beginning of this method
1453//destructor attempts to reset the tree, which is deleted by the loader
1454//need to be redesign
8c2b3fd7 1455 if(GetDigitsArray()) delete GetDigitsArray();
1456 SetDigitsArray(0x0);
f8cf550c 1457}
0a61bf9d 1458//__________________________________________________________________
88cb7938 1459
0a61bf9d 1460void AliTPC::Hits2SDigits()
1461{
1462
1463 //-----------------------------------------------------------
1464 // summable digits - 16 bit "ADC", no noise, no saturation
1465 //-----------------------------------------------------------
1466
506e60a6 1467 if (!fTPCParam->IsGeoRead()){
1468 //
1469 // read transformation matrices for gGeoManager
1470 //
1471 fTPCParam->ReadGeoMatrices();
1472 }
1473
85a5290f 1474 fLoader->LoadHits("read");
1475 fLoader->LoadSDigits("recreate");
1476 AliRunLoader* runLoader = fLoader->GetRunLoader();
0a61bf9d 1477
85a5290f 1478 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1479 runLoader->GetEvent(iEvent);
1480 SetTreeAddress();
1481 SetActiveSectors();
1482 Hits2SDigits2(iEvent);
1483 }
f546a4b4 1484
85a5290f 1485 fLoader->UnloadHits();
1486 fLoader->UnloadSDigits();
23effa8e 1487 if (fDebugStreamer) {
1488 delete fDebugStreamer;
1489 fDebugStreamer=0;
1490 }
0a61bf9d 1491}
fe4da5cc 1492//_____________________________________________________________________________
88cb7938 1493
8c555625 1494void AliTPC::Hits2DigitsSector(Int_t isec)
fe4da5cc 1495{
8c555625 1496 //-------------------------------------------------------------------
fe4da5cc 1497 // TPC conversion from hits to digits.
8c555625 1498 //-------------------------------------------------------------------
1499
1500 //-----------------------------------------------------------------
1501 // Origin: Marek Kowalski IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
1502 //-----------------------------------------------------------------
1503
fe4da5cc 1504 //-------------------------------------------------------
8c555625 1505 // Get the access to the track hits
fe4da5cc 1506 //-------------------------------------------------------
8c555625 1507
afc42102 1508 // check if the parameters are set - important if one calls this method
1509 // directly, not from the Hits2Digits
1510
1511 if(fDefaults == 0) SetDefaults();
cc80f89e 1512
cb5b8b21 1513 TTree *tH = fLoader->TreeH(); // pointer to the hits tree
8c2b3fd7 1514 if (tH == 0x0) {
1515 AliFatal("Can not find TreeH in folder");
1516 return;
1517 }
88cb7938 1518
73042f01 1519 Stat_t ntracks = tH->GetEntries();
8c555625 1520
8c555625 1521
8c555625 1522
8c555625 1523 TObjArray **row;
fe4da5cc 1524
8c2b3fd7 1525 Int_t nrows =fTPCParam->GetNRow(isec);
8c555625 1526
8c2b3fd7 1527 row= new TObjArray* [nrows+2]; // 2 extra rows for cross talk
fe4da5cc 1528
8c2b3fd7 1529 MakeSector(isec,nrows,tH,ntracks,row);
8c555625 1530
8c2b3fd7 1531 //--------------------------------------------------------
1532 // Digitize this sector, row by row
1533 // row[i] is the pointer to the TObjArray of TVectors,
1534 // each one containing electrons accepted on this
1535 // row, assigned into tracks
1536 //--------------------------------------------------------
8c555625 1537
8c2b3fd7 1538 Int_t i;
8c555625 1539
8c2b3fd7 1540 if (fDigitsArray->GetTree()==0) {
1541 AliFatal("Tree not set in fDigitsArray");
1542 }
8c555625 1543
8c2b3fd7 1544 for (i=0;i<nrows;i++){
1545
1546 AliDigits * dig = fDigitsArray->CreateRow(isec,i);
8c555625 1547
8c2b3fd7 1548 DigitizeRow(i,isec,row);
8c555625 1549
8c2b3fd7 1550 fDigitsArray->StoreRow(isec,i);
8c555625 1551
8c2b3fd7 1552 Int_t ndig = dig->GetDigitSize();
88cb7938 1553
8c2b3fd7 1554 AliDebug(10,
1555 Form("*** Sector, row, compressed digits %d %d %d ***\n",
1556 isec,i,ndig));
792bb11c 1557
8c2b3fd7 1558 fDigitsArray->ClearRow(isec,i);
8c555625 1559
cc80f89e 1560
8c2b3fd7 1561 } // end of the sector digitization
8c555625 1562
8c2b3fd7 1563 for(i=0;i<nrows+2;i++){
1564 row[i]->Delete();
1565 delete row[i];
1566 }
cc80f89e 1567
8c2b3fd7 1568 delete [] row; // delete the array of pointers to TObjArray-s
8c555625 1569
8c555625 1570
cc80f89e 1571} // end of Hits2DigitsSector
8c555625 1572
8c555625 1573
8c555625 1574//_____________________________________________________________________________
cc80f89e 1575void AliTPC::DigitizeRow(Int_t irow,Int_t isec,TObjArray **rows)
8c555625 1576{
1577 //-----------------------------------------------------------
1578 // Single row digitization, coupling from the neighbouring
1579 // rows taken into account
1580 //-----------------------------------------------------------
1581
1582 //-----------------------------------------------------------------
1583 // Origin: Marek Kowalski IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
cc80f89e 1584 // Modified: Marian Ivanov GSI Darmstadt, m.ivanov@gsi.de
8c555625 1585 //-----------------------------------------------------------------
1586
8c555625 1587 Float_t zerosup = fTPCParam->GetZeroSup();
28c55807 1588 AliTPCCalPad * gainTPC = AliTPCcalibDB::Instance()->GetDedxGainFactor();
1589 AliTPCCalPad * noiseTPC = AliTPCcalibDB::Instance()->GetPadNoise();
1590 AliTPCCalROC * gainROC = gainTPC->GetCalROC(isec); // pad gains per given sector
1591 AliTPCCalROC * noiseROC = noiseTPC->GetCalROC(isec); // noise per given sector
1592
8c2b3fd7 1593
cc80f89e 1594 fCurrentIndex[1]= isec;
8c555625 1595
8c555625 1596
73042f01 1597 Int_t nofPads = fTPCParam->GetNPads(isec,irow);
1598 Int_t nofTbins = fTPCParam->GetMaxTBin();
1599 Int_t indexRange[4];
8c555625 1600 //
1601 // Integrated signal for this row
1602 // and a single track signal
cc80f89e 1603 //
de61d5d5 1604
e8d02863 1605 TMatrixF *m1 = new TMatrixF(0,nofPads,0,nofTbins); // integrated
1606 TMatrixF *m2 = new TMatrixF(0,nofPads,0,nofTbins); // single
8c555625 1607 //
e8d02863 1608 TMatrixF &total = *m1;
8c555625 1609
1610 // Array of pointers to the label-signal list
1611
73042f01 1612 Int_t nofDigits = nofPads*nofTbins; // number of digits for this row
1613 Float_t **pList = new Float_t* [nofDigits];
8c555625 1614
1615 Int_t lp;
cc80f89e 1616 Int_t i1;
73042f01 1617 for(lp=0;lp<nofDigits;lp++)pList[lp]=0; // set all pointers to NULL
8c555625 1618 //
cc80f89e 1619 //calculate signal
1620 //
b584c7dd 1621 Int_t row1=irow;
1622 Int_t row2=irow+2;
cc80f89e 1623 for (Int_t row= row1;row<=row2;row++){
1624 Int_t nTracks= rows[row]->GetEntries();
1625 for (i1=0;i1<nTracks;i1++){
1626 fCurrentIndex[2]= row;
b584c7dd 1627 fCurrentIndex[3]=irow+1;
1628 if (row==irow+1){
cc80f89e 1629 m2->Zero(); // clear single track signal matrix
73042f01 1630 Float_t trackLabel = GetSignal(rows[row],i1,m2,m1,indexRange);
1631 GetList(trackLabel,nofPads,m2,indexRange,pList);
cc80f89e 1632 }
73042f01 1633 else GetSignal(rows[row],i1,0,m1,indexRange);
cc80f89e 1634 }
8c555625 1635 }
cc80f89e 1636
8c555625 1637 Int_t tracks[3];
8c555625 1638
cc80f89e 1639 AliDigits *dig = fDigitsArray->GetRow(isec,irow);
de61d5d5 1640 Int_t gi=-1;
1641 Float_t fzerosup = zerosup+0.5;
1642 for(Int_t it=0;it<nofTbins;it++){
de61d5d5 1643 for(Int_t ip=0;ip<nofPads;ip++){
1644 gi++;
8c2b3fd7 1645 Float_t q=total(ip,it);
28c55807 1646 if(fDigitsSwitch == 0){
1647 Float_t gain = gainROC->GetValue(irow,ip); // get gain for given - pad-row pad
1648 Float_t noisePad = noiseROC->GetValue(irow,ip);
1649 //
1650 q*=gain;
1651 q+=GetNoise()*noisePad;
de61d5d5 1652 if(q <=fzerosup) continue; // do not fill zeros
68771f83 1653 q = TMath::Nint(q);
e93a083a 1654 if(q >= fTPCParam->GetADCSat()) q = fTPCParam->GetADCSat() - 1; // saturation
cc80f89e 1655
f8cf550c 1656 }
1657
1658 else {
8c2b3fd7 1659 if(q <= 0.) continue; // do not fill zeros
1660 if(q>2000.) q=2000.;
1661 q *= 16.;
1662 q = TMath::Nint(q);
f8cf550c 1663 }
8c555625 1664
1665 //
1666 // "real" signal or electronic noise (list = -1)?
1667 //
1668
1669 for(Int_t j1=0;j1<3;j1++){
407ff276 1670 tracks[j1] = (pList[gi]) ?(Int_t)(*(pList[gi]+j1)) : -2;
8c555625 1671 }
1672
cc80f89e 1673//Begin_Html
1674/*
1675 <A NAME="AliDigits"></A>
1676 using of AliDigits object
1677*/
1678//End_Html
1679 dig->SetDigitFast((Short_t)q,it,ip);
8c2b3fd7 1680 if (fDigitsArray->IsSimulated()) {
1681 ((AliSimDigits*)dig)->SetTrackIDFast(tracks[0],it,ip,0);
1682 ((AliSimDigits*)dig)->SetTrackIDFast(tracks[1],it,ip,1);
1683 ((AliSimDigits*)dig)->SetTrackIDFast(tracks[2],it,ip,2);
1684 }
8c555625 1685
1686 } // end of loop over time buckets
1687 } // end of lop over pads
58c45a5e 1688 //
1689 // test
1690 //
1691 //
8c555625 1692
58c45a5e 1693 // glitch filters if normal simulated digits
1694 //
1695 if(!fDigitsSwitch) ((AliSimDigits*)dig)->GlitchFilter();
8c555625 1696 //
1697 // This row has been digitized, delete nonused stuff
1698 //
1699
73042f01 1700 for(lp=0;lp<nofDigits;lp++){
8c555625 1701 if(pList[lp]) delete [] pList[lp];
1702 }
1703
1704 delete [] pList;
1705
1706 delete m1;
1707 delete m2;
8c555625 1708
1709} // end of DigitizeRow
cc80f89e 1710
8c555625 1711//_____________________________________________________________________________
cc80f89e 1712
de61d5d5 1713Float_t AliTPC::GetSignal(TObjArray *p1, Int_t ntr,
e8d02863 1714 TMatrixF *m1, TMatrixF *m2,Int_t *indexRange)
8c555625 1715{
1716
1717 //---------------------------------------------------------------
1718 // Calculates 2-D signal (pad,time) for a single track,
1719 // returns a pointer to the signal matrix and the track label
1720 // No digitization is performed at this level!!!
1721 //---------------------------------------------------------------
1722
1723 //-----------------------------------------------------------------
1724 // Origin: Marek Kowalski IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
cc80f89e 1725 // Modified: Marian Ivanov
8c555625 1726 //-----------------------------------------------------------------
1727
8c2b3fd7 1728 TVector *tv;
de61d5d5 1729
8c2b3fd7 1730 tv = (TVector*)p1->At(ntr); // pointer to a track
1731 TVector &v = *tv;
8c555625 1732
1733 Float_t label = v(0);
506e60a6 1734 Int_t centralPad = (fTPCParam->GetNPads(fCurrentIndex[1],fCurrentIndex[3]-1))/2;
8c555625 1735
e61fd20d 1736 Int_t nElectrons = (tv->GetNrows()-1)/5;
73042f01 1737 indexRange[0]=9999; // min pad
1738 indexRange[1]=-1; // max pad
1739 indexRange[2]=9999; //min time
1740 indexRange[3]=-1; // max time
8c555625 1741
e8d02863 1742 TMatrixF &signal = *m1;
1743 TMatrixF &total = *m2;
8c555625 1744 //
f546a4b4 1745 // Get LHC clock phase
1746 //
1747 TParameter<float> *ph;
1748 if(fDigitsSwitch){// s-digits
1749 ph = (TParameter<float>*)fLoader->TreeS()->GetUserInfo()->FindObject("lhcphase0");
1750 }
1751 else{ // normal digits
1752 ph = (TParameter<float>*)fLoader->TreeD()->GetUserInfo()->FindObject("lhcphase0");
1753 }
8c555625 1754 // Loop over all electrons
1755 //
8c555625 1756 for(Int_t nel=0; nel<nElectrons; nel++){
e61fd20d 1757 Int_t idx=nel*5;
cc80f89e 1758 Float_t aval = v(idx+4);
1759 Float_t eltoadcfac=aval*fTPCParam->GetTotalNormFac();
e61fd20d 1760 Float_t xyz[4]={v(idx+1),v(idx+2),v(idx+3),v(idx+5)};
f546a4b4 1761 Int_t n = ((AliTPCParamSR*)fTPCParam)->CalcResponseFast(xyz,fCurrentIndex,
1762 fCurrentIndex[3],ph->GetVal());
de61d5d5 1763
1764 Int_t *index = fTPCParam->GetResBin(0);
1765 Float_t *weight = & (fTPCParam->GetResWeight(0));
1766
1767 if (n>0) for (Int_t i =0; i<n; i++){
e30e6a9c 1768 Int_t pad=index[1]+centralPad; //in digit coordinates central pad has coordinate 0
8c2b3fd7 1769
1770 if (pad>=0){
1771 Int_t time=index[2];
1772 Float_t qweight = *(weight)*eltoadcfac;
1773
1774 if (m1!=0) signal(pad,time)+=qweight;
1775 total(pad,time)+=qweight;
1776 if (indexRange[0]>pad) indexRange[0]=pad;
1777 if (indexRange[1]<pad) indexRange[1]=pad;
1778 if (indexRange[2]>time) indexRange[2]=time;
1779 if (indexRange[3]<time) indexRange[3]=time;
1780
1781 index+=3;
1782 weight++;
1783
1784 }
cc80f89e 1785 }
8c555625 1786 } // end of loop over electrons
cc80f89e 1787
8c555625 1788 return label; // returns track label when finished
1789}
1790
1791//_____________________________________________________________________________
e8d02863 1792void AliTPC::GetList(Float_t label,Int_t np,TMatrixF *m,
de61d5d5 1793 Int_t *indexRange, Float_t **pList)
8c555625 1794{
1795 //----------------------------------------------------------------------
1796 // Updates the list of tracks contributing to digits for a given row
1797 //----------------------------------------------------------------------
1798
1799 //-----------------------------------------------------------------
1800 // Origin: Marek Kowalski IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
1801 //-----------------------------------------------------------------
1802
e8d02863 1803 TMatrixF &signal = *m;
8c555625 1804
1805 // lop over nonzero digits
1806
73042f01 1807 for(Int_t it=indexRange[2];it<indexRange[3]+1;it++){
1808 for(Int_t ip=indexRange[0];ip<indexRange[1]+1;ip++){
8c555625 1809
1810
8c2b3fd7 1811 // accept only the contribution larger than 500 electrons (1/2 s_noise)
921bf71a 1812
8c2b3fd7 1813 if(signal(ip,it)<0.5) continue;
921bf71a 1814
8c2b3fd7 1815 Int_t globalIndex = it*np+ip; // globalIndex starts from 0!
8c555625 1816
8c2b3fd7 1817 if(!pList[globalIndex]){
8c555625 1818
8c2b3fd7 1819 //
1820 // Create new list (6 elements - 3 signals and 3 labels),
1821 //
8c555625 1822
8c2b3fd7 1823 pList[globalIndex] = new Float_t [6];
8c555625 1824
8c2b3fd7 1825 // set list to -1
1826
1827 *pList[globalIndex] = -1.;
1828 *(pList[globalIndex]+1) = -1.;
1829 *(pList[globalIndex]+2) = -1.;
1830 *(pList[globalIndex]+3) = -1.;
1831 *(pList[globalIndex]+4) = -1.;
1832 *(pList[globalIndex]+5) = -1.;
1833
1834 *pList[globalIndex] = label;
1835 *(pList[globalIndex]+3) = signal(ip,it);
1836 }
1837 else {
8c555625 1838
8c2b3fd7 1839 // check the signal magnitude
8c555625 1840
8c2b3fd7 1841 Float_t highest = *(pList[globalIndex]+3);
1842 Float_t middle = *(pList[globalIndex]+4);
1843 Float_t lowest = *(pList[globalIndex]+5);
1844
1845 //
1846 // compare the new signal with already existing list
1847 //
1848
1849 if(signal(ip,it)<lowest) continue; // neglect this track
8c555625 1850
8c2b3fd7 1851 //
8c555625 1852
8c2b3fd7 1853 if (signal(ip,it)>highest){
1854 *(pList[globalIndex]+5) = middle;
1855 *(pList[globalIndex]+4) = highest;
1856 *(pList[globalIndex]+3) = signal(ip,it);
1857
1858 *(pList[globalIndex]+2) = *(pList[globalIndex]+1);
1859 *(pList[globalIndex]+1) = *pList[globalIndex];
1860 *pList[globalIndex] = label;
1861 }
1862 else if (signal(ip,it)>middle){
1863 *(pList[globalIndex]+5) = middle;
1864 *(pList[globalIndex]+4) = signal(ip,it);
1865
1866 *(pList[globalIndex]+2) = *(pList[globalIndex]+1);
1867 *(pList[globalIndex]+1) = label;
1868 }
1869 else{
1870 *(pList[globalIndex]+5) = signal(ip,it);
1871 *(pList[globalIndex]+2) = label;
1872 }
1873 }
1874
8c555625 1875 } // end of loop over pads
1876 } // end of loop over time bins
1877
8c555625 1878}//end of GetList
1879//___________________________________________________________________
1880void AliTPC::MakeSector(Int_t isec,Int_t nrows,TTree *TH,
1881 Stat_t ntracks,TObjArray **row)
1882{
1883
1884 //-----------------------------------------------------------------
1885 // Prepares the sector digitization, creates the vectors of
1886 // tracks for each row of this sector. The track vector
1887 // contains the track label and the position of electrons.
1888 //-----------------------------------------------------------------
1889
79fa3dc8 1890 //
1891 // The trasport of the electrons through TPC drift volume
1892 // Drift (drift velocity + velocity map(not yet implemented)))
1893 // Application of the random processes (diffusion, gas gain)
1894 // Systematic effects (ExB effect in drift volume + ROCs)
1895 //
1896 // Algorithm:
1897 // Loop over primary electrons:
1898 // Creation of the secondary electrons
1899 // Loop over electrons (primary+ secondaries)
1900 // Global coordinate frame:
1901 // 1. Skip electrons if attached
1902 // 2. ExB effect in drift volume
1903 // a.) Simulation calib->GetExB()->CorrectInverse(dxyz0,dxyz1);
1904 // b.) Reconstruction - calib->GetExB()->CorrectInverse(dxyz0,dxyz1);
1905 // 3. Generation of gas gain (Random - Exponential distribution)
1906 // 4. TransportElectron function (diffusion)
1907 //
1908 // 5. Conversion to the local coordinate frame pad-row, pad, timebin
1909 // 6. Apply Time0 shift - AliTPCCalPad class
1910 // a.) Plus sign in simulation
1911 // b.) Minus sign in reconstruction
1912 // end of loop
1913 //
8c555625 1914 //-----------------------------------------------------------------
1915 // Origin: Marek Kowalski IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
79fa3dc8 1916 // Origin: Marian Ivanov, marian.ivanov@cern.ch
8c555625 1917 //-----------------------------------------------------------------
79fa3dc8 1918 AliTPCcalibDB* const calib=AliTPCcalibDB::Instance();
7bab76cf 1919 if (gAlice){ // Set correctly the magnetic field in the ExB calculation
9ac65c7c 1920 if (!calib->GetExB()){
1921 AliMagF * field = ((AliMagF*)TGeoGlobalMagField::Instance()->GetField());
1922 if (field) {
1923 calib->SetExBField(field);
1924 }
7bab76cf 1925 }
1926 }
8c555625 1927
cc80f89e 1928 Float_t gasgain = fTPCParam->GetGasGain();
8261c673 1929 gasgain = gasgain/fGainFactor;
8c555625 1930 Int_t i;
e61fd20d 1931 Float_t xyz[5];
8c555625 1932
1933 AliTPChit *tpcHit; // pointer to a sigle TPC hit
39c8eb58 1934 //MI change
1935 TBranch * branch=0;
792bb11c 1936 if (fHitType>1) branch = TH->GetBranch("TPC2");
39c8eb58 1937 else branch = TH->GetBranch("TPC");
1938
8c555625 1939
1940 //----------------------------------------------
1941 // Create TObjArray-s, one for each row,
8c2b3fd7 1942 // each TObjArray will store the TVectors
1943 // of electrons, one TVectors per each track.
8c555625 1944 //----------------------------------------------
1945
b584c7dd 1946 Int_t *nofElectrons = new Int_t [nrows+2]; // electron counter for each row
8c2b3fd7 1947 TVector **tracks = new TVector* [nrows+2]; //pointers to the track vectors
de61d5d5 1948
b584c7dd 1949 for(i=0; i<nrows+2; i++){
8c555625 1950 row[i] = new TObjArray;
f74bb6f5 1951 nofElectrons[i]=0;
1952 tracks[i]=0;
8c555625 1953 }
8c555625 1954
37831078 1955
1956
8c555625 1957 //--------------------------------------------------------------------
1958 // Loop over tracks, the "track" contains the full history
1959 //--------------------------------------------------------------------
8c2b3fd7 1960
8c555625 1961 Int_t previousTrack,currentTrack;
1962 previousTrack = -1; // nothing to store so far!
1963
1964 for(Int_t track=0;track<ntracks;track++){
39c8eb58 1965 Bool_t isInSector=kTRUE;
8c555625 1966 ResetHits();
792bb11c 1967 isInSector = TrackInVolume(isec,track);
39c8eb58 1968 if (!isInSector) continue;
1969 //MI change
1970 branch->GetEntry(track); // get next track
8c2b3fd7 1971
39c8eb58 1972 //M.I. changes
8c555625 1973
39c8eb58 1974 tpcHit = (AliTPChit*)FirstHit(-1);
8c555625 1975
1976 //--------------------------------------------------------------
1977 // Loop over hits
1978 //--------------------------------------------------------------
1979
8c555625 1980
39c8eb58 1981 while(tpcHit){
8c555625 1982
1983 Int_t sector=tpcHit->fSector; // sector number
39c8eb58 1984 if(sector != isec){
1985 tpcHit = (AliTPChit*) NextHit();
1986 continue;
1987 }
8c555625 1988
daf14717 1989 // Remove hits which arrive before the TPC opening gate signal
a1ec4d07 1990 if(((fTPCParam->GetZLength(isec)-TMath::Abs(tpcHit->Z()))
daf14717 1991 /fTPCParam->GetDriftV()+tpcHit->Time())<fTPCParam->GetGateDelay()) {
1992 tpcHit = (AliTPChit*) NextHit();
1993 continue;
1994 }
1995
8c2b3fd7 1996 currentTrack = tpcHit->Track(); // track number
39c8eb58 1997
8c2b3fd7 1998 if(currentTrack != previousTrack){
8c555625 1999
8c2b3fd7 2000 // store already filled fTrack
8c555625 2001
8c2b3fd7 2002 for(i=0;i<nrows+2;i++){
2003 if(previousTrack != -1){
2004 if(nofElectrons[i]>0){
2005 TVector &v = *tracks[i];
2006 v(0) = previousTrack;
e61fd20d 2007 tracks[i]->ResizeTo(5*nofElectrons[i]+1); // shrink if necessary
8c2b3fd7 2008 row[i]->Add(tracks[i]);
2009 }
2010 else {
2011 delete tracks[i]; // delete empty TVector
2012 tracks[i]=0;
2013 }
2014 }
2015
2016 nofElectrons[i]=0;
e61fd20d 2017 tracks[i] = new TVector(601); // TVectors for the next fTrack
8c2b3fd7 2018
2019 } // end of loop over rows
8c555625 2020
8c2b3fd7 2021 previousTrack=currentTrack; // update track label
2022 }
8c555625 2023
8c2b3fd7 2024 Int_t qI = (Int_t) (tpcHit->fQ); // energy loss (number of electrons)
8c555625 2025
8c2b3fd7 2026 //---------------------------------------------------
2027 // Calculate the electron attachment probability
2028 //---------------------------------------------------
8c555625 2029
cc80f89e 2030
a1ec4d07 2031 Float_t time = 1.e6*(fTPCParam->GetZLength(isec)-TMath::Abs(tpcHit->Z()))
8c2b3fd7 2032 /fTPCParam->GetDriftV();
2033 // in microseconds!
2034 Float_t attProb = fTPCParam->GetAttCoef()*
2035 fTPCParam->GetOxyCont()*time; // fraction!
8c555625 2036
8c2b3fd7 2037 //-----------------------------------------------
2038 // Loop over electrons
2039 //-----------------------------------------------
2040 Int_t index[3];
2041 index[1]=isec;
2042 for(Int_t nel=0;nel<qI;nel++){
2043 // skip if electron lost due to the attachment
2044 if((gRandom->Rndm(0)) < attProb) continue; // electron lost!
79fa3dc8 2045
2046 //
2047 // ExB effect
2048 //
2049 Double_t dxyz0[3],dxyz1[3];
2050 dxyz0[0]=tpcHit->X();
2051 dxyz0[1]=tpcHit->Y();
2052 dxyz0[2]=tpcHit->Z();
2053 if (calib->GetExB()){
2054 calib->GetExB()->CorrectInverse(dxyz0,dxyz1);
2055 }else{
2056 AliError("Not valid ExB calibration");
2057 dxyz1[0]=tpcHit->X();
2058 dxyz1[1]=tpcHit->Y();
2059 dxyz1[2]=tpcHit->Z();
2060 }
2061 xyz[0]=dxyz1[0];
2062 xyz[1]=dxyz1[1];
2063 xyz[2]=dxyz1[2];
2064 //
2065 //
8c2b3fd7 2066 //
2067 // protection for the nonphysical avalanche size (10**6 maximum)
2068 //
2069 Double_t rn=TMath::Max(gRandom->Rndm(0),1.93e-22);
2070 xyz[3]= (Float_t) (-gasgain*TMath::Log(rn));
2071 index[0]=1;
79fa3dc8 2072
8c2b3fd7 2073 TransportElectron(xyz,index);
2074 Int_t rowNumber;
23effa8e 2075 Int_t padrow = fTPCParam->GetPadRow(xyz,index);
79fa3dc8 2076 //
2077 // Add Time0 correction due unisochronity
2078 // xyz[0] - pad row coordinate
2079 // xyz[1] - pad coordinate
2080 // xyz[2] - is in now time bin coordinate system
2081 Float_t correction =0;
2082 if (calib->GetPadTime0()){
23effa8e 2083 if (!calib->GetPadTime0()->GetCalROC(isec)) continue;
2084 Int_t npads = fTPCParam->GetNPads(isec,padrow);
17807087 2085 // Int_t pad = TMath::Nint(xyz[1]+fTPCParam->GetNPads(isec,TMath::Nint(xyz[0]))*0.5);
23effa8e 2086 // pad numbering from -npads/2 .. npads/2-1
2087 Int_t pad = TMath::Nint(xyz[1]+npads/2);
17807087 2088 if (pad<0) pad=0;
2089 if (pad>=npads) pad=npads-1;
23effa8e 2090 correction = calib->GetPadTime0()->GetCalROC(isec)->GetValue(padrow,pad);
2091 // printf("%d\t%d\t%d\t%f\n",isec,padrow,pad,correction);
2092 if (fDebugStreamer){
2093 (*fDebugStreamer)<<"Time0"<<
2094 "isec="<<isec<<
2095 "padrow="<<padrow<<
2096 "pad="<<pad<<
2097 "x0="<<xyz[0]<<
2098 "x1="<<xyz[1]<<
2099 "x2="<<xyz[2]<<
2100 "hit.="<<tpcHit<<
2101 "cor="<<correction<<
2102 "\n";
2103 }
79fa3dc8 2104 }
2105 xyz[2]+=correction;
c2940e77 2106 xyz[2]+=fTPCParam->GetNTBinsL1(); // adding Level 1 time bin offset
79fa3dc8 2107 //
e61fd20d 2108 // Electron track time (for pileup simulation)
c2940e77 2109 xyz[2]+=tpcHit->Time()/fTPCParam->GetTSample(); // adding time of flight
2110 xyz[4] =0;
2111
2112 //
8c2b3fd7 2113 // row 0 - cross talk from the innermost row
2114 // row fNRow+1 cross talk from the outermost row
2115 rowNumber = index[2]+1;
2116 //transform position to local digit coordinates
2117 //relative to nearest pad row
2118 if ((rowNumber<0)||rowNumber>fTPCParam->GetNRow(isec)+1) continue;
4555e0d2 2119 /* Float_t x1,y1;
8c2b3fd7 2120 if (isec <fTPCParam->GetNInnerSector()) {
2121 x1 = xyz[1]*fTPCParam->GetInnerPadPitchWidth();
2122 y1 = fTPCParam->GetYInner(rowNumber);
2123 }
2124 else{
2125 x1=xyz[1]*fTPCParam->GetOuterPadPitchWidth();
2126 y1 = fTPCParam->GetYOuter(rowNumber);
2127 }
2128 // gain inefficiency at the wires edges - linear
2129 x1=TMath::Abs(x1);
2130 y1-=1.;
4555e0d2 2131 if(x1>y1) xyz[3]*=TMath::Max(1.e-6,(y1-x1+1.)); */
8c2b3fd7 2132
2133 nofElectrons[rowNumber]++;
2134 //----------------------------------
2135 // Expand vector if necessary
2136 //----------------------------------
2137 if(nofElectrons[rowNumber]>120){
2138 Int_t range = tracks[rowNumber]->GetNrows();
e61fd20d 2139 if((nofElectrons[rowNumber])>(range-1)/5){
8c2b3fd7 2140
e61fd20d 2141 tracks[rowNumber]->ResizeTo(range+500); // Add 100 electrons
fe4da5cc 2142 }
8c2b3fd7 2143 }
2144
2145 TVector &v = *tracks[rowNumber];
e61fd20d 2146 Int_t idx = 5*nofElectrons[rowNumber]-4;
8c2b3fd7 2147 Real_t * position = &(((TVector&)v)(idx)); //make code faster
e61fd20d 2148 memcpy(position,xyz,5*sizeof(Float_t));
8c2b3fd7 2149
2150 } // end of loop over electrons
39c8eb58 2151
8c2b3fd7 2152 tpcHit = (AliTPChit*)NextHit();
2153
2154 } // end of loop over hits
2155 } // end of loop over tracks
8c555625 2156
2157 //
2158 // store remaining track (the last one) if not empty
2159 //
8c2b3fd7 2160
2161 for(i=0;i<nrows+2;i++){
2162 if(nofElectrons[i]>0){
2163 TVector &v = *tracks[i];
2164 v(0) = previousTrack;
e61fd20d 2165 tracks[i]->ResizeTo(5*nofElectrons[i]+1); // shrink if necessary
8c2b3fd7 2166 row[i]->Add(tracks[i]);
2167 }
2168 else{
2169 delete tracks[i];
2170 tracks[i]=0;
2171 }
2172 }
2173
2174 delete [] tracks;
2175 delete [] nofElectrons;
8c555625 2176
cc80f89e 2177} // end of MakeSector
8c555625 2178
fe4da5cc 2179
2180//_____________________________________________________________________________
2181void AliTPC::Init()
2182{
2183 //
2184 // Initialise TPC detector after definition of geometry
2185 //
8c2b3fd7 2186 AliDebug(1,"*********************************************");
fe4da5cc 2187}
2188
fe4da5cc 2189//_____________________________________________________________________________
2190void AliTPC::ResetDigits()
2191{
2192 //
2193 // Reset number of digits and the digits array for this detector
fe4da5cc 2194 //
2195 fNdigits = 0;
cc80f89e 2196 if (fDigits) fDigits->Clear();
fe4da5cc 2197}
2198
fe4da5cc 2199
fe4da5cc 2200
2201//_____________________________________________________________________________
2202void AliTPC::SetSens(Int_t sens)
2203{
8c555625 2204
2205 //-------------------------------------------------------------
2206 // Activates/deactivates the sensitive strips at the center of
2207 // the pad row -- this is for the space-point resolution calculations
2208 //-------------------------------------------------------------
2209
2210 //-----------------------------------------------------------------
2211 // Origin: Marek Kowalski IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
2212 //-----------------------------------------------------------------
2213
fe4da5cc 2214 fSens = sens;
2215}
2b06d5c3 2216
4b0fdcad 2217
73042f01 2218void AliTPC::SetSide(Float_t side=0.)
4b0fdcad 2219{
73042f01 2220 // choice of the TPC side
2221
4b0fdcad 2222 fSide = side;
2223
2224}
cc80f89e 2225//_____________________________________________________________________________
2226
2227void AliTPC::TransportElectron(Float_t *xyz, Int_t *index)
2228{
2229 //
2230 // electron transport taking into account:
2231 // 1. diffusion,
2232 // 2.ExB at the wires
2233 // 3. nonisochronity
2234 //
2235 // xyz and index must be already transformed to system 1
2236 //
2237
23effa8e 2238 fTPCParam->Transform1to2(xyz,index); // mis-alignment applied in this step
cc80f89e 2239
2240 //add diffusion
2241 Float_t driftl=xyz[2];
2242 if(driftl<0.01) driftl=0.01;
2243 driftl=TMath::Sqrt(driftl);
73042f01 2244 Float_t sigT = driftl*(fTPCParam->GetDiffT());
2245 Float_t sigL = driftl*(fTPCParam->GetDiffL());
2246 xyz[0]=gRandom->Gaus(xyz[0],sigT);
2247 xyz[1]=gRandom->Gaus(xyz[1],sigT);
2248 xyz[2]=gRandom->Gaus(xyz[2],sigL);
cc80f89e 2249
2250 // ExB
2251
2252 if (fTPCParam->GetMWPCReadout()==kTRUE){
b584c7dd 2253 Float_t dx = fTPCParam->Transform2to2NearestWire(xyz,index);
cc80f89e 2254 xyz[1]+=dx*(fTPCParam->GetOmegaTau());
2255 }
79fa3dc8 2256 //add nonisochronity (not implemented yet)
2257
2258
1283eee5 2259}
fe4da5cc 2260
fe4da5cc 2261ClassImp(AliTPChit)
179c6296 2262 //______________________________________________________________________
2263 AliTPChit::AliTPChit()
2264 :AliHit(),
2265 fSector(0),
2266 fPadRow(0),
2267 fQ(0),
2268 fTime(0)
2269{
2270 //
2271 // default
2272 //
2273
2274}
fe4da5cc 2275//_____________________________________________________________________________
179c6296 2276AliTPChit::AliTPChit(Int_t shunt, Int_t track, Int_t *vol, Float_t *hits)
2277 :AliHit(shunt,track),
2278 fSector(0),
2279 fPadRow(0),
2280 fQ(0),
2281 fTime(0)
fe4da5cc 2282{
2283 //
2284 // Creates a TPC hit object
2285 //
2286 fSector = vol[0];
2287 fPadRow = vol[1];
2288 fX = hits[0];
2289 fY = hits[1];
2290 fZ = hits[2];
2291 fQ = hits[3];
e61fd20d 2292 fTime = hits[4];
fe4da5cc 2293}
2294
39c8eb58 2295//________________________________________________________________________
792bb11c 2296// Additional code because of the AliTPCTrackHitsV2
39c8eb58 2297
a6005888 2298void AliTPC::MakeBranch(Option_t *option)
39c8eb58 2299{
2300 //
2301 // Create a new branch in the current Root Tree
2302 // The branch of fHits is automatically split
2303 // MI change 14.09.2000
8c2b3fd7 2304 AliDebug(1,"");
792bb11c 2305 if (fHitType<2) return;
39c8eb58 2306 char branchname[10];
2307 sprintf(branchname,"%s2",GetName());
2308 //
2309 // Get the pointer to the header
5cf7bbad 2310 const char *cH = strstr(option,"H");
39c8eb58 2311 //
cb5b8b21 2312 if (fTrackHits && fLoader->TreeH() && cH && fHitType&4) {
8c2b3fd7 2313 AliDebug(1,"Making branch for Type 4 Hits");
cb5b8b21 2314 fLoader->TreeH()->Branch(branchname,"AliTPCTrackHitsV2",&fTrackHits,fBufferSize,99);
8c2b3fd7 2315 }
792bb11c 2316
cb5b8b21 2317// if (fTrackHitsOld && fLoader->TreeH() && cH && fHitType&2) {
be5ffbfe 2318// AliDebug(1,"Making branch for Type 2 Hits");
2319// AliObjectBranch * branch = new AliObjectBranch(branchname,"AliTPCTrackHits",&fTrackHitsOld,
cb5b8b21 2320// fLoader->TreeH(),fBufferSize,99);
2321// fLoader->TreeH()->GetListOfBranches()->Add(branch);
be5ffbfe 2322// }
39c8eb58 2323}
2324
2325void AliTPC::SetTreeAddress()
2326{
8c2b3fd7 2327 //Sets tree address for hits
2328 if (fHitType<=1) {
2329 if (fHits == 0x0 ) fHits = new TClonesArray("AliTPChit", 176);//skowron 20.06.03
2330 AliDetector::SetTreeAddress();
2331 }
792bb11c 2332 if (fHitType>1) SetTreeAddress2();
39c8eb58 2333}
2334
2335void AliTPC::SetTreeAddress2()
2336{
2337 //
2338 // Set branch address for the TrackHits Tree
2339 //
8c2b3fd7 2340 AliDebug(1,"");
88cb7938 2341
39c8eb58 2342 TBranch *branch;
2343 char branchname[20];
2344 sprintf(branchname,"%s2",GetName());
2345 //
2346 // Branch address for hit tree
cb5b8b21 2347 TTree *treeH = fLoader->TreeH();
792bb11c 2348 if ((treeH)&&(fHitType&4)) {
39c8eb58 2349 branch = treeH->GetBranch(branchname);
8c2b3fd7 2350 if (branch) {
2351 branch->SetAddress(&fTrackHits);
2352 AliDebug(1,"fHitType&4 Setting");
2353 }
f2a509af 2354 else
8c2b3fd7 2355 AliDebug(1,"fHitType&4 Failed (can not find branch)");
f2a509af 2356
39c8eb58 2357 }
be5ffbfe 2358 // if ((treeH)&&(fHitType&2)) {
2359// branch = treeH->GetBranch(branchname);
2360// if (branch) {
2361// branch->SetAddress(&fTrackHitsOld);
2362// AliDebug(1,"fHitType&2 Setting");
2363// }
2364// else
2365// AliDebug(1,"fHitType&2 Failed (can not find branch)");
2366// }
39c8eb58 2367}
2368
2369void AliTPC::FinishPrimary()
2370{
792bb11c 2371 if (fTrackHits &&fHitType&4) fTrackHits->FlushHitStack();
be5ffbfe 2372 // if (fTrackHitsOld && fHitType&2) fTrackHitsOld->FlushHitStack();
39c8eb58 2373}
2374
2375
2376void AliTPC::AddHit2(Int_t track, Int_t *vol, Float_t *hits)
2377{
2378 //
a6005888 2379 // add hit to the list
2380
39c8eb58 2381 Int_t rtrack;
2382 if (fIshunt) {
5d12ce38 2383 int primary = gAlice->GetMCApp()->GetPrimary(track);
2384 gAlice->GetMCApp()->Particle(primary)->SetBit(kKeepBit);
39c8eb58 2385 rtrack=primary;
2386 } else {
2387 rtrack=track;
5d12ce38 2388 gAlice->GetMCApp()->FlagTrack(track);
39c8eb58 2389 }
792bb11c 2390 if (fTrackHits && fHitType&4)
39c8eb58 2391 fTrackHits->AddHitKartez(vol[0],rtrack, hits[0],
e61fd20d 2392 hits[1],hits[2],(Int_t)hits[3],hits[4]);
be5ffbfe 2393 // if (fTrackHitsOld &&fHitType&2 )
2394// fTrackHitsOld->AddHitKartez(vol[0],rtrack, hits[0],
2395// hits[1],hits[2],(Int_t)hits[3]);
792bb11c 2396
39c8eb58 2397}
2398
2399void AliTPC::ResetHits()
88cb7938 2400{
39c8eb58 2401 if (fHitType&1) AliDetector::ResetHits();
792bb11c 2402 if (fHitType>1) ResetHits2();
39c8eb58 2403}
2404
2405void AliTPC::ResetHits2()
2406{
2407 //
2408 //reset hits
792bb11c 2409 if (fTrackHits && fHitType&4) fTrackHits->Clear();
be5ffbfe 2410 // if (fTrackHitsOld && fHitType&2) fTrackHitsOld->Clear();
792bb11c 2411
39c8eb58 2412}
2413
2414AliHit* AliTPC::FirstHit(Int_t track)
2415{
792bb11c 2416 if (fHitType>1) return FirstHit2(track);
39c8eb58 2417 return AliDetector::FirstHit(track);
2418}
2419AliHit* AliTPC::NextHit()
2420{
9bdd974b 2421 //
2422 // gets next hit
2423 //
792bb11c 2424 if (fHitType>1) return NextHit2();
2425
39c8eb58 2426 return AliDetector::NextHit();
2427}
2428
2429AliHit* AliTPC::FirstHit2(Int_t track)
2430{
2431 //
2432 // Initialise the hit iterator
2433 // Return the address of the first hit for track
2434 // If track>=0 the track is read from disk
2435 // while if track<0 the first hit of the current
2436 // track is returned
2437 //
2438 if(track>=0) {
3e2e3ece 2439 gAlice->GetMCApp()->ResetHits();
cb5b8b21 2440 fLoader->TreeH()->GetEvent(track);
39c8eb58 2441 }
2442 //
792bb11c 2443 if (fTrackHits && fHitType&4) {
39c8eb58 2444 fTrackHits->First();
2445 return fTrackHits->GetHit();
2446 }
be5ffbfe 2447 // if (fTrackHitsOld && fHitType&2) {
2448// fTrackHitsOld->First();
2449// return fTrackHitsOld->GetHit();
2450// }
792bb11c 2451
39c8eb58 2452 else return 0;
2453}
2454
2455AliHit* AliTPC::NextHit2()
2456{
2457 //
2458 //Return the next hit for the current track
2459
792bb11c 2460
be5ffbfe 2461// if (fTrackHitsOld && fHitType&2) {
2462// fTrackHitsOld->Next();
2463// return fTrackHitsOld->GetHit();
2464// }
39c8eb58 2465 if (fTrackHits) {
2466 fTrackHits->Next();
2467 return fTrackHits->GetHit();
2468 }
2469 else
2470 return 0;
2471}
2472
39c8eb58 2473void AliTPC::RemapTrackHitIDs(Int_t *map)
2474{
9bdd974b 2475 //
2476 // remapping
2477 //
39c8eb58 2478 if (!fTrackHits) return;
39c8eb58 2479
be5ffbfe 2480// if (fTrackHitsOld && fHitType&2){
2481// AliObjectArray * arr = fTrackHitsOld->fTrackHitsInfo;
2482// for (UInt_t i=0;i<arr->GetSize();i++){
2483// AliTrackHitsInfo * info = (AliTrackHitsInfo *)(arr->At(i));
2484// info->fTrackID = map[info->fTrackID];
2485// }
2486// }
2487// if (fTrackHitsOld && fHitType&4){
2488 if (fTrackHits && fHitType&4){
792bb11c 2489 TClonesArray * arr = fTrackHits->GetArray();;
2490 for (Int_t i=0;i<arr->GetEntriesFast();i++){
2491 AliTrackHitsParamV2 * info = (AliTrackHitsParamV2 *)(arr->At(i));
b9671574 2492 info->SetTrackID(map[info->GetTrackID()]);
792bb11c 2493 }
2494 }
39c8eb58 2495}
2496
792bb11c 2497Bool_t AliTPC::TrackInVolume(Int_t id,Int_t track)
2498{
2499 //return bool information - is track in given volume
2500 //load only part of the track information
2501 //return true if current track is in volume
2502 //
2503 // return kTRUE;
be5ffbfe 2504 // if (fTrackHitsOld && fHitType&2) {
cb5b8b21 2505// TBranch * br = fLoader->TreeH()->GetBranch("fTrackHitsInfo");
be5ffbfe 2506// br->GetEvent(track);
2507// AliObjectArray * ar = fTrackHitsOld->fTrackHitsInfo;
2508// for (UInt_t j=0;j<ar->GetSize();j++){
2509// if ( ((AliTrackHitsInfo*)ar->At(j))->fVolumeID==id) return kTRUE;
2510// }
2511// }
792bb11c 2512
2513 if (fTrackHits && fHitType&4) {
cb5b8b21 2514 TBranch * br1 = fLoader->TreeH()->GetBranch("fVolumes");
2515 TBranch * br2 = fLoader->TreeH()->GetBranch("fNVolumes");
792bb11c 2516 br2->GetEvent(track);
2517 br1->GetEvent(track);
2518 Int_t *volumes = fTrackHits->GetVolumes();
2519 Int_t nvolumes = fTrackHits->GetNVolumes();
2520 if (!volumes && nvolumes>0) {
8c2b3fd7 2521 AliWarning(Form("Problematic track\t%d\t%d",track,nvolumes));
792bb11c 2522 return kFALSE;
2523 }
2524 for (Int_t j=0;j<nvolumes; j++)
2525 if (volumes[j]==id) return kTRUE;
2526 }
2527
2528 if (fHitType&1) {
cb5b8b21 2529 TBranch * br = fLoader->TreeH()->GetBranch("fSector");
792bb11c 2530 br->GetEvent(track);
2531 for (Int_t j=0;j<fHits->GetEntriesFast();j++){
2532 if ( ((AliTPChit*)fHits->At(j))->fSector==id) return kTRUE;
2533 }
2534 }
2535 return kFALSE;
2536
2537}
39c8eb58 2538
39c8eb58 2539
88cb7938 2540AliLoader* AliTPC::MakeLoader(const char* topfoldername)
2541{
8c2b3fd7 2542 //Makes TPC loader
2543 fLoader = new AliTPCLoader(GetName(),topfoldername);
2544 return fLoader;
88cb7938 2545}
2546
42157e55 2547////////////////////////////////////////////////////////////////////////
2548AliTPCParam* AliTPC::LoadTPCParam(TFile *file) {
2549//
2550// load TPC paarmeters from a given file or create new if the object
2551// is not found there
88cb7938 2552// 12/05/2003 This method should be moved to the AliTPCLoader
2553// and one has to decide where to store the TPC parameters
2554// M.Kowalski
42157e55 2555 char paramName[50];
2556 sprintf(paramName,"75x40_100x60_150x60");
2557 AliTPCParam *paramTPC=(AliTPCParam*)file->Get(paramName);
2558 if (paramTPC) {
8c2b3fd7 2559 AliDebugClass(1,Form("TPC parameters %s found.",paramName));
42157e55 2560 } else {
8c2b3fd7 2561 AliWarningClass("TPC parameters not found. Create new (they may be incorrect)");
db2fdcfb 2562 //paramTPC = new AliTPCParamSR;
2563 paramTPC = AliTPCcalibDB::Instance()->GetParameters();
2564 if (!paramTPC->IsGeoRead()){
2565 //
2566 // read transformation matrices for gGeoManager
2567 //
2568 paramTPC->ReadGeoMatrices();
2569 }
2570
42157e55 2571 }
2572 return paramTPC;
2573
2574// the older version of parameters can be accessed with this code.
2575// In some cases, we have old parameters saved in the file but
2576// digits were created with new parameters, it can be distinguish
2577// by the name of TPC TreeD. The code here is just for the case
2578// we would need to compare with old data, uncomment it if needed.
2579//
2580// char paramName[50];
2581// sprintf(paramName,"75x40_100x60");
2582// AliTPCParam *paramTPC=(AliTPCParam*)in->Get(paramName);
2583// if (paramTPC) {
2584// cout<<"TPC parameters "<<paramName<<" found."<<endl;
2585// } else {
2586// sprintf(paramName,"75x40_100x60_150x60");
2587// paramTPC=(AliTPCParam*)in->Get(paramName);
2588// if (paramTPC) {
2589// cout<<"TPC parameters "<<paramName<<" found."<<endl;
2590// } else {
2591// cerr<<"TPC parameters not found. Create new (they may be incorrect)."
2592// <<endl;
2593// paramTPC = new AliTPCParamSR;
2594// }
2595// }
2596// return paramTPC;
2597
2598}
2599
85a5290f 2600