]> git.uio.no Git - u/mrichter/AliRoot.git/blob - PHOS/AliPHOSSDigitizer.cxx
Do not write Tree's but do an AutoSave() to fit in the framework
[u/mrichter/AliRoot.git] / PHOS / AliPHOSSDigitizer.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, 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 /* $Id$ */
17
18 //_________________________________________________________________________
19 // This is a TTask that makes SDigits out of Hits
20 // The name of the TTask is also the title of the branch that will contain 
21 // the created SDigits
22 // The title of the TTAsk is the name of the file that contains the hits from
23 // which the SDigits are created
24 // A Summable Digits is the sum of all hits originating 
25 // from one primary in one active cell
26 // A threshold for assignment of the primary to SDigit is applied 
27 // SDigits are written to TreeS, branch "PHOS"
28 // AliPHOSSDigitizer with all current parameters is written 
29 // to TreeS branch "AliPHOSSDigitizer".
30 // Both branches have the same title. If necessary one can produce 
31 // another set of SDigits with different parameters. Two versions
32 // can be distunguished using titles of the branches.
33 // User case:
34 //  root [0] AliPHOSSDigitizer * s = new AliPHOSSDigitizer("galice.root")
35 //  Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
36 //  root [1] s->ExecuteTask()
37 //             // Makes SDigitis for all events stored in galice.root
38 //  root [2] s->SetPedestalParameter(0.001)
39 //             // One can change parameters of digitization
40 // root [3] s->SetSDigitsBranch("Pedestal 0.001")
41 //             // and write them into the new branch
42 // root [4] s->ExecuteTask("deb all tim")
43 //             // available parameters:
44 //             deb - print # of produced SDigitis
45 //             deb all  - print # and list of produced SDigits
46 //             tim - print benchmarking information
47 //
48 //*-- Author :  Dmitri Peressounko (SUBATECH & KI) 
49 //////////////////////////////////////////////////////////////////////////////
50
51
52 // --- ROOT system ---
53 #include "TFile.h"
54 #include "TTask.h"
55 #include "TTree.h"
56 #include "TSystem.h"
57 #include "TROOT.h"
58 #include "TFolder.h"
59 #include "TBenchmark.h"
60 // --- Standard library ---
61 #include <iomanip.h>
62
63 // --- AliRoot header files ---
64 #include "AliRun.h"
65 #include "AliPHOSDigit.h"
66 #include "AliPHOSGeometry.h"
67 #include "AliPHOSGetter.h"
68 #include "AliPHOSHit.h"
69 #include "AliPHOSSDigitizer.h"
70
71
72 ClassImp(AliPHOSSDigitizer)
73
74            
75 //____________________________________________________________________________ 
76   AliPHOSSDigitizer::AliPHOSSDigitizer():TTask("","") 
77 {
78   // ctor
79   fA             = 0;
80   fB             = 10000000.;
81   fPrimThreshold = 0.01 ;
82   fSDigitsInRun  = 0 ; 
83 }
84
85 //____________________________________________________________________________ 
86 AliPHOSSDigitizer::AliPHOSSDigitizer(const char * headerFile, const char * sDigitsTitle):TTask(sDigitsTitle, headerFile)
87 {
88   // ctor
89   fA             = 0;
90   fB             = 10000000.;
91   fPrimThreshold = 0.01 ;
92   fSDigitsInRun  = 0 ; 
93   Init();
94 }
95
96
97 //____________________________________________________________________________ 
98 void AliPHOSSDigitizer::Init()
99 {
100   // Initialization: open root-file, allocate arrays for hits and sdigits,
101   // attach task SDigitizer to the list of PHOS tasks
102   // 
103   // Initialization can not be done in the default constructor
104   //============================================================= YS
105   //  The initialisation is now done by AliPHOSGetter
106   
107   if( strcmp(GetTitle(), "") == 0 )
108     SetTitle("galice.root") ;
109   
110   AliPHOSGetter * gime = AliPHOSGetter::GetInstance(GetTitle(), GetName()) ;     
111   if ( gime == 0 ) {
112     cerr << "ERROR: AliPHOSSDigitizer::Init -> Could not obtain the Getter object !" << endl ; 
113     return ;
114   } 
115   
116   gime->PostSDigits( GetName(), GetTitle() ) ; 
117
118   TString sdname(GetName() );
119   sdname.Append(":") ;
120   sdname.Append(GetTitle() ) ;
121   SetName(sdname) ;
122   gime->PostSDigitizer(this) ;
123      
124 }
125
126 //____________________________________________________________________________
127 void AliPHOSSDigitizer::Exec(Option_t *option) 
128
129   // Collects all hits in the same active volume into digit
130   if( strcmp(GetName(), "") == 0 )
131     Init() ;
132   
133   if (strstr(option, "print") ) {
134     Print("") ; 
135     return ; 
136   }
137
138   if(strstr(option,"tim"))
139     gBenchmark->Start("PHOSSDigitizer");
140
141   //Check, if this branch already exits
142   gAlice->GetEvent(0) ;
143
144   TString sdname(GetName()) ;
145   sdname.Remove(sdname.Index(GetTitle())-1) ;
146   
147   if(gAlice->TreeS() ) {
148     TObjArray * lob = static_cast<TObjArray*>(gAlice->TreeS()->GetListOfBranches()) ;
149     TIter next(lob) ; 
150     TBranch * branch = 0 ;  
151     Bool_t phosfound = kFALSE, sdigitizerfound = kFALSE ; 
152     
153     while ( (branch = (static_cast<TBranch*>(next()))) && (!phosfound || !sdigitizerfound) ) {
154       if ( (strcmp(branch->GetName(), "PHOS")==0) && (strcmp(branch->GetTitle(), GetName())==0) ) 
155         phosfound = kTRUE ;
156       
157       else if ( (strcmp(branch->GetName(), "AliPHOSSDigitizer")==0) && (strcmp(branch->GetTitle(), GetName())==0) ) 
158         sdigitizerfound = kTRUE ; 
159     }
160     
161     if ( phosfound || sdigitizerfound ) {
162       cerr << "WARNING: AliPHOSSDigitizer::Exec -> SDigits and/or SDigitizer branch with name " << GetName() 
163            << " already exits" << endl ;
164       return ; 
165     }   
166   }  
167
168   AliPHOSGetter * gime = AliPHOSGetter::GetInstance() ; 
169
170   Int_t nevents = (Int_t) gAlice->TreeE()->GetEntries() ; 
171   
172   Int_t ievent ;
173   for(ievent = 0; ievent < nevents; ievent++){
174     gime->Event(ievent,"H") ;
175     const TClonesArray * hits = gime->Hits() ;
176     TClonesArray * sdigits = gime->SDigits(sdname.Data()) ;
177     sdigits->Clear();
178     Int_t nSdigits = 0 ;
179     
180     
181     //Now make SDigits from hits, for PHOS it is the same, so just copy    
182
183 //******************** CHECK HERE
184 //YS DOES NOT UNDERSTAND THE NEED FOR THE FOLLOWING LOOP
185 //     Int_t itrack ;
186 //     for (itrack=0; itrack < gAlice->GetNtrack(); itrack++){
187 //       //=========== Get the PHOS branch from Hits Tree for the Primary track itrack
188 //       gime->Track(itrack) ;
189       Int_t i;
190       for ( i = 0 ; i < hits->GetEntries() ; i++ ) {
191         AliPHOSHit * hit = dynamic_cast<AliPHOSHit *>(hits->At(i)) ;
192         // Assign primary number only if contribution is significant
193         
194         if( hit->GetEnergy() > fPrimThreshold)
195           new((*sdigits)[nSdigits]) AliPHOSDigit(hit->GetPrimary(),hit->GetId(),
196                                                   Digitize(hit->GetEnergy()), hit->GetTime()) ;
197         else
198           new((*sdigits)[nSdigits]) AliPHOSDigit( -1              , hit->GetId(), 
199                                                    Digitize(hit->GetEnergy()), hit->GetTime()) ;
200         nSdigits++ ;    
201         
202       }
203 //    } // loop over tracks
204     
205     sdigits->Sort() ;
206     
207     nSdigits = sdigits->GetEntriesFast() ;
208     fSDigitsInRun += nSdigits ;  
209     sdigits->Expand(nSdigits) ;
210 //    Int_t i ;
211     for (i = 0 ; i < nSdigits ; i++) { 
212       AliPHOSDigit * digit = dynamic_cast<AliPHOSDigit *>(sdigits->At(i)) ; 
213       digit->SetIndexInList(i) ;     
214     }
215
216     if(gAlice->TreeS() == 0)
217       gAlice->MakeTree("S") ;
218     
219     //Make (if necessary) branches    
220     char * file =0;
221     if(gSystem->Getenv("CONFIG_SPLIT_FILE")){ //generating file name
222       file = new char[strlen(gAlice->GetBaseFile())+20] ;
223       sprintf(file,"%s/PHOS.SDigits.root",gAlice->GetBaseFile()) ;
224     }
225     
226     TDirectory *cwd = gDirectory;
227     
228     //First list of sdigits
229     Int_t bufferSize = 32000 ;    
230     TBranch * sdigitsBranch = gAlice->TreeS()->Branch("PHOS",&sdigits,bufferSize);
231     sdigitsBranch->SetTitle(sdname);
232     if (file) {
233       sdigitsBranch->SetFile(file);
234       TIter next( sdigitsBranch->GetListOfBranches());
235       TBranch * subbr;
236       while ((subbr=static_cast<TBranch*>(next()))) {
237         subbr->SetFile(file);
238       }   
239       cwd->cd();
240     } 
241       
242     //Next - SDigitizer
243     Int_t splitlevel = 0 ;
244     AliPHOSSDigitizer * sd = this ;
245     TBranch * sdigitizerBranch = gAlice->TreeS()->Branch("AliPHOSSDigitizer","AliPHOSSDigitizer",
246                                                &sd,bufferSize,splitlevel); 
247     sdigitizerBranch->SetTitle(sdname);
248     if (file) {
249       sdigitizerBranch->SetFile(file);
250       TIter next( sdigitizerBranch->GetListOfBranches());
251       TBranch * subbr ;
252       while ((subbr=static_cast<TBranch*>(next()))) {
253         subbr->SetFile(file);
254       }   
255       cwd->cd();
256       delete file;
257     }
258
259     sdigitsBranch->Fill() ; 
260     gAlice->TreeS()->AutoSave() ;
261     
262     if(strstr(option,"deb"))
263       PrintSDigits(option) ;
264     
265   }
266   
267   if(strstr(option,"tim")){
268     gBenchmark->Stop("PHOSSDigitizer");
269     cout << "AliPHOSSDigitizer:" << endl ;
270     cout << "   took " << gBenchmark->GetCpuTime("PHOSSDigitizer") << " seconds for SDigitizing " 
271          <<  gBenchmark->GetCpuTime("PHOSSDigitizer")/nevents << " seconds per event " << endl ;
272     cout << endl ;
273   }
274   
275   
276 }
277 //__________________________________________________________________
278 void AliPHOSSDigitizer::SetSDigitsBranch(const char * title )
279 {
280   // Setting title to branch SDigits 
281
282   TString stitle(title) ;
283
284   // check if branch with title already exists
285   TBranch * sdigitsBranch    = 
286     static_cast<TBranch*>(gAlice->TreeS()->GetListOfBranches()->FindObject("PHOS")) ; 
287   TBranch * sdigitizerBranch =  
288     static_cast<TBranch*>(gAlice->TreeS()->GetListOfBranches()->FindObject("AliPHOSSDigitizer")) ;
289   const char * sdigitsTitle    = sdigitsBranch ->GetTitle() ;  
290   const char * sdigitizerTitle = sdigitizerBranch ->GetTitle() ;
291   if ( stitle.CompareTo(sdigitsTitle)==0 || stitle.CompareTo(sdigitizerTitle)==0 ){
292     cerr << "ERROR: AliPHOSSdigitizer::SetSDigitsBranch -> Cannot overwrite existing branch with title " << title << endl ;
293     return ;
294   }
295   
296   cout << "AliPHOSSdigitizer::SetSDigitsBranch -> Changing SDigits file from " << GetName() << " to " << title << endl ;
297
298   SetName(title) ; 
299     
300   // Post to the WhiteBoard
301   AliPHOSGetter * gime = AliPHOSGetter::GetInstance() ; 
302   gime->PostSDigits( title, GetTitle()) ; 
303 }
304
305 //__________________________________________________________________
306 void AliPHOSSDigitizer::Print(Option_t* option)const
307 {
308   // Prints parameters of SDigitizer
309   cout << "------------------- "<< GetName() << " -------------" << endl ;
310   cout << "   Writing SDigits to branch with title  " << GetName() << endl ;
311   cout << "   with digitization parameters  A = " << fA << endl ;
312   cout << "                                 B = " << fB << endl ;
313   cout << "   Threshold for Primary assignment= " << fPrimThreshold << endl ; 
314   cout << "---------------------------------------------------"<<endl ;
315   
316 }
317 //__________________________________________________________________
318 Bool_t AliPHOSSDigitizer::operator==( AliPHOSSDigitizer const &sd )const
319 {
320   // Equal operator.
321   // SDititizers are equal if their pedestal, slope and threshold are equal
322
323   if( (fA==sd.fA)&&(fB==sd.fB)&&(fPrimThreshold==sd.fPrimThreshold))
324     return kTRUE ;
325   else
326     return kFALSE ;
327 }
328 //__________________________________________________________________
329 //__________________________________________________________________
330 void AliPHOSSDigitizer::PrintSDigits(Option_t * option)
331 {
332   // Prints list of digits produced in the current pass of AliPHOSDigitizer
333
334
335   AliPHOSGetter * gime = AliPHOSGetter::GetInstance() ; 
336   TString sdname(GetName()) ;
337   sdname.Remove(sdname.Index(GetTitle())-1) ;
338   const TClonesArray * sdigits = gime->SDigits(sdname.Data()) ; 
339
340   cout << "AliPHOSSDigitiser: event " << gAlice->GetEvNumber() << endl ;
341   cout << "      Number of entries in SDigits list " << sdigits->GetEntriesFast() << endl ;
342   cout << endl ;
343   if(strstr(option,"all")||strstr(option,"EMC")){
344     
345     //loop over digits
346     AliPHOSDigit * digit;
347     cout << "EMC sdigits " << endl ;
348     cout << "Digit Id    Amplitude     Index     Nprim  Primaries list " <<  endl;      
349     Int_t maxEmc = gime->PHOSGeometry()->GetNModules()*gime->PHOSGeometry()->GetNCristalsInModule() ;
350     Int_t index ;
351     for (index = 0 ; (index < sdigits->GetEntriesFast()) && 
352          (((AliPHOSDigit * )  sdigits->At(index))->GetId() <= maxEmc) ; index++) {
353       digit = (AliPHOSDigit * )  sdigits->At(index) ;
354       if(digit->GetNprimary() == 0) continue;
355       cout << setw(6)  <<  digit->GetId() << "   "  <<  setw(10)  <<  digit->GetAmp() <<   "    "  
356            << setw(6)  <<  digit->GetIndexInList() << "    "   
357            << setw(5)  <<  digit->GetNprimary() <<"    ";
358       
359       Int_t iprimary;
360       for (iprimary=0; iprimary<digit->GetNprimary(); iprimary++)
361         cout << setw(5)  <<  digit->GetPrimary(iprimary+1) << "    ";
362       cout << endl;      
363     }    
364     cout << endl;
365   }
366
367   if(strstr(option,"all")||strstr(option,"CPV")){
368     
369     //loop over CPV digits
370     AliPHOSDigit * digit;
371     cout << "CPV sdigits " << endl ;
372     cout << "Digit Id    Amplitude     Index     Nprim  Primaries list " <<  endl;      
373     Int_t maxEmc = gime->PHOSGeometry()->GetNModules()*gime->PHOSGeometry()->GetNCristalsInModule() ;
374     Int_t index ;
375     for (index = 0 ; index < sdigits->GetEntriesFast(); index++) {
376       digit = (AliPHOSDigit * )  sdigits->At(index) ;
377       if(digit->GetId() > maxEmc){
378         cout << setw(6)  <<  digit->GetId() << "   "  <<        setw(10)  <<  digit->GetAmp() <<   "    "  
379              << setw(6)  <<  digit->GetIndexInList() << "    "   
380              << setw(5)  <<  digit->GetNprimary() <<"    ";
381         
382         Int_t iprimary;
383         for (iprimary=0; iprimary<digit->GetNprimary(); iprimary++)
384           cout << setw(5)  <<  digit->GetPrimary(iprimary+1) << "    ";
385         cout << endl;    
386       }    
387     }
388   }
389
390 }
391
392 //____________________________________________________________________________ 
393 void AliPHOSSDigitizer::UseHitsFrom(const char * filename)
394 {
395   SetTitle(filename) ; 
396   Init() ; 
397 }