0c2d89c2dd6d427350d65956217b9aa1172dc1bc
[u/mrichter/AliRoot.git] / PWG2 / FORWARD / analysis2 / AliFMDCorrELossFit.cxx
1 // Object holding the Energy loss fit 'correction' 
2 // 
3 // These are generated from Monte-Carlo or real ESDs. 
4 //
5 #include "AliFMDCorrELossFit.h"
6 #include "AliForwardUtil.h"
7 #include <TF1.h>
8 #include <TBrowser.h>
9 #include <TVirtualPad.h>
10 #include <THStack.h>
11 #include <TH1D.h>
12 #include <AliLog.h>
13 #include <TMath.h>
14 #include <TList.h>
15 #include <iostream>
16 #include <iomanip>
17
18 Double_t AliFMDCorrELossFit::ELossFit::fgMaxRelError = .12;
19 Double_t AliFMDCorrELossFit::ELossFit::fgLeastWeight = 1e-5;
20 Double_t AliFMDCorrELossFit::ELossFit::fgMaxChi2nu   = 20;
21
22 //____________________________________________________________________
23 AliFMDCorrELossFit::ELossFit::ELossFit()
24   : fN(0),
25     fNu(0),
26     fChi2(0),
27     fC(0),
28     fDelta(0),
29     fXi(0),
30     fSigma(0),
31     fSigmaN(0),
32     fA(0),
33     fEC(0),
34     fEDelta(0),
35     fEXi(0),
36     fESigma(0),
37     fESigmaN(0),
38     fEA(0),
39     fQuality(0), 
40     fDet(0), 
41     fRing('\0'),
42     fBin(0)
43 {
44   //
45   // Default constructor 
46   // 
47   //
48 }
49 //____________________________________________________________________
50 AliFMDCorrELossFit::ELossFit::ELossFit(Int_t quality, const TF1& f)
51   : fN(f.GetNpar() > AliForwardUtil::ELossFitter::kN ? 
52        Int_t(f.GetParameter(AliForwardUtil::ELossFitter::kN)) : 
53        1),
54     fNu(f.GetNDF()),
55     fChi2(f.GetChisquare()),
56     fC(f.GetParameter(AliForwardUtil::ELossFitter::kC)),
57     fDelta(f.GetParameter(AliForwardUtil::ELossFitter::kDelta)),
58     fXi(f.GetParameter(AliForwardUtil::ELossFitter::kXi)),
59     fSigma(f.GetParameter(AliForwardUtil::ELossFitter::kSigma)),
60     fSigmaN(f.GetParameter(AliForwardUtil::ELossFitter::kSigmaN)),
61     fA(0),
62     fEC(f.GetParError(AliForwardUtil::ELossFitter::kC)),
63     fEDelta(f.GetParError(AliForwardUtil::ELossFitter::kDelta)),
64     fEXi(f.GetParError(AliForwardUtil::ELossFitter::kXi)),
65     fESigma(f.GetParError(AliForwardUtil::ELossFitter::kSigma)),
66     fESigmaN(f.GetParError(AliForwardUtil::ELossFitter::kSigmaN)),
67     fEA(0),
68     fQuality(quality),
69     fDet(0), 
70     fRing('\0'),
71     fBin(0)
72 {
73   // 
74   // Construct from a function
75   // 
76   // Parameters:
77   //    quality Quality flag
78   //    f       Function
79   //
80   if (fN <= 0) return;
81   fA  = new Double_t[fN];
82   fEA = new Double_t[fN];
83   for (Int_t i = 0; i < fN-1; i++) { 
84     fA[i]  = f.GetParameter(AliForwardUtil::ELossFitter::kA+i);
85     fEA[i] = f.GetParError(AliForwardUtil::ELossFitter::kA+i);
86   }
87   fA[fN-1]  = -9999;
88   fEA[fN-1] = -9999;
89 }
90
91 //____________________________________________________________________
92 AliFMDCorrELossFit::ELossFit::ELossFit(Int_t     quality,UShort_t  n, 
93                                        Double_t  chi2,   UShort_t  nu, 
94                                        Double_t  c,      Double_t  ec, 
95                                        Double_t  delta,  Double_t  edelta, 
96                                        Double_t  xi,     Double_t  exi,
97                                        Double_t  sigma,  Double_t  esigma, 
98                                        Double_t  sigman, Double_t  esigman, 
99                                        const Double_t* a,const Double_t* ea)
100   : fN(n),
101     fNu(nu),
102     fChi2(chi2),
103     fC(c),
104     fDelta(delta),
105     fXi(xi),
106     fSigma(sigma),
107     fSigmaN(sigman),
108     fA(0),
109     fEC(ec),
110     fEDelta(edelta),
111     fEXi(exi),
112     fESigma(esigma),
113     fESigmaN(esigman),
114     fEA(0),
115     fQuality(quality),
116     fDet(0), 
117     fRing('\0'),
118     fBin(0)
119 {
120   // 
121   // Constructor with full parameter set
122   // 
123   // Parameters:
124   //    quality   Quality flag
125   //    n         @f$ N@f$ - Number of fitted peaks
126   //    chi2      @f$ \chi^2 @f$
127   //    nu        @f$ \nu @f$ - number degrees of freedom
128   //    c         @f$ C@f$ - scale constant
129   //    ec        @f$ \delta C@f$ - error on @f$ C@f$ 
130   //    delta     @f$ \Delta@f$ - Most probable value             
131   //    edelta    @f$ \delta\Delta@f$ - error on @f$\Delta@f$ 
132   //    xi        @f$ \xi@f$ - width  
133   //    exi       @f$ \delta\xi@f$ - error on @f$\xi@f$ 
134   //    sigma     @f$ \sigma@f$ - Width of Gaussian                
135   //    esigma    @f$ \delta\sigma@f$ - error on @f$\sigma@f$ 
136   //    sigman    @f$ \sigma_n@f$ - Noise width           
137   //    esigman   @f$ \delta\sigma_n@f$ - error on @f$\sigma_n@f$ 
138   //    a         Array of @f$ N-1@f$ weights @f$ a_i@f$ for 
139   //                  @f$ i=2,\ldots@f$ 
140   //    ea        Array of @f$ N-1@f$ error on the weights @f$ a_i@f$ for 
141   //                  @f$ i=2,\ldots@f$ 
142   //
143   if (fN <= 0) return;
144   fA  = new Double_t[fN];
145   fEA = new Double_t[fN];
146   for (Int_t i = 0; i < fN-1; i++) { 
147     fA[i]  = a[i];
148     fEA[i] = ea[i];
149   }
150   fA[fN-1]  = -9999;
151   fEA[fN-1] = -9999;
152 }
153 //____________________________________________________________________
154 AliFMDCorrELossFit::ELossFit::ELossFit(const ELossFit& o)
155   : TObject(o), 
156     fN(o.fN),
157     fNu(o.fNu),
158     fChi2(o.fChi2),
159     fC(o.fC),
160     fDelta(o.fDelta),
161     fXi(o.fXi),
162     fSigma(o.fSigma),
163     fSigmaN(o.fSigmaN),
164     fA(0),
165     fEC(o.fEC),
166     fEDelta(o.fEDelta),
167     fEXi(o.fEXi),
168     fESigma(o.fESigma),
169     fESigmaN(o.fESigmaN),
170     fEA(0),
171     fQuality(o.fQuality),
172     fDet(o.fDet), 
173     fRing(o.fRing),
174     fBin(o.fBin)
175 {
176   // 
177   // Copy constructor 
178   // 
179   // Parameters:
180   //    o Object to copy from 
181   //
182   if (fN <= 0) return;
183   fA  = new Double_t[fN];
184   fEA = new Double_t[fN];
185   for (Int_t i = 0; i < fN-1; i++) { 
186     fA[i]  = o.fA[i];
187     fEA[i] = o.fEA[i];
188   }
189   fA[fN-1]  = -9999;
190   fEA[fN-1] = -9999;
191 }
192
193 //____________________________________________________________________
194 AliFMDCorrELossFit::ELossFit&
195 AliFMDCorrELossFit::ELossFit::operator=(const ELossFit& o)
196 {
197   // 
198   // Assignment operator 
199   // 
200   // Parameters:
201   //    o Object to assign from 
202   // 
203   // Return:
204   //    Reference to this object 
205   //
206   fN       = o.fN;
207   fNu      = o.fNu;
208   fChi2    = o.fChi2;
209   fC       = o.fC;
210   fDelta   = o.fDelta;
211   fXi      = o.fXi;
212   fSigma   = o.fSigma;
213   fSigmaN  = o.fSigmaN;
214   fEC      = o.fEC;
215   fEDelta  = o.fEDelta;
216   fEXi     = o.fEXi;
217   fESigma  = o.fESigma;
218   fESigmaN = o.fESigmaN;
219   fQuality = o.fQuality;
220   fDet     = o.fDet; 
221   fRing    = o.fRing;
222   fBin     = o.fBin;
223   if (fA)  delete [] fA;
224   if (fEA) delete [] fEA; 
225   fA  = 0;
226   fEA = 0;
227
228   if (fN <= 0) return *this;
229   fA  = new Double_t[fN];
230   fEA = new Double_t[fN];
231   for (Int_t i = 0; i < fN; i++) { 
232     fA[i]  = o.fA[i];
233     fEA[i] = o.fEA[i];
234   }
235
236   return *this;
237 }
238
239 //____________________________________________________________________
240 AliFMDCorrELossFit::ELossFit::~ELossFit()
241 {
242   if (fA)  delete[] fA;
243   if (fEA) delete[] fEA;
244 }
245
246
247 //____________________________________________________________________
248 Int_t
249 AliFMDCorrELossFit::ELossFit::FindMaxWeight(Double_t maxRelError, 
250                                             Double_t leastWeight, 
251                                             UShort_t  maxN) const
252 {
253   // 
254   // Find the maximum weight to use.  The maximum weight is the
255   // largest i for which 
256   // 
257   // - @f$ i \leq \max{N}@f$ 
258   // - @f$ a_i > \min{a}@f$ 
259   // - @f$ \delta a_i/a_i > \delta_{max}@f$ 
260   // 
261   // Parameters:
262   //    maxRelError @f$ \min{a}@f$ 
263   //    leastWeight @f$ \delta_{max}@f$ 
264   //    maxN        @f$ \max{N}@f$      
265   // 
266   // Return:
267   //    The largest index @f$ i@f$ for which the above
268   // conditions hold.  Will never return less than 1. 
269   //
270   Int_t n = TMath::Min(maxN, UShort_t(fN-1));
271   Int_t m = 1;
272   // fN is one larger than we have data 
273   for (Int_t i = 0; i < n-1; i++, m++) { 
274     if (fA[i] < leastWeight)  break;
275     if (fEA[i] / fA[i] > maxRelError) break;
276   }
277   return m;
278 }
279
280 //____________________________________________________________________
281 Double_t 
282 AliFMDCorrELossFit::ELossFit::Evaluate(Double_t x, 
283                                        UShort_t maxN) const
284 {
285   // 
286   // Evaluate 
287   // @f[ 
288   //  f_N(x;\Delta,\xi,\sigma') = 
289   //     \sum_{i=1}^{n} a_i f(x;\Delta_i,\xi_i,\sigma_i')
290   // @f] 
291   //
292   // (see AliForwardUtil::NLandauGaus) for the maximum @f$ N @f$
293   // that fulfills the requirements 
294   // 
295   // Parameters:
296   //    x           Where to evaluate 
297   //    maxN      @f$ \max{N}@f$    
298   // 
299   // Return:
300   //    @f$ f_N(x;\Delta,\xi,\sigma')@f$ 
301   //
302   return AliForwardUtil::NLandauGaus(x, fDelta, fXi, fSigma, fSigmaN, 
303                                      TMath::Min(maxN, UShort_t(fN)), fA);
304 }
305
306 //____________________________________________________________________
307 Double_t 
308 AliFMDCorrELossFit::ELossFit::EvaluateWeighted(Double_t x, 
309                                                UShort_t maxN) const
310 {                                                                       
311   // 
312   // Evaluate 
313   // @f[ 
314   //   f_W(x;\Delta,\xi,\sigma') = 
315   //   \frac{\sum_{i=1}^{n} i a_i f_i(x;\Delta,\xi,\sigma')}{
316   //     f_N(x;\Delta,\xi,\sigma')} = 
317   //   \frac{\sum_{i=1}^{n} i a_i f(x;\Delta_i,\xi_i,\sigma_i')}{
318   //     \sum_{i=1}^{n} a_i f(x;\Delta_i,\xi_i,\sigma_i')}
319   // @f] 
320   // where @f$ n@f$ fulfills the requirements (see FindMaxWeight). 
321   //
322   // If the denominator is zero, then 1 is returned. 
323   //
324   // See also AliForwardUtil::ILandauGaus and AliForwardUtil::NLandauGaus
325   // for more information on the evaluated functions. 
326   // 
327   // Parameters:
328   //    x           Where to evaluate 
329   //    maxN      @f$ \max{N}@f$      
330   // 
331   // Return:
332   //    @f$ f_W(x;\Delta,\xi,\sigma')@f$.  
333   //
334   UShort_t n   = TMath::Min(maxN, UShort_t(fN-1));
335   Double_t num = 0;
336   Double_t den = 0;
337   for (Int_t i = 1; i <= n; i++) {
338     Double_t a = (i == 1 ? 1 : fA[i-1]);
339     if (fA[i-1] < 0) break;
340     Double_t f = AliForwardUtil::ILandauGaus(x,fDelta,fXi,fSigma,fSigmaN,i);
341     num += i * a * f;
342     den += a * f;
343   }
344   if (den <= 0) return 1;
345   return num / den;
346 }
347
348
349 #define OUTPAR(N,V,E)                   \
350   std::setprecision(9) <<               \
351   std::setw(12) << N << ": " <<         \
352   std::setw(14) << V << " +/- " <<      \
353   std::setw(14) << E << " (" <<         \
354   std::setprecision(-1) <<              \
355   std::setw(5)  << 100*(V>0?E/V:1) << "%)\n"
356   
357
358 //____________________________________________________________________
359 Int_t
360 AliFMDCorrELossFit::ELossFit::Compare(const TObject* o) const
361 {
362   // 
363   // Compare to another ELossFit object. 
364   // 
365   // - +1, if this quality is better (larger) than other objects quality
366   // - -1, if this quality is worse (smaller) than other objects quality
367   // - +1, if this @f$|\chi^2/\nu-1|@f$ is smaller than the same for other
368   // - -1, if this @f$|\chi^2/\nu-1|@f$ is larger than the same for other
369   // - 0 otherwise 
370   // 
371   // Parameters:
372   //    o Other object to compare to 
373   //
374   const ELossFit* other = static_cast<const ELossFit*>(o);
375   if (this->fQuality < other->fQuality) return -1;
376   if (this->fQuality > other->fQuality) return +1;
377   Double_t chi2nu  = (fNu == 0 ? 1000 : fChi2 / fNu);
378   Double_t oChi2nu = (other->fNu == 0 ? 1000 : other->fChi2 / other->fNu);
379   if (TMath::Abs(chi2nu-1) < TMath::Abs(oChi2nu-1)) return +1;
380   if (TMath::Abs(chi2nu-1) > TMath::Abs(oChi2nu-1)) return -1;
381   return 0;
382 }
383
384 //____________________________________________________________________
385 void
386 AliFMDCorrELossFit::ELossFit::Print(Option_t*) const
387 {
388   // 
389   // Information to standard output 
390   // 
391   // Parameters:
392   //    option Not used 
393   //
394   std::cout << GetName() << ":\n"
395             << " chi^2/nu = " << fChi2 << "/" << fNu << " = " 
396             << (fNu == 0 ? 999 : fChi2 / fNu) << "\n"
397             << "     Quality:   " << fQuality << "\n" 
398             << "  NParticles:   " << fN << "  (" << FindMaxWeight() << ")\n"
399             << OUTPAR("Delta", fDelta, fEDelta) 
400             << OUTPAR("xi", fXi, fEXi)
401             << OUTPAR("sigma", fSigma, fESigma)
402             << OUTPAR("sigma_n", fSigmaN, fESigmaN);
403   for (Int_t i = 0; i < fN-1; i++) 
404     std::cout << OUTPAR(Form("a%d", i+2), fA[i], fEA[i]);
405   std::cout << std::flush;
406 }
407
408 //____________________________________________________________________
409 const Char_t*
410 AliFMDCorrELossFit::ELossFit::GetName() const 
411 {
412   // 
413   // Get the name of this object 
414   // 
415   // 
416   // Return:
417   //    
418   //
419   return Form("FMD%d%c_etabin%03d", fDet, fRing, fBin);
420 }
421
422 //____________________________________________________________________
423 void
424 AliFMDCorrELossFit::ELossFit::Browse(TBrowser* b)
425 {
426   // 
427   // Browse this object 
428   // 
429   // Parameters:
430   //    b Browser
431   //
432   Draw(b ? b->GetDrawOption() : "comp");
433   gPad->SetLogy();
434   gPad->Update();
435 }
436      
437 //____________________________________________________________________
438 void
439 AliFMDCorrELossFit::ELossFit::Draw(Option_t* option)
440 {
441   // 
442   // Draw this fit 
443   // 
444   // Parameters:
445   //    option Options 
446   //  - COMP  Draw components too 
447   //
448   TString opt(option);
449   opt.ToUpper();
450   bool comp = false;
451   if (opt.Contains("COMP")) { 
452     opt.ReplaceAll("COMP","");
453     comp = true;
454   }
455   if (!opt.Contains("SAME")) { 
456     gPad->Clear();
457   }
458
459   TObjArray cleanup;
460   TF1* tot = AliForwardUtil::MakeNLandauGaus(1, 
461                                              fDelta, fXi, 
462                                              fSigma, fSigmaN, 
463                                              fN,     fA, 
464                                              0.01,   10);
465   tot->SetLineColor(kBlack);
466   tot->SetLineWidth(2);
467   tot->SetLineStyle(1);
468   tot->SetTitle(GetName());
469   Double_t max = tot->GetMaximum();
470
471   if (!opt.Contains("SAME")) {
472     TH1* frame = new TH1F(GetName(), 
473                           Form("FMD%d%c, eta bin %d",fDet,fRing,fBin),
474                           100, 0, 10);
475     frame->SetMinimum(max/10000);
476     frame->SetMaximum(max*1.4);
477     frame->SetXTitle("#Delta / #Delta_{mip}");
478     frame->Draw();
479     opt.Append(" SAME");
480   }
481   tot->DrawCopy(opt.Data());
482   cleanup.Add(tot);
483
484   if (!comp) { 
485     gPad->cd();
486     return;
487   }
488
489   Double_t min = max;
490   opt.Append(" same");
491   Int_t maxW = FindMaxWeight();
492   for (Int_t i=1; i <= fN; i++) { 
493     TF1* f = AliForwardUtil::MakeILandauGaus((i == 1 ? 1 : fA[i-2]), 
494                                              fDelta, fXi, 
495                                              fSigma, fSigmaN, 
496                                              i,      0.01, 10);
497     f->SetLineWidth(2);
498     f->SetLineStyle(i > maxW ? 2 : 1);
499     min = TMath::Min(f->GetMaximum(), min);
500     f->DrawCopy(opt.Data());
501     cleanup.Add(f);
502   }
503   min /= 100;
504   tot->GetHistogram()->SetMaximum(max);
505   tot->GetHistogram()->SetMinimum(min);
506   tot->GetHistogram()->GetYaxis()->SetRangeUser(min, max);
507
508   gPad->cd();
509 }
510
511                                                  
512 //____________________________________________________________________
513 #define CHECKPAR(V,E,T) ((V > 0) && (E / V < T))
514
515 //____________________________________________________________________
516 Double_t
517 AliFMDCorrELossFit::ELossFit::GetLowerBound(Double_t f) const
518 {
519   // 
520   // Return 
521   //    Delta * f
522   return f * fDelta;
523 }
524 //____________________________________________________________________
525 Double_t
526 AliFMDCorrELossFit::ELossFit::GetLowerBound(Double_t f, 
527                                             Bool_t includeSigma) const
528 {
529   // 
530   // Return 
531   //    Delta - f * (xi + sigma)
532   return fDelta - f * (fXi + (includeSigma ? fSigma : 0));
533 }
534
535 //____________________________________________________________________
536 void 
537 AliFMDCorrELossFit::ELossFit::CalculateQuality(Double_t maxChi2nu, 
538                                                Double_t maxRelError, 
539                                                Double_t leastWeight)
540 {
541   // 
542   // Calculate the quality 
543   //
544   Int_t qual = 0;
545   if (fNu > 0 && fChi2 / fNu < maxChi2nu) qual += 4;;
546   if (CHECKPAR(fDelta,  fEDelta,  maxRelError)) qual++;
547   if (CHECKPAR(fXi,     fEXi,     maxRelError)) qual++;
548   if (CHECKPAR(fSigma,  fESigma,  maxRelError)) qual++;
549   if (CHECKPAR(fSigmaN, fESigmaN, maxRelError)) qual++;
550   qual += FindMaxWeight(1.5*maxRelError, leastWeight, fN);
551   fQuality = qual;
552 }
553
554 //____________________________________________________________________
555 AliFMDCorrELossFit::AliFMDCorrELossFit()
556   : TObject(), 
557     fRings(), 
558     fEtaAxis(0,0,0), 
559     fLowCut(0)
560 {
561   // 
562   // Default constructor 
563   //
564   fRings.SetOwner(kTRUE);
565   fEtaAxis.SetTitle("#eta");
566   fEtaAxis.SetName("etaAxis");
567   fRings.SetName("rings");
568 }
569
570 //____________________________________________________________________
571 AliFMDCorrELossFit::AliFMDCorrELossFit(const AliFMDCorrELossFit& o)
572   : TObject(o), 
573     fRings(o.fRings),
574     fEtaAxis(o.fEtaAxis.GetNbins(),o.fEtaAxis.GetXmin(),o.fEtaAxis.GetXmax()), 
575     fLowCut(0)
576 {
577   // 
578   // Copy constructor 
579   // 
580   // Parameters:
581   //    o Object to copy from 
582   //
583   fEtaAxis.SetTitle("#eta");
584   fEtaAxis.SetName("etaAxis");
585 }
586
587 //____________________________________________________________________
588 AliFMDCorrELossFit::~AliFMDCorrELossFit()
589 {
590   // 
591   // Destructor 
592   //
593   fRings.Clear();
594 }
595
596 //____________________________________________________________________
597 AliFMDCorrELossFit&
598 AliFMDCorrELossFit::operator=(const AliFMDCorrELossFit& o)
599 {
600   // 
601   // Assignment operator 
602   // 
603   // Parameters:
604   //    o Object to assign from 
605   // 
606   // Return:
607   //    Reference to this object 
608   //
609   fRings = o.fRings;
610   fLowCut = o.fLowCut;
611   SetEtaAxis(o.fEtaAxis.GetNbins(), o.fEtaAxis.GetXmin(), o.fEtaAxis.GetXmax());
612
613   return *this;
614 }
615 //____________________________________________________________________
616 Int_t
617 AliFMDCorrELossFit::FindEtaBin(Double_t eta) const
618 {
619   // 
620   // Find the eta bin corresponding to the given eta 
621   // 
622   // Parameters:
623   //    eta  Eta value 
624   // 
625   // Return:
626   //    Bin (in @f$[1,N_{bins}]@f$) corresponding to the given
627   // eta, or 0 if out of range.
628   //
629   if (TMath::Abs(fEtaAxis.GetXmin() - fEtaAxis.GetXmax()) < 1e-6 
630       || fEtaAxis.GetNbins() == 0) {
631     AliWarning("No eta axis defined");
632     return -1;
633   }
634   Int_t bin = const_cast<TAxis&>(fEtaAxis).FindBin(eta);
635   if (bin <= 0 || bin > fEtaAxis.GetNbins()) return 0;
636   return bin;
637 }
638
639 //____________________________________________________________________
640 Bool_t
641 AliFMDCorrELossFit::SetFit(UShort_t d, Char_t r, Int_t etaBin, ELossFit* fit)
642 {
643   // 
644   // Set the fit parameters from a function 
645   // 
646   // Parameters:
647   //    d       Detector
648   //    r       Ring 
649   //    etaBin  Eta (bin number, 1->nBins)
650   //    f       ELoss fit result - note, the object will take ownership
651   //  
652   TObjArray* ringArray = GetOrMakeRingArray(d, r);
653   if (!ringArray) { 
654     AliError(Form("Failed to make ring array for FMD%d%c", d, r));
655     return kFALSE;
656   }
657   if (etaBin <= 0 || etaBin >= fEtaAxis.GetNbins()+1) { 
658     AliError(Form("bin=%d is out of range [%d,%d]", 
659                   etaBin, 1, fEtaAxis.GetNbins()));
660     return kFALSE;
661   }
662   // AliInfo(Form("Adding fit %p at %3d", fit, etaBin));
663   ringArray->AddAtAndExpand(fit, etaBin);
664   return kTRUE;
665 }
666
667 //____________________________________________________________________
668 Bool_t
669 AliFMDCorrELossFit::SetFit(UShort_t d, Char_t r, Double_t eta, ELossFit* fit)
670 {
671   // 
672   // Set the fit parameters from a function 
673   // 
674   // Parameters:
675   //    d    Detector
676   //    r    Ring 
677   //    eta  Eta 
678   //    f    ELoss fit result - note, the object will take ownership
679   //  
680   Int_t bin = FindEtaBin(eta);
681   if (bin <= 0) { 
682     AliError(Form("eta=%f is out of range [%f,%f]", 
683                   eta, fEtaAxis.GetXmin(), fEtaAxis.GetXmax()));
684     return kFALSE;
685   }
686
687   return SetFit(d, r, bin, fit);
688 }
689 //____________________________________________________________________
690 Bool_t
691 AliFMDCorrELossFit::SetFit(UShort_t  d,      Char_t    r, 
692                            Double_t  eta, 
693                            Int_t     quality,UShort_t  n, 
694                            Double_t  chi2,   UShort_t  nu, 
695                            Double_t  c,      Double_t  ec, 
696                            Double_t  delta,  Double_t  edelta, 
697                            Double_t  xi,     Double_t  exi,
698                            Double_t  sigma,  Double_t  esigma, 
699                            Double_t  sigman, Double_t  esigman, 
700                            Double_t* a,      Double_t* ea)
701 {
702   // 
703   // Set the fit parameters from a function 
704   // 
705   // Parameters:
706   //    d         Detector number
707   //    r         Ring identifier 
708   //    eta       Eta value
709   //    quality   Quality flag
710   //    n         @f$ N@f$ - Number of fitted peaks
711   //    chi2      @f$ \chi^2 @f$
712   //    nu        @f$ \nu @f$ - number degrees of freedom
713   //    c         @f$ C@f$ - scale constant
714   //    ec        @f$ \delta C@f$ - error on @f$ C@f$ 
715   //    delta     @f$ \Delta@f$ - most probable value
716   //    edelta    @f$ \delta\Delta@f$ - error on @f$\Delta@f$ 
717   //    xi        @f$ \xi@f$ - Landau width               
718   //    exi       @f$ \delta\xi@f$ - error on @f$\xi@f$ 
719   //    sigma     @f$ \sigma@f$ - Gaussian width
720   //    esigma    @f$ \delta\sigma@f$ - error on @f$\sigma@f$ 
721   //    sigman    @f$ \sigma_n@f$ - Noise width           
722   //    esigman   @f$ \delta\sigma_n@f$ - error on @f$\sigma_n@f$ 
723   //    a         Array of @f$ N-1@f$ weights @f$ a_i@f$ for 
724   //                  @f$ i=2,\ldots@f$ 
725   //    ea        Array of @f$ N-1@f$ errors on weights @f$ a_i@f$ for 
726   //                  @f$ i=2,\ldots@f$ 
727   //
728   ELossFit* e = new ELossFit(quality, n, 
729                              chi2,    nu,
730                              c,       ec,
731                              delta,   edelta,
732                              xi,      exi,
733                              sigma,   esigma,
734                              sigman,  esigman,
735                              a,       ea);
736   if (!SetFit(d, r, eta, e)) { 
737     delete e;
738     return kFALSE;
739   }
740   return kTRUE;
741 }
742 //____________________________________________________________________
743 Bool_t
744 AliFMDCorrELossFit::SetFit(UShort_t  d, Char_t r, Double_t eta, 
745                            Int_t quality, const TF1& f)
746 {
747   // 
748   // Set the fit parameters from a function 
749   // 
750   // Parameters:
751   //    d        Detector
752   //    r        Ring 
753   //    eta      Eta 
754   //    quality  Quality flag
755   //    f        Function from fit 
756   //  
757   ELossFit* e = new ELossFit(quality, f);
758   if (!SetFit(d, r, eta, e)) { 
759     delete e;
760     return kFALSE;
761   }
762   return kTRUE;
763 }
764 //____________________________________________________________________
765 AliFMDCorrELossFit::ELossFit*
766 AliFMDCorrELossFit::GetFit(UShort_t  d, Char_t r, Int_t etabin) const
767 {
768   // 
769   // Get the fit corresponding to the specified parameters 
770   // 
771   // Parameters:
772   //    d      Detector 
773   //    r      Ring 
774   //    etabin Eta bin (1 based)
775   // 
776   // Return:
777   //    Fit parameters or null in case of problems 
778   //
779   TObjArray* ringArray = GetRingArray(d, r);
780   if (!ringArray)                                   return 0;
781   if (etabin <= 0 || etabin >= fEtaAxis.GetNbins()) return 0;
782   if      (etabin > ringArray->GetEntriesFast())    return 0;
783   else if (etabin >= ringArray->GetEntriesFast())   etabin--;
784   else if (!ringArray->At(etabin))                  etabin++;
785   return static_cast<ELossFit*>(ringArray->At(etabin));
786 }
787 //____________________________________________________________________
788 AliFMDCorrELossFit::ELossFit*
789 AliFMDCorrELossFit::GetFit(UShort_t  d, Char_t r, Double_t eta) const
790 {
791   // 
792   // Find the fit corresponding to the specified parameters 
793   // 
794   // Parameters:
795   //    d   Detector 
796   //    r   Ring 
797   //    eta Eta value 
798   // 
799   // Return:
800   //    Fit parameters or null in case of problems 
801   //
802   Int_t etabin = FindEtaBin(eta);
803   return GetFit(d, r, etabin);
804 }
805
806 //____________________________________________________________________
807 AliFMDCorrELossFit::ELossFit*
808 AliFMDCorrELossFit::FindFit(UShort_t  d, Char_t r, Int_t etabin) const
809 {
810   // 
811   // Find the fit corresponding to the specified parameters 
812   // 
813   // Parameters:
814   //    d      Detector 
815   //    r      Ring 
816   //    etabin Eta bin (1 based)
817   // 
818   // Return:
819   //    Fit parameters or null in case of problems 
820   //
821   TObjArray* ringArray = GetRingArray(d, r);
822   if (!ringArray) { 
823     AliError(Form("Failed to make ring array for FMD%d%c", d, r));
824     return 0;
825   }
826   if (etabin <= 0 || etabin >= fEtaAxis.GetNbins()) { 
827     // AliError(Form("Eta bin=%3d out of bounds [%d,%d] for FMD%d%c", 
828     //               etabin, 1, fEtaAxis.GetNbins(), d, r));
829     return 0;
830   }
831   if (etabin > ringArray->GetEntriesFast()) { 
832     // AliError(Form("Eta bin=%3d out of bounds [%d,%d] for FMD%d%c", 
833     //               etabin, 1, ringArray->GetEntriesFast(), d, r));
834     return 0;
835   }
836   else if (etabin >= ringArray->GetEntriesFast()) { 
837     // AliWarning(Form("Eta bin=%3d out of bounds by +1 [%d,%d] for FMD%d%c, " 
838     //              "trying %3d", etabin, 1, ringArray->GetEntriesFast(), d, r,
839     //              etabin-1));
840     etabin--;
841   }
842   else if (!ringArray->At(etabin)) { 
843     // AliWarning(Form("Eta bin=%d has no fit for FMD%d%c, trying %03d", 
844     //              etabin, d, r, etabin+1));
845     etabin++;
846   }
847   return static_cast<ELossFit*>(ringArray->At(etabin));
848 }
849 //____________________________________________________________________
850 AliFMDCorrELossFit::ELossFit*
851 AliFMDCorrELossFit::FindFit(UShort_t  d, Char_t r, Double_t eta) const
852 {
853   // 
854   // Find the fit corresponding to the specified parameters 
855   // 
856   // Parameters:
857   //    d   Detector 
858   //    r   Ring 
859   //    eta Eta value 
860   // 
861   // Return:
862   //    Fit parameters or null in case of problems 
863   //
864   Int_t etabin = FindEtaBin(eta);
865   return FindFit(d, r, etabin);
866 }
867 //____________________________________________________________________
868 TObjArray*
869 AliFMDCorrELossFit::GetRingArray(UShort_t d, Char_t r) const
870 {
871   // 
872   // Get the ring array corresponding to the specified ring
873   // 
874   // Parameters:
875   //    d Detector 
876   //    r Ring 
877   // 
878   // Return:
879   //    Pointer to ring array, or null in case of problems
880   //
881   Int_t idx = -1;
882   switch (d) { 
883   case 1:   idx = 0; break;
884   case 2:   idx = (r == 'i' || r == 'I') ? 1 : 2; break;
885   case 3:   idx = (r == 'o' || r == 'I') ? 3 : 4; break;
886   }
887   if (idx < 0 || idx >= fRings.GetEntriesFast()) return 0;
888   return static_cast<TObjArray*>(fRings.At(idx));
889 }
890 //____________________________________________________________________
891 TObjArray*
892 AliFMDCorrELossFit::GetOrMakeRingArray(UShort_t d, Char_t r)
893 {
894   // 
895   // Get the ring array corresponding to the specified ring
896   // 
897   // Parameters:
898   //    d Detector 
899   //    r Ring 
900   // 
901   // Return:
902   //    Pointer to ring array, or newly created container 
903   //
904   Int_t idx = -1;
905   switch (d) { 
906   case 1:   idx = 0; break;
907   case 2:   idx = (r == 'i' || r == 'I') ? 1 : 2; break;
908   case 3:   idx = (r == 'o' || r == 'I') ? 3 : 4; break;
909   }
910   if (idx < 0) return 0;
911   if (idx >= fRings.GetEntriesFast() || !fRings.At(idx)) {
912     TObjArray* a = new TObjArray(0);
913     // TOrdCollection* a = new TOrdCollection(fEtaAxis.GetNbins());
914     a->SetName(Form("FMD%d%c", d, r));
915     a->SetOwner();
916     fRings.AddAtAndExpand(a, idx);
917   }
918   return static_cast<TObjArray*>(fRings.At(idx));
919 }
920
921 //____________________________________________________________________
922 Double_t
923 AliFMDCorrELossFit::GetLowerBound(UShort_t  d, Char_t r, Int_t etabin,
924                                   Double_t f) const
925 {
926   ELossFit* fit = GetFit(d, r, etabin);
927   if (!fit) return -1024;
928   return fit->GetLowerBound(f);
929 }
930 //____________________________________________________________________
931 Double_t
932 AliFMDCorrELossFit::GetLowerBound(UShort_t  d, Char_t r, Double_t eta,
933                                   Double_t f) const
934 {
935   Int_t bin = FindEtaBin(eta);
936   if (bin <= 0) return -1024;
937   return GetLowerBound(d, r, Int_t(bin), f);
938 }
939 //____________________________________________________________________
940 Double_t
941 AliFMDCorrELossFit::GetLowerBound(UShort_t  d, Char_t r, Int_t etabin,
942                                   Double_t f, Bool_t showErrors, 
943                                   Bool_t includeSigma) const
944 {
945   ELossFit* fit = GetFit(d, r, etabin);
946   if (!fit) { 
947     if (showErrors) {
948       AliWarning(Form("No fit for FMD%d%c @ etabin=%d", d, r, etabin));
949     }
950     return -1024;
951   }
952   return fit->GetLowerBound(f, includeSigma);
953 }
954
955 //____________________________________________________________________
956 Double_t
957 AliFMDCorrELossFit::GetLowerBound(UShort_t  d, Char_t r, Double_t eta,
958                                   Double_t f, Bool_t showErrors, 
959                                   Bool_t includeSigma) const
960 {
961   Int_t bin = FindEtaBin(eta);
962   if (bin <= 0) { 
963     if (showErrors)
964       AliError(Form("eta=%f out of bounds for FMD%d%c", eta, d, r));
965     return -1024;
966   }
967   return GetLowerBound(d, r, bin, f, showErrors, includeSigma);
968 }
969
970 //____________________________________________________________________
971 namespace { 
972   TH1D* MakeHist(const TAxis& axis, const char* name, const char* title, 
973                  Int_t color)
974   {
975     TH1D* h = new TH1D(Form("%s_%s", name, title), 
976                        Form("%s %s", name, title), 
977                        axis.GetNbins(), axis.GetXmin(), axis.GetXmax());
978     h->SetDirectory(0);
979     h->SetMarkerStyle(20);
980     h->SetMarkerColor(color);
981     h->SetMarkerSize(0.5);
982     h->SetFillColor(color);
983     h->SetFillStyle(3001);
984     h->SetLineColor(color);
985     return h;
986   }
987 }
988
989   
990
991 #define IDX2RING(I) (i == 0 || i == 1 || i == 3 ? 'I' : 'O')
992 #define IDX2DET(I)  (i == 0 ? 1 : (i == 1 || i == 2 ? 2 : 3))
993 //____________________________________________________________________
994 TList*
995 AliFMDCorrELossFit::GetStacks(Bool_t err, Bool_t rel, UShort_t maxN) const
996 {
997   // Get a list of THStacks 
998   Int_t nRings = fRings.GetEntriesFast();
999   Int_t nPad   = 6+maxN-1; // 7 regular params, and maxN-1 weights
1000
1001   enum { 
1002     kChi2nu = 0, 
1003     kC      = 1, 
1004     kDelta  = 2, 
1005     kXi     = 3, 
1006     kSigma  = 4, 
1007     kN      = 5 
1008   };
1009   
1010   THStack* sChi2nu;
1011   THStack* sC;
1012   THStack* sDelta;
1013   THStack* sXi;
1014   THStack* sSigma;
1015   // THStack* sigman;
1016   THStack* n;
1017   TList* stacks = new TList;
1018   stacks->AddAt(sChi2nu= new THStack("chi2",   "#chi^{2}/#nu"), kChi2nu);
1019   stacks->AddAt(sC     = new THStack("c",       "C"),           kC);
1020   stacks->AddAt(sDelta = new THStack("delta",  "#Delta_{mp}"),  kDelta);
1021   stacks->AddAt(sXi    = new THStack("xi",     "#xi"),          kXi);
1022   stacks->AddAt(sSigma = new THStack("sigma",  "#sigma"),       kSigma);
1023   //stacks->AddAt(sigman= new THStack("sigman", "#sigma_{n}"),   5);
1024   stacks->AddAt(n     = new THStack("n",      "N"),            kN);
1025   for (Int_t i = 1; i <= maxN; i++) {
1026     stacks->AddAt(new THStack(Form("a_%02d", i+1), Form("a_{%d}", i+1)), kN+i);
1027   }
1028   
1029   TArrayD min(nPad);
1030   TArrayD max(nPad);
1031   min.Reset(100000);
1032   max.Reset(-100000);
1033
1034   for (Int_t i = 0; i < nRings; i++) { 
1035     if (!fRings.At(i)) continue;
1036     TObjArray* a     = static_cast<TObjArray*>(fRings.At(i));
1037     Int_t      nFits = a->GetEntriesFast();
1038     Int_t      color = AliForwardUtil::RingColor(IDX2DET(i), IDX2RING(i));
1039
1040     TH1D* hChi    = MakeHist(fEtaAxis,a->GetName(), "chi2",   color);
1041     TH1D* hC      = MakeHist(fEtaAxis,a->GetName(), "c",      color);
1042     TH1D* hDelta  = MakeHist(fEtaAxis,a->GetName(), "delta",  color);
1043     TH1D* hXi     = MakeHist(fEtaAxis,a->GetName(), "xi",     color);
1044     TH1D* hSigma  = MakeHist(fEtaAxis,a->GetName(), "sigma",  color);
1045     // TH1D* hSigmaN = MakeHist(fEtaAxis,a->GetName(), "sigman", color);
1046     TH1D* hN      = MakeHist(fEtaAxis,a->GetName(), "n",      color);
1047     TH1D* hA[maxN];
1048     const char* ho = (rel || !err ? "hist" : "e");
1049     sChi2nu->Add(hChi,   "hist"); // 0
1050     sC     ->Add(hC,     ho);     // 1
1051     sDelta ->Add(hDelta, ho);     // 2
1052     sXi    ->Add(hXi,    ho);     // 3
1053     sSigma ->Add(hSigma, ho);     // 4
1054     // sigman->Add(hSigmaN,ho);     // 5
1055     n     ->Add(hN,     "hist"); // 5
1056     hChi->SetFillColor(color);
1057     hChi->SetFillStyle(3001);
1058     hN->SetFillColor(color);
1059     hN->SetFillStyle(3001);
1060
1061     for (Int_t k = 1; k <= maxN; k++) { 
1062       hA[k-1] = MakeHist(fEtaAxis,a->GetName(), Form("a%02d",k+1), color);
1063       static_cast<THStack*>(stacks->At(kN+k))->Add(hA[k-1], ho);
1064     }
1065                           
1066     for (Int_t j = 0; j < nFits; j++) {
1067       ELossFit* f = static_cast<ELossFit*>(a->At(j));
1068       if (!f) continue;
1069
1070       Int_t     b      = f->fBin;
1071       Int_t     nW     = f->FindMaxWeight();
1072       Double_t  vChi2nu = (f->fNu <= 0 ? 0 : f->fChi2 / f->fNu);
1073       Double_t  vC      = (rel ? (f->fC     >0 ?f->fEC     /f->fC      :0) 
1074                            : f->fC);
1075       Double_t  vDelta  = (rel ? (f->fDelta >0 ?f->fEDelta /f->fDelta  :0) 
1076                            : f->fDelta);
1077       Double_t  vXi     = (rel ? (f->fXi    >0 ?f->fEXi    /f->fXi     :0) 
1078                            : f->fXi);
1079       Double_t  vSigma  = (rel ? (f->fSigma >0 ?f->fESigma /f->fSigma  :0) 
1080                            : f->fSigma);
1081       // Double_t  sigman = (rel ? (f->fSigmaN>0 ?f->fESigmaN/f->fSigmaN :0) 
1082       //                     : f->SigmaN); 
1083       hChi   ->SetBinContent(b, vChi2nu);
1084       hN     ->SetBinContent(b, nW);
1085       hC     ->SetBinContent(b, vC);
1086       hDelta ->SetBinContent(b, vDelta);
1087       hXi    ->SetBinContent(b, vXi);
1088       hSigma ->SetBinContent(b, vSigma);
1089       if (vChi2nu > 1e-12) {
1090         min[kChi2nu] = TMath::Min(min[kChi2nu], vChi2nu);
1091         max[kChi2nu] = TMath::Max(max[kChi2nu], vChi2nu);
1092       }
1093       if (vC > 1e-12) {
1094         min[kC] = TMath::Min(min[kC], vC);
1095         max[kC] = TMath::Max(max[kC], vC);
1096       }
1097       if (vDelta > 1e-12) {
1098         min[kDelta] = TMath::Min(min[kDelta], vDelta);
1099         max[kDelta] = TMath::Max(max[kDelta], vDelta);
1100       }
1101       if (vXi > 1e-12) {
1102         min[kXi] = TMath::Min(min[kXi], vXi);
1103         max[kXi] = TMath::Max(max[kXi], vXi);
1104       }
1105       if (vSigma > 1e-12) {
1106         min[kSigma] = TMath::Min(min[kSigma], vSigma);
1107         max[kSigma] = TMath::Max(max[kSigma], vSigma);
1108       }
1109       if (nW > 1e-12) { 
1110         min[kN] = TMath::Min(min[kN], Double_t(nW));
1111         max[kN] = TMath::Max(max[kN], Double_t(nW));
1112       }
1113       // hSigmaN->SetBinContent(b,  sigman);
1114       if (!rel) {
1115         hC     ->SetBinError(b, f->fEC);
1116         hDelta ->SetBinError(b, f->fEDelta);
1117         hXi    ->SetBinError(b, f->fEXi);
1118         hSigma ->SetBinError(b, f->fESigma);
1119         // hSigmaN->SetBinError(b, f->fESigmaN);
1120       }
1121       for (Int_t k = 0; k < f->fN-1 && k < maxN; k++) { 
1122         Double_t vA = (rel ? (f->fA[k] > 0 ? f->fEA[k] / f->fA[k] : 0) 
1123                        : f->fA[k]);
1124         hA[k]->SetBinContent(b, vA);
1125         if (!rel) hA[k]->SetBinError(b, f->fEA[k]);
1126         if (vA > 1e-12) {
1127           min[kN+1+k] = TMath::Min(min[kN+1+k], vA);
1128           max[kN+1+k] = TMath::Max(max[kN+1+k], vA);
1129         }
1130       }
1131     }
1132   }
1133   return stacks;
1134 }
1135
1136 //____________________________________________________________________
1137 void
1138 AliFMDCorrELossFit::Draw(Option_t* option)
1139 {
1140   // 
1141   // Draw this object 
1142   // 
1143   // Parameters:
1144   //    option Options.  Possible values are 
1145   //  - err Plot error bars 
1146   //
1147   TString opt(Form("nostack %s", option));
1148   opt.ToLower();
1149   Bool_t  rel = (opt.Contains("relative"));
1150   Bool_t  err = (opt.Contains("error"));
1151   if (rel) opt.ReplaceAll("relative","");
1152   if (err) opt.ReplaceAll("error","");
1153
1154   UShort_t maxN   = 0;
1155   Int_t nRings = fRings.GetEntriesFast();
1156   for (Int_t i = 0; i < nRings; i++) { 
1157     if (!fRings.At(i)) continue;
1158     TObjArray* a     = static_cast<TObjArray*>(fRings.At(i));
1159     Int_t      nFits = a->GetEntriesFast();
1160
1161     for (Int_t j = 0; j < nFits; j++) {
1162       ELossFit* fit = static_cast<ELossFit*>(a->At(j));
1163       if (!fit) continue;
1164       maxN          = TMath::Max(maxN, UShort_t(fit->fN));
1165     }
1166   }
1167   // AliInfo(Form("Maximum N is %d", maxN));
1168   Int_t nPad = 6+maxN-1; // 7 regular params, and maxN-1 weights
1169   TVirtualPad* pad = gPad;
1170   pad->Divide(2, (nPad+1)/2, 0.1, 0, 0);
1171
1172   TList* stacks = GetStacks(err, rel, maxN);
1173
1174   Int_t nPad2 = (nPad+1) / 2;
1175   for (Int_t i = 0; i < nPad; i++) {
1176     Int_t iPad = 1 + i/nPad2 + 2 * (i % nPad2);
1177     TVirtualPad* p = pad->cd(iPad);
1178     p->SetLeftMargin(.15);
1179     p->SetFillColor(0);
1180     p->SetFillStyle(0);
1181     p->SetGridx();
1182     p->SetGridy();
1183     if (rel && i != 0 && i != 6 && i != 5 && i != 4) p->SetLogy();
1184
1185
1186     THStack* stack = static_cast<THStack*>(stacks->At(i));
1187
1188     // Double_t powMax = TMath::Log10(max[i]);
1189     // Double_t powMin = min[i] <= 0 ? powMax : TMath::Log10(min[i]);
1190     // if (powMax-powMin > 2. && min[i] != 0) p->SetLogy();
1191
1192     // stack->SetMinimum(min[i]);
1193     // stack->SetMaximum(max[i]);
1194     stack->Draw(opt.Data());
1195
1196     TString tit(stack->GetTitle());
1197     if (rel && i != 0 && i != 6)
1198       tit = Form("#delta %s/%s", tit.Data(), tit.Data());
1199     TH1*   hist  = stack->GetHistogram();
1200     TAxis* yaxis = hist->GetYaxis();
1201     yaxis->SetTitle(tit.Data());
1202     yaxis->SetTitleSize(0.15);
1203     yaxis->SetLabelSize(0.08);
1204     yaxis->SetTitleOffset(0.35);
1205     yaxis->SetTitleFont(132);
1206     yaxis->SetLabelFont(132);
1207     yaxis->SetNdivisions(5);
1208
1209
1210     TAxis* xaxis = stack->GetHistogram()->GetXaxis();
1211     xaxis->SetTitle("#eta");
1212     xaxis->SetTitleSize(0.15);
1213     xaxis->SetLabelSize(0.08);
1214     xaxis->SetTitleOffset(0.35);
1215     xaxis->SetTitleFont(132);
1216     xaxis->SetLabelFont(132);
1217     xaxis->SetNdivisions(10);
1218
1219     stack->Draw(opt.Data());
1220   }
1221   pad->cd();      
1222 }
1223
1224 //____________________________________________________________________
1225 void
1226 AliFMDCorrELossFit::Print(Option_t* option) const
1227 {
1228   // 
1229   // Print this object.  
1230   // 
1231   // Parameters:
1232   //    option Options 
1233   //   - R   Print recursive  
1234   //
1235   //
1236   TString opt(option);
1237   opt.ToUpper();
1238   Int_t nRings = fRings.GetEntriesFast();
1239   bool recurse = opt.Contains("R");
1240
1241   std::cout << "Low cut in fit range: " << fLowCut << "\n"
1242             << "Eta axis:             " << fEtaAxis.GetNbins() 
1243             << " bins, range [" << fEtaAxis.GetXmin() << "," 
1244             << fEtaAxis.GetXmax() << "]" << std::endl;
1245   
1246   for (Int_t i = 0; i < nRings; i++) { 
1247     if (!fRings.At(i)) continue;
1248     TObjArray* a     = static_cast<TObjArray*>(fRings.At(i));
1249     Int_t      nFits = a->GetEntriesFast();
1250
1251     std::cout << a->GetName() << " [" << nFits << " entries]" 
1252               << (recurse ? ":\n" : "\t");
1253     Int_t min = fEtaAxis.GetNbins()+1;
1254     Int_t max = 0;
1255     for (Int_t j = 0; j < nFits; j++) {
1256       if (!a->At(j)) continue;
1257       
1258       min = TMath::Min(j, min);
1259       max = TMath::Max(j, max);
1260
1261       if (recurse) {
1262         std::cout << "Bin # " << j << "\t";
1263         ELossFit* fit = static_cast<ELossFit*>(a->At(j));
1264         fit->Print(option);
1265       }
1266     }
1267     if (!recurse) 
1268       std::cout << " bin range: " << std::setw(3) << min 
1269                 << "-" << std::setw(3) << max << " " << std::setw(3) 
1270                 << (max-min+1) << " bins" << std::endl;
1271   }
1272 }
1273
1274 //____________________________________________________________________
1275 void
1276 AliFMDCorrELossFit::Browse(TBrowser* b)
1277 {
1278   // 
1279   // Browse this object 
1280   // 
1281   // Parameters:
1282   //    b 
1283   //
1284   b->Add(&fRings);
1285   b->Add(&fEtaAxis);
1286 }
1287
1288
1289
1290 //____________________________________________________________________
1291 //
1292 // EOF
1293 //