0a1afaa610f1afbeebdb2a89b2d809eb36334525
[u/mrichter/AliRoot.git] / PHOS / testsimglobal.C
1 #include "AliPHOSGetter.h"
2 #include "TH2.h"
3 #include "TH1.h"
4 #include "TFile.h"
5 #include "TTree.h"
6 #include "TBranch.h"
7 #include "TClonesArray.h"
8 #include "TCanvas.h"
9 #include "TSystem.h"
10 #include "AliPHOSHit.h"
11 #include "TF1.h"
12 #include "TFormula.h"
13 #include "TFolder.h"
14 #include "TStopwatch.h"
15 #include "TObjArray.h"
16 #include "AliPHOSGeometry.h"
17 #include "AliPHOSDigit.h"
18 #include "AliPHOSSDigitizer.h"
19 #include "AliPHOSDigitizer.h"
20 #include "AliPHOSClusterizer.h"
21 #include "AliPHOSClusterizerv1.h"
22 #include "AliPHOSTrackSegmentMaker.h"
23 #include "AliPHOSTrackSegmentMakerv1.h"
24 #include "AliPHOSPID.h"
25 #include "AliPHOSPIDv1.h"
26
27
28 void testsimglobal(Int_t nevent = 1, const char *config="testconfig.C")
29 {
30   //
31   // Simple macro to run aliroot in a batch mode
32   //
33   cerr<<" __________________________________________________________________ "<<endl;
34   cerr<< " " <<endl;
35   cerr<<"             MESS ==> Beginning of the simulation examination."<<endl;
36   cerr<<" __________________________________________________________________ "<<endl;
37   // Definition of the alarm bounds
38      const Float_t maxAlaHitsM = 12.79  ;  // total multiplicity
39      const Float_t maxAlaTotEn = 19.34  ;  // total energy
40      const Float_t maxAlaTotEnB = 18.35 ;  // per block multiplicity
41      const Float_t maxAlaHitsMB = 11.1  ;  // per block energy
42    
43        // boolean which will test if there is an alarm
44    Int_t boolerror = 0;
45    
46    // define the array in which the events that have not reach th EMCA will be put.
47   
48    
49
50    TStopwatch timer;
51    timer.Start();
52
53    // Get the number of events generated in the simulation
54
55   AliPHOSGetter * gime = AliPHOSGetter::GetInstance("testPHOS.root") ; 
56   Int_t maxevent = gime->MaxEvent() ; 
57   gime->Event(0,"Q");
58  
59   // Examine the alarms
60
61    TObjArray * alahm = (TObjArray*)(gime->Alarms()->FindObject("HitsM"));
62    Float_t ratiohm = 100.0*(Float_t)alahm->GetEntries()/(Float_t)maxevent ;  
63  
64   
65    TObjArray * alaet = (TObjArray*)(gime->Alarms()->FindObject("TotEn"));
66    Float_t ratioet = 100.0*(Float_t)alaet->GetEntries()/(Float_t)maxevent ; 
67   
68
69    // Define the alarms per block and examine them
70    char namemul[80], namen[80];
71    TObjArray* alahmb[5];
72    TObjArray* alaenb[5];
73    Float_t ratiohmb[5], ratioenb[5];
74    for (Int_t i = 0 ; i < 5 ; i++)
75      {
76         sprintf(namemul,"HitsMB%d",i+1);
77         sprintf(namen,"TotEnB%d",i+1);
78         alahmb[i] = (TObjArray*) (gime->Alarms()->FindObject(namemul));
79         ratiohmb[i] = 100.0*(Float_t)alahmb[i]->GetEntries()/(Float_t)maxevent;        
80         alaenb[i] = (TObjArray*)(gime->Alarms()->FindObject(namen));
81         ratioenb[i] = 100.0*(Float_t)alaenb[i]->GetEntries()/(Float_t)maxevent;      
82         if (ratiohmb[i]>maxAlaHitsMB){
83           boolerror = 1; 
84           cerr<<" _____________________________________________________________ "<<endl;
85           cerr<< " " <<endl;
86           cerr << "             MESS ==> Examination detected an error in "<<namemul << endl ;    
87           cerr<<" _____________________________________________________________ "<<endl;
88         }
89         if (ratioenb[i]>maxAlaTotEnB) {
90           boolerror = 1;
91           cerr<<" _____________________________________________________________ "<<endl;
92           cerr<< " " <<endl;
93           cerr << "             MESS ==> Examination detected an error in "<<namen << endl ;
94           cerr<<" _____________________________________________________________ "<<endl;
95         }
96             
97      }
98  
99
100   timer.Stop();
101   timer.Print();
102
103       
104   if (ratiohm>maxAlaHitsM){
105      boolerror = 1;
106      cerr<<" _____________________________________________________________ "<<endl;
107      cerr<< " " <<endl;
108      cerr << "             MESS ==> Examination detected an error in HitsM." << endl ;
109      cerr<<" _____________________________________________________________ "<<endl;
110   }
111
112   if (ratioet>maxAlaTotEn){
113      boolerror = 1;
114      cerr<<" _____________________________________________________________ "<<endl;
115      cerr<< " " <<endl;
116      cerr << "             MESS ==> Examination detected an error in TotEn." << endl ;
117      cerr<<" _____________________________________________________________ "<<endl;
118   }
119
120  
121   // Condition that will launch the general loop that builds the histograms in order to allow a further analysis.
122   
123     if ( boolerror == 1 ) {
124   cerr<<" _____________________________________________________________ "<<endl;
125   cerr<< " " <<endl;
126   cerr << "             MESS ==> Examination sets up the file that will be sent to PHOS director (30s). " << endl ;
127   cerr<<" _____________________________________________________________ "<<endl;    
128       Int_t index = 0 ; 
129       Int_t nhits = 0 ; 
130     
131       TH1F * his = new TH1F("Total Multiplicity", "Total Multiplicity in PHOS", 200, 0., 200.) ;
132       TH1F * hisnrg = new TH1F("Total Energy", "Total energy in PHOS", 200, 0., 12.) ;
133     
134      // name and define the different histograms per block involved in the analysis
135       TClonesArray hisba("TH1F") ;
136       TClonesArray hisbanrg("TH1F") ;  
137       for (Int_t i = 0 ; i < 6 ; i++)
138        {
139          char name[80], title[80] ; 
140          sprintf(name,"multiplicity for block %d",i) ; 
141          sprintf(title,"multiplicity per blocks, block %d",i) ;
142          new(hisba[i]) TH1F(name, title, 100, 0., 100.) ; 
143       
144
145          char namenrg[80], titlenrg[80] ; 
146          sprintf(namenrg,"total energy for block %d",i) ; 
147          sprintf(titlenrg,"total energy per block, block %d",i) ;           
148          new(hisbanrg[i]) TH1F(namenrg, titlenrg, 200, 0., 12.) ; 
149        }
150       // define the global histograms, the hit pointer and give the means to get the actual block reached by the photon
151       AliPHOSHit * hit;
152       TH1F * hist ; 
153       TH1F * histnrg;
154       TH2F * hbiz = new TH2F ("hbiz","hbiz",200.,0.,200.,200.,0.,12.);
155       AliPHOSGeometry * geom = AliPHOSGeometry::GetInstance() ; 
156
157
158       // the very general loop
159       for (index = 0 ; index < maxevent ; index++)
160        {
161          //get the number of the event
162          gime->Event(index,"H") ; 
163          // get the number of cells reached during this event and fill the total multiplicity histogram
164          Int_t n = gime->Hits()->GetEntries() ;
165          nhits += n ; 
166          his->Fill(n) ;  
167          // Get the data per block      
168          TClonesArray * hita = (TClonesArray *) gime -> Hits();
169          TIter next(hita);
170          Float_t Et = 0.;
171          Int_t id = 0, block = 0;
172          Int_t nhit[6], rid[4];
173          Float_t etblock[6]; 
174          for ( Int_t i = 0; i < 6 ; i++) { 
175            nhit[i] = 0 ; 
176            etblock[i] = 0 ;
177          }
178      
179          while ( (hit = (AliPHOSHit *) next()) ) {
180                  id = hit->GetId();
181                  if  (geom->IsInEMC(id) ) { 
182                    Et += ( hit -> GetEnergy());
183                    geom->AbsToRelNumbering(id,rid) ;
184                    block = rid[0];
185                    nhit[block]++  ;
186                    etblock[block] +=  ( hit -> GetEnergy());
187                  }
188                 }
189
190
191          //Fill all the histograms but total multiplicity, already done
192          hist = static_cast<TH1F*>(hisba.At(block)) ; 
193          hist->Fill(nhit[block]) ; 
194          histnrg = static_cast<TH1F*>(hisbanrg.At(block)) ;
195          histnrg->Fill(etblock[block]);
196          hisnrg -> Fill(Et);
197          hbiz->Fill(n,Et);               
198        }   
199         
200       nhits /= maxevent ; 
201       cerr << "av = " << nhits << endl ;
202       TFile * file = gROOT -> GetFile("testPHOS.root");
203       file -> Write();
204       his->Draw() ;
205       hisnrg->Draw() ;
206       hbiz->Draw();
207
208       //Put the histograms in the root file
209       for (Int_t i = 0 ; i < 6 ; i++)   { 
210          hist = static_cast<TH1F*>(hisba.At(i)) ;
211          histnrg = static_cast<TH1F*>(hisbanrg.At(i));
212      
213         
214          // cout << hist << endl << i << endl ; 
215          hist->Draw();
216          histnrg->Draw(); 
217          
218    
219        }
220       
221       file->Write() ;
222       hisba.Delete() ;
223       hisbanrg.Delete() ; 
224       file->Close();
225         gSystem->Exec("uuencode $ALICE_ROOT/PHOS/testPHOS.root testPHOS.root | mail -s 'PHOS INSTALLATION ERROR  ' almajor@caramail.com"); 
226    
227     }
228   cerr<<" _____________________________________________________________ "<<endl;
229   cerr<< " " <<endl;
230   cerr << "             MESS ==> Examination ended successfully. " << endl ;
231   cerr<<" _____________________________________________________________ "<<endl;   
232
233
234
235   const Float_t maxSDigits = 62.89 ;
236   const Float_t widSDigits = TMath::Sqrt(maxSDigits) ;
237   const Float_t maxDigits = 3489.41 ;
238   const Float_t widDigits = TMath::Sqrt(maxDigits) ;
239   const Float_t maxRecPoints = 222.83 ;
240   const Float_t widRecPoints = TMath::Sqrt(maxRecPoints) ;
241   const Float_t maxTrackSegments = 1 ;
242   const Float_t maxRecParticles = 1 ;
243   TString reconame = "test suite" ;
244   if (boolerror == 1){
245   cerr<<" ___________________________________________________________________ "<<endl;
246   cerr<<" "<<endl;
247   cerr<<"             MESS ==> Please find the error and try the whole simulation again. "<<endl;
248   cerr<<" ___________________________________________________________________ "<<endl;
249   }
250
251   else {
252
253   cerr<<" ___________________________________________________________________ "<<endl;
254   cerr<<" "<<endl;
255   cerr<<"             MESS ==> Beginning of the PHOS reconstruction. "<<endl;
256   cerr<<" ___________________________________________________________________ "<<endl;
257
258   //SDigits process
259    AliPHOSSDigitizer *sd = new AliPHOSSDigitizer("testPHOS.root",reconame.Data());
260    AliPHOSGetter * gime = AliPHOSGetter::GetInstance() ;
261    sd->ExecuteTask("deb"); 
262    Float_t nSDigits =  (Float_t) (gime->SDigitizer()->GetSDigitsInRun()) / gime->MaxEvent(); 
263    if ( nSDigits < maxSDigits-widSDigits || nSDigits > maxSDigits+widSDigits ) {
264     boolerror = 1 ;  
265     cerr<<"__________________________________________________________________"<<endl;
266     cerr<<" "<<endl;
267     cerr<<"             MESS ==> Error detected in the SDigits process. Sending error file to PHOS director."<<endl;
268     cerr<<"__________________________________________________________________"<<endl;
269     gSystem->Exec("uuencode $ALICE_ROOT/PHOS/testPHOS.root testPHOS.root | mail -s 'PHOS INSTALLATION ERROR' almajor@caramail.com");
270     }
271
272   cerr<<"__________________________________________________________________"<<endl;
273   cerr<<" "<<endl;
274   cerr<<"             MESS ==> SDigits process ended successfully."<<endl;
275   cerr<<"__________________________________________________________________"<<endl;
276
277   //Digits process
278      if (boolerror == 0){
279         AliPHOSDigitizer *d = new AliPHOSDigitizer("testPHOS.root",reconame.Data());
280         d->ExecuteTask("deb"); 
281          Float_t nDigits = (Float_t)(gime->Digitizer()->GetDigitsInRun()) / gime->MaxEvent();
282   
283               if ( nDigits < maxDigits-widDigits || nDigits > maxDigits+widDigits ) {
284                   boolerror =1 ;
285                   cerr<<"__________________________________________________________________"<<endl;
286                   cerr<<" "<<endl;
287                   cerr<<"             MESS ==> Error detected in the Digits process. Sending error file to PHOS director."<<endl;
288                   cerr<<"__________________________________________________________________"<<endl;
289    
290     gSystem->Exec("uuencode $ALICE_ROOT/PHOS/testPHOS.root testPHOS.root | mail -s 'PHOS INSTALLATION ERROR' almajor@caramail.com");
291               }
292   
293     cerr<<"__________________________________________________________________"<<endl;
294     cerr<<" "<<endl;
295     cerr<<"             MESS ==> Digits process ended successfully."<<endl;
296     cerr<<"__________________________________________________________________"<<endl;
297   }
298   
299 //RecPoints process
300  if (boolerror == 0){
301  AliPHOSClusterizer * cluster = new  AliPHOSClusterizerv1("testPHOS.root", reconame.Data());
302  
303   cluster->ExecuteTask("deb");
304   Float_t nRecPoints =  (Float_t) (gime->Clusterizer(reconame.Data())->GetRecPointsInRun()) / gime->MaxEvent();
305  
306    if ( nRecPoints < maxRecPoints-widRecPoints || nRecPoints > maxRecPoints+widRecPoints ) {
307     boolerror = 1;
308     cerr<<"__________________________________________________________________"<<endl;
309     cerr<<" "<<endl;
310     cerr<<"             MESS ==> Error detected in the Clusterizing process. Sending error file to PHOS director."<<endl;
311     cerr<<"__________________________________________________________________"<<endl;
312   
313       gSystem->Exec("uuencode $ALICE_ROOT/PHOS/testPHOS.root testPHOS.root | mail -s 'PHOS INSTALLATION ERROR' almajor@caramail.com");
314     }
315    
316   cerr<<"__________________________________________________________________"<<endl;
317   cerr<<" "<<endl;
318   cerr<<"             MESS ==> Cluster process ended successfully."<<endl;
319   cerr<<"__________________________________________________________________"<<endl;
320  }
321
322 //TrackSegments process
323  if(boolerror == 0){ 
324   AliPHOSTrackSegmentMaker * tracks = new AliPHOSTrackSegmentMakerv1("testPHOS.root",reconame.Data());
325   tracks->ExecuteTask("deb");
326   Float_t nTrackSegments =  (Float_t) (gime->TrackSegmentMaker(reconame.Data())->GetTrackSegmentsInRun()) / gime->MaxEvent();
327  
328    if ( nTrackSegments < maxTrackSegments-0.25 || nTrackSegments > maxTrackSegments+0.25 ) {
329     boolerror = 1;
330     cerr<<"__________________________________________________________________"<<endl;
331     cerr<<" "<<endl;
332     cerr<<"             MESS ==> Error detected in the TrackSegments process. Sending error file to PHOS director."<<endl;
333     cerr<<"__________________________________________________________________"<<endl;
334     
335       gSystem->Exec("uuencode $ALICE_ROOT/PHOS/testPHOS.root testPHOS.root | mail -s 'PHOS INSTALLATION ERROR' almajor@caramail.com");
336     }
337    
338   cerr<<"__________________________________________________________________"<<endl;
339   cerr<<" "<<endl;
340   cerr<<"             MESS ==> TrackSegments process ended successfully."<<endl;
341   cerr<<"__________________________________________________________________"<<endl;
342  }
343 //RecParticles process
344  if(boolerror == 0){
345   AliPHOSPID * pid = new AliPHOSPIDv1("testPHOS.root",reconame.Data()); pid->ExecuteTask("deb");
346  
347   Float_t nRecParticles =  (Float_t) (gime->PID(reconame.Data())->GetRecParticlesInRun())/gime->MaxEvent();
348  
349  
350    if ( nRecParticles < maxRecParticles-0.25 || nRecParticles > maxRecParticles+0.25 ) {
351     boolerror = 1;
352     cerr<<"__________________________________________________________________"<<endl;
353     cerr<<" "<<endl;
354     cerr<<"             MESS ==> Error detected in the RecParticles process. Sending error file to PHOS director.Stop reconstruction."<<endl;
355     cerr<<"__________________________________________________________________"<<endl;
356     
357     gSystem->Exec("uuencode $ALICE_ROOT/PHOS/testPHOS.root testPHOS.root | mail -s 'PHOS INSTALLATION ERROR' almajor@caramail.com");
358     }
359  
360   cerr<<"__________________________________________________________________"<<endl;
361   cerr<<" "<<endl;
362   cerr<<"             MESS ==> RecParticles process ended successfully."<<endl;
363   cerr<<"__________________________________________________________________"<<endl;
364  }
365  if(boolerror == 0){
366   cerr<<"_____________________________________________________________________"<<endl;
367   cerr<<" "<<endl;
368   cerr<<"             MESS ==> reconstruction ended successfully."<<endl;
369   cerr<<"_____________________________________________________________________"<<endl;
370  }
371  else {
372  cerr<<"_____________________________________________________________________"<<endl;
373   cerr<<" "<<endl;
374   cerr<<"             MESS ==> Due to an error, reconstruction has not been completed.Please try again when the error has been found."<<endl;
375   cerr<<"_____________________________________________________________________"<<endl;
376   }
377  }
378 }