7d54fa91f333744cf2d473cb806fad9b20c4d50d
[u/mrichter/AliRoot.git] / TPC / AliTPCDigitizer.cxx
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
16
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
34 ClassImp(AliTPCDigitizer)
35
36 //___________________________________________
37 AliTPCDigitizer::AliTPCDigitizer() :AliDigitizer()
38 {
39 // Default ctor - don't use it
40 }
41
42 //___________________________________________
43 AliTPCDigitizer::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 //------------------------------------------------------------------------
53 AliTPCDigitizer::~AliTPCDigitizer()
54 {
55 // Destructor
56 }
57
58
59
60 //------------------------------------------------------------------------
61 Bool_t AliTPCDigitizer::Init()
62 {
63 // Initialization 
64     
65  return kTRUE;
66 }
67
68
69 //------------------------------------------------------------------------
70 void 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();
84   pTPC->GenerNoise(500000); //create teble with noise
85   printf("noise %f \n",  param->GetNoise()*param->GetNoiseNormFac());
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
94   //create digits array for given sectors
95   AliSimDigits ** digarr = new AliSimDigits*[nInputs]; 
96   for (Int_t i1=0;i1<nInputs; i1++){
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();
111   //if (tree->GetBranch("Segment") ) tree->GetBranch("Segment")->SetAddress(&digrow);
112   //else
113   tree->Branch("Segment","AliSimDigits",&digrow);
114   //
115
116   param->SetZeroSup(2);
117
118   Int_t zerosup = param->GetZeroSup();
119   //Loop over segments of the TPC
120     
121   for (Int_t n=0; n<nentries; n++) {
122   //    for (Int_t n=0; n<300; n++) {
123            
124     for (Int_t i=0;i<nInputs; i++){ 
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
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           
169           for (Int_t tr=0;tr<3;tr++) {
170             //             Int_t lab = digarr[i]->GetTrackIDFast(rows,col,tr);
171             Int_t lab = ptr[i][tr*nElems];
172             if ( lab > 1) {
173               label[labptr]=lab+masks[i];
174               labptr++;
175             }      
176           }
177           pdig[i]++;
178           ptr[i]++;
179           
180         }
181         q/=16.;  //conversion factor
182         //      Float_t noise  = gRandom->Gaus(0,param->GetNoise()*param->GetNoiseNormFac());  
183         Float_t noise  = pTPC->GetNoise();
184         q+=noise;
185         q=TMath::Nint(q);
186         if (q > zerosup){ 
187           
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++;
202     }
203     
204     digrow->CompresBuffer(1,zerosup);
205     digrow->CompresTrackBuffer(1);
206     tree->Fill();
207     if (fDebug>0) cerr<<sec<<"\t"<<row<<"\n";  
208   }
209   fManager->GetTreeDTPC()->Write(0,TObject::kOverwrite);
210   delete digrow;     
211   for (Int_t i1=0;i1<nInputs; i1++) delete digarr[i1];
212   delete []masks;
213   delete digarr;  
214 }
215