X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=PWG2%2FFORWARD%2Fanalysis2%2FAliFMDEnergyFitter.cxx;h=8c7eb80e004ca4aa33c8106c21492eecc8d33776;hb=f4a1b795965b59e87ae7c5b8b62485b592600a52;hp=caedd41a879d72ed2c939818218c445b5acf8255;hpb=95985fc421ff3b83b2567d6a0ea340cd7ecb18e0;p=u%2Fmrichter%2FAliRoot.git diff --git a/PWG2/FORWARD/analysis2/AliFMDEnergyFitter.cxx b/PWG2/FORWARD/analysis2/AliFMDEnergyFitter.cxx index caedd41a879..8c7eb80e004 100644 --- a/PWG2/FORWARD/analysis2/AliFMDEnergyFitter.cxx +++ b/PWG2/FORWARD/analysis2/AliFMDEnergyFitter.cxx @@ -16,7 +16,6 @@ #include #include #include -#include ClassImp(AliFMDEnergyFitter) #if 0 @@ -31,13 +30,14 @@ namespace { AliFMDEnergyFitter::AliFMDEnergyFitter() : TNamed(), fRingHistos(), - fLowCut(0.3), - fNParticles(3), - fMinEntries(100), + fLowCut(0.4), + fNParticles(5), + fMinEntries(1000), fFitRangeBinWidth(4), fDoFits(false), fDoMakeObject(false), fEtaAxis(), + fCentralityAxis(), fMaxE(10), fNEbins(300), fUseIncreasingBins(true), @@ -55,13 +55,14 @@ AliFMDEnergyFitter::AliFMDEnergyFitter() AliFMDEnergyFitter::AliFMDEnergyFitter(const char* title) : TNamed("fmdEnergyFitter", title), fRingHistos(), - fLowCut(0.3), - fNParticles(3), - fMinEntries(100), + 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), @@ -78,6 +79,8 @@ AliFMDEnergyFitter::AliFMDEnergyFitter(const char* title) // 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')); @@ -96,6 +99,7 @@ AliFMDEnergyFitter::AliFMDEnergyFitter(const AliFMDEnergyFitter& o) fDoFits(o.fDoFits), fDoMakeObject(o.fDoMakeObject), fEtaAxis(o.fEtaAxis), + fCentralityAxis(o.fCentralityAxis), fMaxE(o.fMaxE), fNEbins(o.fNEbins), fUseIncreasingBins(o.fUseIncreasingBins), @@ -145,7 +149,17 @@ AliFMDEnergyFitter::operator=(const AliFMDEnergyFitter& o) fFitRangeBinWidth= o.fFitRangeBinWidth; fDoFits = o.fDoFits; fDoMakeObject = o.fDoMakeObject; - fEtaAxis.Set(o.fEtaAxis.GetNbins(),o.fEtaAxis.GetXmin(),o.fEtaAxis.GetXmax()); + 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; @@ -199,10 +213,16 @@ AliFMDEnergyFitter::Init(const TAxis& eAxis) if (fEtaAxis.GetNbins() == 0 || 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(next()))) - o->Init(fEtaAxis, fMaxE, fNEbins, fUseIncreasingBins); + o->Init(fEtaAxis, fCentralityAxis, fMaxE, fNEbins, fUseIncreasingBins); } //____________________________________________________________________ void @@ -238,10 +258,18 @@ AliFMDEnergyFitter::SetEtaAxis(Int_t nBins, Double_t etaMin, Double_t etaMax) // 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) { // @@ -249,11 +277,15 @@ AliFMDEnergyFitter::Accumulate(const AliESDFMD& input, // // 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++) { @@ -276,7 +308,7 @@ AliFMDEnergyFitter::Accumulate(const AliESDFMD& input, 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 @@ -334,13 +366,6 @@ AliFMDEnergyFitter::Fit(const TList* dir) if (!fDoMakeObject) return; MakeCorrectionsObject(d); - d->ls(); - TDirectory* savdir = gDirectory; - TFile* tmp = TFile::Open("elossfits.root", "RECREATE"); - d->Write(); - tmp->Write(); - tmp->Close(); - savdir->cd(); } //____________________________________________________________________ @@ -370,6 +395,7 @@ AliFMDEnergyFitter::MakeCorrectionsObject(TList* d) 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()); } @@ -425,7 +451,7 @@ AliFMDEnergyFitter::Print(Option_t*) const char ind[gROOT->GetDirLevel()+1]; for (Int_t i = 0; i < gROOT->GetDirLevel(); i++) ind[i] = ' '; ind[gROOT->GetDirLevel()] = '\0'; - std::cout << ind << "AliFMDEnergyFitter: " << GetName() << '\n' + std::cout << ind << ClassName() << ": " << GetName() << '\n' << ind << " Low cut: " << fLowCut << " E/E_mip\n" << ind << " Max(particles): " << fNParticles << '\n' << ind << " Min(entries): " << fMinEntries << '\n' @@ -555,14 +581,16 @@ AliFMDEnergyFitter::RingHistos::~RingHistos() //____________________________________________________________________ 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 + // ieta Eta bin (0-based) + // icent Centrality bin (1-based) // mult Signal // if (empty) { @@ -571,6 +599,9 @@ AliFMDEnergyFitter::RingHistos::Fill(Bool_t empty, Int_t ieta, Double_t mult) } 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(fEtaEDists.At(ieta)); @@ -729,6 +760,7 @@ AliFMDEnergyFitter::RingHistos::MakeTotal(const char* name, //____________________________________________________________________ void AliFMDEnergyFitter::RingHistos::Init(const TAxis& eAxis, + const TAxis& /* cAxis */, Double_t maxDE, Int_t nDEbins, Bool_t useIncrBin) { @@ -749,13 +781,67 @@ AliFMDEnergyFitter::RingHistos::Init(const TAxis& eAxis, 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* @@ -1005,7 +1091,9 @@ AliFMDEnergyFitter::RingHistos::FitHist(TH1* dist, 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 @@ -1013,15 +1101,15 @@ AliFMDEnergyFitter::RingHistos::FitHist(TH1* dist, // 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--) { good[i] = 0; - TF1* ff = static_cast(f.fFunctions.At(i)); + TF1* ff = static_cast(f.GetFunctions().At(i)); // ff->SetLineColor(Color()); ff->SetRange(0, maxRange); dist->GetListOfFunctions()->Add(new TF1(*ff)); - if (CheckResult(static_cast(f.fFitResults.At(i)), + if (CheckResult(static_cast(f.GetFitResults().At(i)), relErrorCut, chi2nuCut)) { good[i] = ff; ff->SetLineWidth(2); @@ -1029,20 +1117,20 @@ AliFMDEnergyFitter::RingHistos::FitHist(TH1* dist, } } // If all else fails, use the 1 particle fit - TF1* ret = static_cast(f.fFunctions.At(0)); + TF1* ret = static_cast(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.fFitResults.At(i)->Print(); + f.GetFitResults().At(i)->Print(); } ret = good[i]; break; } // Give a warning if we're using fall-back - if (ret == f.fFunctions.At(0)) { + 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) @@ -1133,7 +1221,7 @@ AliFMDEnergyFitter::RingHistos::CheckResult(TFitResult* r, //__________________________________________________________________ void -AliFMDEnergyFitter::RingHistos::FindBestFits(TList* d, +AliFMDEnergyFitter::RingHistos::FindBestFits(const TList* d, AliFMDCorrELossFit& obj, const TAxis& eta, Double_t relErrorCut, @@ -1189,7 +1277,7 @@ AliFMDEnergyFitter::RingHistos::FindBestFits(TList* d, //__________________________________________________________________ AliFMDCorrELossFit::ELossFit* -AliFMDEnergyFitter::RingHistos::FindBestFit(TH1* dist, +AliFMDEnergyFitter::RingHistos::FindBestFit(const TH1* dist, Double_t relErrorCut, Double_t chi2nuCut, Double_t minWeightCut)