L1 Gating signal is now used to remove the TPC hits which arrive before the gating...
[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)
de61d5d5 82//_____________________________________________________________________________
fe4da5cc 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
188//_____________________________________________________________________________
fe4da5cc 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//_____________________________________________________________________________
fe4da5cc 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//_____________________________________________________________________________
fe4da5cc 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;
0421c3d1 783
784 fLoader->LoadDigits();
785 TTree* digits = fLoader->TreeD();
786 if (!digits) {
8c2b3fd7 787 AliError("No digits tree");
0421c3d1 788 return;
789 }
790
791 AliSimDigits digarr;
792 AliSimDigits* digrow = &digarr;
793 digits->GetBranch("Segment")->SetAddress(&digrow);
794
795 const char* fileName = "AliTPCDDL.dat";
796 AliTPCBuffer* buffer = new AliTPCBuffer(fileName);
797 //Verbose level
798 // 0: Silent
799 // 1: cout messages
800 // 2: txt files with digits
801 //BE CAREFUL, verbose level 2 MUST be used only for debugging and
802 //it is highly suggested to use this mode only for debugging digits files
803 //reasonably small, because otherwise the size of the txt files can reach
804 //quickly several MB wasting time and disk space.
805 buffer->SetVerbose(0);
806
807 Int_t nEntries = Int_t(digits->GetEntries());
808 Int_t previousSector = -1;
809 Int_t subSector = 0;
810 for (Int_t i = 0; i < nEntries; i++) {
811 digits->GetEntry(i);
812 Int_t sector, row;
813 fTPCParam->AdjustSectorRow(digarr.GetID(), sector, row);
814 if(previousSector != sector) {
815 subSector = 0;
816 previousSector = sector;
817 }
818
819 if (sector < 36) { //inner sector [0;35]
820 if (row != 30) {
821 //the whole row is written into the output file
822 buffer->WriteRowBinary(kThreshold, digrow, 0, 0, 0,
823 sector, subSector, row);
824 } else {
825 //only the pads in the range [37;48] are written into the output file
826 buffer->WriteRowBinary(kThreshold, digrow, 37, 48, 1,
827 sector, subSector, row);
828 subSector = 1;
829 //only the pads outside the range [37;48] are written into the output file
830 buffer->WriteRowBinary(kThreshold, digrow, 37, 48, 2,
831 sector, subSector, row);
832 }//end else
833
834 } else { //outer sector [36;71]
835 if (row == 54) subSector = 2;
836 if ((row != 27) && (row != 76)) {
837 buffer->WriteRowBinary(kThreshold, digrow, 0, 0, 0,
838 sector, subSector, row);
839 } else if (row == 27) {
8c2b3fd7 840 //only the pads outside the range [43;46] are written into the output file
841 buffer->WriteRowBinary(kThreshold, digrow, 43, 46, 2,
0421c3d1 842 sector, subSector, row);
8c2b3fd7 843 subSector = 1;
844 //only the pads in the range [43;46] are written into the output file
845 buffer->WriteRowBinary(kThreshold, digrow, 43, 46, 1,
0421c3d1 846 sector, subSector, row);
847 } else if (row == 76) {
8c2b3fd7 848 //only the pads outside the range [33;88] are written into the output file
849 buffer->WriteRowBinary(kThreshold, digrow, 33, 88, 2,
850 sector, subSector, row);
851 subSector = 3;
852 //only the pads in the range [33;88] are written into the output file
853 buffer->WriteRowBinary(kThreshold, digrow, 33, 88, 1,
0421c3d1 854 sector, subSector, row);
855 }
856 }//end else
857 }//end for
858
859 delete buffer;
860 fLoader->UnloadDigits();
861
862 AliTPCDDLRawData rawWriter;
863 rawWriter.SetVerbose(0);
864
865 rawWriter.RawData(fileName);
866 gSystem->Unlink(fileName);
867
0421c3d1 868}
869
870
0a61bf9d 871
85a5290f 872//______________________________________________________________________
9bdd974b 873AliDigitizer* AliTPC::CreateDigitizer(AliRunDigitizer* manager) const
85a5290f 874{
875 return new AliTPCDigitizer(manager);
876}
0a61bf9d 877//__
176aff27 878void AliTPC::SDigits2Digits2(Int_t /*eventnumber*/)
0a61bf9d 879{
880 //create digits from summable digits
407ff276 881 GenerNoise(500000); //create teble with noise
0a61bf9d 882
883 //conect tree with sSDigits
88cb7938 884 TTree *t = fLoader->TreeS();
885
8c2b3fd7 886 if (t == 0x0) {
887 fLoader->LoadSDigits("READ");
888 t = fLoader->TreeS();
889 if (t == 0x0) {
890 AliError("Can not get input TreeS");
891 return;
892 }
893 }
88cb7938 894
895 if (fLoader->TreeD() == 0x0) fLoader->MakeTree("D");
896
0a61bf9d 897 AliSimDigits digarr, *dummy=&digarr;
88cb7938 898 TBranch* sdb = t->GetBranch("Segment");
8c2b3fd7 899 if (sdb == 0x0) {
900 AliError("Can not find branch with segments in TreeS.");
901 return;
902 }
88cb7938 903
904 sdb->SetAddress(&dummy);
905
0a61bf9d 906 Stat_t nentries = t->GetEntries();
907
5f16d237 908 // set zero suppression
909
910 fTPCParam->SetZeroSup(2);
911
912 // get zero suppression
913
914 Int_t zerosup = fTPCParam->GetZeroSup();
915
916 //make tree with digits
917
0a61bf9d 918 AliTPCDigitsArray *arr = new AliTPCDigitsArray;
919 arr->SetClass("AliSimDigits");
920 arr->Setup(fTPCParam);
88cb7938 921 arr->MakeTree(fLoader->TreeD());
0a61bf9d 922
88cb7938 923 AliTPCParam * par = fTPCParam;
5f16d237 924
0a61bf9d 925 //Loop over segments of the TPC
5f16d237 926
0a61bf9d 927 for (Int_t n=0; n<nentries; n++) {
928 t->GetEvent(n);
929 Int_t sec, row;
930 if (!par->AdjustSectorRow(digarr.GetID(),sec,row)) {
8c2b3fd7 931 AliWarning(Form("Invalid segment ID ! %d",digarr.GetID()));
0a61bf9d 932 continue;
933 }
8c2b3fd7 934 if (!IsSectorActive(sec)) continue;
935
0a61bf9d 936 AliSimDigits * digrow =(AliSimDigits*) arr->CreateRow(sec,row);
937 Int_t nrows = digrow->GetNRows();
938 Int_t ncols = digrow->GetNCols();
939
940 digrow->ExpandBuffer();
941 digarr.ExpandBuffer();
942 digrow->ExpandTrackBuffer();
943 digarr.ExpandTrackBuffer();
944
5f16d237 945
407ff276 946 Short_t * pamp0 = digarr.GetDigits();
947 Int_t * ptracks0 = digarr.GetTracks();
948 Short_t * pamp1 = digrow->GetDigits();
949 Int_t * ptracks1 = digrow->GetTracks();
950 Int_t nelems =nrows*ncols;
951 Int_t saturation = fTPCParam->GetADCSat();
952 //use internal structure of the AliDigits - for speed reason
953 //if you cahnge implementation
954 //of the Alidigits - it must be rewriten -
955 for (Int_t i= 0; i<nelems; i++){
407ff276 956 Float_t q = TMath::Nint(Float_t(*pamp0)/16.+GetNoise());
957 if (q>zerosup){
958 if (q>saturation) q=saturation;
959 *pamp1=(Short_t)q;
8c2b3fd7 960
407ff276 961 ptracks1[0]=ptracks0[0];
962 ptracks1[nelems]=ptracks0[nelems];
963 ptracks1[2*nelems]=ptracks0[2*nelems];
964 }
965 pamp0++;
966 pamp1++;
967 ptracks0++;
968 ptracks1++;
969 }
5f16d237 970
5f16d237 971 arr->StoreRow(sec,row);
972 arr->ClearRow(sec,row);
5f16d237 973 }
0a61bf9d 974
0a61bf9d 975
5f16d237 976 //write results
88cb7938 977 fLoader->WriteDigits("OVERWRITE");
5f16d237 978
88cb7938 979 delete arr;
afc42102 980}
981//__________________________________________________________________
982void AliTPC::SetDefaults(){
9bdd974b 983 //
984 // setting the defaults
985 //
afc42102 986
afc42102 987 // Set response functions
988
88cb7938 989 //
024a7e64 990 AliRunLoader* rl = (AliRunLoader*)fLoader->GetEventFolder()->FindObject(AliRunLoader::GetRunLoaderName());
88cb7938 991 rl->CdGAFile();
2ab0c725 992 AliTPCParamSR *param=(AliTPCParamSR*)gDirectory->Get("75x40_100x60");
7a09f434 993 if(param){
8c2b3fd7 994 AliInfo("You are using 2 pad-length geom hits with 3 pad-lenght geom digits...");
7a09f434 995 delete param;
996 param = new AliTPCParamSR();
997 }
998 else {
999 param=(AliTPCParamSR*)gDirectory->Get("75x40_100x60_150x60");
1000 }
1001 if(!param){
8c2b3fd7 1002 AliFatal("No TPC parameters found");
7a09f434 1003 }
1004
1005
2ab0c725 1006 AliTPCPRF2D * prfinner = new AliTPCPRF2D;
f03e3423 1007 AliTPCPRF2D * prfouter1 = new AliTPCPRF2D;
1008 AliTPCPRF2D * prfouter2 = new AliTPCPRF2D;
2ab0c725 1009 AliTPCRF1D * rf = new AliTPCRF1D(kTRUE);
2ab0c725 1010 rf->SetGauss(param->GetZSigma(),param->GetZWidth(),1.);
1011 rf->SetOffset(3*param->GetZSigma());
1012 rf->Update();
afc42102 1013
1014 TDirectory *savedir=gDirectory;
2ab0c725 1015 TFile *f=TFile::Open("$ALICE_ROOT/TPC/AliTPCprf2d.root");
8c2b3fd7 1016 if (!f->IsOpen())
1017 AliFatal("Can't open $ALICE_ROOT/TPC/AliTPCprf2d.root !");
2a336e15 1018
1019 TString s;
2ab0c725 1020 prfinner->Read("prf_07504_Gati_056068_d02");
2a336e15 1021 //PH Set different names
1022 s=prfinner->GetGRF()->GetName();
1023 s+="in";
1024 prfinner->GetGRF()->SetName(s.Data());
1025
f03e3423 1026 prfouter1->Read("prf_10006_Gati_047051_d03");
2a336e15 1027 s=prfouter1->GetGRF()->GetName();
1028 s+="out1";
1029 prfouter1->GetGRF()->SetName(s.Data());
1030
f03e3423 1031 prfouter2->Read("prf_15006_Gati_047051_d03");
2a336e15 1032 s=prfouter2->GetGRF()->GetName();
1033 s+="out2";
1034 prfouter2->GetGRF()->SetName(s.Data());
1035
2ab0c725 1036 f->Close();
afc42102 1037 savedir->cd();
1038
2ab0c725 1039 param->SetInnerPRF(prfinner);
f03e3423 1040 param->SetOuter1PRF(prfouter1);
1041 param->SetOuter2PRF(prfouter2);
2ab0c725 1042 param->SetTimeRF(rf);
1043
afc42102 1044 // set fTPCParam
1045
2ab0c725 1046 SetParam(param);
1047
afc42102 1048
1049 fDefaults = 1;
1050
1051}
1052//__________________________________________________________________
85a5290f 1053void AliTPC::Hits2Digits()
1054{
9bdd974b 1055 //
1056 // creates digits from hits
1057 //
506e60a6 1058 if (!fTPCParam->IsGeoRead()){
1059 //
1060 // read transformation matrices for gGeoManager
1061 //
1062 fTPCParam->ReadGeoMatrices();
1063 }
9bdd974b 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
506e60a6 1209 if (!fTPCParam->IsGeoRead()){
1210 //
1211 // read transformation matrices for gGeoManager
1212 //
1213 fTPCParam->ReadGeoMatrices();
1214 }
1215
85a5290f 1216 fLoader->LoadHits("read");
1217 fLoader->LoadSDigits("recreate");
1218 AliRunLoader* runLoader = fLoader->GetRunLoader();
0a61bf9d 1219
85a5290f 1220 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1221 runLoader->GetEvent(iEvent);
1222 SetTreeAddress();
1223 SetActiveSectors();
1224 Hits2SDigits2(iEvent);
1225 }
0a61bf9d 1226
85a5290f 1227 fLoader->UnloadHits();
1228 fLoader->UnloadSDigits();
0a61bf9d 1229}
fe4da5cc 1230//_____________________________________________________________________________
88cb7938 1231
8c555625 1232void AliTPC::Hits2DigitsSector(Int_t isec)
fe4da5cc 1233{
8c555625 1234 //-------------------------------------------------------------------
fe4da5cc 1235 // TPC conversion from hits to digits.
8c555625 1236 //-------------------------------------------------------------------
1237
1238 //-----------------------------------------------------------------
1239 // Origin: Marek Kowalski IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
1240 //-----------------------------------------------------------------
1241
fe4da5cc 1242 //-------------------------------------------------------
8c555625 1243 // Get the access to the track hits
fe4da5cc 1244 //-------------------------------------------------------
8c555625 1245
afc42102 1246 // check if the parameters are set - important if one calls this method
1247 // directly, not from the Hits2Digits
1248
1249 if(fDefaults == 0) SetDefaults();
cc80f89e 1250
88cb7938 1251 TTree *tH = TreeH(); // pointer to the hits tree
8c2b3fd7 1252 if (tH == 0x0) {
1253 AliFatal("Can not find TreeH in folder");
1254 return;
1255 }
88cb7938 1256
73042f01 1257 Stat_t ntracks = tH->GetEntries();
8c555625 1258
1259 if( ntracks > 0){
1260
1261 //-------------------------------------------
1262 // Only if there are any tracks...
1263 //-------------------------------------------
1264
8c555625 1265 TObjArray **row;
fe4da5cc 1266
8c2b3fd7 1267 Int_t nrows =fTPCParam->GetNRow(isec);
8c555625 1268
8c2b3fd7 1269 row= new TObjArray* [nrows+2]; // 2 extra rows for cross talk
fe4da5cc 1270
8c2b3fd7 1271 MakeSector(isec,nrows,tH,ntracks,row);
8c555625 1272
8c2b3fd7 1273 //--------------------------------------------------------
1274 // Digitize this sector, row by row
1275 // row[i] is the pointer to the TObjArray of TVectors,
1276 // each one containing electrons accepted on this
1277 // row, assigned into tracks
1278 //--------------------------------------------------------
8c555625 1279
8c2b3fd7 1280 Int_t i;
8c555625 1281
8c2b3fd7 1282 if (fDigitsArray->GetTree()==0) {
1283 AliFatal("Tree not set in fDigitsArray");
1284 }
8c555625 1285
8c2b3fd7 1286 for (i=0;i<nrows;i++){
1287
1288 AliDigits * dig = fDigitsArray->CreateRow(isec,i);
8c555625 1289
8c2b3fd7 1290 DigitizeRow(i,isec,row);
8c555625 1291
8c2b3fd7 1292 fDigitsArray->StoreRow(isec,i);
8c555625 1293
8c2b3fd7 1294 Int_t ndig = dig->GetDigitSize();
88cb7938 1295
8c2b3fd7 1296 AliDebug(10,
1297 Form("*** Sector, row, compressed digits %d %d %d ***\n",
1298 isec,i,ndig));
792bb11c 1299
8c2b3fd7 1300 fDigitsArray->ClearRow(isec,i);
8c555625 1301
cc80f89e 1302
8c2b3fd7 1303 } // end of the sector digitization
8c555625 1304
8c2b3fd7 1305 for(i=0;i<nrows+2;i++){
1306 row[i]->Delete();
1307 delete row[i];
1308 }
cc80f89e 1309
8c2b3fd7 1310 delete [] row; // delete the array of pointers to TObjArray-s
8c555625 1311
1312 } // ntracks >0
8c555625 1313
cc80f89e 1314} // end of Hits2DigitsSector
8c555625 1315
8c555625 1316
8c555625 1317//_____________________________________________________________________________
cc80f89e 1318void AliTPC::DigitizeRow(Int_t irow,Int_t isec,TObjArray **rows)
8c555625 1319{
1320 //-----------------------------------------------------------
1321 // Single row digitization, coupling from the neighbouring
1322 // rows taken into account
1323 //-----------------------------------------------------------
1324
1325 //-----------------------------------------------------------------
1326 // Origin: Marek Kowalski IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
cc80f89e 1327 // Modified: Marian Ivanov GSI Darmstadt, m.ivanov@gsi.de
8c555625 1328 //-----------------------------------------------------------------
1329
8c555625 1330 Float_t zerosup = fTPCParam->GetZeroSup();
8c2b3fd7 1331
cc80f89e 1332 fCurrentIndex[1]= isec;
8c555625 1333
8c555625 1334
73042f01 1335 Int_t nofPads = fTPCParam->GetNPads(isec,irow);
1336 Int_t nofTbins = fTPCParam->GetMaxTBin();
1337 Int_t indexRange[4];
8c555625 1338 //
1339 // Integrated signal for this row
1340 // and a single track signal
cc80f89e 1341 //
de61d5d5 1342
e8d02863 1343 TMatrixF *m1 = new TMatrixF(0,nofPads,0,nofTbins); // integrated
1344 TMatrixF *m2 = new TMatrixF(0,nofPads,0,nofTbins); // single
8c555625 1345 //
e8d02863 1346 TMatrixF &total = *m1;
8c555625 1347
1348 // Array of pointers to the label-signal list
1349
73042f01 1350 Int_t nofDigits = nofPads*nofTbins; // number of digits for this row
1351 Float_t **pList = new Float_t* [nofDigits];
8c555625 1352
1353 Int_t lp;
cc80f89e 1354 Int_t i1;
73042f01 1355 for(lp=0;lp<nofDigits;lp++)pList[lp]=0; // set all pointers to NULL
8c555625 1356 //
cc80f89e 1357 //calculate signal
1358 //
b584c7dd 1359 Int_t row1=irow;
1360 Int_t row2=irow+2;
cc80f89e 1361 for (Int_t row= row1;row<=row2;row++){
1362 Int_t nTracks= rows[row]->GetEntries();
1363 for (i1=0;i1<nTracks;i1++){
1364 fCurrentIndex[2]= row;
b584c7dd 1365 fCurrentIndex[3]=irow+1;
1366 if (row==irow+1){
cc80f89e 1367 m2->Zero(); // clear single track signal matrix
73042f01 1368 Float_t trackLabel = GetSignal(rows[row],i1,m2,m1,indexRange);
1369 GetList(trackLabel,nofPads,m2,indexRange,pList);
cc80f89e 1370 }
73042f01 1371 else GetSignal(rows[row],i1,0,m1,indexRange);
cc80f89e 1372 }
8c555625 1373 }
cc80f89e 1374
8c555625 1375 Int_t tracks[3];
8c555625 1376
cc80f89e 1377 AliDigits *dig = fDigitsArray->GetRow(isec,irow);
de61d5d5 1378 Int_t gi=-1;
1379 Float_t fzerosup = zerosup+0.5;
1380 for(Int_t it=0;it<nofTbins;it++){
de61d5d5 1381 for(Int_t ip=0;ip<nofPads;ip++){
1382 gi++;
8c2b3fd7 1383 Float_t q=total(ip,it);
f8cf550c 1384 if(fDigitsSwitch == 0){
407ff276 1385 q+=GetNoise();
de61d5d5 1386 if(q <=fzerosup) continue; // do not fill zeros
68771f83 1387 q = TMath::Nint(q);
f8cf550c 1388 if(q > fTPCParam->GetADCSat()) q = fTPCParam->GetADCSat(); // saturation
cc80f89e 1389
f8cf550c 1390 }
1391
1392 else {
8c2b3fd7 1393 if(q <= 0.) continue; // do not fill zeros
1394 if(q>2000.) q=2000.;
1395 q *= 16.;
1396 q = TMath::Nint(q);
f8cf550c 1397 }
8c555625 1398
1399 //
1400 // "real" signal or electronic noise (list = -1)?
1401 //
1402
1403 for(Int_t j1=0;j1<3;j1++){
407ff276 1404 tracks[j1] = (pList[gi]) ?(Int_t)(*(pList[gi]+j1)) : -2;
8c555625 1405 }
1406
cc80f89e 1407//Begin_Html
1408/*
1409 <A NAME="AliDigits"></A>
1410 using of AliDigits object
1411*/
1412//End_Html
1413 dig->SetDigitFast((Short_t)q,it,ip);
8c2b3fd7 1414 if (fDigitsArray->IsSimulated()) {
1415 ((AliSimDigits*)dig)->SetTrackIDFast(tracks[0],it,ip,0);
1416 ((AliSimDigits*)dig)->SetTrackIDFast(tracks[1],it,ip,1);
1417 ((AliSimDigits*)dig)->SetTrackIDFast(tracks[2],it,ip,2);
1418 }
fe4da5cc 1419
8c555625 1420 } // end of loop over time buckets
1421 } // end of lop over pads
1422
1423 //
1424 // This row has been digitized, delete nonused stuff
1425 //
1426
73042f01 1427 for(lp=0;lp<nofDigits;lp++){
8c555625 1428 if(pList[lp]) delete [] pList[lp];
1429 }
1430
1431 delete [] pList;
1432
1433 delete m1;
1434 delete m2;
8c555625 1435
1436} // end of DigitizeRow
cc80f89e 1437
8c555625 1438//_____________________________________________________________________________
cc80f89e 1439
de61d5d5 1440Float_t AliTPC::GetSignal(TObjArray *p1, Int_t ntr,
e8d02863 1441 TMatrixF *m1, TMatrixF *m2,Int_t *indexRange)
8c555625 1442{
1443
1444 //---------------------------------------------------------------
1445 // Calculates 2-D signal (pad,time) for a single track,
1446 // returns a pointer to the signal matrix and the track label
1447 // No digitization is performed at this level!!!
1448 //---------------------------------------------------------------
1449
1450 //-----------------------------------------------------------------
1451 // Origin: Marek Kowalski IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
cc80f89e 1452 // Modified: Marian Ivanov
8c555625 1453 //-----------------------------------------------------------------
1454
8c2b3fd7 1455 TVector *tv;
de61d5d5 1456
8c2b3fd7 1457 tv = (TVector*)p1->At(ntr); // pointer to a track
1458 TVector &v = *tv;
8c555625 1459
1460 Float_t label = v(0);
506e60a6 1461 Int_t centralPad = (fTPCParam->GetNPads(fCurrentIndex[1],fCurrentIndex[3]-1))/2;
8c555625 1462
e61fd20d 1463 Int_t nElectrons = (tv->GetNrows()-1)/5;
73042f01 1464 indexRange[0]=9999; // min pad
1465 indexRange[1]=-1; // max pad
1466 indexRange[2]=9999; //min time
1467 indexRange[3]=-1; // max time
8c555625 1468
e8d02863 1469 TMatrixF &signal = *m1;
1470 TMatrixF &total = *m2;
8c555625 1471 //
1472 // Loop over all electrons
1473 //
8c555625 1474 for(Int_t nel=0; nel<nElectrons; nel++){
e61fd20d 1475 Int_t idx=nel*5;
cc80f89e 1476 Float_t aval = v(idx+4);
1477 Float_t eltoadcfac=aval*fTPCParam->GetTotalNormFac();
e61fd20d 1478 Float_t xyz[4]={v(idx+1),v(idx+2),v(idx+3),v(idx+5)};
de61d5d5 1479 Int_t n = ((AliTPCParamSR*)fTPCParam)->CalcResponseFast(xyz,fCurrentIndex,fCurrentIndex[3]);
1480
1481 Int_t *index = fTPCParam->GetResBin(0);
1482 Float_t *weight = & (fTPCParam->GetResWeight(0));
1483
1484 if (n>0) for (Int_t i =0; i<n; i++){
8c2b3fd7 1485 Int_t pad=index[1]+centralPad; //in digit coordinates central pad has coordinate 0
1486
1487 if (pad>=0){
1488 Int_t time=index[2];
1489 Float_t qweight = *(weight)*eltoadcfac;
1490
1491 if (m1!=0) signal(pad,time)+=qweight;
1492 total(pad,time)+=qweight;
1493 if (indexRange[0]>pad) indexRange[0]=pad;
1494 if (indexRange[1]<pad) indexRange[1]=pad;
1495 if (indexRange[2]>time) indexRange[2]=time;
1496 if (indexRange[3]<time) indexRange[3]=time;
1497
1498 index+=3;
1499 weight++;
1500
1501 }
cc80f89e 1502 }
8c555625 1503 } // end of loop over electrons
cc80f89e 1504
8c555625 1505 return label; // returns track label when finished
1506}
1507
1508//_____________________________________________________________________________
e8d02863 1509void AliTPC::GetList(Float_t label,Int_t np,TMatrixF *m,
de61d5d5 1510 Int_t *indexRange, Float_t **pList)
8c555625 1511{
1512 //----------------------------------------------------------------------
1513 // Updates the list of tracks contributing to digits for a given row
1514 //----------------------------------------------------------------------
1515
1516 //-----------------------------------------------------------------
1517 // Origin: Marek Kowalski IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
1518 //-----------------------------------------------------------------
1519
e8d02863 1520 TMatrixF &signal = *m;
8c555625 1521
1522 // lop over nonzero digits
1523
73042f01 1524 for(Int_t it=indexRange[2];it<indexRange[3]+1;it++){
1525 for(Int_t ip=indexRange[0];ip<indexRange[1]+1;ip++){
8c555625 1526
1527
8c2b3fd7 1528 // accept only the contribution larger than 500 electrons (1/2 s_noise)
921bf71a 1529
8c2b3fd7 1530 if(signal(ip,it)<0.5) continue;
921bf71a 1531
8c2b3fd7 1532 Int_t globalIndex = it*np+ip; // globalIndex starts from 0!
8c555625 1533
8c2b3fd7 1534 if(!pList[globalIndex]){
8c555625 1535
8c2b3fd7 1536 //
1537 // Create new list (6 elements - 3 signals and 3 labels),
1538 //
8c555625 1539
8c2b3fd7 1540 pList[globalIndex] = new Float_t [6];
8c555625 1541
8c2b3fd7 1542 // set list to -1
1543
1544 *pList[globalIndex] = -1.;
1545 *(pList[globalIndex]+1) = -1.;
1546 *(pList[globalIndex]+2) = -1.;
1547 *(pList[globalIndex]+3) = -1.;
1548 *(pList[globalIndex]+4) = -1.;
1549 *(pList[globalIndex]+5) = -1.;
1550
1551 *pList[globalIndex] = label;
1552 *(pList[globalIndex]+3) = signal(ip,it);
1553 }
1554 else {
8c555625 1555
8c2b3fd7 1556 // check the signal magnitude
8c555625 1557
8c2b3fd7 1558 Float_t highest = *(pList[globalIndex]+3);
1559 Float_t middle = *(pList[globalIndex]+4);
1560 Float_t lowest = *(pList[globalIndex]+5);
1561
1562 //
1563 // compare the new signal with already existing list
1564 //
1565
1566 if(signal(ip,it)<lowest) continue; // neglect this track
8c555625 1567
8c2b3fd7 1568 //
8c555625 1569
8c2b3fd7 1570 if (signal(ip,it)>highest){
1571 *(pList[globalIndex]+5) = middle;
1572 *(pList[globalIndex]+4) = highest;
1573 *(pList[globalIndex]+3) = signal(ip,it);
1574
1575 *(pList[globalIndex]+2) = *(pList[globalIndex]+1);
1576 *(pList[globalIndex]+1) = *pList[globalIndex];
1577 *pList[globalIndex] = label;
1578 }
1579 else if (signal(ip,it)>middle){
1580 *(pList[globalIndex]+5) = middle;
1581 *(pList[globalIndex]+4) = signal(ip,it);
1582
1583 *(pList[globalIndex]+2) = *(pList[globalIndex]+1);
1584 *(pList[globalIndex]+1) = label;
1585 }
1586 else{
1587 *(pList[globalIndex]+5) = signal(ip,it);
1588 *(pList[globalIndex]+2) = label;
1589 }
1590 }
1591
8c555625 1592 } // end of loop over pads
1593 } // end of loop over time bins
1594
8c555625 1595}//end of GetList
1596//___________________________________________________________________
1597void AliTPC::MakeSector(Int_t isec,Int_t nrows,TTree *TH,
1598 Stat_t ntracks,TObjArray **row)
1599{
1600
1601 //-----------------------------------------------------------------
1602 // Prepares the sector digitization, creates the vectors of
1603 // tracks for each row of this sector. The track vector
1604 // contains the track label and the position of electrons.
1605 //-----------------------------------------------------------------
1606
1607 //-----------------------------------------------------------------
1608 // Origin: Marek Kowalski IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
1609 //-----------------------------------------------------------------
1610
cc80f89e 1611 Float_t gasgain = fTPCParam->GetGasGain();
8c555625 1612 Int_t i;
e61fd20d 1613 Float_t xyz[5];
8c555625 1614
1615 AliTPChit *tpcHit; // pointer to a sigle TPC hit
39c8eb58 1616 //MI change
1617 TBranch * branch=0;
792bb11c 1618 if (fHitType>1) branch = TH->GetBranch("TPC2");
39c8eb58 1619 else branch = TH->GetBranch("TPC");
1620
8c555625 1621
1622 //----------------------------------------------
1623 // Create TObjArray-s, one for each row,
8c2b3fd7 1624 // each TObjArray will store the TVectors
1625 // of electrons, one TVectors per each track.
8c555625 1626 //----------------------------------------------
fe4da5cc 1627
b584c7dd 1628 Int_t *nofElectrons = new Int_t [nrows+2]; // electron counter for each row
8c2b3fd7 1629 TVector **tracks = new TVector* [nrows+2]; //pointers to the track vectors
de61d5d5 1630
b584c7dd 1631 for(i=0; i<nrows+2; i++){
8c555625 1632 row[i] = new TObjArray;
f74bb6f5 1633 nofElectrons[i]=0;
1634 tracks[i]=0;
8c555625 1635 }
8c555625 1636
37831078 1637
1638
8c555625 1639 //--------------------------------------------------------------------
1640 // Loop over tracks, the "track" contains the full history
1641 //--------------------------------------------------------------------
8c2b3fd7 1642
8c555625 1643 Int_t previousTrack,currentTrack;
1644 previousTrack = -1; // nothing to store so far!
1645
1646 for(Int_t track=0;track<ntracks;track++){
39c8eb58 1647 Bool_t isInSector=kTRUE;
8c555625 1648 ResetHits();
792bb11c 1649 isInSector = TrackInVolume(isec,track);
39c8eb58 1650 if (!isInSector) continue;
1651 //MI change
1652 branch->GetEntry(track); // get next track
8c2b3fd7 1653
39c8eb58 1654 //M.I. changes
8c555625 1655
39c8eb58 1656 tpcHit = (AliTPChit*)FirstHit(-1);
8c555625 1657
1658 //--------------------------------------------------------------
1659 // Loop over hits
1660 //--------------------------------------------------------------
1661
8c555625 1662
39c8eb58 1663 while(tpcHit){
fe4da5cc 1664
8c555625 1665 Int_t sector=tpcHit->fSector; // sector number
39c8eb58 1666 if(sector != isec){
1667 tpcHit = (AliTPChit*) NextHit();
1668 continue;
1669 }
8c555625 1670
daf14717 1671 // Remove hits which arrive before the TPC opening gate signal
1672 if(((fTPCParam->GetZLength()-TMath::Abs(tpcHit->Z()))
1673 /fTPCParam->GetDriftV()+tpcHit->Time())<fTPCParam->GetGateDelay()) {
1674 tpcHit = (AliTPChit*) NextHit();
1675 continue;
1676 }
1677
8c2b3fd7 1678 currentTrack = tpcHit->Track(); // track number
39c8eb58 1679
8c2b3fd7 1680 if(currentTrack != previousTrack){
8c555625 1681
8c2b3fd7 1682 // store already filled fTrack
8c555625 1683
8c2b3fd7 1684 for(i=0;i<nrows+2;i++){
1685 if(previousTrack != -1){
1686 if(nofElectrons[i]>0){
1687 TVector &v = *tracks[i];
1688 v(0) = previousTrack;
e61fd20d 1689 tracks[i]->ResizeTo(5*nofElectrons[i]+1); // shrink if necessary
8c2b3fd7 1690 row[i]->Add(tracks[i]);
1691 }
1692 else {
1693 delete tracks[i]; // delete empty TVector
1694 tracks[i]=0;
1695 }
1696 }
1697
1698 nofElectrons[i]=0;
e61fd20d 1699 tracks[i] = new TVector(601); // TVectors for the next fTrack
8c2b3fd7 1700
1701 } // end of loop over rows
8c555625 1702
8c2b3fd7 1703 previousTrack=currentTrack; // update track label
1704 }
8c555625 1705
8c2b3fd7 1706 Int_t qI = (Int_t) (tpcHit->fQ); // energy loss (number of electrons)
8c555625 1707
8c2b3fd7 1708 //---------------------------------------------------
1709 // Calculate the electron attachment probability
1710 //---------------------------------------------------
8c555625 1711
cc80f89e 1712
8c2b3fd7 1713 Float_t time = 1.e6*(fTPCParam->GetZLength()-TMath::Abs(tpcHit->Z()))
1714 /fTPCParam->GetDriftV();
1715 // in microseconds!
1716 Float_t attProb = fTPCParam->GetAttCoef()*
1717 fTPCParam->GetOxyCont()*time; // fraction!
8c555625 1718
8c2b3fd7 1719 //-----------------------------------------------
1720 // Loop over electrons
1721 //-----------------------------------------------
1722 Int_t index[3];
1723 index[1]=isec;
1724 for(Int_t nel=0;nel<qI;nel++){
1725 // skip if electron lost due to the attachment
1726 if((gRandom->Rndm(0)) < attProb) continue; // electron lost!
1727 xyz[0]=tpcHit->X();
1728 xyz[1]=tpcHit->Y();
1729 xyz[2]=tpcHit->Z();
1730 //
1731 // protection for the nonphysical avalanche size (10**6 maximum)
1732 //
1733 Double_t rn=TMath::Max(gRandom->Rndm(0),1.93e-22);
1734 xyz[3]= (Float_t) (-gasgain*TMath::Log(rn));
1735 index[0]=1;
1736
1737 TransportElectron(xyz,index);
1738 Int_t rowNumber;
1739 fTPCParam->GetPadRow(xyz,index);
e61fd20d 1740 // Electron track time (for pileup simulation)
1741 xyz[4] = tpcHit->Time()/fTPCParam->GetTSample();
8c2b3fd7 1742 // row 0 - cross talk from the innermost row
1743 // row fNRow+1 cross talk from the outermost row
1744 rowNumber = index[2]+1;
1745 //transform position to local digit coordinates
1746 //relative to nearest pad row
1747 if ((rowNumber<0)||rowNumber>fTPCParam->GetNRow(isec)+1) continue;
1748 Float_t x1,y1;
1749 if (isec <fTPCParam->GetNInnerSector()) {
1750 x1 = xyz[1]*fTPCParam->GetInnerPadPitchWidth();
1751 y1 = fTPCParam->GetYInner(rowNumber);
1752 }
1753 else{
1754 x1=xyz[1]*fTPCParam->GetOuterPadPitchWidth();
1755 y1 = fTPCParam->GetYOuter(rowNumber);
1756 }
1757 // gain inefficiency at the wires edges - linear
1758 x1=TMath::Abs(x1);
1759 y1-=1.;
1760 if(x1>y1) xyz[3]*=TMath::Max(1.e-6,(y1-x1+1.));
1761
1762 nofElectrons[rowNumber]++;
1763 //----------------------------------
1764 // Expand vector if necessary
1765 //----------------------------------
1766 if(nofElectrons[rowNumber]>120){
1767 Int_t range = tracks[rowNumber]->GetNrows();
e61fd20d 1768 if((nofElectrons[rowNumber])>(range-1)/5){
8c2b3fd7 1769
e61fd20d 1770 tracks[rowNumber]->ResizeTo(range+500); // Add 100 electrons
fe4da5cc 1771 }
8c2b3fd7 1772 }
1773
1774 TVector &v = *tracks[rowNumber];
e61fd20d 1775 Int_t idx = 5*nofElectrons[rowNumber]-4;
8c2b3fd7 1776 Real_t * position = &(((TVector&)v)(idx)); //make code faster
e61fd20d 1777 memcpy(position,xyz,5*sizeof(Float_t));
8c2b3fd7 1778
1779 } // end of loop over electrons
39c8eb58 1780
8c2b3fd7 1781 tpcHit = (AliTPChit*)NextHit();
1782
1783 } // end of loop over hits
1784 } // end of loop over tracks
8c555625 1785
1786 //
1787 // store remaining track (the last one) if not empty
1788 //
8c2b3fd7 1789
1790 for(i=0;i<nrows+2;i++){
1791 if(nofElectrons[i]>0){
1792 TVector &v = *tracks[i];
1793 v(0) = previousTrack;
e61fd20d 1794 tracks[i]->ResizeTo(5*nofElectrons[i]+1); // shrink if necessary
8c2b3fd7 1795 row[i]->Add(tracks[i]);
1796 }
1797 else{
1798 delete tracks[i];
1799 tracks[i]=0;
1800 }
1801 }
1802
1803 delete [] tracks;
1804 delete [] nofElectrons;
8c555625 1805
cc80f89e 1806} // end of MakeSector
8c555625 1807
fe4da5cc 1808
1809//_____________________________________________________________________________
1810void AliTPC::Init()
1811{
1812 //
1813 // Initialise TPC detector after definition of geometry
1814 //
8c2b3fd7 1815 AliDebug(1,"*********************************************");
fe4da5cc 1816}
1817
1818//_____________________________________________________________________________
88cb7938 1819void AliTPC::MakeBranch(Option_t* option)
fe4da5cc 1820{
1821 //
1822 // Create Tree branches for the TPC.
1823 //
8c2b3fd7 1824 AliDebug(1,"");
fe4da5cc 1825 Int_t buffersize = 4000;
1826 char branchname[10];
1827 sprintf(branchname,"%s",GetName());
88cb7938 1828
1829 const char *h = strstr(option,"H");
8c2b3fd7 1830
88cb7938 1831 if ( h && (fHitType<=1) && (fHits == 0x0)) fHits = new TClonesArray("AliTPChit", 176);//skowron 20.06.03
1832
1833 AliDetector::MakeBranch(option);
8c2b3fd7 1834
5cf7bbad 1835 const char *d = strstr(option,"D");
8c2b3fd7 1836
1837 if (fDigits && fLoader->TreeD() && d) {
1838 MakeBranchInTree(gAlice->TreeD(), branchname, &fDigits, buffersize, 0);
1839 }
fe4da5cc 1840
88cb7938 1841 if (fHitType>1) MakeBranch2(option,0); // MI change 14.09.2000
fe4da5cc 1842}
1843
1844//_____________________________________________________________________________
1845void AliTPC::ResetDigits()
1846{
1847 //
1848 // Reset number of digits and the digits array for this detector
fe4da5cc 1849 //
1850 fNdigits = 0;
cc80f89e 1851 if (fDigits) fDigits->Clear();
fe4da5cc 1852}
1853
fe4da5cc 1854
fe4da5cc 1855
1856//_____________________________________________________________________________
1857void AliTPC::SetSens(Int_t sens)
1858{
8c555625 1859
1860 //-------------------------------------------------------------
1861 // Activates/deactivates the sensitive strips at the center of
1862 // the pad row -- this is for the space-point resolution calculations
1863 //-------------------------------------------------------------
1864
1865 //-----------------------------------------------------------------
1866 // Origin: Marek Kowalski IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
1867 //-----------------------------------------------------------------
1868
fe4da5cc 1869 fSens = sens;
1870}
2b06d5c3 1871
4b0fdcad 1872
73042f01 1873void AliTPC::SetSide(Float_t side=0.)
4b0fdcad 1874{
73042f01 1875 // choice of the TPC side
1876
4b0fdcad 1877 fSide = side;
1878
1879}
fe4da5cc 1880//_____________________________________________________________________________
cc80f89e 1881
1882void AliTPC::TransportElectron(Float_t *xyz, Int_t *index)
1883{
1884 //
1885 // electron transport taking into account:
1886 // 1. diffusion,
1887 // 2.ExB at the wires
1888 // 3. nonisochronity
1889 //
1890 // xyz and index must be already transformed to system 1
1891 //
1892
1893 fTPCParam->Transform1to2(xyz,index);
1894
1895 //add diffusion
1896 Float_t driftl=xyz[2];
1897 if(driftl<0.01) driftl=0.01;
1898 driftl=TMath::Sqrt(driftl);
73042f01 1899 Float_t sigT = driftl*(fTPCParam->GetDiffT());
1900 Float_t sigL = driftl*(fTPCParam->GetDiffL());
1901 xyz[0]=gRandom->Gaus(xyz[0],sigT);
1902 xyz[1]=gRandom->Gaus(xyz[1],sigT);
1903 xyz[2]=gRandom->Gaus(xyz[2],sigL);
cc80f89e 1904
1905 // ExB
1906
1907 if (fTPCParam->GetMWPCReadout()==kTRUE){
b584c7dd 1908 Float_t dx = fTPCParam->Transform2to2NearestWire(xyz,index);
cc80f89e 1909 xyz[1]+=dx*(fTPCParam->GetOmegaTau());
1910 }
b584c7dd 1911 //add nonisochronity (not implemented yet)
cc80f89e 1912}
fe4da5cc 1913
fe4da5cc 1914ClassImp(AliTPChit)
1915
1916//_____________________________________________________________________________
1917AliTPChit::AliTPChit(Int_t shunt, Int_t track, Int_t *vol, Float_t *hits):
1918AliHit(shunt,track)
1919{
1920 //
1921 // Creates a TPC hit object
1922 //
1923 fSector = vol[0];
1924 fPadRow = vol[1];
1925 fX = hits[0];
1926 fY = hits[1];
1927 fZ = hits[2];
1928 fQ = hits[3];
e61fd20d 1929 fTime = hits[4];
fe4da5cc 1930}
1931
39c8eb58 1932//________________________________________________________________________
792bb11c 1933// Additional code because of the AliTPCTrackHitsV2
39c8eb58 1934
176aff27 1935void AliTPC::MakeBranch2(Option_t *option,const char */*file*/)
39c8eb58 1936{
1937 //
1938 // Create a new branch in the current Root Tree
1939 // The branch of fHits is automatically split
1940 // MI change 14.09.2000
8c2b3fd7 1941 AliDebug(1,"");
792bb11c 1942 if (fHitType<2) return;
39c8eb58 1943 char branchname[10];
1944 sprintf(branchname,"%s2",GetName());
1945 //
1946 // Get the pointer to the header
5cf7bbad 1947 const char *cH = strstr(option,"H");
39c8eb58 1948 //
8c2b3fd7 1949 if (fTrackHits && TreeH() && cH && fHitType&4) {
1950 AliDebug(1,"Making branch for Type 4 Hits");
88cb7938 1951 TreeH()->Branch(branchname,"AliTPCTrackHitsV2",&fTrackHits,fBufferSize,99);
8c2b3fd7 1952 }
792bb11c 1953
be5ffbfe 1954// if (fTrackHitsOld && TreeH() && cH && fHitType&2) {
1955// AliDebug(1,"Making branch for Type 2 Hits");
1956// AliObjectBranch * branch = new AliObjectBranch(branchname,"AliTPCTrackHits",&fTrackHitsOld,
1957// TreeH(),fBufferSize,99);
1958// TreeH()->GetListOfBranches()->Add(branch);
1959// }
39c8eb58 1960}
1961
1962void AliTPC::SetTreeAddress()
1963{
8c2b3fd7 1964 //Sets tree address for hits
1965 if (fHitType<=1) {
1966 if (fHits == 0x0 ) fHits = new TClonesArray("AliTPChit", 176);//skowron 20.06.03
1967 AliDetector::SetTreeAddress();
1968 }
792bb11c 1969 if (fHitType>1) SetTreeAddress2();
39c8eb58 1970}
1971
1972void AliTPC::SetTreeAddress2()
1973{
1974 //
1975 // Set branch address for the TrackHits Tree
1976 //
8c2b3fd7 1977 AliDebug(1,"");
88cb7938 1978
39c8eb58 1979 TBranch *branch;
1980 char branchname[20];
1981 sprintf(branchname,"%s2",GetName());
1982 //
1983 // Branch address for hit tree
88cb7938 1984 TTree *treeH = TreeH();
792bb11c 1985 if ((treeH)&&(fHitType&4)) {
39c8eb58 1986 branch = treeH->GetBranch(branchname);
8c2b3fd7 1987 if (branch) {
1988 branch->SetAddress(&fTrackHits);
1989 AliDebug(1,"fHitType&4 Setting");
1990 }
f2a509af 1991 else
8c2b3fd7 1992 AliDebug(1,"fHitType&4 Failed (can not find branch)");
f2a509af 1993
39c8eb58 1994 }
be5ffbfe 1995 // if ((treeH)&&(fHitType&2)) {
1996// branch = treeH->GetBranch(branchname);
1997// if (branch) {
1998// branch->SetAddress(&fTrackHitsOld);
1999// AliDebug(1,"fHitType&2 Setting");
2000// }
2001// else
2002// AliDebug(1,"fHitType&2 Failed (can not find branch)");
2003// }
b6e0d3fe 2004 //set address to TREETR
8c2b3fd7 2005
88cb7938 2006 TTree *treeTR = TreeTR();
b6e0d3fe 2007 if (treeTR && fTrackReferences) {
2008 branch = treeTR->GetBranch(GetName());
2009 if (branch) branch->SetAddress(&fTrackReferences);
2010 }
2011
39c8eb58 2012}
2013
2014void AliTPC::FinishPrimary()
2015{
792bb11c 2016 if (fTrackHits &&fHitType&4) fTrackHits->FlushHitStack();
be5ffbfe 2017 // if (fTrackHitsOld && fHitType&2) fTrackHitsOld->FlushHitStack();
39c8eb58 2018}
2019
2020
2021void AliTPC::AddHit2(Int_t track, Int_t *vol, Float_t *hits)
2022{
2023 //
2024 // add hit to the list
39c8eb58 2025 Int_t rtrack;
2026 if (fIshunt) {
5d12ce38 2027 int primary = gAlice->GetMCApp()->GetPrimary(track);
2028 gAlice->GetMCApp()->Particle(primary)->SetBit(kKeepBit);
39c8eb58 2029 rtrack=primary;
2030 } else {
2031 rtrack=track;
5d12ce38 2032 gAlice->GetMCApp()->FlagTrack(track);
39c8eb58 2033 }
792bb11c 2034 if (fTrackHits && fHitType&4)
39c8eb58 2035 fTrackHits->AddHitKartez(vol[0],rtrack, hits[0],
e61fd20d 2036 hits[1],hits[2],(Int_t)hits[3],hits[4]);
be5ffbfe 2037 // if (fTrackHitsOld &&fHitType&2 )
2038// fTrackHitsOld->AddHitKartez(vol[0],rtrack, hits[0],
2039// hits[1],hits[2],(Int_t)hits[3]);
792bb11c 2040
39c8eb58 2041}
2042
2043void AliTPC::ResetHits()
88cb7938 2044{
39c8eb58 2045 if (fHitType&1) AliDetector::ResetHits();
792bb11c 2046 if (fHitType>1) ResetHits2();
39c8eb58 2047}
2048
2049void AliTPC::ResetHits2()
2050{
2051 //
2052 //reset hits
792bb11c 2053 if (fTrackHits && fHitType&4) fTrackHits->Clear();
be5ffbfe 2054 // if (fTrackHitsOld && fHitType&2) fTrackHitsOld->Clear();
792bb11c 2055
39c8eb58 2056}
2057
2058AliHit* AliTPC::FirstHit(Int_t track)
2059{
792bb11c 2060 if (fHitType>1) return FirstHit2(track);
39c8eb58 2061 return AliDetector::FirstHit(track);
2062}
2063AliHit* AliTPC::NextHit()
2064{
9bdd974b 2065 //
2066 // gets next hit
2067 //
792bb11c 2068 if (fHitType>1) return NextHit2();
2069
39c8eb58 2070 return AliDetector::NextHit();
2071}
2072
2073AliHit* AliTPC::FirstHit2(Int_t track)
2074{
2075 //
2076 // Initialise the hit iterator
2077 // Return the address of the first hit for track
2078 // If track>=0 the track is read from disk
2079 // while if track<0 the first hit of the current
2080 // track is returned
2081 //
2082 if(track>=0) {
2083 gAlice->ResetHits();
88cb7938 2084 TreeH()->GetEvent(track);
39c8eb58 2085 }
2086 //
792bb11c 2087 if (fTrackHits && fHitType&4) {
39c8eb58 2088 fTrackHits->First();
2089 return fTrackHits->GetHit();
2090 }
be5ffbfe 2091 // if (fTrackHitsOld && fHitType&2) {
2092// fTrackHitsOld->First();
2093// return fTrackHitsOld->GetHit();
2094// }
792bb11c 2095
39c8eb58 2096 else return 0;
2097}
2098
2099AliHit* AliTPC::NextHit2()
2100{
2101 //
2102 //Return the next hit for the current track
2103
792bb11c 2104
be5ffbfe 2105// if (fTrackHitsOld && fHitType&2) {
2106// fTrackHitsOld->Next();
2107// return fTrackHitsOld->GetHit();
2108// }
39c8eb58 2109 if (fTrackHits) {
2110 fTrackHits->Next();
2111 return fTrackHits->GetHit();
2112 }
2113 else
2114 return 0;
2115}
2116
2117void AliTPC::LoadPoints(Int_t)
2118{
2119 //
2120 Int_t a = 0;
8c2b3fd7 2121
f2e8b846 2122 if(fHitType==1) AliDetector::LoadPoints(a);
2123 else LoadPoints2(a);
39c8eb58 2124}
2125
2126
2127void AliTPC::RemapTrackHitIDs(Int_t *map)
2128{
9bdd974b 2129 //
2130 // remapping
2131 //
39c8eb58 2132 if (!fTrackHits) return;
39c8eb58 2133
be5ffbfe 2134// if (fTrackHitsOld && fHitType&2){
2135// AliObjectArray * arr = fTrackHitsOld->fTrackHitsInfo;
2136// for (UInt_t i=0;i<arr->GetSize();i++){
2137// AliTrackHitsInfo * info = (AliTrackHitsInfo *)(arr->At(i));
2138// info->fTrackID = map[info->fTrackID];
2139// }
2140// }
2141// if (fTrackHitsOld && fHitType&4){
2142 if (fTrackHits && fHitType&4){
792bb11c 2143 TClonesArray * arr = fTrackHits->GetArray();;
2144 for (Int_t i=0;i<arr->GetEntriesFast();i++){
2145 AliTrackHitsParamV2 * info = (AliTrackHitsParamV2 *)(arr->At(i));
2146 info->fTrackID = map[info->fTrackID];
2147 }
2148 }
39c8eb58 2149}
2150
792bb11c 2151Bool_t AliTPC::TrackInVolume(Int_t id,Int_t track)
2152{
2153 //return bool information - is track in given volume
2154 //load only part of the track information
2155 //return true if current track is in volume
2156 //
2157 // return kTRUE;
be5ffbfe 2158 // if (fTrackHitsOld && fHitType&2) {
2159// TBranch * br = TreeH()->GetBranch("fTrackHitsInfo");
2160// br->GetEvent(track);
2161// AliObjectArray * ar = fTrackHitsOld->fTrackHitsInfo;
2162// for (UInt_t j=0;j<ar->GetSize();j++){
2163// if ( ((AliTrackHitsInfo*)ar->At(j))->fVolumeID==id) return kTRUE;
2164// }
2165// }
792bb11c 2166
2167 if (fTrackHits && fHitType&4) {
88cb7938 2168 TBranch * br1 = TreeH()->GetBranch("fVolumes");
2169 TBranch * br2 = TreeH()->GetBranch("fNVolumes");
792bb11c 2170 br2->GetEvent(track);
2171 br1->GetEvent(track);
2172 Int_t *volumes = fTrackHits->GetVolumes();
2173 Int_t nvolumes = fTrackHits->GetNVolumes();
2174 if (!volumes && nvolumes>0) {
8c2b3fd7 2175 AliWarning(Form("Problematic track\t%d\t%d",track,nvolumes));
792bb11c 2176 return kFALSE;
2177 }
2178 for (Int_t j=0;j<nvolumes; j++)
2179 if (volumes[j]==id) return kTRUE;
2180 }
2181
2182 if (fHitType&1) {
88cb7938 2183 TBranch * br = TreeH()->GetBranch("fSector");
792bb11c 2184 br->GetEvent(track);
2185 for (Int_t j=0;j<fHits->GetEntriesFast();j++){
2186 if ( ((AliTPChit*)fHits->At(j))->fSector==id) return kTRUE;
2187 }
2188 }
2189 return kFALSE;
2190
2191}
39c8eb58 2192
2193//_____________________________________________________________________________
2194void AliTPC::LoadPoints2(Int_t)
2195{
2196 //
2197 // Store x, y, z of all hits in memory
2198 //
be5ffbfe 2199 // if (fTrackHits == 0 && fTrackHitsOld==0) return;
2200 if (fTrackHits == 0 ) return;
39c8eb58 2201 //
792bb11c 2202 Int_t nhits =0;
2203 if (fHitType&4) nhits = fTrackHits->GetEntriesFast();
be5ffbfe 2204 // if (fHitType&2) nhits = fTrackHitsOld->GetEntriesFast();
792bb11c 2205
39c8eb58 2206 if (nhits == 0) return;
5d12ce38 2207 Int_t tracks = gAlice->GetMCApp()->GetNtrack();
39c8eb58 2208 if (fPoints == 0) fPoints = new TObjArray(tracks);
2209 AliHit *ahit;
2210 //
2211 Int_t *ntrk=new Int_t[tracks];
2212 Int_t *limi=new Int_t[tracks];
2213 Float_t **coor=new Float_t*[tracks];
2214 for(Int_t i=0;i<tracks;i++) {
2215 ntrk[i]=0;
2216 coor[i]=0;
2217 limi[i]=0;
2218 }
2219 //
2220 AliPoints *points = 0;
2221 Float_t *fp=0;
2222 Int_t trk;
2223 Int_t chunk=nhits/4+1;
2224 //
2225 // Loop over all the hits and store their position
2226 //
2227 ahit = FirstHit2(-1);
39c8eb58 2228 while (ahit){
39c8eb58 2229 trk=ahit->GetTrack();
2230 if(ntrk[trk]==limi[trk]) {
2231 //
2232 // Initialise a new track
2233 fp=new Float_t[3*(limi[trk]+chunk)];
2234 if(coor[trk]) {
2235 memcpy(fp,coor[trk],sizeof(Float_t)*3*limi[trk]);
2236 delete [] coor[trk];
2237 }
2238 limi[trk]+=chunk;
2239 coor[trk] = fp;
2240 } else {
2241 fp = coor[trk];
2242 }
2243 fp[3*ntrk[trk] ] = ahit->X();
2244 fp[3*ntrk[trk]+1] = ahit->Y();
2245 fp[3*ntrk[trk]+2] = ahit->Z();
2246 ntrk[trk]++;
2247 ahit = NextHit2();
2248 }
792bb11c 2249
2250
2251
39c8eb58 2252 //
2253 for(trk=0; trk<tracks; ++trk) {
2254 if(ntrk[trk]) {
2255 points = new AliPoints();
2256 points->SetMarkerColor(GetMarkerColor());
2257 points->SetMarkerSize(GetMarkerSize());
2258 points->SetDetector(this);
2259 points->SetParticle(trk);
2260 points->SetPolyMarker(ntrk[trk],coor[trk],GetMarkerStyle());
2261 fPoints->AddAt(points,trk);
2262 delete [] coor[trk];
2263 coor[trk]=0;
2264 }
2265 }
2266 delete [] coor;
2267 delete [] ntrk;
2268 delete [] limi;
2269}
2270
2271
2272//_____________________________________________________________________________
2273void AliTPC::LoadPoints3(Int_t)
2274{
2275 //
2276 // Store x, y, z of all hits in memory
2277 // - only intersection point with pad row
2278 if (fTrackHits == 0) return;
2279 //
2280 Int_t nhits = fTrackHits->GetEntriesFast();
2281 if (nhits == 0) return;
5d12ce38 2282 Int_t tracks = gAlice->GetMCApp()->GetNtrack();
39c8eb58 2283 if (fPoints == 0) fPoints = new TObjArray(2*tracks);
2284 fPoints->Expand(2*tracks);
2285 AliHit *ahit;
2286 //
2287 Int_t *ntrk=new Int_t[tracks];
2288 Int_t *limi=new Int_t[tracks];
2289 Float_t **coor=new Float_t*[tracks];
2290 for(Int_t i=0;i<tracks;i++) {
2291 ntrk[i]=0;
2292 coor[i]=0;
2293 limi[i]=0;
2294 }
2295 //
2296 AliPoints *points = 0;
2297 Float_t *fp=0;
2298 Int_t trk;
2299 Int_t chunk=nhits/4+1;
2300 //
2301 // Loop over all the hits and store their position
2302 //
2303 ahit = FirstHit2(-1);
39c8eb58 2304
2305 Int_t lastrow = -1;
2306 while (ahit){
39c8eb58 2307 trk=ahit->GetTrack();
2308 Float_t x[3]={ahit->X(),ahit->Y(),ahit->Z()};
2309 Int_t index[3]={1,((AliTPChit*)ahit)->fSector,0};
2310 Int_t currentrow = fTPCParam->GetPadRow(x,index) ;
2311 if (currentrow!=lastrow){
2312 lastrow = currentrow;
2313 //later calculate intersection point
2314 if(ntrk[trk]==limi[trk]) {
2315 //
2316 // Initialise a new track
2317 fp=new Float_t[3*(limi[trk]+chunk)];
2318 if(coor[trk]) {
2319 memcpy(fp,coor[trk],sizeof(Float_t)*3*limi[trk]);
2320 delete [] coor[trk];
2321 }
2322 limi[trk]+=chunk;
2323 coor[trk] = fp;
2324 } else {
2325 fp = coor[trk];
2326 }
2327 fp[3*ntrk[trk] ] = ahit->X();
2328 fp[3*ntrk[trk]+1] = ahit->Y();
2329 fp[3*ntrk[trk]+2] = ahit->Z();
2330 ntrk[trk]++;
2331 }
2332 ahit = NextHit2();
2333 }
2334
2335 //
2336 for(trk=0; trk<tracks; ++trk) {
2337 if(ntrk[trk]) {
2338 points = new AliPoints();
2339 points->SetMarkerColor(GetMarkerColor()+1);
2340 points->SetMarkerStyle(5);
2341 points->SetMarkerSize(0.2);
2342 points->SetDetector(this);
2343 points->SetParticle(trk);
39c8eb58 2344 points->SetPolyMarker(ntrk[trk],coor[trk],30);
2345 fPoints->AddAt(points,tracks+trk);
2346 delete [] coor[trk];
2347 coor[trk]=0;
2348 }
2349 }
2350 delete [] coor;
2351 delete [] ntrk;
2352 delete [] limi;
2353}
2354
2355
2356
88cb7938 2357AliLoader* AliTPC::MakeLoader(const char* topfoldername)
2358{
8c2b3fd7 2359 //Makes TPC loader
2360 fLoader = new AliTPCLoader(GetName(),topfoldername);
2361 return fLoader;
88cb7938 2362}
2363
42157e55 2364////////////////////////////////////////////////////////////////////////
2365AliTPCParam* AliTPC::LoadTPCParam(TFile *file) {
2366//
2367// load TPC paarmeters from a given file or create new if the object
2368// is not found there
88cb7938 2369// 12/05/2003 This method should be moved to the AliTPCLoader
2370// and one has to decide where to store the TPC parameters
2371// M.Kowalski
42157e55 2372 char paramName[50];
2373 sprintf(paramName,"75x40_100x60_150x60");
2374 AliTPCParam *paramTPC=(AliTPCParam*)file->Get(paramName);
2375 if (paramTPC) {
8c2b3fd7 2376 AliDebugClass(1,Form("TPC parameters %s found.",paramName));
42157e55 2377 } else {
8c2b3fd7 2378 AliWarningClass("TPC parameters not found. Create new (they may be incorrect)");
42157e55 2379 paramTPC = new AliTPCParamSR;
2380 }
2381 return paramTPC;
2382
2383// the older version of parameters can be accessed with this code.
2384// In some cases, we have old parameters saved in the file but
2385// digits were created with new parameters, it can be distinguish
2386// by the name of TPC TreeD. The code here is just for the case
2387// we would need to compare with old data, uncomment it if needed.
2388//
2389// char paramName[50];
2390// sprintf(paramName,"75x40_100x60");
2391// AliTPCParam *paramTPC=(AliTPCParam*)in->Get(paramName);
2392// if (paramTPC) {
2393// cout<<"TPC parameters "<<paramName<<" found."<<endl;
2394// } else {
2395// sprintf(paramName,"75x40_100x60_150x60");
2396// paramTPC=(AliTPCParam*)in->Get(paramName);
2397// if (paramTPC) {
2398// cout<<"TPC parameters "<<paramName<<" found."<<endl;
2399// } else {
2400// cerr<<"TPC parameters not found. Create new (they may be incorrect)."
2401// <<endl;
2402// paramTPC = new AliTPCParamSR;
2403// }
2404// }
2405// return paramTPC;
2406
2407}
2408
85a5290f 2409