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