char * to TString
[u/mrichter/AliRoot.git] / ESDCheck / AliFMDQATask.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 // An analysis task to check the FMD data in simulated data
20 //
21 //*-- Hans Hjersing Dalsgaard 
22 //////////////////////////////////////////////////////////////////////////////
23
24 #include <TCanvas.h> 
25 #include <TChain.h>
26 #include <TF1.h> 
27 #include <TFile.h> 
28 #include <TH1D.h> 
29 #include <TROOT.h>
30 #include <TString.h> 
31
32 #include "AliFMDQATask.h" 
33 #include "AliESD.h" 
34 #include "AliLog.h"
35
36 //______________________________________________________________________________
37 AliFMDQATask::AliFMDQATask(const char *name) : 
38   AliAnalysisTask(name,""),  
39   fChain(0),
40   fESD(0), 
41   fhFMD1i(0),
42   fhFMD2i(0), 
43   fhFMD2o(0), 
44   fhFMD3i(0), 
45   fhFMD3o(0) 
46 {
47   // Constructor.
48   // Input slot #0 works with an Ntuple
49   DefineInput(0, TChain::Class());
50   // Output slot #0 writes into a TH1 container
51   DefineOutput(0,  TObjArray::Class()) ; 
52 }
53
54 //______________________________________________________________________________
55 AliFMDQATask::~AliFMDQATask()
56 {
57   // dtor
58
59   fOutputContainer->Clear() ; 
60   delete fOutputContainer ; 
61
62   delete fhFMD1i ;
63   delete fhFMD2i ; 
64   delete fhFMD2o ; 
65   delete fhFMD3i ; 
66   delete fhFMD3o ;
67
68 }
69
70 //______________________________________________________________________________
71 void AliFMDQATask::ConnectInputData(const Option_t*)
72 {
73   // Initialisation of branch container and histograms 
74     
75   AliInfo(Form("*** Initialization of %s", GetName())) ; 
76   
77   // Get input data
78   fChain = dynamic_cast<TChain *>(GetInputData(0)) ;
79   if (!fChain) {
80     AliError(Form("Input 0 for %s not found\n", GetName()));
81     return ;
82   }
83   
84   // One should first check if the branch address was taken by some other task
85   char ** address = (char **)GetBranchAddress(0, "ESD");
86   if (address) {
87     fESD = (AliESD*)(*address);
88   } else {
89     fESD = new AliESD();
90     SetBranchAddress(0, "ESD", &fESD);
91   }
92 }
93
94 //________________________________________________________________________
95 void AliFMDQATask::CreateOutputObjects()
96 {  
97   // create histograms 
98
99   OpenFile(0) ; 
100
101   fhFMD1i = new TH1D("FMD1i", "FMD1i", 100, -0.5, 3);
102   fhFMD2i = new TH1D("FMD2i", "FMD2i", 100, -0.5, 3);
103   fhFMD2o = new TH1D("FMD2o", "FMD2o", 100, -0.5, 3);
104   fhFMD3i = new TH1D("FMD3i", "FMD3i", 100, -0.5, 3);
105   fhFMD3o = new TH1D("FMD3o", "FMD3o", 100, -0.5, 3);
106
107   // create output container
108   
109   fOutputContainer = new TObjArray(5) ; 
110   fOutputContainer->SetName(GetName()) ; 
111
112   fOutputContainer->AddAt(fhFMD1i,             0) ; 
113   fOutputContainer->AddAt(fhFMD2i,             1) ; 
114   fOutputContainer->AddAt(fhFMD2o,             2) ; 
115   fOutputContainer->AddAt(fhFMD3i,             3) ; 
116   fOutputContainer->AddAt(fhFMD3o,             4) ; 
117 }
118
119 //______________________________________________________________________________
120 void AliFMDQATask::Exec(Option_t *) 
121 {
122   // Processing of one event
123     
124   Long64_t entry = fChain->GetReadEntry() ;
125   
126   if (!fESD) {
127     AliError("fESD is not connected to the input!") ; 
128     return ; 
129   }
130   
131   TFile * currentFile = (dynamic_cast<TChain *>(fChain))->GetFile() ; 
132   if ( !((entry-1)%100) ) 
133     AliInfo(Form("%s ----> Processing event # %lld", currentFile->GetName(), entry)) ; 
134             
135   // ************************  FMD *************************************
136  AliESDFMD * fmd = fESD->GetFMDData() ;
137  if ( ! fmd ) {
138    AliError("No FMD found in ESD") ; 
139    return ; 
140  }
141  fmd->CheckNeedUShort(currentFile);
142   Int_t nFMD1 = 0, nFMD2i = 0, nFMD2o = 0, nFMD3i = 0, nFMD3o = 0 ;
143   
144   UShort_t detector = 1 ;
145   for(detector = 1 ; detector <= fmd->MaxDetectors() ; detector++) {
146     Char_t ring = 'O' ;
147     UShort_t sector ; 
148     for(sector = 0 ;sector < fmd->MaxSectors() ; sector++) {
149       UShort_t strip ;
150       for(strip = 0 ; strip < fmd->MaxStrips(); strip++) {
151         if(fmd->Multiplicity(detector, ring, sector, strip) != AliESDFMD::kInvalidMult)
152           RingSelector(detector, ring, fmd->Multiplicity(detector, ring, sector, strip)) ;
153         if( (fmd->Multiplicity(detector, ring, sector, strip) == AliESDFMD::kInvalidMult) && detector == 2 )
154           nFMD2o++ ;
155         if( (fmd->Multiplicity(detector, ring, sector, strip)==AliESDFMD::kInvalidMult) && detector == 3 )
156           nFMD3o++ ;
157       }
158     }
159     ring='I';
160     for(sector = 0; sector < fmd->MaxSectors(); sector++) {
161       UShort_t strip ;
162       for(strip = 0 ; strip < fmd->MaxStrips() ; strip++) {
163         if(fmd->Multiplicity(detector, ring, sector, strip) != AliESDFMD::kInvalidMult)
164           RingSelector(detector, ring, fmd->Multiplicity(detector, ring, sector, strip));
165         if( (fmd->Multiplicity(detector, ring, sector, strip) == AliESDFMD::kInvalidMult) && detector == 1 )
166           nFMD1++;
167         if( (fmd->Multiplicity(detector, ring, sector, strip) == AliESDFMD::kInvalidMult) && detector == 2 )
168           nFMD2i++;
169         if( (fmd->Multiplicity(detector, ring, sector, strip) == AliESDFMD::kInvalidMult) && detector == 3 )
170           nFMD3i++;
171       }     
172     }
173   }
174
175   if(nFMD1>100+10240)
176     AliWarning(Form("number of missing strips in FMD1i too high in event number %lld in file", entry, fChain->GetCurrentFile()->GetName())) ;
177   if(nFMD2i>100+10240)
178     AliWarning(Form("number of missing strips in FMD2i too high in event number %lld in file", entry, fChain->GetCurrentFile()->GetName())) ;
179   if(nFMD2o>100+10240)
180     AliWarning(Form("number of missing strips in FMD2o too high in event number %lld in file", entry, fChain->GetCurrentFile()->GetName())) ;
181   if(nFMD3i>100+10240)
182     AliWarning(Form("number of missing strips in FMD3i too high in event number %lld in file", entry, fChain->GetCurrentFile()->GetName())) ;
183   if(nFMD3o>100+10240)
184     AliWarning(Form("number of missing strips in FMD3o too high in event number %lld in file", entry, fChain->GetCurrentFile()->GetName())) ;
185   
186   PostData(0, fOutputContainer);
187 }
188
189 //______________________________________________________________________________
190 void AliFMDQATask::Terminate(Option_t *)
191 {
192   // Processing when the event loop is ended
193  
194   fOutputContainer = (TObjArray*)GetOutputData(0);
195   if ( ! fOutputContainer ) 
196     return ; 
197
198   fhFMD1i = (TH1D*)fOutputContainer->At(0);
199   fhFMD2i = (TH1D*)fOutputContainer->At(1);
200   fhFMD2o = (TH1D*)fOutputContainer->At(2);
201   fhFMD3i = (TH1D*)fOutputContainer->At(3);
202   fhFMD3o = (TH1D*)fOutputContainer->At(4);
203
204   TCanvas * cFMD1 = new TCanvas("cFMD1", "FMD ESD Test", 400, 10, 600, 700);
205   cFMD1->Divide(3, 2) ; 
206
207   cFMD1->cd(1) ;; 
208   fhFMD1i->Draw() ; 
209   
210   cFMD1->cd(2) ;; 
211   fhFMD2i->Draw() ; 
212
213   cFMD1->cd(3) ;; 
214   fhFMD2o->Draw() ; 
215
216   cFMD1->cd(4) ;; 
217   fhFMD3i->Draw() ; 
218
219   cFMD1->cd(5) ;; 
220   fhFMD3o->Draw() ; 
221   
222   cFMD1->Print("FMD.eps") ;
223
224   Bool_t rv1i = TestHisto(fhFMD1i) ;
225   Bool_t rv2i = TestHisto(fhFMD2i) ;
226   Bool_t rv2o = TestHisto(fhFMD2o) ;
227   Bool_t rv3i = TestHisto(fhFMD2i) ;
228   Bool_t rv3o = TestHisto(fhFMD2o) ;
229   
230   if (  !(rv1i * rv2i * rv2o * rv3i * rv3o) )
231     AliWarning("Possible problem in file !!! Check output!") ;
232
233
234
235   char line[1024] ; 
236   sprintf(line, ".!tar -zcvf %s.tar.gz *.eps", GetName()) ; 
237   gROOT->ProcessLine(line);
238   sprintf(line, ".!rm -fR *.eps"); 
239   gROOT->ProcessLine(line);
240  
241   AliInfo(Form("!!! All the eps files are in %s.tar.gz !!! \n", GetName())) ;
242 }
243
244 //______________________________________________________________________________
245 void AliFMDQATask::RingSelector(const UShort_t detector, const Char_t ring, const Float_t mult) const 
246
247   // fill the histograms for each ring in each detector layer
248
249   if(ring == 'I' && detector == 1)
250     fhFMD1i->Fill(mult) ;
251   if(ring == 'I' && detector == 2)
252     fhFMD2i->Fill(mult) ;
253   if(ring == 'O' && detector == 2)
254     fhFMD2o ->Fill(mult) ;
255   if(ring == 'I' && detector == 3)
256     fhFMD3i ->Fill(mult) ;
257   if(ring == 'O' && detector == 3)
258     fhFMD3o ->Fill(mult) ;
259 }
260 //______________________________________________________________________________
261 Bool_t AliFMDQATask::TestHisto(TH1D * hTest) const 
262 {  
263   // analyses the histogram with a Landau function
264
265   Float_t chiMax = 3, chiLow=0.5 ;
266   Float_t chiSq ;
267   Float_t mpv ;
268   Int_t   ndf ;
269   
270   FitAll(hTest, chiSq, ndf, mpv, chiMax, chiLow);
271
272   if( (chiSq > chiMax || chiSq < chiLow) || mpv < 0.6 || mpv > 1 ) {
273     hTest->Rebin(2) ;
274     FitAll(hTest, chiSq, ndf, mpv, chiMax, chiLow) ;
275   }
276
277   Bool_t   ret   = kFALSE ;
278   Char_t * test  = "not OK";
279   Char_t * test2 = "not OK";
280   
281   if(chiSq < chiMax && chiSq > chiLow)
282     test = "OK" ;
283   if(mpv > 0.6 && mpv < 1)
284     test2 = "OK" ;
285  
286   if(test == "OK" && test2 == "OK")
287     ret = kTRUE;
288   
289   if(test == "not OK" || test2 == "not OK") {
290     AliWarning("Bad fit results") ; 
291     printf("Detector : %s\n", hTest->GetName()) ;
292     printf("Landau fit Chi Square / NDF = %f / %d which is %s\n", chiSq*ndf, ndf, test) ; 
293     printf("Landau fit MPV is: %f which is %s\n", mpv, test2) ;
294   }
295   return ret;  
296 }
297 //______________________________________________________________________________
298 void AliFMDQATask::FitAll(TH1D* hTest, Float_t &chiSq, Int_t &ndf, Float_t &mpv, Float_t chiMax, Float_t chiLow ) const 
299 {
300   // fit a histogram with a Landau distribution and returns chi square
301
302   Float_t fitMax = hTest->GetXaxis()->GetXmax() ;
303
304   TH1D hTmp = *hTest ;
305   hTmp.SetAxisRange(0.4,fitMax) ;
306   Int_t   maxBin = hTmp.GetMaximumBin();
307   Float_t max    = hTmp.GetBinCenter(maxBin);
308  
309   hTest->Fit("landau", "QOI", "", max-0.3, fitMax) ;
310   TF1* fitfunc = hTest->GetFunction("landau") ;
311   chiSq = fitfunc->GetChisquare() / fitfunc->GetNDF() ;
312   mpv   = fitfunc->GetParameter(1) ;
313   ndf   = fitfunc->GetNDF() ;
314
315   if( ( chiSq > chiMax || chiSq < chiLow ) || ( mpv < 0.6 || mpv > 1 ) ) {
316     hTest->Fit("landau", "QOI", "", max-0.2, fitMax) ;
317     fitfunc = hTest->GetFunction("landau") ;
318     chiSq   = fitfunc->GetChisquare() / fitfunc->GetNDF() ;
319     mpv     = fitfunc->GetParameter(1) ;
320     ndf     = fitfunc->GetNDF() ;
321   }
322   if( ( chiSq >chiMax || chiSq < chiLow ) || ( mpv < 0.6 || mpv > 1 ) ) {
323     hTest->Fit("landau", "QOI", "", max-0.1, fitMax) ;
324     fitfunc = hTest->GetFunction("landau") ;
325     chiSq   = fitfunc->GetChisquare() / fitfunc->GetNDF() ;
326     mpv     = fitfunc->GetParameter(1) ;
327     ndf     = fitfunc->GetNDF() ;
328   }
329   if( ( chiSq > chiMax || chiSq <chiLow ) || ( mpv < 0.6 || mpv > 1 ) ) {
330     hTest->Fit("landau", "QOI", "", max, fitMax) ;
331     fitfunc = hTest->GetFunction("landau") ;
332     chiSq   = fitfunc->GetChisquare() / fitfunc->GetNDF() ;
333     mpv     = fitfunc->GetParameter(1) ;
334     ndf     = fitfunc->GetNDF(); 
335   }
336   if( ( chiSq > chiMax || chiSq < chiLow ) || ( mpv < 0.6 ||mpv > 1 ) ) {
337     hTest->Fit("landau", "QOI", "", max+0.1, fitMax) ;
338     fitfunc = hTest->GetFunction("landau") ;
339     chiSq   = fitfunc->GetChisquare() / fitfunc->GetNDF() ;
340     mpv     = fitfunc->GetParameter(1) ;
341     ndf     = fitfunc->GetNDF() ;
342   } 
343 }