From a4097895839c04c83ecd6ae8be564e1f731c3721 Mon Sep 17 00:00:00 2001 From: mfloris Date: Fri, 14 Oct 2011 13:41:01 +0000 Subject: [PATCH] Added more functions and fixed blast wave --- PWG2/SPECTRA/Fit/AliBWFunc.cxx | 196 ++++++++++++++++++++++++++++----- PWG2/SPECTRA/Fit/AliBWFunc.h | 44 ++++++-- PWG2/SPECTRA/Fit/FitParticle.C | 128 +++++++++++++++------ 3 files changed, 298 insertions(+), 70 deletions(-) diff --git a/PWG2/SPECTRA/Fit/AliBWFunc.cxx b/PWG2/SPECTRA/Fit/AliBWFunc.cxx index 0063de4841b..00c43409ff7 100644 --- a/PWG2/SPECTRA/Fit/AliBWFunc.cxx +++ b/PWG2/SPECTRA/Fit/AliBWFunc.cxx @@ -43,6 +43,7 @@ TF1 * AliBWFunc::GetHistoFunc(TH1 * h, const char * name) { // Regardless of the variable type, this returns a function made // from the histo * a multiplicative normalization. + // This uses a bad hack... fLastFunc = new TF1 (name, StaticHistoFunc, 0.0, 10, 2); fLastFunc->SetParameter(0,1); @@ -52,20 +53,36 @@ TF1 * AliBWFunc::GetHistoFunc(TH1 * h, const char * name) { return fLastFunc; + +} +TF1 * AliBWFunc::GetGraphFunc(TGraph * g, const char * name) { + + // Regardless of the variable type, this returns a function made + // from the graph * a multiplicative normalization. + // This uses a bad hack... + + fLastFunc = new TF1 (name, StaticHistoFunc, 0.0, 10, 2); + fLastFunc->SetParameter(0,1); + fLastFunc->FixParameter(1,Double_t(Long64_t(g))); + fLastFunc->SetParNames("norm", "pointer to histo"); + fLastFunc->SetLineWidth(fLineWidth); + return fLastFunc; + + } TF1 * AliBWFunc::GetBGBW(Double_t mass, Double_t beta, Double_t T, - Double_t norm, const char * name){ + Double_t n, Double_t norm, const char * name){ // Boltzmann-Gibbs blast wave switch (fVarType) { case kdNdpt: - return GetBGBWdNdptTimesPt(mass,beta,T,norm,name); + return GetBGBWdNdptTimesPt(mass,beta,T,n,norm,name); break; case kOneOverPtdNdpt: - return GetBGBWdNdpt(mass,beta,T,norm,name); + return GetBGBWdNdpt(mass,beta,T,n,norm,name); break; case kOneOverMtdNdmt: AliFatal("Not implemented"); @@ -144,6 +161,52 @@ TF1 * AliBWFunc::GetMTExp(Double_t mass, Double_t T, Double_t norm, const char * } + +TF1 * AliBWFunc::GetBoseEinstein(Double_t mass, Double_t T, Double_t norm, const char * name){ + + // Bose einstein + switch (fVarType) { + case kdNdpt: + return GetBoseEinsteindNdptTimesPt(mass,T,norm,name); + break; + case kOneOverPtdNdpt: + return GetBoseEinsteindNdpt(mass,T,norm,name); + break; + case kOneOverMtdNdmt: + AliFatal("Not implemented"); + break; + default: + AliFatal("Not implemented"); + } + + return 0; + + +} + +TF1 * AliBWFunc::GetFermiDirac(Double_t mass, Double_t T, Double_t norm, const char * name){ + + // Simple exponential in 1/mt*MT + switch (fVarType) { + case kdNdpt: + return GetFermiDiracdNdptTimesPt(mass,T,norm,name); + break; + case kOneOverPtdNdpt: + return GetFermiDiracdNdpt(mass,T,norm,name); + break; + case kOneOverMtdNdmt: + AliFatal("Not implemented"); + break; + default: + AliFatal("Not implemented"); + } + + return 0; + + +} + + TF1 * AliBWFunc::GetPTExp(Double_t T, Double_t norm, const char * name){ // Simple exponential in 1/mt*MT @@ -269,18 +332,26 @@ Double_t AliBWFunc::StaticHistoFunc(const double * x, const double* p){ double norm = p[0]; - TH1 * h = (TH1*) Long64_t(p[1]); + TObject * h = (TObject*) Long64_t(p[1]); // Int_t bin = h->FindBin(x[0]); // double value = h->GetBinContent(bin); - if (h->FindBin(x[0]) > h->GetNbinsX()) return 0; // static TH1 * oldptr = 0; // static TSpline3 * spl = 0; // if (h!=oldptr) { // FIXME: recheck static pointers - TSpline3 * spl = new TSpline3(h); + TSpline3 * spl = 0; + if(h->InheritsFrom("TH1")) { + if ( ((TH1*)h)->FindBin(x[0]) > ((TH1*)h)->GetNbinsX()) return 0; + spl= new TSpline3((TH1*)h); + } + else if(h->InheritsFrom("TGraph")) spl= new TSpline3("fGraph",(TGraph*)h); + else { + Printf("AliBWFunc::StaticHistoFunc: Unsupported type"); + return 0; + } // } double value = spl->Eval(x[0]); delete spl; @@ -329,21 +400,36 @@ Double_t AliBWFunc::StaticUA1Func(const double * x, const double* p) { Double_t AliBWFunc::IntegrandBG(const double * x, const double* p){ // integrand for boltzman-gibbs blast wave + // x[0] -> r (radius) + // p[0] -> mass + // p[1] -> pT (transverse momentum) + // p[2] -> beta_max (surface velocity) + // p[3] -> T (freezout temperature) + // p[4] -> n (velocity profile) + double x0 = x[0]; - double mass = p[0]; - double pT = p[1]; - double beta = p[2]; - double temp = p[3]; - + double mass = p[0]; + double pT = p[1]; + double beta_max = p[2]; + double temp = p[3]; + Double_t n = p[4]; + + // Keep beta within reasonable limits + Double_t beta = beta_max * TMath::Power(x0, n); + if (beta > 0.9999999999999999) beta = 0.9999999999999999; + double mT = TMath::Sqrt(mass*mass+pT*pT); - double rho0 = TMath::ATanH(beta*x0); + double rho0 = TMath::ATanH(beta); double arg00 = pT*TMath::SinH(rho0)/temp; + if (arg00 > 700.) arg00 = 700.; // avoid FPE double arg01 = mT*TMath::CosH(rho0)/temp; double f0 = x0*mT*TMath::BesselI0(arg00)*TMath::BesselK1(arg01); + // printf("r=%f, pt=%f, beta_max=%f, temp=%f, n=%f, mt=%f, beta=%f, rho=%f, argI0=%f, argK1=%f\n", x0, pT, beta_max, temp, n, mT, beta, rho0, arg00, arg01); + return f0; } @@ -356,17 +442,20 @@ Double_t AliBWFunc::StaticBGdNdPt(const double * x, const double* p) { double pT = x[0];; - double mass = p[0]; - double beta = p[1]; + double mass = p[0]; + double beta = p[1]; double temp = p[2]; + double n = p[3]; + double norm = p[4]; static TF1 * fIntBG = 0; if(!fIntBG) - fIntBG = new TF1 ("fIntBG", IntegrandBG, 0, 1, 4); + fIntBG = new TF1 ("fIntBG", IntegrandBG, 0, 1, 5); - fIntBG->SetParameters(mass, pT, beta, temp); + fIntBG->SetParameters(mass, pT, beta, temp,n); double result = fIntBG->Integral(0,1); - return result*p[3];//*1e30;; + // printf ("[%4.4f], Int :%f\n", pT, result); + return result*norm;//*1e30;; } @@ -377,13 +466,14 @@ Double_t AliBWFunc::StaticBGdNdPtTimesPt(const double * x, const double* p) { TF1 * AliBWFunc::GetBGBWdNdpt(Double_t mass, Double_t beta, Double_t temp, - Double_t norm, const char * name){ + Double_t n, Double_t norm, const char * name){ // BGBW 1/pt dNdpt - fLastFunc = new TF1 (name, StaticBGdNdPt, 0.0, 10, 4); - fLastFunc->SetParameters(mass,beta,temp,norm); - fLastFunc->SetParNames("mass", "#beta", "T", "norm"); + fLastFunc = new TF1 (name, StaticBGdNdPt, 0.0, 10, 5); + fLastFunc->SetParameters(mass,beta,temp,n,norm); + fLastFunc->FixParameter(0,mass); + fLastFunc->SetParNames("mass", "#beta", "T", "n", "norm"); fLastFunc->SetLineWidth(fLineWidth); return fLastFunc; @@ -484,14 +574,15 @@ TF1 * AliBWFunc::GetTsallisBWdNdpt(Double_t mass, Double_t beta, Double_t temp, // Times Pt funcs // Boltzmann-Gibbs Blast Wave -TF1 * AliBWFunc::GetBGBWdNdptTimesPt(Double_t mass, Double_t beta, Double_t temp, +TF1 * AliBWFunc::GetBGBWdNdptTimesPt(Double_t mass, Double_t beta, Double_t temp, Double_t n, Double_t norm, const char * name){ // BGBW, dNdpt - fLastFunc = new TF1 (name, StaticBGdNdPtTimesPt, 0.0, 10, 4); - fLastFunc->SetParameters(mass,beta,temp,norm); - fLastFunc->SetParNames("mass", "#beta", "temp", "norm"); + fLastFunc = new TF1 (name, StaticBGdNdPtTimesPt, 0.0, 10, 5); + fLastFunc->SetParameters(mass,beta,temp,n,norm); + fLastFunc->FixParameter(0,mass); + fLastFunc->SetParNames("mass", "#beta", "temp", "n", "norm"); fLastFunc->SetLineWidth(fLineWidth); return fLastFunc; @@ -531,6 +622,37 @@ TF1 * AliBWFunc::GetMTExpdNdptTimesPt(Double_t mass, Double_t temp, Double_t nor } +TF1 * AliBWFunc::GetBoseEinsteindNdptTimesPt(Double_t mass, Double_t temp, Double_t norm, const char * name){ + + // Bose einstein distribution as a function of dNdpt + char formula[500]; + snprintf(formula,500,"[0]*x*1./(exp(sqrt(x**2+%f**2)/[1])-1)", mass); + fLastFunc=new TF1(name,formula,0,10); + fLastFunc->SetParameters(norm, temp); + fLastFunc->SetParLimits(1, 0.01, 10); + fLastFunc->SetParNames("norm", "T"); + fLastFunc->SetLineWidth(fLineWidth); + return fLastFunc; + + +} + +TF1 * AliBWFunc::GetFermiDiracdNdptTimesPt(Double_t mass, Double_t temp, Double_t norm, const char * name){ + + // Bose einstein distribution as a function of dNdpt + char formula[500]; + snprintf(formula,500,"[0]*x*1./(exp(sqrt(x**2+%f**2)/[1])+1)", mass); + fLastFunc=new TF1(name,formula,0,10); + fLastFunc->SetParameters(norm, temp); + fLastFunc->SetParLimits(1, 0.01, 10); + fLastFunc->SetParNames("norm", "T"); + fLastFunc->SetLineWidth(fLineWidth); + return fLastFunc; + + +} + + TF1 * AliBWFunc::GetPTExpdNdptTimesPt(Double_t temp, Double_t norm, const char * name){ @@ -748,6 +870,30 @@ TF1 * AliBWFunc::GetMTExpdNdpt(Double_t mass, Double_t temp, Double_t norm, cons return fLastFunc; } +TF1 * AliBWFunc::GetBoseEinsteindNdpt(Double_t mass, Double_t temp, Double_t norm, const char * name){ + // bose einstein + char formula[500]; + snprintf(formula,500,"[0]*1./(exp(sqrt(x**2+%f**2)/[1])-1)", mass); + fLastFunc=new TF1(name,formula,0,10); + fLastFunc->SetParameters(norm, temp); + fLastFunc->SetParLimits(1, 0.01, 10); + fLastFunc->SetParNames("norm", "T"); + fLastFunc->SetLineWidth(fLineWidth); + return fLastFunc; +} + +TF1 * AliBWFunc::GetFermiDiracdNdpt(Double_t mass, Double_t temp, Double_t norm, const char * name){ + // bose einstein + char formula[500]; + snprintf(formula,500,"[0]*1./(exp(sqrt(x**2+%f**2)/[1])+1)", mass); + fLastFunc=new TF1(name,formula,0,10); + fLastFunc->SetParameters(norm, temp); + fLastFunc->SetParLimits(1, 0.01, 10); + fLastFunc->SetParNames("norm", "T"); + fLastFunc->SetLineWidth(fLineWidth); + return fLastFunc; +} + // // Simple tsallis (a la CMS) // TF1 * AliBWFunc::GetTsallisdNdpt(Double_t mass, Double_t temp, Double_t q, Double_t norm, const char * name){ diff --git a/PWG2/SPECTRA/Fit/AliBWFunc.h b/PWG2/SPECTRA/Fit/AliBWFunc.h index 27a2f7aadfc..d774687fad9 100644 --- a/PWG2/SPECTRA/Fit/AliBWFunc.h +++ b/PWG2/SPECTRA/Fit/AliBWFunc.h @@ -23,7 +23,7 @@ class TF1; class TH1; - +class TGraph; #endif @@ -38,9 +38,10 @@ public: AliBWFunc(); ~AliBWFunc(); - // Boltzmann-Gibbs blast wave - TF1 * GetBGBW(Double_t mass, Double_t beta, Double_t T, + // Boltzmann-Gibbs Blast Wave + TF1 * GetBGBW(Double_t mass, Double_t beta, Double_t T, Double_t n, Double_t norm, const char * name = "fBGBW"); + // Boltzmann TF1 * GetBoltzmann(Double_t mass, Double_t T, Double_t norm, const char * name ="fBoltzmann"); @@ -50,10 +51,10 @@ public: Double_t norm, Double_t ymax = 0.5, const char * name = "fTsallisBW"); // Simple exponential in 1/mt*dNdmt - TF1 * GetMTExp(Double_t mass, Double_t T, Double_t norm, const char * name ="fExp"); + TF1 * GetMTExp(Double_t mass, Double_t T, Double_t norm, const char * name ="fMtExp"); // Simple exponential in 1/pt*dNdpt - TF1 * GetPTExp(Double_t T, Double_t norm, const char * name ="fExp"); + TF1 * GetPTExp(Double_t T, Double_t norm, const char * name ="fPtExp"); // Tsallis (no BW, a la CMS) TF1 * GetTsallis(Double_t mass, Double_t T, Double_t q, Double_t norm, const char * name="fTsallis") @@ -68,8 +69,18 @@ public: // Function derived from a histo TF1 * GetHistoFunc(TH1 * h, const char * name = "fHisto"); + // Function derived from a graph + TF1 * GetGraphFunc(TGraph * h, const char * name = "fHisto"); + // Power law TF1 * GetPowerLaw(Double_t pt0, Double_t n, Double_t norm, const char * name="fPowerLaw"); + + + // Bose-Einstein + TF1 * GetBoseEinstein(Double_t mass, Double_t T, Double_t norm, const char * name="fBoseEinstein"); + + // Fermi-Dirac + TF1 * GetFermiDirac(Double_t mass, Double_t T, Double_t norm, const char * name="fFermiDirac"); void SetVarType(VarType_t tp) {fVarType=tp;} @@ -78,9 +89,10 @@ protected: // dNdpt here means 1/pt dN/dpt - // Boltzmann-Gibbs Blast Wave - TF1 * GetBGBWdNdpt(Double_t mass, Double_t beta, Double_t T, - Double_t norm, const char * name = "fBGBW"); + + // Boltzmann-Gibbs blast wave + TF1 * GetBGBWdNdpt(Double_t mass, Double_t beta, Double_t temp, + Double_t n, Double_t norm, const char * name= "fBGBW"); // Tsallis blast wave TF1 * GetTsallisBWdNdpt(Double_t mass, Double_t beta, Double_t T, Double_t q, @@ -89,6 +101,13 @@ protected: // Simple exponential in 1/mt*MT TF1 * GetMTExpdNdpt(Double_t mass, Double_t T, Double_t norm, const char * name ="fExp"); + // Bose-Einstein + TF1 * GetBoseEinsteindNdpt(Double_t mass, Double_t T, Double_t norm, const char * name="fBoseEinstein"); + + // Fermi-Dirac + TF1 * GetFermiDiracdNdpt(Double_t mass, Double_t T, Double_t norm, const char * name="fFermiDirac"); + + // Tsallis (no BW, a la CMS) TF1 * GetTsallisdNdpt(Double_t mass, Double_t T, Double_t q, Double_t norm, const char * name="fTsallis"); @@ -104,7 +123,7 @@ protected: // TimesPt means dNdpt // Boltzmann-Gibbs Blast Wave - TF1 * GetBGBWdNdptTimesPt(Double_t mass, Double_t beta, Double_t T, + TF1 * GetBGBWdNdptTimesPt(Double_t mass, Double_t beta, Double_t T, Double_t n, Double_t norm, const char * name = "fBGBWTimesPt"); // Tsallis blast wave @@ -122,6 +141,13 @@ protected: // Simple exponential in 1/mt*dNdmT TF1 * GetMTExpdNdptTimesPt(Double_t mass, Double_t T, Double_t norm, const char * name ="fMtExpTimesPt"); + // Bose-Einstein + TF1 * GetBoseEinsteindNdptTimesPt(Double_t mass, Double_t T, Double_t norm, const char * name="fBoseEinstein"); + + // Fermi-Dirac + TF1 * GetFermiDiracdNdptTimesPt(Double_t mass, Double_t T, Double_t norm, const char * name="fFermiDirac"); + + // Simple exponential in 1/mp*dNdpT TF1 * GetPTExpdNdptTimesPt(Double_t T, Double_t norm, const char * name ="fPtExpTimesPt"); diff --git a/PWG2/SPECTRA/Fit/FitParticle.C b/PWG2/SPECTRA/Fit/FitParticle.C index d7ff6ad6212..236854e4f1b 100644 --- a/PWG2/SPECTRA/Fit/FitParticle.C +++ b/PWG2/SPECTRA/Fit/FitParticle.C @@ -14,18 +14,19 @@ #include "TSystem.h" #include "THashList.h" #include "TMinuit.h" +#include "TLatex.h" using namespace std; #endif -enum {kFitExpPt, kFitLevi, fFitExpMt, kFitBoltzmann, kFitBlastWave}; +enum {kFitExpPt, kFitLevi, fFitExpMt, kFitBoltzmann, kFitBlastWave, kFitBoseEinstein, kFitFermiDirac}; -void FitParticle(const char * file, const char * histo, const char * partName, const char * listname=0, Float_t min = 0, Float_t max =3, Float_t scaleHisto = -1., Int_t fitFunc = kFitLevi, Int_t vartype = AliBWFunc::kdNdpt) { +void FitParticle(const char * file, const char * histo, const char * partName, const char * listname=0, Float_t min = 0, Float_t max =3, Float_t scaleHisto = -1., Int_t fitFunc = kFitLevi, Int_t vartype = AliBWFunc::kdNdpt, const char * fileOut = 0, Bool_t wait = 0) { // Generic Macro to fit any particle using the PWG2/SPECTRA/Fit macros // // You have to provide: - // - file: file name + // - file: file name. If not set, uses gFile // - histo: histogram name name (assumend to be in the main folder of the file) // - listname: if different from 0, histo is looked in this tlist // - partName: it is used to get the mass of the particle from @@ -78,28 +79,69 @@ void FitParticle(const char * file, const char * histo, const char * partName, Double_t mass = TDatabasePDG::Instance()->GetParticle(partName)->Mass(); cout << "Fitting ["<GetLevi (mass, 0.4, 20,3); - func->SetParLimits(1,1,100); + func->SetParLimits(1,0.0001,20000); + table = new AliLatexTable(9,"c|cccccccc"); + table->InsertCustomRow("Part & Func & Yield & FD Slope & $\\Chi^2$/NDF & Min X & Fit Min & Fit Max & Frac Below "); normPar = 0; + slopePar = 2; } if (fitFunc == kFitExpPt) { func = fm->GetPTExp(0.2, 20); + table = new AliLatexTable(9,"c|cccccccc"); + table->InsertCustomRow("Part & Func & Yield & pT Slope & $\\Chi^2$/NDF & Min X & Fit Min & Fit Max & Frac Below "); + normPar = 0; + slopePar = 1; } if (fitFunc == fFitExpMt) { func = fm->GetMTExp(mass,0.2, 20); + table = new AliLatexTable(9,"c|cccccccc"); + table->InsertCustomRow("Part & Func & Yield & mT Slope & $\\Chi^2$/NDF & Min X & Fit Min & Fit Max & Frac Below "); + normPar = 0; + slopePar = 1; } if (fitFunc == kFitBoltzmann) { func = fm->GetBoltzmann(mass, 0.2, 20); + table = new AliLatexTable(9,"c|cccccccc"); + table->InsertCustomRow("Part & Func & Yield & Bz Slope & $\\Chi^2$/NDF & Min X & Fit Min & Fit Max & Frac Below "); + normPar = 0; + slopePar = 1; + } if (fitFunc == kFitBlastWave) { - func = fm->GetBGBW(mass,0.6,0.3, 20); + func = fm->GetBGBW(mass,0.6,0.3, 1, 1e5);// beta, T, n, norm + table = new AliLatexTable(11,"cc|ccccccccccc"); + table->InsertCustomRow("Part & Func & Yield & T & beta & n & $\\Chi^2$/NDF & Min X & Fit Min & Fit Max & Frac Below "); + func->SetParLimits(1, 0.1, 0.99); + func->SetParLimits(2, 0.01, 1); + func->SetParLimits(3, 0.01, 2); + normPar = 4; + slopePar = 2; + } + if (fitFunc == kFitBoseEinstein) { + func = fm->GetBoseEinstein(mass,0.3,20); + table = new AliLatexTable(9,"c|cccccccc"); + table->InsertCustomRow("Part & Func & Yield & BE Slope & $\\Chi^2$/NDF & Min X & Fit Min & Fit Max & Frac Below "); + normPar = 0; + slopePar = 1; + } + if (fitFunc == kFitFermiDirac) { + func = fm->GetFermiDirac(mass,0.3,20); + table = new AliLatexTable(9,"c|cccccccc"); + table->InsertCustomRow("Part & Func & Yield & FD Slope & $\\Chi^2$/NDF & Min X & Fit Min & Fit Max & Frac Below "); + normPar = 0; + slopePar = 1; } - - TFile * f = new TFile(file); + + TFile * f = file ? new TFile(file) : gFile; TH1 * h = 0; if(listname){ TList * l = (TList*) gDirectory->Get(listname); @@ -132,45 +174,59 @@ void FitParticle(const char * file, const char * histo, const char * partName, func->Draw("same"); cout << "4" << endl; - // Print results nicely - AliLatexTable table(9,"c|ccccccc"); - table.InsertCustomRow("Part & Integral & T Slope & n & $\\Chi^2$/NDF & Min X & Frac Above & \\langle p_{t} \\rangle & \\langle p_{t}^{2} \\rangle"); // populate table - // Float_t yield = func->Integral(0,100); - // Float_t yieldE = func->IntegralError(0,100); + Double_t yield=0,yieldE=0; cout << "Y" << endl; AliBWTools::GetYield(h,func,yield,yieldE); cout << "YE" << endl; + TLatex * l = new TLatex(2,(h->GetMaximum()+h->GetMinimum())/2,Form("%3.3f #pm %3.3f (%s)", yield,yieldE, func->GetName())); + l->Draw(); + + if(wait) c->WaitPrimitive(); + if(fileOut) { + TFile * fout = new TFile(fileOut, "update"); + c->SetName(Form("c_%s_%s_fit_%2.2f_%2.2f",histo, func->GetName(), min, max)); + c->Write(); + fout->Purge(); // remove old canvases + fout->Close(); + } + // Float_t yield = func->Integral(0.45,1.05); // Float_t yieldE = func->IntegralError(0.45,1.05); - Double_t tslope = func->GetParameter(2); - Double_t tslopeE = func->GetParError(2); - - table.SetNextCol(partName); - table.SetNextCol(yield,yieldE,-4); - table.SetNextCol(tslope,tslopeE,-4); - table.SetNextCol(func->GetParameter(1),func->GetParError(1),-4); - table.SetNextCol(Form("%2.2f/%d",func->GetChisquare(),func->GetNDF())); - cout << "5" << endl; + cout << "Slope Par: "<< slopePar << endl; + + Double_t tslope = func->GetParameter(slopePar); + Double_t tslopeE = func->GetParError(slopePar); + + table->SetNextCol(Form("%s (%s)", partName, h->GetName())); + table->SetNextCol(func->GetName()); + table->SetNextCol(yield,yieldE,-4); + table->SetNextCol(tslope,tslopeE,-4); + if(fitFunc == kFitBlastWave) { + table->SetNextCol(func->GetParameter(1),func->GetParError(1),-4); + table->SetNextCol(func->GetParameter(3),func->GetParError(3),-4); + } + // table->SetNextCol(func->GetParameter(1),func->GetParError(1),-4); + table->SetNextCol(Form("%2.2f/%d",func->GetChisquare(),func->GetNDF())); Float_t lowestPoint = TMath::Max(AliBWTools::GetLowestNotEmptyBinEdge(h),min); - cout << "6" << endl; - //Float_t yieldAbove = 0; - Float_t yieldAbove = func->Integral(lowestPoint,100); - table.SetNextCol(lowestPoint,-3); - table.SetNextCol(yieldAbove/yield,-3); - Float_t mean=0, meane=0; - Float_t mean2=0, mean2e=0; - cout << "6" << endl; - AliBWTools::GetMean (func, mean, meane , 0.,100., normPar); - AliBWTools::GetMeanSquare(func, mean2, mean2e, 0.,100., normPar); - cout << "8" << endl; - table.SetNextCol(mean, meane ,-4); - table.SetNextCol(mean2, mean2e,-4); + // Float_t yieldAbove = func->Integral(lowestPoint,100); + Float_t yieldBelow = func->Integral(0,lowestPoint); + table->SetNextCol(lowestPoint,-3); + table->SetNextCol(min,-2); + table->SetNextCol(max,-2); + table->SetNextCol(yieldBelow/yield,-3); + + // Float_t mean=0, meane=0; + // Float_t mean2=0, mean2e=0; + // AliBWTools::GetMean (func, mean, meane , 0.,100., normPar); + // AliBWTools::GetMeanSquare(func, mean2, mean2e, 0.,100., normPar); + // table->SetNextCol(mean, meane ,-4); + // table->SetNextCol(mean2, mean2e,-4); // fMean2->IntegralError(0,100)/func->Integral(0,100),-7); - table.InsertRow(); - table.PrintTable("ASCII"); + table->InsertRow(); + table->PrintTable("ASCII"); -- 2.39.3