77c21172e5c487adc122366e4995b944a26fe785
[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 $Log$
18 Revision 1.8  2002/10/23 07:17:33  alibrary
19 Introducing Riostream.h
20
21 Revision 1.7  2002/10/14 14:57:42  hristov
22 Merging the VirtualMC branch to the main development branch (HEAD)
23
24 Revision 1.3.6.1  2002/06/10 15:26:11  hristov
25 Merged with v3-08-02
26
27 Revision 1.6  2002/04/06 14:41:04  kowal2
28 Added #include<stdlib.h> and log
29
30 */
31
32
33
34
35 #include <stdlib.h>
36 #include <TTree.h> 
37 #include <TObjArray.h>
38 #include <TFile.h>
39 #include <TDirectory.h>
40 #include <Riostream.h>
41
42 #include "AliTPCDigitizer.h"
43
44 #include "AliTPC.h"
45 #include "AliTPCParam.h"
46 #include "AliTPCParamSR.h" 
47 #include "AliRun.h"
48 #include "AliPDG.h"
49 #include "AliRunDigitizer.h"
50 #include "AliSimDigits.h"
51
52
53 ClassImp(AliTPCDigitizer)
54
55 //___________________________________________
56 AliTPCDigitizer::AliTPCDigitizer() :AliDigitizer()
57 {
58 // Default ctor - don't use it
59   fDebug =0;
60 }
61
62 //___________________________________________
63 AliTPCDigitizer::AliTPCDigitizer(AliRunDigitizer* manager) 
64     :AliDigitizer(manager)
65 {
66 // ctor which should be used
67   fDebug =0;
68   if (GetDebug()>2) 
69     cerr<<"AliTPCDigitizer::AliTPCDigitizer"
70         <<"(AliRunDigitizer* manager) was processed"<<endl;
71 }
72
73 //------------------------------------------------------------------------
74 AliTPCDigitizer::~AliTPCDigitizer()
75 {
76 // Destructor
77 }
78
79
80
81 //------------------------------------------------------------------------
82 Bool_t AliTPCDigitizer::Init()
83 {
84 // Initialization 
85     
86  return kTRUE;
87 }
88
89
90 //------------------------------------------------------------------------
91 void AliTPCDigitizer::Exec(Option_t* option)
92 {
93   ExecFast(option);  
94 }
95 //------------------------------------------------------------------------
96 void AliTPCDigitizer::ExecFast(Option_t* option)
97 {
98   
99   // merge input tree's with summable digits
100   //output stored in TreeTPCD
101   char s[100]; 
102   char ss[100];
103   TString optionString = option;
104   if (optionString.Data() == "deb") {
105     cout<<"AliTPCDigitizer::Exec: called with option deb "<<endl;
106     fDebug = 3;
107   }
108   //get detector and geometry
109   AliTPC *pTPC  = (AliTPC *) gAlice->GetModule("TPC");
110   AliTPCParam * param = pTPC->GetParam();
111   
112   sprintf(s,param->GetTitle());
113   sprintf(ss,"75x40_100x60");
114   if(strcmp(s,ss)==0){
115     printf("2 pad-length geom hits with 3 pad-lenght geom digits...\n");
116     delete param;
117     param=new AliTPCParamSR();
118   }
119   else{
120    sprintf(ss,"75x40_100x60_150x60");
121    if(strcmp(s,ss)!=0) {
122      printf("No TPC parameters found...\n");
123      exit(2); 
124    }
125   }
126   
127   pTPC->GenerNoise(500000); //create teble with noise
128   //
129   Int_t nInputs = fManager->GetNinputs();
130   Int_t * masks = new Int_t[nInputs];
131   for (Int_t i=0; i<nInputs;i++)
132     masks[i]= fManager->GetMask(i);
133   Short_t **pdig= new Short_t*[nInputs];   //pointers to the expanded digits array
134   Int_t **ptr=  new Int_t*[nInputs];       //pointers to teh expanded tracks array
135
136   //create digits array for given sectors
137   // make indexes
138   AliSimDigits ** digarr = new AliSimDigits*[nInputs]; 
139   for (Int_t i1=0;i1<nInputs; i1++){
140     digarr[i1]=0;
141     //    intree[i1]
142     TTree * treear =  fManager->GetInputTreeTPCS(i1);
143     if (!treear) {
144       cerr<<"AliTPCDigitizer: Input tree with SDigits not found in"
145           <<" input "<< i1<<endl;
146       return;
147     }
148     treear->GetBranch("Segment")->SetAddress(&digarr[i1]);
149     if (treear->GetIndex()==0) 
150       treear->BuildIndex("fSegmentID","fSegmentID");
151   }
152   Stat_t nentries = fManager->GetInputTreeTPCS(0)->GetEntries();
153   
154
155   //create branch's in TPC treeD
156   AliSimDigits * digrow = new AliSimDigits;
157   TTree * tree  = fManager->GetTreeDTPC();
158   tree->Branch("Segment","AliSimDigits",&digrow);
159   //
160
161   param->SetZeroSup(2);
162
163   Int_t zerosup = param->GetZeroSup(); 
164   //
165   //Loop over segments of the TPC
166     
167   for (Int_t n=0; n<nentries; n++) {
168   //    for (Int_t n=0; n<300; n++) {
169     fManager->GetInputTreeTPCS(0)->GetEvent(n);      
170     digarr[0]->ExpandBuffer();
171     digarr[0]->ExpandTrackBuffer();
172            
173     for (Int_t i=1;i<nInputs; i++){ 
174       fManager->GetInputTreeTPCS(i)->GetEntryWithIndex(digarr[0]->GetID(),digarr[0]->GetID());
175       if ((digarr[0]->GetID()-digarr[i]->GetID())>0) 
176         printf("problem - not corresponding segment in background event\n");
177       
178       digarr[i]->ExpandBuffer();
179       digarr[i]->ExpandTrackBuffer();
180     }   
181     Int_t sec, row;
182     if (!param->AdjustSectorRow(digarr[0]->GetID(),sec,row)) {
183       cerr<<"AliTPC warning: invalid segment ID ! "<<digarr[0]->GetID()<<endl;
184       continue;
185     }
186
187     digrow->SetID(digarr[0]->GetID());
188
189     Int_t nrows = digarr[0]->GetNRows();
190     Int_t ncols = digarr[0]->GetNCols();
191     digrow->Allocate(nrows,ncols);
192     digrow->AllocateTrack(3);
193
194     Float_t q=0;
195     Int_t label[1000]; //stack for 300 events 
196     Int_t labptr = 0;
197
198     Int_t nElems = nrows*ncols;     
199  
200     for (Int_t i=0;i<nInputs; i++){ 
201       pdig[i] = digarr[i]->GetDigits();
202       ptr[i]  = digarr[i]->GetTracks();
203     }
204     Short_t *pdig1= digrow->GetDigits();
205     Int_t   *ptr1= digrow->GetTracks() ;
206
207     
208
209     //    for (Int_t rows=0;rows<nrows; rows++){
210     //  for (Int_t col=0;col<ncols; col++){
211     for (Int_t elem=0;elem<nElems; elem++){    
212       //for (Int_t elem=nElems;elem<nElems; elem++){
213
214         q=0;
215         labptr=0;
216         // looop over digits 
217         for (Int_t i=0;i<nInputs; i++){ 
218           //          q  += digarr[i]->GetDigitFast(rows,col);
219           q  += *(pdig[i]);
220           
221           for (Int_t tr=0;tr<3;tr++) {
222             //             Int_t lab = digarr[i]->GetTrackIDFast(rows,col,tr);
223             Int_t lab = ptr[i][tr*nElems];
224             if ( (lab > 1) && *(pdig[i])>zerosup) {
225               label[labptr]=lab+masks[i];
226               labptr++;
227             }      
228           }
229           pdig[i]++;
230           ptr[i]++;
231           
232         }
233         q/=16.;  //conversion factor
234         //      Float_t noise  = gRandom->Gaus(0,param->GetNoise()*param->GetNoiseNormFac());  
235         Float_t noise  = pTPC->GetNoise();
236         q+=noise;
237         q=TMath::Nint(q);
238         if (q > zerosup){ 
239           
240           if(q > param->GetADCSat()) q = (Short_t)(param->GetADCSat());
241           //digrow->SetDigitFast((Short_t)q,rows,col);  
242           *pdig1 =Short_t(q);
243           for (Int_t tr=0;tr<3;tr++){
244             if (tr<labptr) 
245               // ((AliSimDigits*)digrow)->SetTrackIDFast(label[tr],rows,col,tr);
246               ptr1[tr*nElems] = label[tr];
247             //else
248               //            ((AliSimDigits*)digrow)->SetTrackIDFast(-1,rows,col,tr);          
249             //  ptr1[tr*nElems] = 1;
250           }
251         }
252         pdig1++;
253         ptr1++;
254     }
255     
256     digrow->CompresBuffer(1,zerosup);
257     digrow->CompresTrackBuffer(1);
258     tree->Fill();
259     if (fDebug>0) cerr<<sec<<"\t"<<row<<"\n";  
260   } 
261   fManager->GetTreeDTPC()->Write(0,TObject::kOverwrite);
262   delete digrow;     
263   for (Int_t i1=0;i1<nInputs; i1++) delete digarr[i1];
264   delete []masks;
265   delete digarr;  
266 }
267
268
269
270 //------------------------------------------------------------------------
271 void AliTPCDigitizer::ExecSave(Option_t* option)
272 {
273   
274   // merge input tree's with summable digits
275   //output stored in TreeTPCD
276
277   TString optionString = option;
278   if (optionString.Data() == "deb") {
279     cout<<"AliTPCDigitizer::Exec: called with option deb "<<endl;
280     fDebug = 3;
281   }
282   //get detector and geometry 
283   printf("TPC merging -1  -Tree %s\t%p\n",fManager->GetInputTreeH(0)->GetName(),fManager->GetInputTreeH(0)->GetListOfBranches()->At(3));
284   AliTPC *pTPC  = (AliTPC *) gAlice->GetModule("TPC");
285   AliTPCParam * param = pTPC->GetParam();
286   pTPC->GenerNoise(500000); //create teble with noise
287   printf("noise %f \n",  param->GetNoise()*param->GetNoiseNormFac());
288   //
289   Int_t nInputs = fManager->GetNinputs();
290   Int_t * masks = new Int_t[nInputs];
291   for (Int_t i=0; i<nInputs;i++)
292     masks[i]= fManager->GetMask(i);
293   //  Short_t **pdig= new Short_t*[nInputs];   //pointers to the expanded digits array
294   //Int_t **ptr=  new Int_t*[nInputs];       //pointers to teh expanded tracks array
295
296   //create digits array for given sectors
297   // make indexes
298    printf("TPC merging -2  -Tree %s\t%p\n",fManager->GetInputTreeH(0)->GetName(),fManager->GetInputTreeH(0)->GetListOfBranches()->At(3));
299   AliSimDigits ** digarr = new AliSimDigits*[nInputs]; 
300   for (Int_t i1=0;i1<nInputs; i1++){
301     digarr[i1]=0;
302     //    intree[i1]
303     TTree * treear =  fManager->GetInputTreeTPCS(i1);
304     printf("TPC merging -2.7  -Tree %s\t%p\n",fManager->GetInputTreeH(0)->GetName(),fManager->GetInputTreeH(0)->GetListOfBranches()->At(3));
305     TBranch * br = treear->GetBranch("fSegmentID");
306     if (br) br->GetFile()->cd();
307     printf("TPC merging -2.75  -Tree %s\t%p\n",fManager->GetInputTreeH(0)->GetName(),fManager->GetInputTreeH(0)->GetListOfBranches()->At(3));
308     //treear->BuildIndex("fSegmentID","fSegmentID");
309     if (!treear) {      
310       cerr<<" TPC -  not existing input = \n"<<i1<<" ";      
311     } 
312     printf("TPC merging -2.8  -Tree %s\t%p\n",fManager->GetInputTreeH(0)->GetName(),fManager->GetInputTreeH(0)->GetListOfBranches()->At(3));
313     treear->GetBranch("Segment")->SetAddress(&digarr[i1]);
314      printf("TPC merging -2.9  -Tree %s\t%p\n",fManager->GetInputTreeH(0)->GetName(),fManager->GetInputTreeH(0)->GetListOfBranches()->At(3));
315   }
316   Stat_t nentries = fManager->GetInputTreeTPCS(0)->GetEntries();
317   
318
319   //create branch's in TPC treeD
320   AliSimDigits * digrow = new AliSimDigits;
321   TTree * tree  = fManager->GetTreeDTPC();
322   //if (tree->GetBranch("Segment") ) tree->GetBranch("Segment")->SetAddress(&digrow);
323   //else
324   tree->Branch("Segment","AliSimDigits",&digrow);
325   //
326  printf("TPC merging -3  -Tree %s\t%p\n",fManager->GetInputTreeH(0)->GetName(),fManager->GetInputTreeH(0)->GetListOfBranches()->At(3));
327   param->SetZeroSup(2);
328
329   Int_t zerosup = param->GetZeroSup();
330   //Loop over segments of the TPC
331     
332   for (Int_t n=0; n<nentries; n++) {
333   //    for (Int_t n=0; n<300; n++) {
334     fManager->GetInputTreeTPCS(0)->GetEvent(n);      
335     digarr[0]->ExpandBuffer();
336     digarr[0]->ExpandTrackBuffer();
337            
338     for (Int_t i=1;i<nInputs; i++){ 
339       fManager->GetInputTreeTPCS(i)->GetEntryWithIndex(digarr[0]->GetID(),digarr[0]->GetID());      
340       //fManager->GetInputTreeTPCS(i)->GetEntryWithIndex(digarr[0]->GetID(),1);      
341       digarr[i]->ExpandBuffer();
342       digarr[i]->ExpandTrackBuffer();
343       if ((digarr[0]->GetID()-digarr[i]->GetID())>0) 
344         printf("problem\n");
345     
346     }   
347     
348     Int_t sec, row;
349     if (!param->AdjustSectorRow(digarr[0]->GetID(),sec,row)) {
350       cerr<<"AliTPC warning: invalid segment ID ! "<<digarr[0]->GetID()<<endl;
351       continue;
352     }
353
354     digrow->SetID(digarr[0]->GetID());
355
356     Int_t nrows = digarr[0]->GetNRows();
357     Int_t ncols = digarr[0]->GetNCols();
358     digrow->Allocate(nrows,ncols);
359     digrow->AllocateTrack(3);
360
361     Float_t q=0;
362     Int_t label[1000]; //stack for 300 events 
363     Int_t labptr = 0;
364
365     
366
367     for (Int_t rows=0;rows<nrows; rows++){
368       for (Int_t col=0;col<ncols; col++){
369     
370         q=0;
371         labptr=0;
372         // looop over digits 
373         for (Int_t i=0;i<nInputs; i++){ 
374           q  += digarr[i]->GetDigitFast(rows,col);
375           //q  += *(pdig[i]);
376           
377           for (Int_t tr=0;tr<3;tr++) {
378             Int_t lab = digarr[i]->GetTrackIDFast(rows,col,tr);
379             //Int_t lab = ptr[i][tr*nElems];
380             if ( (lab > 1) ) {
381               label[labptr]=lab+masks[i];
382               labptr++;
383             }      
384           }
385           // pdig[i]++;
386           //ptr[i]++;
387           
388         }
389         q/=16.;  //conversion factor
390         //      Float_t noise  = gRandom->Gaus(0,param->GetNoise()*param->GetNoiseNormFac());  
391         Float_t noise  = pTPC->GetNoise();
392         q+=noise;
393         q=TMath::Nint(q);
394         if (q > zerosup){ 
395           
396           if(q > param->GetADCSat()) q = (Short_t)(param->GetADCSat());
397           digrow->SetDigitFast((Short_t)q,rows,col);  
398           // *pdig1 =Short_t(q);
399           for (Int_t tr=0;tr<3;tr++){
400             if (tr<labptr) 
401               ((AliSimDigits*)digrow)->SetTrackIDFast(label[tr],rows,col,tr);
402             //ptr1[tr*nElems] = label[tr];
403             //else
404               //            ((AliSimDigits*)digrow)->SetTrackIDFast(-1,rows,col,tr);          
405             //  ptr1[tr*nElems] = 1;
406           }
407         }
408         //pdig1++;
409         //ptr1++;
410     }
411     }
412     
413     digrow->CompresBuffer(1,zerosup);
414     digrow->CompresTrackBuffer(1);
415     tree->Fill();
416     if (fDebug>0) cerr<<sec<<"\t"<<row<<"\n";  
417   } 
418   printf("end TPC merging - end -Tree %s\t%p\n",fManager->GetInputTreeH(0)->GetName(),fManager->GetInputTreeH(0)->GetListOfBranches()->At(3));
419   fManager->GetTreeDTPC()->Write(0,TObject::kOverwrite);
420   delete digrow;     
421   for (Int_t i1=0;i1<nInputs; i1++) delete digarr[i1];
422   delete []masks;
423   delete digarr;  
424 }