Better protection.
[u/mrichter/AliRoot.git] / TPC / AliTPCDigitizer.cxx
CommitLineData
3c038d07 1/**************************************************************************
2 * Copyright(c) 1998-2000, 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
407ff276 16
3c038d07 17
18#include <TTree.h>
19#include <TObjArray.h>
20#include <TFile.h>
21#include <TDirectory.h>
22#include <iostream.h>
23
24#include "AliTPCDigitizer.h"
25
26#include "AliTPC.h"
27#include "AliTPCParam.h"
28#include "AliRun.h"
29#include "AliPDG.h"
30#include "AliRunDigitizer.h"
31#include "AliSimDigits.h"
32
33
34ClassImp(AliTPCDigitizer)
35
36//___________________________________________
37AliTPCDigitizer::AliTPCDigitizer() :AliDigitizer()
38{
39// Default ctor - don't use it
40}
41
42//___________________________________________
43AliTPCDigitizer::AliTPCDigitizer(AliRunDigitizer* manager)
44 :AliDigitizer(manager)
45{
46// ctor which should be used
47 if (GetDebug()>2)
48 cerr<<"AliTPCDigitizer::AliTPCDigitizer"
49 <<"(AliRunDigitizer* manager) was processed"<<endl;
50}
51
52//------------------------------------------------------------------------
53AliTPCDigitizer::~AliTPCDigitizer()
54{
55// Destructor
56}
57
58
59
60//------------------------------------------------------------------------
61Bool_t AliTPCDigitizer::Init()
62{
63// Initialization
64
65 return kTRUE;
66}
67
68
69//------------------------------------------------------------------------
70void AliTPCDigitizer::Exec(Option_t* option)
71{
72
73 // merge input tree's with summable digits
74 //output stored in TreeTPCD
75
76 TString optionString = option;
77 if (optionString.Data() == "deb") {
78 cout<<"AliTPCDigitizer::Exec: called with option deb "<<endl;
79 fDebug = 3;
80 }
81 //get detector and geometry
82 AliTPC *pTPC = (AliTPC *) gAlice->GetModule("TPC");
83 AliTPCParam * param = pTPC->GetParam();
407ff276 84 pTPC->GenerNoise(500000); //create teble with noise
3c038d07 85 printf("noise %f \n", param->GetNoise()*param->GetNoiseNormFac());
407ff276 86 //
87 Int_t nInputs = fManager->GetNinputs();
88 Int_t * masks = new Int_t[nInputs];
89 for (Int_t i=0; i<nInputs;i++)
90 masks[i]= fManager->GetMask(i);
91 Short_t **pdig= new Short_t*[nInputs]; //pointers to the expanded digits array
92 Int_t **ptr= new Int_t*[nInputs]; //pointers to teh expanded tracks array
93
3c038d07 94 //create digits array for given sectors
407ff276 95 AliSimDigits ** digarr = new AliSimDigits*[nInputs];
96 for (Int_t i1=0;i1<nInputs; i1++){
3c038d07 97 digarr[i1]=0;
98 // intree[i1]
99 TTree * treear = fManager->GetInputTreeTPCS(i1);
100 if (!treear) {
101 cerr<<" TPC - not existing input = \n"<<i1<<" ";
102 }
103 treear->GetBranch("Segment")->SetAddress(&digarr[i1]);
104 }
105 Stat_t nentries = fManager->GetInputTreeTPCS(0)->GetEntries();
106
107
108 //create branch's in TPC treeD
109 AliSimDigits * digrow = new AliSimDigits;
110 TTree * tree = fManager->GetTreeDTPC();
407ff276 111 //if (tree->GetBranch("Segment") ) tree->GetBranch("Segment")->SetAddress(&digrow);
112 //else
113 tree->Branch("Segment","AliSimDigits",&digrow);
3c038d07 114 //
115
407ff276 116 param->SetZeroSup(2);
117
3c038d07 118 Int_t zerosup = param->GetZeroSup();
119 //Loop over segments of the TPC
120
121 for (Int_t n=0; n<nentries; n++) {
407ff276 122 // for (Int_t n=0; n<300; n++) {
123
124 for (Int_t i=0;i<nInputs; i++){
3c038d07 125 fManager->GetInputTreeTPCS(i)->GetEvent(n);
126 digarr[i]->ExpandBuffer();
127 digarr[i]->ExpandTrackBuffer();
128 }
129 if ((digarr[0]->GetID()-digarr[1]->GetID())>0)
130 printf("problem\n");
131 Int_t sec, row;
132 if (!param->AdjustSectorRow(digarr[0]->GetID(),sec,row)) {
133 cerr<<"AliTPC warning: invalid segment ID ! "<<digarr[0]->GetID()<<endl;
134 continue;
135 }
136
137 digrow->SetID(digarr[0]->GetID());
138
139 Int_t nrows = digarr[0]->GetNRows();
140 Int_t ncols = digarr[0]->GetNCols();
141 digrow->Allocate(nrows,ncols);
142 digrow->AllocateTrack(3);
143
407ff276 144 Float_t q=0;
145 Int_t label[1000]; //stack for 300 events
146 Int_t labptr = 0;
147
148
149 for (Int_t i=0;i<nInputs; i++){
150 pdig[i] = digarr[i]->GetDigits();
151 ptr[i] = digarr[i]->GetTracks();
152 }
153 Short_t *pdig1= digrow->GetDigits();
154 Int_t *ptr1= digrow->GetTracks() ;
155
156 Int_t nElems = nrows*ncols;
157
158 // for (Int_t rows=0;rows<nrows; rows++){
159 // for (Int_t col=0;col<ncols; col++){
160 for (Int_t elem=0;elem<nElems; elem++){
161
162 q=0;
163 labptr=0;
164 // looop over digits
165 for (Int_t i=0;i<nInputs; i++){
166 // q += digarr[i]->GetDigitFast(rows,col);
167 q += *(pdig[i]);
168
3c038d07 169 for (Int_t tr=0;tr<3;tr++) {
407ff276 170 // Int_t lab = digarr[i]->GetTrackIDFast(rows,col,tr);
171 Int_t lab = ptr[i][tr*nElems];
3c038d07 172 if ( lab > 1) {
407ff276 173 label[labptr]=lab+masks[i];
3c038d07 174 labptr++;
407ff276 175 }
3c038d07 176 }
407ff276 177 pdig[i]++;
178 ptr[i]++;
179
3c038d07 180 }
3c038d07 181 q/=16.; //conversion factor
407ff276 182 // Float_t noise = gRandom->Gaus(0,param->GetNoise()*param->GetNoiseNormFac());
183 Float_t noise = pTPC->GetNoise();
3c038d07 184 q+=noise;
185 q=TMath::Nint(q);
407ff276 186 if (q > zerosup){
3c038d07 187
407ff276 188 if(q > param->GetADCSat()) q = (Short_t)(param->GetADCSat());
189 //digrow->SetDigitFast((Short_t)q,rows,col);
190 *pdig1 =Short_t(q);
191 for (Int_t tr=0;tr<3;tr++){
192 if (tr<labptr)
193 // ((AliSimDigits*)digrow)->SetTrackIDFast(label[tr],rows,col,tr);
194 ptr1[tr*nElems] = label[tr];
195 //else
196 // ((AliSimDigits*)digrow)->SetTrackIDFast(-1,rows,col,tr);
197 // ptr1[tr*nElems] = 1;
198 }
199 }
200 pdig1++;
201 ptr1++;
3c038d07 202 }
407ff276 203
3c038d07 204 digrow->CompresBuffer(1,zerosup);
205 digrow->CompresTrackBuffer(1);
206 tree->Fill();
207 if (fDebug>0) cerr<<sec<<"\t"<<row<<"\n";
208 }
407ff276 209 fManager->GetTreeDTPC()->Write(0,TObject::kOverwrite);
3c038d07 210 delete digrow;
407ff276 211 for (Int_t i1=0;i1<nInputs; i1++) delete digarr[i1];
212 delete []masks;
3c038d07 213 delete digarr;
214}
215