// Class to do the energy correction of FMD ESD data
//
#include "AliFMDEnergyFitter.h"
+#include "AliForwardCorrectionManager.h"
#include <AliESDFMD.h>
#include <TAxis.h>
#include <TList.h>
#include <TH1.h>
#include <TF1.h>
#include <TMath.h>
-#include "AliFMDAnaParameters.h"
#include <AliLog.h>
#include <TClonesArray.h>
#include <TFitResult.h>
#include <THStack.h>
+#include <TROOT.h>
+#include <iostream>
+#include <iomanip>
ClassImp(AliFMDEnergyFitter)
#if 0
; // This is for Emacs
#endif
-
-#define N_A(N) (4+N-2)
-#define N2_A(N) (4+(N-2)*3)
-#define N2_D(N) (4+(N-2)*3+1)
-#define N2_X(N) (4+(N-2)*3+2)
-
-//____________________________________________________________________
namespace {
- Double_t
- NLandau(Double_t* xp, Double_t* pp)
- {
- Double_t e = xp[0];
- Double_t constant = pp[0];
- Double_t mpv = pp[1];
- Double_t fwhm = pp[2];
- Int_t n = Int_t(pp[3]);
- Double_t result = 0;
- for (Int_t i = 1; i <= n; i++) {
- Double_t mpvi = i*(mpv+fwhm*TMath::Log(i));
- Double_t fwhmi = i*fwhm;
- Double_t ai = (i == 1 ? 1 : pp[N_A(i)]);
- result += ai * TMath::Landau(e,mpvi,fwhmi,kTRUE);
- }
- result *= constant;
- return result;
- }
-
- Double_t
- NLandau2(Double_t* xp, Double_t* pp)
- {
- Double_t e = xp[0];
- Double_t constant = pp[0];
- Double_t mpv = pp[1];
- Double_t fwhm = pp[2];
- Int_t n = Int_t(pp[3]);
- Double_t result = TMath::Landau(e,mpv,fwhm,kTRUE);
- for (Int_t i = 2; i <= n; i++) {
- Double_t ai = pp[N2_A(i)];
- Double_t mpvi = pp[N2_D(i)];
- Double_t fwhmi = pp[N2_X(i)];
- result += ai * TMath::Landau(e,mpvi,fwhmi,kTRUE);
- }
- result *= constant;
- return result;
- }
-
- Double_t
- TripleLandau(Double_t *x, Double_t *par)
- {
- Double_t energy = x[0];
- Double_t constant = par[0];
- Double_t mpv = par[1];
- Double_t fwhm = par[2];
- Double_t alpha = par[3];
- Double_t beta = par[4];
- Double_t mpv2 = 2*mpv+2*fwhm*TMath::Log(2);
- Double_t mpv3 = 3*mpv+3*fwhm*TMath::Log(3);
-
- Double_t f = constant * (TMath::Landau(energy,mpv,fwhm,kTRUE)+
- alpha * TMath::Landau(energy,mpv2,2*fwhm,kTRUE)+
- beta * TMath::Landau(energy,mpv3,3*fwhm,kTRUE));
-
- return f;
- }
- Double_t
- DoubleLandau(Double_t *x, Double_t *par)
- {
- Double_t energy = x[0];
- Double_t constant = par[0];
- Double_t mpv = par[1];
- Double_t fwhm = par[2];
- Double_t alpha = par[3];
- Double_t mpv2 = 2*mpv+2*fwhm*TMath::Log(2);
-
- Double_t f = constant * (TMath::Landau(energy,mpv,fwhm,kTRUE)+
- alpha * TMath::Landau(energy,mpv2,2*fwhm,kTRUE));
-
- return f;
- }
- Double_t
- SingleLandau(Double_t *x, Double_t *par)
- {
- Double_t energy = x[0];
- Double_t constant = par[0];
- Double_t mpv = par[1];
- Double_t fwhm = par[2];
-
- Double_t f = constant * TMath::Landau(energy,mpv,fwhm,kTRUE);
-
- return f;
- }
+ const char* fgkEDistFormat = "%s_etabin%03d";
}
+
//____________________________________________________________________
AliFMDEnergyFitter::AliFMDEnergyFitter()
: TNamed(),
fRingHistos(),
- fLowCut(0.3),
- fNLandau(3),
- fMinEntries(100),
- fBinsToSubtract(4),
+ fLowCut(0.4),
+ fNParticles(5),
+ fMinEntries(1000),
+ fFitRangeBinWidth(4),
fDoFits(false),
+ fDoMakeObject(false),
fEtaAxis(),
+ fCentralityAxis(),
fMaxE(10),
fNEbins(300),
fUseIncreasingBins(true),
+ fMaxRelParError(.25),
+ fMaxChi2PerNDF(20),
+ fMinWeight(1e-7),
fDebug(0)
-{}
+{
+ //
+ // Default Constructor - do not use
+ //
+}
//____________________________________________________________________
AliFMDEnergyFitter::AliFMDEnergyFitter(const char* title)
: TNamed("fmdEnergyFitter", title),
fRingHistos(),
- fLowCut(0.3),
- fNLandau(3),
- fMinEntries(100),
- fBinsToSubtract(4),
+ fLowCut(0.4),
+ fNParticles(5),
+ fMinEntries(1000),
+ fFitRangeBinWidth(4),
fDoFits(false),
+ fDoMakeObject(false),
fEtaAxis(0,0,0),
+ fCentralityAxis(0,0,0),
fMaxE(10),
fNEbins(300),
fUseIncreasingBins(true),
+ fMaxRelParError(.25),
+ fMaxChi2PerNDF(20),
+ fMinWeight(1e-7),
fDebug(3)
{
+ //
+ // Constructor
+ //
+ // Parameters:
+ // title Title of object - not significant
+ //
fEtaAxis.SetName("etaAxis");
fEtaAxis.SetTitle("#eta");
+ fCentralityAxis.SetName("centralityAxis");
+ fCentralityAxis.SetTitle("Centrality [%]");
fRingHistos.Add(new RingHistos(1, 'I'));
fRingHistos.Add(new RingHistos(2, 'I'));
fRingHistos.Add(new RingHistos(2, 'O'));
: TNamed(o),
fRingHistos(),
fLowCut(o.fLowCut),
- fNLandau(o.fNLandau),
+ fNParticles(o.fNParticles),
fMinEntries(o.fMinEntries),
- fBinsToSubtract(o.fBinsToSubtract),
+ fFitRangeBinWidth(o.fFitRangeBinWidth),
fDoFits(o.fDoFits),
+ fDoMakeObject(o.fDoMakeObject),
fEtaAxis(o.fEtaAxis),
+ fCentralityAxis(o.fCentralityAxis),
fMaxE(o.fMaxE),
fNEbins(o.fNEbins),
fUseIncreasingBins(o.fUseIncreasingBins),
+ fMaxRelParError(o.fMaxRelParError),
+ fMaxChi2PerNDF(o.fMaxChi2PerNDF),
+ fMinWeight(o.fMinWeight),
fDebug(o.fDebug)
{
+ //
+ // Copy constructor
+ //
+ // Parameters:
+ // o Object to copy from
+ //
TIter next(&o.fRingHistos);
TObject* obj = 0;
while ((obj = next())) fRingHistos.Add(obj);
//____________________________________________________________________
AliFMDEnergyFitter::~AliFMDEnergyFitter()
{
+ //
+ // Destructor
+ //
fRingHistos.Delete();
}
AliFMDEnergyFitter&
AliFMDEnergyFitter::operator=(const AliFMDEnergyFitter& o)
{
+ //
+ // Assignment operator
+ //
+ // Parameters:
+ // o Object to assign from
+ //
+ // Return:
+ // Reference to this
+ //
TNamed::operator=(o);
fLowCut = o.fLowCut;
- fNLandau = o.fNLandau;
+ fNParticles = o.fNParticles;
fMinEntries = o.fMinEntries;
- fBinsToSubtract= o.fBinsToSubtract;
+ fFitRangeBinWidth= o.fFitRangeBinWidth;
fDoFits = o.fDoFits;
- fEtaAxis.Set(o.fEtaAxis.GetNbins(),o.fEtaAxis.GetXmin(),o.fEtaAxis.GetXmax());
+ fDoMakeObject = o.fDoMakeObject;
+ fEtaAxis.Set(o.fEtaAxis.GetNbins(),
+ o.fEtaAxis.GetXmin(),
+ o.fEtaAxis.GetXmax());
+ if (o.fCentralityAxis.GetXbins()) {
+ const TArrayD* bins = o.fCentralityAxis.GetXbins();
+ fCentralityAxis.Set(bins->GetSize()-1,bins->GetArray());
+ }
+ else
+ fCentralityAxis.Set(o.fCentralityAxis.GetNbins(),
+ o.fCentralityAxis.GetXmin(),
+ o.fCentralityAxis.GetXmax());
fDebug = o.fDebug;
+ fMaxRelParError= o.fMaxRelParError;
+ fMaxChi2PerNDF = o.fMaxChi2PerNDF;
+ fMinWeight = o.fMinWeight;
fRingHistos.Delete();
TIter next(&o.fRingHistos);
AliFMDEnergyFitter::RingHistos*
AliFMDEnergyFitter::GetRingHistos(UShort_t d, Char_t r) const
{
+ //
+ // Get the ring histogram container
+ //
+ // Parameters:
+ // d Detector
+ // r Ring
+ //
+ // Return:
+ // Ring histogram container
+ //
Int_t idx = -1;
switch (d) {
case 1: idx = 0; break;
void
AliFMDEnergyFitter::Init(const TAxis& eAxis)
{
+ //
+ // Initialise the task
+ //
+ // Parameters:
+ // etaAxis The eta axis to use. Note, that if the eta axis
+ // has already been set (using SetEtaAxis), then this parameter will be
+ // ignored
+ //
if (fEtaAxis.GetNbins() == 0 ||
- fEtaAxis.GetXmin() == fEtaAxis.GetXmax())
+ TMath::Abs(fEtaAxis.GetXmax() - fEtaAxis.GetXmin()) < 1e-7)
SetEtaAxis(eAxis);
+ if (fCentralityAxis.GetNbins() == 0) {
+ UShort_t n = 12;
+ Double_t bins[] = { 0., 5., 10., 15., 20., 30.,
+ 40., 50., 60., 70., 80., 100. };
+ SetCentralityAxis(n, bins);
+ }
TIter next(&fRingHistos);
RingHistos* o = 0;
while ((o = static_cast<RingHistos*>(next())))
- o->Init(fEtaAxis, fMaxE, fNEbins, fUseIncreasingBins);
+ o->Init(fEtaAxis, fCentralityAxis, fMaxE, fNEbins, fUseIncreasingBins);
}
//____________________________________________________________________
void
AliFMDEnergyFitter::SetEtaAxis(const TAxis& eAxis)
{
+ //
+ // Set the eta axis to use. This will force the code to use this
+ // eta axis definition - irrespective of whatever axis is passed to
+ // the Init member function. Therefore, this member function can be
+ // used to force another eta axis than one found in the correction
+ // objects.
+ //
+ // Parameters:
+ // etaAxis Eta axis to use
+ //
SetEtaAxis(eAxis.GetNbins(),eAxis.GetXmin(),eAxis.GetXmax());
}
//____________________________________________________________________
void
AliFMDEnergyFitter::SetEtaAxis(Int_t nBins, Double_t etaMin, Double_t etaMax)
{
+ //
+ // Set the eta axis to use. This will force the code to use this
+ // eta axis definition - irrespective of whatever axis is passed to
+ // the Init member function. Therefore, this member function can be
+ // used to force another eta axis than one found in the correction
+ // objects.
+ //
+ // Parameters:
+ // nBins Number of bins
+ // etaMin Minimum of the eta axis
+ // etaMax Maximum of the eta axis
+ //
fEtaAxis.Set(nBins,etaMin,etaMax);
}
+//____________________________________________________________________
+void
+AliFMDEnergyFitter::SetCentralityAxis(UShort_t n, Double_t* bins)
+{
+ fCentralityAxis.Set(n-1, bins);
+}
+
//____________________________________________________________________
Bool_t
-AliFMDEnergyFitter::Accumulate(const AliESDFMD& input,
+AliFMDEnergyFitter::Accumulate(const AliESDFMD& input,
+ Double_t cent,
Bool_t empty)
{
+ //
+ // Fitter the input AliESDFMD object
+ //
+ // Parameters:
+ // input Input
+ // cent Centrality
+ // empty Whether the event is 'empty'
+ //
+ // Return:
+ // True on success, false otherwise
+ //
+ Int_t icent = fCentralityAxis.FindBin(cent);
+ if (icent < 1 || icent > fCentralityAxis.GetNbins()) icent = 0;
+
for(UShort_t d = 1; d <= 3; d++) {
Int_t nRings = (d == 1 ? 1 : 2);
for (UShort_t q = 0; q < nRings; q++) {
Int_t ieta = fEtaAxis.FindBin(eta1);
if (ieta < 1 || ieta > fEtaAxis.GetNbins()) continue;
- histos->Fill(empty, ieta-1, mult);
+ histos->Fill(empty, ieta-1, icent, mult);
} // for strip
} // for sector
//____________________________________________________________________
void
-AliFMDEnergyFitter::Fit(TList* dir)
+AliFMDEnergyFitter::Fit(const TList* dir)
{
+ //
+ // Scale the histograms to the total number of events
+ //
+ // Parameters:
+ // dir Where the histograms are
+ //
if (!fDoFits) return;
TList* d = static_cast<TList*>(dir->FindObject(GetName()));
// +1 for chi^2
// +3 for 1 landau
// +1 for N
- // +fNLandau-1 for weights
- Int_t nStack = 1+3+1+fNLandau-1;
+ // +fNParticles-1 for weights
+ Int_t nStack = kN+fNParticles;
THStack* stack[nStack];
- stack[0] = new THStack("chi2", "#chi^{2}/#nu");
- stack[1] = new THStack("c", "constant");
- stack[2] = new THStack("mpv", "#Delta_{p}");
- stack[3] = new THStack("w", "w");
- stack[4] = new THStack("n", "# of Landaus");
- for (Int_t i = 2; i <= fNLandau; i++)
- stack[i-1+4] = new THStack(Form("a%d", i),
- Form("Weight of %d signal", i));
+ stack[0] = new THStack("chi2", "#chi^{2}/#nu");
+ stack[kC +1] = new THStack("c", "Constant");
+ stack[kDelta +1] = new THStack("delta", "#Delta_{p}");
+ stack[kXi +1] = new THStack("xi", "#xi");
+ stack[kSigma +1] = new THStack("sigma", "#sigma");
+ stack[kSigmaN+1] = new THStack("sigman", "#sigma_{n}");
+ stack[kN +1] = new THStack("n", "# Particles");
+ for (Int_t i = 2; i <= fNParticles; i++)
+ stack[kN+i] = new THStack(Form("a%d", i), Form("a_{%d}", i));
for (Int_t i = 0; i < nStack; i++)
d->Add(stack[i]);
TIter next(&fRingHistos);
RingHistos* o = 0;
while ((o = static_cast<RingHistos*>(next()))) {
- TObjArray* l = o->Fit(d, fEtaAxis, fLowCut, fNLandau,
- fMinEntries, fBinsToSubtract);
+ TObjArray* l = o->Fit(d, fEtaAxis, fLowCut, fNParticles,
+ fMinEntries, fFitRangeBinWidth,
+ fMaxRelParError, fMaxChi2PerNDF);
if (!l) continue;
for (Int_t i = 0; i < l->GetEntriesFast(); i++) {
stack[i % nStack]->Add(static_cast<TH1*>(l->At(i)));
}
}
+
+ if (!fDoMakeObject) return;
+
+ MakeCorrectionsObject(d);
+}
+
+//____________________________________________________________________
+void
+AliFMDEnergyFitter::MakeCorrectionsObject(TList* d)
+{
+ //
+ // Generate the corrections object
+ //
+ // Parameters:
+ // dir List to analyse
+ //
+ AliForwardCorrectionManager& mgr = AliForwardCorrectionManager::Instance();
+
+ AliFMDCorrELossFit* obj = new AliFMDCorrELossFit;
+ obj->SetEtaAxis(fEtaAxis);
+ obj->SetLowCut(fLowCut);
+
+ TIter next(&fRingHistos);
+ RingHistos* o = 0;
+ while ((o = static_cast<RingHistos*>(next()))) {
+ o->FindBestFits(d, *obj, fEtaAxis, fMaxRelParError,
+ fMaxChi2PerNDF, fMinWeight);
+ }
+
+ TString oName(mgr.GetObjectName(AliForwardCorrectionManager::kELossFits));
+ TString fileName(mgr.GetFilePath(AliForwardCorrectionManager::kELossFits));
+ AliInfo(Form("Object %s created in output - should be extracted and copied "
+ "to %s", oName.Data(), fileName.Data()));
+ d->Add(new TNamed("filename", fileName.Data()));
+ d->Add(obj, oName.Data());
}
//____________________________________________________________________
void
AliFMDEnergyFitter::DefineOutput(TList* dir)
{
+ //
+ // Define the output histograms. These are put in a sub list of the
+ // passed list. The histograms are merged before the parent task calls
+ // AliAnalysisTaskSE::Terminate
+ //
+ // Parameters:
+ // dir Directory to add to
+ //
TList* d = new TList;
d->SetName(GetName());
dir->Add(d);
void
AliFMDEnergyFitter::SetDebug(Int_t dbg)
{
+ //
+ // Set the debug level. The higher the value the more output
+ //
+ // Parameters:
+ // dbg Debug level
+ //
fDebug = dbg;
TIter next(&fRingHistos);
RingHistos* o = 0;
while ((o = static_cast<RingHistos*>(next())))
o->fDebug = dbg;
}
+
+//____________________________________________________________________
+void
+AliFMDEnergyFitter::Print(Option_t*) const
+{
+ //
+ // Print information
+ //
+ // Parameters:
+ // option Not used
+ //
+ char ind[gROOT->GetDirLevel()+1];
+ for (Int_t i = 0; i < gROOT->GetDirLevel(); i++) ind[i] = ' ';
+ ind[gROOT->GetDirLevel()] = '\0';
+ std::cout << ind << ClassName() << ": " << GetName() << '\n'
+ << ind << " Low cut: " << fLowCut << " E/E_mip\n"
+ << ind << " Max(particles): " << fNParticles << '\n'
+ << ind << " Min(entries): " << fMinEntries << '\n'
+ << ind << " Fit range: "
+ << fFitRangeBinWidth << " bins\n"
+ << ind << " Make fits: "
+ << (fDoFits ? "yes\n" : "no\n")
+ << ind << " Make object: "
+ << (fDoMakeObject ? "yes\n" : "no\n")
+ << ind << " Max E/E_mip: " << fMaxE << '\n'
+ << ind << " N bins: " << fNEbins << '\n'
+ << ind << " Increasing bins: "
+ << (fUseIncreasingBins ?"yes\n" : "no\n")
+ << ind << " max(delta p/p): " << fMaxRelParError << '\n'
+ << ind << " max(chi^2/nu): " << fMaxChi2PerNDF << '\n'
+ << ind << " min(a_i): " << fMinWeight << std::endl;
+}
//====================================================================
AliFMDEnergyFitter::RingHistos::RingHistos()
fEmpty(0),
fEtaEDists(),
fList(0),
+ fFits("AliFMDCorrELossFit::ELossFit"),
fDebug(0)
-{}
+{
+ //
+ // Default CTOR
+ //
+}
//____________________________________________________________________
AliFMDEnergyFitter::RingHistos::RingHistos(UShort_t d, Char_t r)
fEmpty(0),
fEtaEDists(),
fList(0),
+ fFits("AliFMDCorrELossFit::ELossFit"),
fDebug(0)
{
+ //
+ // Constructor
+ //
+ // Parameters:
+ // d detector
+ // r ring
+ //
fEtaEDists.SetName("EDists");
}
//____________________________________________________________________
fEmpty(o.fEmpty),
fEtaEDists(),
fList(0),
+ fFits("AliFMDCorrELossFit::ELossFit"),
fDebug(0)
{
+ //
+ // Copy constructor
+ //
+ // Parameters:
+ // o Object to copy from
+ //
+ fFits.Clear();
TIter next(&o.fEtaEDists);
TObject* obj = 0;
while ((obj = next())) fEtaEDists.Add(obj->Clone());
AliFMDEnergyFitter::RingHistos&
AliFMDEnergyFitter::RingHistos::operator=(const RingHistos& o)
{
+ //
+ // Assignment operator
+ //
+ // Parameters:
+ // o Object to assign from
+ //
+ // Return:
+ // Reference to this
+ //
AliForwardUtil::RingHistos::operator=(o);
if (fEDist) delete fEDist;
if (fEmpty) delete fEmpty;
fEtaEDists.Delete();
if (fList) fList->Delete();
+ fFits.Clear();
fEDist = static_cast<TH1D*>(o.fEDist->Clone());
fEmpty = static_cast<TH1D*>(o.fEmpty->Clone());
//____________________________________________________________________
AliFMDEnergyFitter::RingHistos::~RingHistos()
{
+ //
+ // Destructor
+ //
if (fEDist) delete fEDist;
fEtaEDists.Delete();
}
//____________________________________________________________________
void
-AliFMDEnergyFitter::RingHistos::Fill(Bool_t empty, Int_t ieta, Double_t mult)
+AliFMDEnergyFitter::RingHistos::Fill(Bool_t empty, Int_t ieta,
+ Int_t /* icent */, Double_t mult)
{
+ //
+ // Fill histogram
+ //
+ // Parameters:
+ // empty True if event is empty
+ // ieta Eta bin (0-based)
+ // icent Centrality bin (1-based)
+ // mult Signal
+ //
if (empty) {
fEmpty->Fill(mult);
return;
}
fEDist->Fill(mult);
+ // Here, we should first look up in a centrality array, and then in
+ // that array look up the eta bin - or perhaps we should do it the
+ // other way around?
if (ieta < 0 || ieta >= fEtaEDists.GetEntries()) return;
TH1D* h = static_cast<TH1D*>(fEtaEDists.At(ieta));
AliFMDEnergyFitter::RingHistos::MakeIncreasingAxis(Int_t n, Double_t min,
Double_t max) const
{
- // Service function to define a logarithmic axis.
- // Parameters:
- // n Number of bins
- // min Minimum of axis
- // max Maximum of axis
+ //
+ // Make an axis with increasing bins
+ //
+ // Parameters:
+ // n Number of bins
+ // min Minimum
+ // max Maximum
+ //
+ // Return:
+ // An axis with quadratically increasing bin size
+ //
+
TArrayD bins(n+1);
Double_t dx = (max-min) / n;
bins[0] = min;
//____________________________________________________________________
void
-AliFMDEnergyFitter::RingHistos::Make(Int_t ieta, Double_t emin, Double_t emax,
- Double_t deMax, Int_t nDeBins,
- Bool_t incr)
+AliFMDEnergyFitter::RingHistos::Make(Int_t ieta,
+ Double_t emin,
+ Double_t emax,
+ Double_t deMax,
+ Int_t nDeBins,
+ Bool_t incr)
{
+ //
+ // Make E/E_mip histogram
+ //
+ // Parameters:
+ // ieta Eta bin
+ // eMin Least signal
+ // eMax Largest signal
+ // deMax Maximum energy loss
+ // nDeBins Number energy loss bins
+ // incr Whether to make bins of increasing size
+ //
TH1D* h = 0;
- TString name = Form("%s_etabin%03d", fName.Data(), ieta);
+ TString name = Form(fgkEDistFormat, fName.Data(), ieta);
TString title = Form("#DeltaE/#DeltaE_{mip} for %s %+5.3f#leq#eta<%+5.3f "
"(#eta bin %d)", fName.Data(), emin, emax, ieta);
if (!incr)
fEtaEDists.AddAt(h, ieta-1);
}
//____________________________________________________________________
-TH1D*
-AliFMDEnergyFitter::RingHistos::MakePar(const char* name,
- const char* title,
- const TAxis& eta) const
+TH1D* AliFMDEnergyFitter::RingHistos::MakePar(const char* name,
+ const char* title,
+ const TAxis& eta) const
{
+ //
+ // Make a parameter histogram
+ //
+ // Parameters:
+ // name Name of histogram.
+ // title Title of histogram.
+ // eta Eta axis
+ //
+ // Return:
+ //
+ //
TH1D* h = new TH1D(Form("%s_%s", fName.Data(), name),
Form("%s for %s", title, fName.Data()),
eta.GetNbins(), eta.GetXmin(), eta.GetXmax());
Double_t val,
Double_t err) const
{
+ //
+ // Make a histogram that contains the results of the fit over the full ring
+ //
+ // Parameters:
+ // name Name
+ // title Title
+ // eta Eta axis
+ // low Least bin
+ // high Largest bin
+ // val Value of parameter
+ // err Error on parameter
+ //
+ // Return:
+ // The newly allocated histogram
+ //
Double_t xlow = eta.GetBinLowEdge(low);
Double_t xhigh = eta.GetBinUpEdge(high);
TH1D* h = new TH1D(Form("%s_%s", fName.Data(), name),
//____________________________________________________________________
void
AliFMDEnergyFitter::RingHistos::Init(const TAxis& eAxis,
+ const TAxis& /* cAxis */,
Double_t maxDE, Int_t nDEbins,
Bool_t useIncrBin)
{
+ //
+ // Initialise object
+ //
+ // Parameters:
+ // eAxis Eta axis
+ // maxDE Max energy loss to consider
+ // nDEbins Number of bins
+ // useIncrBin Whether to use an increasing bin size
+ //
fEDist = new TH1D(Form("%s_edist", fName.Data()),
Form("#DeltaE/#DeltaE_{mip} for %s", fName.Data()),
200, 0, 6);
fName.Data()), 200, 0, 6);
fList->Add(fEDist);
fList->Add(fEmpty);
+ // Here, we should make an array of centrality bins ...
+ // In that case, the structure will be
+ //
+ // -+- Ring1 -+- Centrality1 -+- Eta1
+ // | | +- Eta2
+ // ... ... ...
+ // | | +- EtaM
+ // | +- Centrality2 -+- Eta1
+ // ... ... ...
+ // | | +- EtaM
+ // ... ...
+ // | +- CentralityN -+- Eta1
+ // ... ... ...
+ // | +- EtaM
+ // -+- Ring2 -+- Centrality1 -+- Eta1
+ // | | +- Eta2
+ // ... ... ...
+ // | | +- EtaM
+ // | +- Centrality2 -+- Eta1
+ // ... ... ...
+ // | | +- EtaM
+ // ... ...
+ // | +- CentralityN -+- Eta1
+ // ... ... ...
+ // | +- EtaM
+ // ... ... ...
+ // -+- RingO -+- Centrality1 -+- Eta1
+ // | +- Eta2
+ // ... ...
+ // | +- EtaM
+ // +- Centrality2 -+- Eta1
+ // ... ...
+ // | +- EtaM
+ // ...
+ // +- CentralityN -+- Eta1
+ // ... ...
+ // +- EtaM
+ //
// fEtaEDists.Expand(eAxis.GetNbins());
for (Int_t i = 1; i <= eAxis.GetNbins(); i++) {
Double_t min = eAxis.GetBinLowEdge(i);
Double_t max = eAxis.GetBinUpEdge(i);
+ // Or may we should do it here? In that case we'd have the
+ // following structure:
+ // Ring1 -+- Eta1 -+- Centrality1
+ // | +- Centrality2
+ // ... ...
+ // | +- CentralityN
+ // +- Eta2 -+- Centrality1
+ // | +- Centrality2
+ // ... ...
+ // | +- CentralityN
+ // ...
+ // +- EtaM -+- Centrality1
+ // +- Centrality2
+ // ...
+ // +- CentralityN
Make(i, min, max, maxDE, nDEbins, useIncrBin);
}
fList->Add(&fEtaEDists);
+ // fEtaEDists.ls();
}
//____________________________________________________________________
TObjArray*
-AliFMDEnergyFitter::RingHistos::Fit(TList* dir, const TAxis& eta,
- Double_t lowCut, UShort_t nLandau,
- UShort_t minEntries,
- UShort_t minusBins) const
+AliFMDEnergyFitter::RingHistos::Fit(TList* dir,
+ const TAxis& eta,
+ Double_t lowCut,
+ UShort_t nParticles,
+ UShort_t minEntries,
+ UShort_t minusBins,
+ Double_t relErrorCut,
+ Double_t chi2nuCut) const
{
+ //
+ // Fit each histogram to up to @a nParticles particle responses.
+ //
+ // Parameters:
+ // dir Output list
+ // eta Eta axis
+ // lowCut Lower cut
+ // nParticles Max number of convolved landaus to fit
+ // minEntries Minimum number of entries
+ // minusBins Number of bins from peak to subtract to
+ // get the fit range
+ // relErrorCut Cut applied to relative error of parameter.
+ // Note, for multi-particle weights, the cut
+ // is loosend by a factor of 2
+ // chi2nuCut Cut on @f$ \chi^2/\nu@f$ -
+ // the reduced @f$\chi^2@f$
+ //
+
+ // Get our ring list from the output container
TList* l = GetOutputList(dir);
if (!l) return 0;
+ // Get the energy distributions from the output container
TList* dists = static_cast<TList*>(l->FindObject("EDists"));
if (!dists) {
AliWarning(Form("Didn't find %s_EtaEDists in %s",
return 0;
}
- TObjArray* pars = new TObjArray(3+nLandau+1);
+ // Container of the fit results histograms
+ TObjArray* pars = new TObjArray(3+nParticles+1);
pars->SetName("FitResults");
l->Add(pars);
- TH1* hChi2 = 0;
- TH1* hC = 0;
- TH1* hMpv = 0;
- TH1* hW = 0;
- TH1* hS = 0;
- TH1* hN = 0;
- TH1* hA[nLandau-1];
- pars->Add(hChi2 = MakePar("chi2", "#chi^{2}/#nu", eta));
- pars->Add(hC = MakePar("c", "Constant", eta));
- pars->Add(hMpv = MakePar("mpv", "#Delta_{p}", eta));
- pars->Add(hW = MakePar("w", "#xi", eta));
- pars->Add(hS = MakePar("s", "#sigma", eta));
- pars->Add(hN = MakePar("n", "N", eta));
- for (UShort_t i = 1; i < nLandau; i++)
+ // Result objects stored in separate list on the output
+ TH1* hChi2 = 0;
+ TH1* hC = 0;
+ TH1* hDelta = 0;
+ TH1* hXi = 0;
+ TH1* hSigma = 0;
+ TH1* hSigmaN = 0;
+ TH1* hN = 0;
+ TH1* hA[nParticles-1];
+ pars->Add(hChi2 = MakePar("chi2", "#chi^{2}/#nu", eta));
+ pars->Add(hC = MakePar("c", "Constant", eta));
+ pars->Add(hDelta = MakePar("delta", "#Delta_{p}", eta));
+ pars->Add(hXi = MakePar("xi", "#xi", eta));
+ pars->Add(hSigma = MakePar("sigma", "#sigma", eta));
+ pars->Add(hSigmaN = MakePar("sigman", "#sigma_{n}", eta));
+ pars->Add(hN = MakePar("n", "N", eta));
+ for (UShort_t i = 1; i < nParticles; i++)
pars->Add(hA[i-1] = MakePar(Form("a%d",i+1), Form("a_{%d}",i+1), eta));
Int_t high = 0;
for (Int_t i = 0; i < nDists; i++) {
TH1D* dist = static_cast<TH1D*>(dists->At(i));
- if (!dist || dist->GetEntries() <= minEntries) continue;
+ // Ignore empty histograms altoghether
+ if (!dist || dist->GetEntries() <= 0) continue;
+ // Scale to the bin-width
+ dist->Scale(1., "width");
+
+ // Normalize peak to 1
+ Double_t max = dist->GetMaximum();
+ if (max <= 0) continue;
+ dist->Scale(1/max);
+
+ // Check that we have enough entries
+ if (dist->GetEntries() <= minEntries) {
+ AliWarning(Form("Histogram at %3d (%s) has too few entries (%d <= %d)",
+ i, dist->GetName(), Int_t(dist->GetEntries()),
+ minEntries));
+ continue;
+ }
- TF1* res = FitHist(dist,lowCut,nLandau,minusBins);
+ // Now fit
+ TF1* res = FitHist(dist,lowCut,nParticles,minusBins,
+ relErrorCut,chi2nuCut);
if (!res) continue;
-
+ // dist->GetListOfFunctions()->Add(res);
+
+ // Store eta limits
low = TMath::Min(low,i+1);
high = TMath::Max(high,i+1);
+ // Get the reduced chi square
Double_t chi2 = res->GetChisquare();
Int_t ndf = res->GetNDF();
- hChi2->SetBinContent(i+1, ndf > 0 ? chi2 / ndf : 0);
- hC ->SetBinContent(i+1, res->GetParameter(0));
- hMpv->SetBinContent(i+1, res->GetParameter(1));
- hW ->SetBinContent(i+1, res->GetParameter(2));
- hN ->SetBinContent(i+1, res->GetParameter(3));
-
- hC ->SetBinError(i+1, res->GetParError(1));
- hMpv->SetBinError(i+1, res->GetParError(2));
- hW ->SetBinError(i+1, res->GetParError(2));
- // hN ->SetBinError(i, res->GetParError(3));
-
- for (UShort_t j = 0; j < nLandau-1; j++) {
- hA[j]->SetBinContent(i+1, res->GetParameter(4+j));
- hA[j]->SetBinError(i+1, res->GetParError(4+j));
+
+ // Store results of best fit in output histograms
+ res->SetLineWidth(4);
+ hChi2 ->SetBinContent(i+1, ndf > 0 ? chi2 / ndf : 0);
+ hC ->SetBinContent(i+1, res->GetParameter(kC));
+ hDelta ->SetBinContent(i+1, res->GetParameter(kDelta));
+ hXi ->SetBinContent(i+1, res->GetParameter(kXi));
+ hSigma ->SetBinContent(i+1, res->GetParameter(kSigma));
+ hSigmaN ->SetBinContent(i+1, res->GetParameter(kSigmaN));
+ hN ->SetBinContent(i+1, res->GetParameter(kN));
+
+ hC ->SetBinError(i+1, res->GetParError(kC));
+ hDelta ->SetBinError(i+1, res->GetParError(kDelta));
+ hXi ->SetBinError(i+1, res->GetParError(kXi));
+ hSigma ->SetBinError(i+1, res->GetParError(kSigma));
+ hSigmaN->SetBinError(i+1, res->GetParError(kSigmaN));
+ hN ->SetBinError(i+1, res->GetParError(kN));
+
+ for (UShort_t j = 0; j < nParticles-1; j++) {
+ hA[j]->SetBinContent(i+1, res->GetParameter(kA+j));
+ hA[j]->SetBinError(i+1, res->GetParError(kA+j));
}
}
+ // Fit the full-ring histogram
TH1* total = GetOutputHist(l, Form("%s_edist", fName.Data()));
if (total && total->GetEntries() >= minEntries) {
- TF1* res = FitHist(total,lowCut,nLandau,minusBins);
+
+ // Scale to the bin-width
+ total->Scale(1., "width");
+
+ // Normalize peak to 1
+ Double_t max = total->GetMaximum();
+ if (max > 0) total->Scale(1/max);
+
+ TF1* res = FitHist(total,lowCut,nParticles,minusBins,
+ relErrorCut,chi2nuCut);
if (res) {
+ // Make histograms for the result of this fit
Double_t chi2 = res->GetChisquare();
Int_t ndf = res->GetNDF();
- pars->Add(MakeTotal("t_chi2", "#chi^{2}/#nu", eta, low, high,
+ pars->Add(MakeTotal("t_chi2", "#chi^{2}/#nu", eta, low, high,
ndf > 0 ? chi2/ndf : 0, 0));
- pars->Add(MakeTotal("t_c", "Constant", eta, low, high,
- res->GetParameter(0),res->GetParError(0)));
- pars->Add(MakeTotal("t_mpv", "#Delta_{p}", eta, low, high,
- res->GetParameter(1),res->GetParError(1)));
- pars->Add(MakeTotal("t_w", "#xi", eta, low, high,
- res->GetParameter(2),res->GetParError(2)));
- pars->Add(MakeTotal("t_n", "N", eta, low, high,
- res->GetParameter(3),0));
- for (UShort_t j = 1; j < nLandau; j++)
- pars->Add(MakeTotal(Form("a%d_t",j+1),
- Form("a_{%d}",j+1), eta, low, high,
- res->GetParameter(3+j), res->GetParError(3+j)));
+ pars->Add(MakeTotal("t_c", "Constant", eta, low, high,
+ res->GetParameter(kC),
+ res->GetParError(kC)));
+ pars->Add(MakeTotal("t_delta", "#Delta_{p}", eta, low, high,
+ res->GetParameter(kDelta),
+ res->GetParError(kDelta)));
+ pars->Add(MakeTotal("t_xi", "#xi", eta, low, high,
+ res->GetParameter(kXi),
+ res->GetParError(kXi)));
+ pars->Add(MakeTotal("t_sigma", "#sigma", eta, low, high,
+ res->GetParameter(kSigma),
+ res->GetParError(kSigma)));
+ pars->Add(MakeTotal("t_sigman", "#sigma_{n}", eta, low, high,
+ res->GetParameter(kSigmaN),
+ res->GetParError(kSigmaN)));
+ pars->Add(MakeTotal("t_n", "N", eta, low, high,
+ res->GetParameter(kN),0));
+ for (UShort_t j = 0; j < nParticles-1; j++)
+ pars->Add(MakeTotal(Form("t_a%d",j+2),
+ Form("a_{%d}",j+2), eta, low, high,
+ res->GetParameter(kA+j),
+ res->GetParError(kA+j)));
}
}
+ // Clean up list of histogram. Histograms with no entries or
+ // no functions are deleted. We have to do this using the TObjLink
+ // objects stored in the list since ROOT cannot guaranty the validity
+ // of iterators when removing from a list - tsck. Should just implement
+ // TIter::Remove().
TObjLink* lnk = dists->FirstLink();
while (lnk) {
TH1* h = static_cast<TH1*>(lnk->GetObject());
- if (h->GetEntries() <= 0 ||
- h->GetListOfFunctions()->GetEntries() <= 0) {
+ bool remove = false;
+ if (h->GetEntries() <= 0) {
+ // AliWarning(Form("No entries in %s - removing", h->GetName()));
+ remove = true;
+ }
+ else if (h->GetListOfFunctions()->GetEntries() <= 0) {
+ // AliWarning(Form("No fuctions associated with %s - removing",
+ // h->GetName()));
+ remove = true;
+ }
+ if (remove) {
TObjLink* keep = lnk->Next();
dists->Remove(lnk);
lnk = keep;
TF1*
AliFMDEnergyFitter::RingHistos::FitHist(TH1* dist,
Double_t lowCut,
- UShort_t nLandau,
- UShort_t minusBins) const
+ UShort_t nParticles,
+ UShort_t minusBins,
+ Double_t relErrorCut,
+ Double_t chi2nuCut) const
{
+ //
+ // Fit a signal histogram. First, the bin @f$ b_{min}@f$ with
+ // maximum bin content in the range @f$ [E_{min},\infty]@f$ is
+ // found. Then the fit range is set to the bin range
+ // @f$ [b_{min}-\Delta b,b_{min}+2\Delta b]@f$, and a 1
+ // particle signal is fitted to that. The parameters of that fit
+ // is then used as seeds for a fit of the @f$ N@f$ particle response
+ // to the data in the range
+ // @f$ [b_{min}-\Delta b,N(\Delta_1+\xi_1\log(N))+2N\xi@f$
+ //
+ // Parameters:
+ // dist Histogram to fit
+ // lowCut Lower cut @f$ E_{min}@f$ on signal
+ // nParticles Max number @f$ N@f$ of convolved landaus to fit
+ // minusBins Number of bins @f$ \Delta b@f$ from peak to
+ // subtract to get the fit range
+ // relErrorCut Cut applied to relative error of parameter.
+ // Note, for multi-particle weights, the cut
+ // is loosend by a factor of 2
+ // chi2nuCut Cut on @f$ \chi^2/\nu@f$ -
+ // the reduced @f$\chi^2@f$
+ //
+ // Return:
+ // The best fit function
+ //
Double_t maxRange = 10;
+ // Create a fitter object
AliForwardUtil::ELossFitter f(lowCut, maxRange, minusBins);
f.Clear();
+
// If we are only asked to fit a single particle, return this fit,
// no matter what.
- if (nLandau == 1) {
+ if (nParticles == 1) {
TF1* r = f.Fit1Particle(dist, 0);
if (!r) return 0;
- return new TF1(*r);
+ TF1* ret = new TF1(*r);
+ dist->GetListOfFunctions()->Add(ret);
+ return ret;
}
// Fit from 2 upto n particles
- for (Int_t i = 2; i <= nLandau; i++) f.FitNParticle(dist, i, 0);
+ for (Int_t i = 2; i <= nParticles; i++) f.FitNParticle(dist, i, 0);
// Now, we need to select the best fit
- Int_t nFits = f.fFitResults.GetEntriesFast();
+ Int_t nFits = f.GetFitResults().GetEntriesFast();
TF1* good[nFits];
for (Int_t i = nFits-1; i >= 0; i--) {
- if (CheckResult(static_cast<TFitResult*>(f.fFitResults.At(i))))
- good[i] = static_cast<TF1*>(f.fFunctions.At(i));
+ good[i] = 0;
+ TF1* ff = static_cast<TF1*>(f.GetFunctions().At(i));
+ // ff->SetLineColor(Color());
+ ff->SetRange(0, maxRange);
+ dist->GetListOfFunctions()->Add(new TF1(*ff));
+ if (CheckResult(static_cast<TFitResult*>(f.GetFitResults().At(i)),
+ relErrorCut, chi2nuCut)) {
+ good[i] = ff;
+ ff->SetLineWidth(2);
+ // f.fFitResults.At(i)->Print("V");
+ }
}
// If all else fails, use the 1 particle fit
- TF1* ret = static_cast<TF1*>(f.fFunctions.At(0));
+ TF1* ret = static_cast<TF1*>(f.GetFunctions().At(0));
+
+ // Find the fit with the most valid particles
for (Int_t i = nFits-1; i > 0; i--) {
if (!good[i]) continue;
+ if (fDebug > 1) {
+ AliInfo(Form("Choosing fit with n=%d", i+1));
+ f.GetFitResults().At(i)->Print();
+ }
ret = good[i];
break;
}
- return new TF1(*ret);
+ // Give a warning if we're using fall-back
+ if (ret == f.GetFunctions().At(0)) {
+ AliWarning("Choosing fall-back 1 particle fit");
+ }
+ // Copy our result and return (the functions are owned by the fitter)
+ TF1* fret = new TF1(*ret);
+ return fret;
}
//____________________________________________________________________
Bool_t
-AliFMDEnergyFitter::RingHistos::CheckResult(TFitResult* r) const
+AliFMDEnergyFitter::RingHistos::CheckResult(TFitResult* r,
+ Double_t relErrorCut,
+ Double_t chi2nuCut) const
{
+ //
+ // Check the result of the fit. Returns true if
+ // - @f$ \chi^2/\nu < \max{\chi^2/\nu}@f$
+ // - @f$ \Delta p_i/p_i < \delta_e@f$ for all parameters. Note,
+ // for multi-particle fits, this requirement is relaxed by a
+ // factor of 2
+ // - @f$ a_{n} > 10^{-7}@f$ when fitting to an @f$ n@f$
+ // particle response
+ //
+ // Parameters:
+ // r Result to check
+ // relErrorCut Cut @f$ \delta_e@f$ applied to relative error
+ // of parameter.
+ // chi2nuCut Cut @f$ \max{\chi^2/\nu}@f$
+ //
+ // Return:
+ // true if fit is good.
+ //
+ if (fDebug > 10) r->Print();
+ TString n = r->GetName();
+ n.ReplaceAll("TFitResult-", "");
Double_t chi2 = r->Chi2();
Int_t ndf = r->Ndf();
// Double_t prob = r.Prob();
- if (ndf <= 0 || chi2 / ndf > 5) {
- if (fDebug > 2)
- AliWarning(Form("Fit %s has chi^2/ndf=%f/%d=%f",
- r->GetName(), chi2, ndf, (ndf<0 ? 0 : chi2/ndf)));
- return kFALSE;
+ Bool_t ret = kTRUE;
+
+ // Check that the reduced chi square isn't larger than 5
+ if (ndf <= 0 || chi2 / ndf > chi2nuCut) {
+ if (fDebug > 2) {
+ AliWarning(Form("%s: chi^2/ndf=%12.5f/%3d=%12.5f>%12.5f",
+ n.Data(), chi2, ndf, (ndf<0 ? 0 : chi2/ndf),
+ chi2nuCut)); }
+ ret = kFALSE;
}
+ // Check each parameter
UShort_t nPar = r->NPar();
for (UShort_t i = 0; i < nPar; i++) {
- if (i == 3) continue;
+ if (i == kN) continue; // Skip the number parameter
- Double_t v = r->Parameter(i);
- Double_t e = r->ParError(i);
+ // Get value and error and check value
+ Double_t v = r->Parameter(i);
+ Double_t e = r->ParError(i);
if (v == 0) continue;
- if (v == 0 || e / v > 0.2) {
- if (fDebug > 2)
- AliWarning(Form("Fit %s has Delta %s/%s=%f/%f=%f%%",
- r->GetName(), r->ParName(i).c_str(),
- r->ParName(i).c_str(), e, v, 100*(v == 0 ? 0 : e/v)));
- return kFALSE;
+
+ // Calculate the relative error and check it
+ Double_t re = e / v;
+ Double_t cut = relErrorCut * (i < kN ? 1 : 2);
+ if (re > cut) {
+ if (fDebug > 2) {
+ AliWarning(Form("%s: Delta %s/%s=%9.5f/%9.5f=%5.1f%%>%5.1f%%",
+ n.Data(), r->ParName(i).c_str(),
+ r->ParName(i).c_str(), e, v,
+ 100*(v == 0 ? 0 : e/v),
+ 100*(cut))); }
+ ret = kFALSE;
}
}
- if (nPar > 5) {
+
+ // Check if we have scale parameters
+ if (nPar > kN) {
+
+ // Check that the last particle has a significant contribution
Double_t lastScale = r->Parameter(nPar-1);
if (lastScale <= 1e-7) {
- if (fDebug)
- AliWarning(Form("Last scale %s is too small %f<1e-7",
- r->ParName(nPar-1).c_str(), lastScale));
- return kFALSE;
+ if (fDebug) {
+ AliWarning(Form("%s: %s=%9.6f<1e-7",
+ n.Data(), r->ParName(nPar-1).c_str(), lastScale)); }
+ ret = kFALSE;
}
}
- return kTRUE;
+ return ret;
+}
+
+
+//__________________________________________________________________
+void
+AliFMDEnergyFitter::RingHistos::FindBestFits(const TList* d,
+ AliFMDCorrELossFit& obj,
+ const TAxis& eta,
+ Double_t relErrorCut,
+ Double_t chi2nuCut,
+ Double_t minWeightCut)
+{
+ //
+ // Find the best fits
+ //
+ // Parameters:
+ // d Parent list
+ // obj Object to add fits to
+ // eta Eta axis
+ // relErrorCut Cut applied to relative error of parameter.
+ // Note, for multi-particle weights, the cut
+ // is loosend by a factor of 2
+ // chi2nuCut Cut on @f$ \chi^2/\nu@f$ -
+ // the reduced @f$\chi^2@f$
+ // minWeightCut Least valid @f$ a_i@f$
+ //
+
+ // Get our ring list from the output container
+ TList* l = GetOutputList(d);
+ if (!l) return;
+
+ // Get the energy distributions from the output container
+ TList* dists = static_cast<TList*>(l->FindObject("EDists"));
+ if (!dists) {
+ AliWarning(Form("Didn't find %s_EtaEDists in %s",
+ fName.Data(), l->GetName()));
+ l->ls();
+ return;
+ }
+ Int_t nBin = eta.GetNbins();
+
+ for (Int_t b = 1; b <= nBin; b++) {
+ TString n(Form(fgkEDistFormat, fName.Data(), b));
+ TH1D* dist = static_cast<TH1D*>(dists->FindObject(n));
+ // Ignore empty histograms altoghether
+ if (!dist || dist->GetEntries() <= 0) continue;
+
+ AliFMDCorrELossFit::ELossFit* best = FindBestFit(dist,
+ relErrorCut,
+ chi2nuCut,
+ minWeightCut);
+ best->fDet = fDet;
+ best->fRing = fRing;
+ best->fBin = b; //
+ // Double_t eta = fAxis->GetBinCenter(b);
+ obj.SetFit(fDet, fRing, b, new AliFMDCorrELossFit::ELossFit(*best));
+ }
+}
+
+//__________________________________________________________________
+AliFMDCorrELossFit::ELossFit*
+AliFMDEnergyFitter::RingHistos::FindBestFit(const TH1* dist,
+ Double_t relErrorCut,
+ Double_t chi2nuCut,
+ Double_t minWeightCut)
+{
+ //
+ // Find the best fit
+ //
+ // Parameters:
+ // dist Histogram
+ // relErrorCut Cut applied to relative error of parameter.
+ // Note, for multi-particle weights, the cut
+ // is loosend by a factor of 2
+ // chi2nuCut Cut on @f$ \chi^2/\nu@f$ -
+ // the reduced @f$\chi^2@f$
+ // minWeightCut Least valid @f$ a_i@f$
+ //
+ // Return:
+ // Best fit
+ //
+ TList* funcs = dist->GetListOfFunctions();
+ TIter next(funcs);
+ TF1* func = 0;
+ fFits.Clear();
+ Int_t i = 0;
+ // Info("FindBestFit", "%s", dist->GetName());
+ while ((func = static_cast<TF1*>(next()))) {
+ AliFMDCorrELossFit::ELossFit* fit =
+ new(fFits[i++]) AliFMDCorrELossFit::ELossFit(0,*func);
+ fit->CalculateQuality(chi2nuCut, relErrorCut, minWeightCut);
+ }
+ fFits.Sort(false);
+ // fFits.Print();
+ return static_cast<AliFMDCorrELossFit::ELossFit*>(fFits.At(0));
}
void
AliFMDEnergyFitter::RingHistos::Output(TList* dir)
{
+ //
+ // Define outputs
+ //
+ // Parameters:
+ // dir
+ //
fList = DefineOutputList(dir);
}