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