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