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