AliPhysics  ff1d528 (ff1d528)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliHFVnVsMassFitter.cxx
Go to the documentation of this file.
1 #include "AliHFVnVsMassFitter.h"
2 
3 #include <TROOT.h>
4 #include <TMath.h>
5 #include <TF1.h>
6 #include <TGraphErrors.h>
7 #include <TColor.h>
8 #include <TLegend.h>
9 #include <TList.h>
10 #include <TCanvas.h>
11 #include <TStyle.h>
12 #include <TVirtualPad.h>
13 #include <TDatabasePDG.h>
14 #include <TPaveText.h>
15 #include "Fit/BinData.h"
16 #include "HFitInterface.h"
17 #include <vector>
18 
19 #include "AliVertexingHFUtils.h"
20 
24 
25 //________________________________________________________________
27  :TObject()
28  ,fMassHisto(0x0)
29  ,fVnVsMassHisto(0x0)
30  ,fMassSgnFuncType(kGaus)
31  ,fMassBkgFuncType(kExpo)
32  ,fVnBkgFuncType(kLin)
33  ,fMassFuncFromPrefit(0x0)
34  ,fMassBkgFunc(0x0)
35  ,fMassSgnFunc(0x0)
36  ,fMassTotFunc(0x0)
37  ,fVnBkgFuncSb(0x0)
38  ,fVnBkgFunc(0x0)
39  ,fVnTotFunc(0x0)
40  ,fMassFitter(0x0)
41  ,fMassMin(1.69)
42  ,fMassMax(2.05)
43  ,fVn(0.)
44  ,fVnUncertainty(0.)
45  ,fSigma(0.)
46  ,fSigmaUncertainty(0.)
47  ,fMean(0.)
48  ,fMeanUncertainty(0.)
49  ,fRawYield(0.)
50  ,fRawYieldUncertainty(0.)
51  ,fChiSquare(0.)
52  ,fNDF(0)
53  ,fProb(0.)
54  ,fNSigmaForSB(3.)
55  ,fSigmaInit(0.012)
56  ,fMeanInit(1.870)
57  ,fMeanFixedFromMassFit(kFALSE)
58  ,fSigmaFixedFromMassFit(kFALSE)
59  ,fMassParticle(1.870)
60  ,fNParsMassSgn(3)
61  ,fNParsMassBkg(2)
62  ,fNParsVnBkg(2)
63  ,fSigmaFixed(0)
64  ,fMeanFixed(0)
65  ,fPolDegreeBkg(3)
66  ,fReflections(kFALSE)
67  ,fNParsRfl(0)
68  ,fRflOverSig(0.)
69  ,fFixRflOverSig(kFALSE)
70  ,fHistoTemplRfl(0x0)
71  ,fHistoTemplRflInit(0x0)
72  ,fMassRflFunc(0x0)
73  ,fMassBkgRflFunc(0x0)
74  ,fRflOpt("1gaus")
75  ,fMinRefl(0.)
76  ,fMaxRefl(0.)
77  ,fSmoothRfl(kFALSE)
78  ,fRawYieldHelp(0.)
79  ,fSecondPeak(kFALSE)
80  ,fMassSecPeakFunc(0x0)
81  ,fNParsSec(0)
82  ,fSecMass(-999.)
83  ,fSecWidth(9999.)
84  ,fFixSecMass(kFALSE)
85  ,fFixSecWidth(kFALSE)
86  ,fDoSecondPeakVn(kFALSE)
87  ,fHarmonic(2) {
88 
89  //default constructor
90 }
91 
92 //________________________________________________________________
93 AliHFVnVsMassFitter::AliHFVnVsMassFitter(TH1F* hMass, TH1F* hvn, Double_t min, Double_t max, Int_t funcMassBkg, Int_t funcMassSgn, Int_t funcvnBkg)
94  :TObject()
95  ,fMassSgnFuncType(funcMassSgn)
96  ,fMassBkgFuncType(funcMassBkg)
97  ,fVnBkgFuncType(funcvnBkg)
98  ,fMassFuncFromPrefit(0x0)
99  ,fMassBkgFunc(0x0)
100  ,fMassSgnFunc(0x0)
101  ,fMassTotFunc(0x0)
102  ,fVnBkgFuncSb(0x0)
103  ,fVnBkgFunc(0x0)
104  ,fVnTotFunc(0x0)
105  ,fMassFitter(0x0)
106  ,fMassMin(min)
107  ,fMassMax(max)
108  ,fVn(0.)
109  ,fVnUncertainty(0.)
110  ,fSigma(0.)
111  ,fSigmaUncertainty(0.)
112  ,fMean(0.)
113  ,fMeanUncertainty(0.)
114  ,fRawYield(0.)
115  ,fRawYieldUncertainty(0.)
116  ,fChiSquare(0.)
117  ,fNDF(0)
118  ,fProb(0)
119  ,fNSigmaForSB(3.)
120  ,fSigmaInit(0.012)
121  ,fMeanInit(1.870)
122  ,fMeanFixedFromMassFit(kFALSE)
123  ,fSigmaFixedFromMassFit(kFALSE)
124  ,fMassParticle(1.870)
125  ,fNParsMassSgn(3)
126  ,fNParsMassBkg(2)
127  ,fNParsVnBkg(2)
128  ,fSigmaFixed(0)
129  ,fMeanFixed(0)
130  ,fPolDegreeBkg(3)
131  ,fReflections(kFALSE)
132  ,fNParsRfl(0)
133  ,fRflOverSig(0.)
134  ,fFixRflOverSig(kFALSE)
135  ,fHistoTemplRfl(0x0)
136  ,fHistoTemplRflInit(0x0)
137  ,fMassRflFunc(0x0)
138  ,fMassBkgRflFunc(0x0)
139  ,fRflOpt("1gaus")
140  ,fMinRefl(0.)
141  ,fMaxRefl(0.)
142  ,fSmoothRfl(kFALSE)
143  ,fRawYieldHelp(0.)
144  ,fSecondPeak(kFALSE)
145  ,fMassSecPeakFunc(0x0)
146  ,fNParsSec(0)
147  ,fSecMass(-999.)
148  ,fSecWidth(9999.)
149  ,fFixSecMass(kFALSE)
150  ,fFixSecWidth(kFALSE)
151  ,fDoSecondPeakVn(kFALSE)
152  ,fHarmonic(2) {
153 
154  //standard constructor
155  fMassHisto = (TH1F*)hMass->Clone("fHistoInvMass");
156  fMassHisto->SetDirectory(0);
157  fVnVsMassHisto = (TH1F*)hvn->Clone(Form("fHistoV%dVsMass",fHarmonic));
158  fVnVsMassHisto->SetDirectory(0);
159 
161 }
162 
163 //________________________________________________________________
165 
166  //destructor
167  if(fMassHisto) delete fMassHisto;
168  if(fVnVsMassHisto) delete fVnVsMassHisto;
170  if(fMassBkgFunc) delete fMassBkgFunc;
171  if(fMassBkgRflFunc) delete fMassBkgRflFunc;
172  if(fMassSgnFunc) delete fMassSgnFunc;
173  if(fMassTotFunc) delete fMassTotFunc;
174  if(fVnBkgFuncSb) delete fVnBkgFuncSb;
175  if(fVnBkgFunc) delete fVnBkgFunc;
176  if(fVnTotFunc) delete fVnTotFunc;
177  if(fMassFitter) delete fMassFitter;
178  if(fHistoTemplRfl) delete fHistoTemplRfl;
180  if(fMassRflFunc) delete fMassRflFunc;
182 }
183 
184 //________________________________________________________________
186 
187  if(!fMassHisto || !fVnVsMassHisto) {AliError("Histograms not set! Exit."); return kFALSE;}
189 
191  Int_t NvnParsSgn = 1;
192  if(fSecondPeak && fDoSecondPeakVn) {NvnParsSgn+=1;}
193  const Int_t nparsvn = nparsmass+fNParsVnBkg+NvnParsSgn;
194 
195  Bool_t massprefit=MassPrefit();
196  if(!massprefit) {AliError("Impossible to perform the mass prefit"); return kFALSE;}
197  Bool_t vnprefit=VnSBPrefit();
198 
199  std::vector<Double_t> initpars;
200  for(Int_t iBkgPar=0; iBkgPar<fNParsMassBkg; iBkgPar++) {
201  initpars.push_back(fMassFuncFromPrefit->GetParameter(iBkgPar));
202  }
203  for(Int_t iSgnPar=0; iSgnPar<fNParsMassSgn; iSgnPar++) {
204  initpars.push_back(fMassFuncFromPrefit->GetParameter(iSgnPar+fNParsMassBkg));
205  }
206  for(Int_t iSecPeakPar=0; iSecPeakPar<fNParsSec; iSecPeakPar++) {
207  initpars.push_back(fMassFuncFromPrefit->GetParameter(iSecPeakPar+fNParsMassBkg+fNParsMassSgn));
208  }
209  for(Int_t iReflPar=0; iReflPar<fNParsRfl; iReflPar++) {
210  initpars.push_back(fMassFuncFromPrefit->GetParameter(iReflPar+fNParsMassBkg+fNParsMassSgn+fNParsSec));
211  }
212  for(Int_t iVnBkgPar=0; iVnBkgPar<fNParsVnBkg; iVnBkgPar++) {
213  if(vnprefit) {initpars.push_back(fVnBkgFuncSb->GetParameter(iVnBkgPar));}
214  else {initpars.push_back(0.05);}
215  }
216  initpars.push_back(0.10); //initial parameter for signal vn
217  if(fSecondPeak && fDoSecondPeakVn) {initpars.push_back(0.10);} //initial parameter for second peak vn
218 
219  fMassTotFunc = new TF1("fMassTotFunc",this,&AliHFVnVsMassFitter::MassFunc,fMassMin,fMassMax,nparsmass,"AliHFVnVsMassFitter","MassFunc");
220  fVnTotFunc = new TF1("fVnTotFunc",this,&AliHFVnVsMassFitter::vnFunc,fMassMin,fMassMax,nparsvn,"AliHFVnVsMassFitter","vnFunc");
221  SetParNames();
222 
223  ROOT::Math::WrappedMultiTF1 wfTotMass(*fMassTotFunc,1);
224  ROOT::Math::WrappedMultiTF1 wfTotVn(*fVnTotFunc,1);
225 
226  // set data options and ranges
227  ROOT::Fit::DataOptions opt;
228  ROOT::Fit::DataRange rangeMass; //same range for two functions
229 
230  rangeMass.SetRange(fMassMin,fMassMax);
231  ROOT::Fit::BinData dataMass(opt,rangeMass);
232  ROOT::Fit::FillData(dataMass, fMassHisto);
233  ROOT::Fit::BinData dataVn(opt,rangeMass);
234  ROOT::Fit::FillData(dataVn, fVnVsMassHisto);
235 
236  //define the 2 chi squares
237  ROOT::Fit::Chi2Function chi2Mass(dataMass, wfTotMass);
238  ROOT::Fit::Chi2Function chi2Vn(dataVn, wfTotVn);
239 
240  //define the global chi square
241  AliHFGlobalChi2 globalChi2(chi2Mass, chi2Vn);
242 
243  //define fitter
244  ROOT::Fit::Fitter fitter;
245  // create before the parameter settings in order to fix or set range on them
246  fitter.Config().SetParamsSettings(nparsvn,initpars.data()); //set initial parameters from prefits
247  if(fMeanFixed==2 || fMeanFixedFromMassFit) {fitter.Config().ParSettings(fNParsMassBkg+1).Fix();}
248  if(fSigmaFixed==2 || fSigmaFixedFromMassFit) {fitter.Config().ParSettings(fNParsMassBkg+2).Fix();}
249  if(fSecondPeak && fFixSecMass) {fitter.Config().ParSettings(fNParsMassBkg+fNParsMassSgn+1).Fix();}
250  if(fSecondPeak && fFixSecWidth) {fitter.Config().ParSettings(fNParsMassBkg+fNParsMassSgn+2).Fix();}
251  if(fReflections && fFixRflOverSig) {fitter.Config().ParSettings(fNParsMassBkg+fNParsMassSgn+fNParsSec).Fix();}
252 
253  fitter.Config().MinimizerOptions().SetPrintLevel(0);
254  fitter.Config().SetMinimizer("Minuit2","Migrad");
255  for(Int_t iPar=0; iPar<nparsvn; iPar++) {fitter.Config().ParSettings(iPar).SetName(fVnTotFunc->GetParName(iPar));}
256  // fit FCN function directly
257  // (specify optionally data size and flag to indicate that is a chi2 fit
258  fitter.FitFCN(nparsvn,globalChi2,0,dataMass.Size()+dataVn.Size(),kFALSE);
259  ROOT::Fit::FitResult result = fitter.Result();
260  result.Print(std::cout);
261 
262  //set parameters in every function
263  fVnBkgFunc = new TF1("fVnBkgFunc",this,&AliHFVnVsMassFitter::vnBkgFunc,fMassMin,fMassMax,fNParsVnBkg,"AliHFVnVsMassFitter","vnBkgFunc");
264  fMassBkgFunc = new TF1("fMassBkgFunc",this,&AliHFVnVsMassFitter::MassBkg,fMassMin,fMassMax,fNParsMassBkg,"AliHFVnVsMassFitter","MassBkg");
265  fMassSgnFunc = new TF1("fMassSgnFunc",this,&AliHFVnVsMassFitter::MassSignal,fMassMin,fMassMax,fNParsMassSgn,"AliHFVnVsMassFitter","MassSignal");
266  if(fReflections) {fMassRflFunc = new TF1("fMassRflFunc",this,&AliHFVnVsMassFitter::MassRfl,fMassMin,fMassMax,fNParsRfl,"AliHFVnVsMassFitter","MassRfl");}
267  if(fReflections) {fMassBkgRflFunc = new TF1("fMassBkgRflFunc",this,&AliHFVnVsMassFitter::MassBkgRfl,fMassMin,fMassMax,fNParsMassBkg+fNParsRfl,"AliHFVnVsMassFitter","MassBkgRfl");}
268  if(fSecondPeak) {fMassSecPeakFunc = new TF1("fMassSecPeakFunc",this,&AliHFVnVsMassFitter::MassSecondPeak,fMassMin,fMassMax,fNParsSec,"AliHFVnVsMassFitter","MassSecondPeak");}
269  for(Int_t iPar=0; iPar<nparsvn; iPar++) {
270  fVnTotFunc->SetParameter(iPar,result.Parameter(iPar));
271  fVnTotFunc->SetParError(iPar,result.ParError(iPar));
272  if(iPar<nparsmass) {
273  fMassTotFunc->SetParameter(iPar,result.Parameter(iPar));
274  fMassTotFunc->SetParError(iPar,result.ParError(iPar));
275  }
276  if(iPar>=nparsmass && iPar<nparsvn-NvnParsSgn) {
277  fVnBkgFunc->SetParameter(iPar-nparsmass,result.Parameter(iPar));
278  fVnBkgFunc->SetParError(iPar-nparsmass,result.ParError(iPar));
279  }
280  if(iPar>=fNParsMassBkg && iPar<fNParsMassBkg+fNParsMassSgn) {
281  fMassSgnFunc->SetParameter(iPar-fNParsMassBkg,result.Parameter(iPar));
282  }
283  if(iPar<fNParsMassBkg) {
284  fMassBkgFunc->SetParameter(iPar,result.Parameter(iPar));
285  fMassBkgFunc->SetParError(iPar,result.ParError(iPar));
286  if(fReflections) {
287  fMassRflFunc->SetParameter(iPar,result.Parameter(iPar));
288  fMassRflFunc->SetParError(iPar,result.ParError(iPar));
289  }
290  }
291  if(fReflections && (iPar>=fNParsMassBkg+fNParsMassSgn+fNParsSec && iPar<fNParsMassBkg+fNParsMassSgn+fNParsSec+fNParsRfl)) {
292  fMassRflFunc->SetParameter(iPar-(fNParsMassBkg+fNParsMassSgn+fNParsSec),result.Parameter(iPar));
293  fMassRflFunc->SetParError(iPar-(fNParsMassBkg+fNParsMassSgn+fNParsSec),result.ParError(iPar));
294  }
295  if(fSecondPeak && (iPar>=fNParsMassBkg+fNParsMassSgn && iPar<fNParsMassBkg+fNParsMassSgn+fNParsSec)) {
296  fMassSecPeakFunc->SetParameter(iPar-(fNParsMassBkg+fNParsMassSgn),result.Parameter(iPar));
297  fMassSecPeakFunc->SetParError(iPar-(fNParsMassBkg+fNParsMassSgn),result.ParError(iPar));
298  }
299  }
300  if(drawFit) {DrawFit();}
301 
302  fVn = fVnTotFunc->GetParameter(fVnTotFunc->GetNpar()-NvnParsSgn);
303  fVnUncertainty = fVnTotFunc->GetParError(fVnTotFunc->GetNpar()-NvnParsSgn);
304  if(fDoSecondPeakVn) {
305  fVnSecPeak = fVnTotFunc->GetParameter(fVnTotFunc->GetNpar()-1);
306  fVnSecPeakUncertainty = fVnTotFunc->GetParError(fVnTotFunc->GetNpar()-1);
307  }
308  fRawYield = fVnTotFunc->GetParameter(fNParsMassBkg)/fMassHisto->GetBinWidth(10);
309  fRawYieldUncertainty = fVnTotFunc->GetParError(fNParsMassBkg)/fMassHisto->GetBinWidth(10);
310  fMean = fVnTotFunc->GetParameter(fNParsMassBkg+1);
311  fMeanUncertainty = fVnTotFunc->GetParError(fNParsMassBkg+1);
312  fSigma = fVnTotFunc->GetParameter(fNParsMassBkg+2);
313  fSigmaUncertainty = fVnTotFunc->GetParError(fNParsMassBkg+2);
314  fChiSquare = result.MinFcnValue();
315  fNDF = result.Ndf();
316  fProb = result.Prob();
317 
318  return kTRUE;
319 }
320 
321 //______________________________________________________________________________
322 void AliHFVnVsMassFitter::DrawHere(TVirtualPad* c){
324 
325  gStyle->SetOptStat(0);
326  gStyle->SetCanvasColor(0);
327  gStyle->SetFrameFillColor(0);
328  c->Divide(1,2);
329 
330  c->cd(1);
331  fMassHisto->SetTitle("");
332  fMassHisto->SetMarkerStyle(20);
333  fMassHisto->SetMarkerSize(1);
334  fMassHisto->SetMarkerColor(kBlack);
335  fMassHisto->SetLineColor(kBlack);
336  fMassHisto->GetYaxis()->SetRangeUser(0.,fMassHisto->GetMaximum()*1.2);
337  fMassHisto->GetXaxis()->SetRangeUser(fMassMin,fMassMax);
338  fMassHisto->Draw("E");
339  if(fMassFuncFromPrefit) {
340  fMassFuncFromPrefit->SetLineColor(kGray+1);
341  fMassFuncFromPrefit->SetLineStyle(7);
343  fMassFuncFromPrefit->Draw("same");
344  }
345  if(fMassBkgFunc) {
346  fMassBkgFunc->SetLineColor(kRed);
347  fMassBkgFunc->SetRange(fMassMin,fMassMax);
348  fMassBkgFunc->Draw("same");
349  }
350  if(fMassRflFunc) {
351  fMassRflFunc->SetLineColor(kGreen+1);
352  fMassRflFunc->SetRange(fMassMin,fMassMax);
353  fMassRflFunc->Draw("same");
354  }
355  if(fMassBkgRflFunc) {
356  fMassBkgRflFunc->SetLineColor(kRed+1);
357  fMassBkgRflFunc->SetLineStyle(7);
358  fMassBkgRflFunc->SetRange(fMassMin,fMassMax);
359  fMassBkgRflFunc->Draw("same");
360  }
361  if(fMassSecPeakFunc) {
362  fMassSecPeakFunc->SetLineColor(kMagenta+1);
363  fMassSecPeakFunc->SetLineStyle(7);
365  fMassSecPeakFunc->Draw("same");
366  }
367  if(fMassTotFunc) {
368  fMassTotFunc->SetLineColor(kBlue);
369  fMassTotFunc->SetRange(fMassMin,fMassMax);
370  fMassTotFunc->Draw("same");
371  }
372  TPaveText* massinfo = new TPaveText(0.45,0.7,1.,0.87,"NDC");
373  massinfo->SetTextFont(42);
374  massinfo->SetTextSize(0.05);
375  massinfo->SetBorderSize(0);
376  massinfo->SetFillStyle(0);
377  massinfo->SetTextColor(kBlue);
378  massinfo->AddText(Form("mean = %.3f #pm %.3f",fVnTotFunc->GetParameter(fNParsMassBkg+1),fVnTotFunc->GetParError(fNParsMassBkg+1)));
379  massinfo->AddText(Form("sigma = %.3f #pm %.3f",fVnTotFunc->GetParameter(fNParsMassBkg+2),fVnTotFunc->GetParError(fNParsMassBkg+2)));
380  if(fMassSgnFuncType==k2Gaus) {
381  massinfo->AddText(Form("sigma2 = %.3f #pm %.3f",fVnTotFunc->GetParameter(fNParsMassBkg+3),fVnTotFunc->GetParError(fNParsMassBkg+3)));
382  }
383  massinfo->Draw("same");
384 
385  c->cd(2);
386  fVnVsMassHisto->SetTitle("");
387  fVnVsMassHisto->SetMarkerStyle(20);
388  fVnVsMassHisto->SetMarkerSize(1);
389  fVnVsMassHisto->SetMarkerColor(kBlack);
390  fVnVsMassHisto->SetLineColor(kBlack);
391  fVnVsMassHisto->GetYaxis()->SetRangeUser(fVnVsMassHisto->GetMinimum()-0.15,fVnVsMassHisto->GetMaximum()+0.20);
392  fVnVsMassHisto->GetXaxis()->SetRangeUser(fMassMin,fMassMax);
393  fVnBkgFuncSb->SetRange(fMassMin,fMassMax);
394  fVnVsMassHisto->Draw("E");
395  if(fVnBkgFuncSb) {
396  fVnBkgFuncSb->SetLineColor(kGray+1);
397  fVnBkgFuncSb->SetLineStyle(7);
398  fVnBkgFuncSb->Draw("same");
399  }
400  if(fVnBkgFunc) {
401  fVnBkgFunc->SetLineColor(kRed);
402  fVnBkgFunc->SetRange(fMassMin,fMassMax);
403  fVnBkgFunc->Draw("same");
404  }
405  if(fVnTotFunc) {
406  fVnTotFunc->SetLineColor(kBlue);
407  fVnTotFunc->SetRange(fMassMin,fMassMax);
408  fVnTotFunc->Draw("same");
409  }
410 
411  TPaveText* vninfo = new TPaveText(-0.45,0.7,1.,0.87,"NDC");
412  vninfo->SetTextFont(42);
413  vninfo->SetTextSize(0.05);
414  vninfo->SetBorderSize(0);
415  vninfo->SetFillStyle(0);
416  Int_t NvnParsSgn = 1;
417  if(fSecondPeak && fDoSecondPeakVn) {NvnParsSgn+=1;}
418  vninfo->AddText(Form("#it{v}_{%d}^{sgn} = %.3f #pm %.3f",fHarmonic,fVnTotFunc->GetParameter(fVnTotFunc->GetNpar()-NvnParsSgn),fVnTotFunc->GetParError(fVnTotFunc->GetNpar()-NvnParsSgn)));
419  if(fSecondPeak && fDoSecondPeakVn) {vninfo->AddText(Form("#it{v}_{%d}^{sec peak} = %.3f #pm %.3f",fHarmonic,fVnTotFunc->GetParameter(fVnTotFunc->GetNpar()-1),fVnTotFunc->GetParError(fVnTotFunc->GetNpar()-1)));}
420  vninfo->AddText(Form("#chi^{2}/#it{ndf} = %.2f/%d",fChiSquare,fNDF));
421  vninfo->Draw("same");
422 
423  c->Update();
424 }
425 
426 //________________________________________________________________
428 
429  //define proper maxs and mins from histos
430  Double_t tmpmin = TMath::Max(fMassHisto->GetBinLowEdge(1),fVnVsMassHisto->GetBinLowEdge(1));
431  fMassMin=TMath::Max(fMassMin,tmpmin);
432  Double_t tmpmax = TMath::Min(fMassHisto->GetBinLowEdge(fMassHisto->GetNbinsX())+fMassHisto->GetBinWidth(fMassHisto->GetNbinsX()),fVnVsMassHisto->GetBinLowEdge(fVnVsMassHisto->GetNbinsX())+fVnVsMassHisto->GetBinWidth(fVnVsMassHisto->GetNbinsX()));
433  fMassMax=TMath::Min(fMassMax,tmpmax);
434 
443  if(fReflections) {
447  }
448  Bool_t status = fMassFitter->MassFitter(kFALSE);
449 
450  if(status) {
452  fMassFuncFromPrefit->SetName("fMassFuncFromPrefit");
453  }
454 
455  return status;
456 }
457 
458 //________________________________________________________________
460 
461  Double_t mean = fMassFitter->GetMean();
463  const Int_t nMassBins = fVnVsMassHisto->GetNbinsX();
464  Double_t SBbins[nMassBins];
465  Int_t nSBbins=0;
466  for(Int_t iBin=0; iBin<nMassBins; iBin++) {
467  Double_t min = fVnVsMassHisto->GetBinLowEdge(iBin+1);
468  Double_t max = fVnVsMassHisto->GetBinLowEdge(iBin+1)+fVnVsMassHisto->GetBinWidth(iBin+1);
469  if(max<mean){
470  if(max<(mean-fNSigmaForSB*sigma)) {SBbins[iBin]=1; nSBbins++;}
471  else {SBbins[iBin]=0;}
472  }
473  if(min>=mean){
474  if(min>(mean+fNSigmaForSB*sigma)) {SBbins[iBin]=1; nSBbins++;}
475  else {SBbins[iBin]=0;}
476  }
477  }
478  TGraphErrors* gVnVsMassSB = new TGraphErrors(nSBbins);
479  for(Int_t iBin=0; iBin<nMassBins; iBin++) {
480  if(SBbins[iBin]==1) {
481  gVnVsMassSB->SetPoint(iBin,fVnVsMassHisto->GetBinCenter(iBin+1),fVnVsMassHisto->GetBinContent(iBin+1));
482  gVnVsMassSB->SetPointError(iBin,fVnVsMassHisto->GetBinWidth(iBin+1)/2,fVnVsMassHisto->GetBinError(iBin+1));
483  }
484  }
485  fVnBkgFuncSb = new TF1("fVnBkgFuncSb",this,&AliHFVnVsMassFitter::vnBkgFunc,fMassMin,fMassMax,fNParsVnBkg,"AliHFVnVsMassFitter","vnBkgFunc");
486  switch(fVnBkgFuncType) {
487  case 1:
488  fVnBkgFuncSb->SetParName(0,"ConstVnBkg");
489  fVnBkgFuncSb->SetParName(1,"SlopeVnBkg");
490  break;
491  case 2:
492  fVnBkgFuncSb->SetParName(0,"ConstVnBkg");
493  fVnBkgFuncSb->SetParName(1,"Coef1VnBkg");
494  fVnBkgFuncSb->SetParName(2,"Coef2VnBkg");
495  break;
496  default:
497  AliError("Error in setting signal par names: check fVnBkgFuncType");
498  break;
499  }
500  gVnVsMassSB->Fit(fVnBkgFuncSb,"","",fMassMin,fMassMax);
501  Bool_t status=kFALSE;
502  if(fVnBkgFuncSb->GetChisquare()<1000) status=kTRUE;
503 
504  delete gVnVsMassSB;
505  return status;
506 }
507 
508 //________________________________________________________________
510 
511  switch(fMassSgnFuncType) {
512  case 0: //single gaus
513  fNParsMassSgn=3;
514  break;
515  case 1: //double gaus
516  fNParsMassSgn=5;
517  break;
518  default:
519  AliError("Error in computing fMassSgnFuncType: check fMassSgnFuncType");
520  break;
521  }
522 
523  switch(fMassBkgFuncType) {
524  case 0: //expo
525  fNParsMassBkg=2;
526  break;
527  case 1: //lin
528  fNParsMassBkg=2;
529  break;
530  case 2: //pol2
531  fNParsMassBkg=3;
532  break;
533  case 3: //no bkg
534  fNParsMassBkg=1;
535  break;
536  case 4: //power law
537  fNParsMassBkg=2;
538  break;
539  case 5: //power expo
540  fNParsMassBkg=3;
541  break;
542  case 6: //high degree pol
544  break;
545  default:
546  AliError("Error in computing fNParsMassBkg: check fMassBkgFuncType");
547  break;
548  }
549 
550  switch(fVnBkgFuncType) {
551  case 1: //lin
552  fNParsVnBkg=2;
553  break;
554  case 2: //pol2
555  fNParsVnBkg=3;
556  break;
557  default:
558  AliError("Error in computing fNParsVnBkg: check fVnBkgFuncType");
559  break;
560  }
561 
562  if(fReflections) fNParsRfl=1;
563  else fNParsRfl=0;
564 
565  if(fSecondPeak) fNParsSec=3;
566  else fNParsSec=0;
567 }
568 
569 //________________________________________________________________
571 
572  switch(fMassSgnFuncType) {
573  case 0: //single gaus
574  fVnTotFunc->SetParName(fNParsMassBkg,"SgnInt");
575  fVnTotFunc->SetParName(fNParsMassBkg+1,"Mean");
576  fVnTotFunc->SetParName(fNParsMassBkg+2,"Sigma");
577  break;
578  case 1: //double gaus
579  fVnTotFunc->SetParName(fNParsMassBkg,"SgnInt");
580  fVnTotFunc->SetParName(fNParsMassBkg+1,"Mean");
581  fVnTotFunc->SetParName(fNParsMassBkg+2,"Sigma1");
582  fVnTotFunc->SetParName(fNParsMassBkg+3,"Frac");
583  fVnTotFunc->SetParName(fNParsMassBkg+4,"Sigma2");
584  break;
585  default:
586  AliError("Error in setting signal par names: check fMassSgnFuncType");
587  break;
588  }
589 
590  switch(fMassBkgFuncType) {
591  case 0: //expo
592  fVnTotFunc->SetParName(0,"BkgInt");
593  fVnTotFunc->SetParName(1,"Slope");
594  break;
595  case 1: //lin
596  fVnTotFunc->SetParName(0,"BkgInt");
597  fVnTotFunc->SetParName(1,"Slope");
598  break;
599  case 2: //pol2
600  fVnTotFunc->SetParName(0,"BkgInt");
601  fVnTotFunc->SetParName(1,"Coef1");
602  fVnTotFunc->SetParName(2,"Coef1");
603  break;
604  case 3: //no bkg
605  fVnTotFunc->SetParName(0,"Const");
606  break;
607  case 4: //power law
608  fVnTotFunc->SetParName(0,"BkgInt");
609  fVnTotFunc->SetParName(1,"Coef1");
610  break;
611  case 5: //power expo
612  fVnTotFunc->SetParName(0,"BkgInt");
613  fVnTotFunc->SetParName(1,"Coef1");
614  fVnTotFunc->SetParName(2,"Coef2");
615  break;
616  case 6: //high degree pol
617  fVnTotFunc->SetParName(0,"BkgInt");
618  for(Int_t iPar=1; iPar<fNParsMassBkg; iPar++) {fVnTotFunc->SetParName(iPar,Form("Coef%d",iPar));}
619  break;
620  default:
621  AliError("Error in setting signal par names: check fMassBkgFuncType");
622  break;
623  }
624 
625  for(Int_t iPar=0; iPar<fNParsVnBkg; iPar++) {fVnTotFunc->SetParName(fNParsMassBkg+fNParsMassSgn+fNParsSec+fNParsRfl+iPar,fVnBkgFuncSb->GetParName(iPar));}
626 
627  if(fReflections) {fVnTotFunc->SetParName(fNParsMassBkg+fNParsMassSgn+fNParsSec,"ReflOverS");}
628 
629  if(fSecondPeak) {
630  fVnTotFunc->SetParName(fNParsMassBkg+fNParsMassSgn,"SecPeakInt");
631  fVnTotFunc->SetParName(fNParsMassBkg+fNParsMassSgn+1,"SecPeakMean");
632  fVnTotFunc->SetParName(fNParsMassBkg+fNParsMassSgn+2,"SecPeakSigma");
633  }
634  fVnTotFunc->SetParName(fNParsMassBkg+fNParsMassSgn+fNParsSec+fNParsVnBkg,Form("v%dSgn",fHarmonic));
635  if(fSecondPeak && fDoSecondPeakVn) {fVnTotFunc->SetParName(fNParsMassBkg+fNParsMassSgn+fNParsSec+fNParsVnBkg+1,Form("v%dSecPeak",fHarmonic));}
636 }
637 
638 //_________________________________________________________________________
639 void AliHFVnVsMassFitter::Signal(Double_t nOfSigma,Double_t &signal,Double_t &errsignal) const {
642 
643  Double_t minMass=fMean-nOfSigma*fSigma;
644  Double_t maxMass=fMean+nOfSigma*fSigma;
645  Signal(minMass,maxMass,signal,errsignal);
646  return;
647 }
648 
649 //_________________________________________________________________________
650 void AliHFVnVsMassFitter::Signal(Double_t min, Double_t max, Double_t &signal,Double_t &errsignal) const {
653  if(!fMassSgnFunc) {signal=-1; errsignal=0; return;}
654 
655  signal=fMassSgnFunc->Integral(min, max)/(Double_t)fMassHisto->GetBinWidth(1);
656  errsignal=(fRawYieldUncertainty/fRawYield)*signal;/*assume relative error is the same as for total integral*/
657 
658  return;
659 }
660 
661 //___________________________________________________________________________
662 void AliHFVnVsMassFitter::Background(Double_t nOfSigma,Double_t &background,Double_t &errbackground) const {
665 
666  Double_t minMass=fMean-nOfSigma*fSigma;
667  Double_t maxMass=fMean+nOfSigma*fSigma;
668  Background(minMass,maxMass,background,errbackground);
669 
670  return;
671 }
672 
673 //___________________________________________________________________________
674 void AliHFVnVsMassFitter::Background(Double_t min, Double_t max, Double_t &background,Double_t &errbackground) const {
677 
678  if(!fMassBkgFunc) {background=-1; errbackground=0; return;}
679 
680  Double_t intB=fMassBkgFunc->GetParameter(0);
681  Double_t intBerr=fMassBkgFunc->GetParError(0);
682  //relative error evaluation: from histo
683 
684  Int_t leftBand=fMassHisto->FindBin(fMean-4*fSigma);
685  Int_t rightBand=fMassHisto->FindBin(fMean+4*fSigma);
686  intB=fMassHisto->Integral(1,leftBand)+fMassHisto->Integral(rightBand,fMassHisto->GetNbinsX());
687  Double_t sum2=0;
688  for(Int_t iBin=1; iBin<=leftBand; iBin++){
689  sum2+=fMassHisto->GetBinError(iBin)*fMassHisto->GetBinError(iBin);
690  }
691  for(Int_t iBin=rightBand; iBin<=fMassHisto->GetNbinsX(); iBin++){
692  sum2+=fMassHisto->GetBinError(iBin)*fMassHisto->GetBinError(iBin);
693  }
694 
695  intBerr=TMath::Sqrt(sum2);
696 
697  background=fMassBkgFunc->Integral(min,max)/(Double_t)fMassHisto->GetBinWidth(1);
698  errbackground=intBerr/intB*background;
699 
700  return;
701 }
702 
703 //__________________________________________________________________________
704 
705 void AliHFVnVsMassFitter::Significance(Double_t nOfSigma,Double_t &significance,Double_t &errsignificance) const {
708 
709  Double_t minMass=fMean-nOfSigma*fSigma;
710  Double_t maxMass=fMean+nOfSigma*fSigma;
711  Significance(minMass, maxMass, significance, errsignificance);
712 
713  return;
714 }
715 
716 //__________________________________________________________________________
717 
718 void AliHFVnVsMassFitter::Significance(Double_t min, Double_t max, Double_t &significance,Double_t &errsignificance) const {
721 
722  Double_t background,errbackground;
723  Background(min,max,background,errbackground);
724 
725  if (fRawYield+background <= 0.){
726  significance=-1;
727  errsignificance=0;
728  return;
729  }
730 
731  AliVertexingHFUtils::ComputeSignificance(fRawYield,fRawYieldUncertainty,background,errbackground,significance,errsignificance);
732 
733  return;
734 }
735 
736 //________________________________________________________________
738 
739  return TMath::Gaus(x,mean,sigma,kTRUE);
740 }
741 
742 //________________________________________________________________
744 
745  return TMath::Exp(x/coeff)/(coeff*(TMath::Exp(fMassMax/coeff)-TMath::Exp(fMassMin/coeff)));
746 }
747 
748 //________________________________________________________________
750 
751  switch(order) {
752  case 0:
753  if(isnorm) {return 1./(fMassMax-fMassMin);}
754  else {return pars[0];}
755  break;
756  case 1:
757  if(isnorm) {return (pars[0]-pars[1]/2*(fMassMax*fMassMax-fMassMin*fMassMin))/(fMassMax-fMassMin)+pars[1]*x;}
758  else {return pars[0]+pars[1]*x;}
759  break;
760  case 2:
761  if(isnorm) {return (pars[0]-pars[1]/2*(fMassMax*fMassMax-fMassMin*fMassMin)-pars[2]/3*(fMassMax*fMassMax*fMassMax-fMassMin*fMassMin*fMassMin))/(fMassMax-fMassMin)+pars[1]*x;}
762  else {return pars[0]+pars[1]*x+pars[2]*x*x;}
763  }
764  return 0;
765 }
766 
767 //________________________________________________________________
769 
770  Double_t mpi = TDatabasePDG::Instance()->GetParticle(211)->Mass();
771  return pars[0]*(pars[1]+1.)/(TMath::Power(fMassMax-mpi,pars[1]+1.)-TMath::Power(fMassMin-mpi,pars[1]+1.))*TMath::Power(x-mpi,pars[1]);
772 }
773 
774 //________________________________________________________________
776 
777  Double_t mpi = TDatabasePDG::Instance()->GetParticle(211)->Mass();
778  return pars[0]*TMath::Sqrt(x - mpi)*TMath::Exp(-1.*pars[1]*(x-mpi));
779 }
780 
781 //________________________________________________________________
783 
784  Double_t total=pars[0];
785  for(Int_t iT=1; iT<=fPolDegreeBkg; iT++){
786  total+=pars[iT]*TMath::Power(x-fMassParticle,iT)/TMath::Factorial(iT);
787  }
788  return total;
789 }
790 
791 //________________________________________________________________
793 
794  switch(fMassSgnFuncType) {
795  case 0:
796  return pars[0]*GetGausPDF(m[0],pars[1],pars[2]);
797  break;
798  case 1:
799  return pars[0]*(pars[3]*GetGausPDF(m[0],pars[1],pars[2])+(1-pars[3])*GetGausPDF(m[0],pars[1],pars[4]));
800  break;
801  }
802  fRawYieldHelp=pars[0]/fMassHisto->GetBinWidth(1);
803 
804  return 0;
805 }
806 
807 //________________________________________________________________
809 
810  switch(fMassBkgFuncType) {
811  case 0: //exponential
812  return pars[0]*GetExpoPDF(m[0],pars[1]);
813  break;
814  case 1: //linear
815  return GetPolPDF(m[0],pars,1,kTRUE);
816  break;
817  case 2: //parabolic
818  return GetPolPDF(m[0],pars,2,kTRUE);
819  break;
820  case 3: //constant
821  return GetPolPDF(m[0],pars,0,kTRUE);
822  break;
823  case 4: //power law
824  return GetPowerFuncPDF(m[0],pars);
825  break;
826  case 5: //power law expo
827  return GetPowerExpoPDF(m[0],pars);
828  break;
829  case 6: //higher order (>=3) polinomial
830  return GetHigherPolFuncPDF(m[0],pars);
831  break;
832  }
833  return 0;
834 }
835 
836 //_________________________________________________________________________
840  if(!fHistoTemplRfl) return 0;
841 
842  Int_t bin =fHistoTemplRfl->FindBin(m[0]);
843  Double_t value=fHistoTemplRfl->GetBinContent(bin);
844  Int_t binmin=fHistoTemplRfl->FindBin(fMassMin*1.00001);
845  Int_t binmax=fHistoTemplRfl->FindBin(fMassMax*0.99999);
846  Double_t norm=fHistoTemplRfl->Integral(binmin,binmax)*fHistoTemplRfl->GetBinWidth(bin);
847  if(TMath::Abs(value)<1.e-14 && fSmoothRfl){// very rough, assume a constant trend, much better would be a pol1 or pol2 over a broader range
848  value+=fHistoTemplRfl->GetBinContent(bin-1)+fHistoTemplRfl->GetBinContent(bin+1);
849  value/=3.;
850  }
851  return pars[0]*value/norm*fRawYieldHelp*fMassHisto->GetBinWidth(1);
852 }
853 
854 //_________________________________________________________________________
856 
857  if(!fHistoTemplRfl) {return MassBkg(m,pars);}
858  else {
859  //bkg mass parameters
860  const Int_t nBkgPars = fNParsMassBkg;
861  Double_t bkgpars[nBkgPars];
862  for(Int_t iPar=0; iPar<fNParsMassBkg; iPar++) {bkgpars[iPar] = pars[iPar];}
863  //reflection parameters
864  Double_t rflpars[1]; //maximum number of parameters for rfl = 1 for the implemented functions
865  for(Int_t iPar=0; iPar<fNParsRfl; iPar++) {rflpars[iPar] = pars[iPar+fNParsMassBkg];}
866  return MassBkg(m,bkgpars)+MassRfl(m,rflpars);
867  }
868 }
869 
870 //_________________________________________________________________________
874 
875  return pars[0]*GetGausPDF(m[0],pars[1],pars[2]);
876 }
877 
878 //________________________________________________________________
880 
881  switch(fVnBkgFuncType) {
882  case 1: //linear
883  return GetPolPDF(m[0],pars,1,kFALSE);
884  break;
885  case 2: //parabolic
886  return GetPolPDF(m[0],pars,2,kFALSE);
887  break;
888  }
889  return 0;
890 }
891 
892 //________________________________________________________________
894 
895  //bkg mass parameters
896  const Int_t nBkgPars = fNParsMassBkg;
897  Double_t bkgpars[nBkgPars];
898  for(Int_t iPar=0; iPar<fNParsMassBkg; iPar++) {bkgpars[iPar] = pars[iPar];}
899  //signal mass parameters
900  Double_t sgnpars[5]; //maximum number of parameters for sgn = 5 for the implemented functions
901  for(Int_t iPar=0; iPar<fNParsMassSgn; iPar++) {sgnpars[iPar] = pars[iPar+fNParsMassBkg];}
902  //second peak parameters
903  Double_t secpeakpars[3]; //maximum number of parameters for second peak = 3 for the implemented functions
904  for(Int_t iPar=0; iPar<fNParsSec; iPar++) {secpeakpars[iPar] = pars[iPar+fNParsMassBkg+fNParsMassSgn];}
905  //reflection parameters
906  Double_t rflpars[1]; //maximum number of parameters for rfl = 1 for the implemented functions
907  for(Int_t iPar=0; iPar<fNParsRfl; iPar++) {rflpars[iPar] = pars[iPar+fNParsMassBkg+fNParsMassSgn+fNParsSec];}
908 
909  Double_t total = MassSignal(m,sgnpars)+MassBkg(m,bkgpars);
910  if(fSecondPeak) {total += MassSecondPeak(m,secpeakpars);}
911  if(fReflections) {total += MassRfl(m,rflpars);}
912 
913  return total;
914 }
915 
916 //________________________________________________________________
918 
919  //bkg mass parameters
920  const Int_t nBkgPars = fNParsMassBkg;
921  Double_t massbkgpars[nBkgPars];
922  for(Int_t iPar=0; iPar<fNParsMassBkg; iPar++) {massbkgpars[iPar] = pars[iPar];}
923  //signal mass parameters
924  Double_t masssgnpars[5]; //maximum number of parameters for mass sgn = 5 for the implemented functions
925  for(Int_t iPar=0; iPar<fNParsMassSgn; iPar++) {masssgnpars[iPar] = pars[iPar+fNParsMassBkg];}
926  //second peak parameters
927  Double_t secpeakpars[3]; //maximum number of parameters for second peak = 3 for the implemented functions
928  for(Int_t iPar=0; iPar<fNParsSec; iPar++) {secpeakpars[iPar] = pars[iPar+fNParsMassBkg+fNParsMassSgn];}
929  //reflection parameters
930  Double_t rflpars[1]; //maximum number of parameters for rfl = 1 for the implemented functions
931  for(Int_t iPar=0; iPar<fNParsRfl; iPar++) {rflpars[iPar] = pars[iPar+fNParsMassBkg+fNParsMassSgn+fNParsSec];}
932  //bkg vn parameters
933  Double_t vnbkgpars[3]; //maximum number of parameters for vn bkg = 3 for the implemented functions
934  for(Int_t iPar=0; iPar<fNParsVnBkg; iPar++) {vnbkgpars[iPar] = pars[iPar+fNParsMassSgn+fNParsMassBkg+fNParsSec+fNParsRfl];}
935  //signal vn parameter
936  Double_t vnSgn = pars[fNParsMassSgn+fNParsMassBkg+fNParsSec+fNParsRfl+fNParsVnBkg];
937  //second peak vn parameter
938  Double_t vnSecPeak = 0;
939  if(fSecondPeak && fDoSecondPeakVn) {vnSecPeak = pars[fNParsMassSgn+fNParsMassBkg+fNParsSec+fNParsRfl+fNParsVnBkg+1];}
940 
941 
942  Double_t vnBkg = vnBkgFunc(m,vnbkgpars);
943  Double_t Sgn = MassSignal(m,masssgnpars);
944  Double_t Bkg = MassBkg(m,massbkgpars);
945  if(fReflections) {Bkg += MassRfl(m,rflpars);}
946  Double_t SecPeak = 0;
947  if(fSecondPeak) {SecPeak += MassSecondPeak(m,secpeakpars);}
948 
949  if(fSecondPeak && fDoSecondPeakVn) {return (vnSgn*Sgn+vnBkg*Bkg+vnSecPeak*SecPeak)/(Sgn+Bkg+SecPeak);}
950  else {return (vnSgn*Sgn+vnBkg*(Bkg+SecPeak))/(Sgn+Bkg+SecPeak);}
951 }
952 
953 //______________________________________________________________________________
957  TCanvas* c0=new TCanvas("c0","",600,800);
958  DrawHere(c0);
959 }
Double_t MassSecondPeak(Double_t *m, Double_t *par)
Int_t fNParsVnBkg
number of parameters in mass bkg fit function
Bool_t fFixSecWidth
flag to fix the position of the 2nd peak
double Double_t
Definition: External.C:58
void SetFixReflOverS(Double_t rovers)
Int_t fNParsRfl
flag use/not use reflections
Bool_t fReflections
degree of polynomial expansion for back fit (option 6 for back)
void IncludeSecondGausPeak(Double_t mass, Bool_t fixm, Double_t width, Bool_t fixw)
Double_t MassRfl(Double_t *m, Double_t *par)
Bool_t SimultaneusFit(Bool_t drawFit=kTRUE)
Double_t GetPowerFuncPDF(Double_t x, Double_t *pars)
Bool_t fMeanFixedFromMassFit
flag to fix peak width from mass prefit
Double_t fMean
uncertainty on mass peak width from simultaneus fit
TH1F * fHistoTemplRfl
switch for fix refl/signal
Double_t fMassParticle
flag to fix peak position from mass prefit
Bool_t fFixSecMass
width of the 2nd peak
Int_t fNParsSec
fit function for second peak
Bool_t fSecondPeak
internal variable for fit with reflections
Bool_t fSigmaFixedFromMassFit
initialization for peak position
Double_t GetHigherPolFuncPDF(Double_t x, Double_t *pars)
void Signal(Double_t nOfSigma, Double_t &signal, Double_t &errsignal) const
Int_t MassFitter(Bool_t draw=kTRUE)
Double_t fSecMass
number of parameters in second peak fit function
TF1 * fMassTotFunc
mass signal fit function (final, after simultaneus fit)
Double_t vnBkgFunc(Double_t *m, Double_t *pars)
Double_t fRflOverSig
fit parameters in reflection fit function
Double_t fMassMax
upper mass limit
AliHFInvMassFitter class for the fit of invariant mass distribution of charm hadrons.
TCanvas * c
Definition: TestFitELoss.C:172
Int_t fMeanFixed
flag to fix peak width
Double_t fVn
lower mass limit
Int_t fMassSgnFuncType
vn vs. mass histogram to fit
void DefineNumberOfParameters()
private methods
Double_t GetSigma() const
void SetPolDegreeForBackgroundFit(Int_t deg)
Double_t fChiSquare
uncertainty raw yield from simultaneus fit
void SetFixGaussianSigma(Double_t sigma)
Double_t fRawYieldUncertainty
raw yield from simultaneus fit
Double_t GetPowerExpoPDF(Double_t x, Double_t *pars)
Int_t fNSigmaForSB
simultaneus fit probability
Double_t * sigma
TF1 * fMassBkgFunc
mass fit function (1st step, from prefit)
Double_t fMaxRefl
minimum for refelction histo
Double_t vnFunc(Double_t *m, Double_t *pars)
Double_t fVnUncertainty
vn of the signal from fit
Double_t fProb
simultaneus fit number of degree of freedom
Double_t MassSignal(Double_t *m, Double_t *pars)
void SetFixGaussianMean(Double_t mean)
int Int_t
Definition: External.C:63
TF1 * fMassSecPeakFunc
switch off/on second peak (for D+->KKpi in Ds)
Double_t GetMean() const
TF1 * fVnBkgFunc
vn bkg fit function (1st step from SB prefit)
Double_t GetExpoPDF(Double_t x, Double_t slope)
Int_t fPolDegreeBkg
flag to fix peak position
Double_t fSecWidth
position of the 2nd peak
TF1 * fVnTotFunc
vn bkg fit function (final, after simultaneus fit)
Int_t fHarmonic
vn uncertainty of second peak from fit
Double_t fSigmaInit
number of sigma for sidebands region (vn bkg prefit)
Int_t fNParsMassBkg
number of parameters in mass signal fit function
Double_t MassBkgRfl(Double_t *m, Double_t *par)
Double_t fMeanUncertainty
mass peak position from simultaneus fit
Double_t GetPolPDF(Double_t x, Double_t *pars, Int_t order, Bool_t isnorm=kTRUE)
TH1F * SetTemplateReflections(const TH1 *h, TString opt, Double_t minRange, Double_t maxRange)
Double_t fSigma
uncertainty on vn of the signal from simultaneus fit
TF1 * fVnBkgFuncSb
mass fit function (final, after simultaneus fit)
Double_t fMeanInit
initialization for peak width
TList * fitter
Definition: DrawAnaELoss.C:26
TH1F * fMassHisto
data members
void DrawHere(TVirtualPad *c)
TH1F * fHistoTemplRflInit
histogram with reflection template
Double_t GetGausPDF(Double_t x, Double_t mean, Double_t sigma)
fit functions
Int_t fNDF
simultaneus fit chi square
static void ComputeSignificance(Double_t signal, Double_t errsignal, Double_t background, Double_t errbackground, Double_t &significance, Double_t &errsignificance)
Significance calculator.
Bool_t fDoSecondPeakVn
vn of second peak from fit
TF1 * fMassRflFunc
initial histogram with reflection template
Double_t MassFunc(Double_t *m, Double_t *pars)
void Background(Double_t nOfSigma, Double_t &background, Double_t &errbackground) const
AliHFInvMassFitter * fMassFitter
vn fit function (final, after simultaneus fit)
TF1 * fMassBkgRflFunc
fit function for reflections
void SetInitialReflOverS(Double_t rovers)
TF1 * fMassFuncFromPrefit
type of vn bkg fit function
TH1F * fVnVsMassHisto
mass histogram to fit
ClassImp(AliAnalysisTaskCRC) AliAnalysisTaskCRC
Int_t fNParsMassSgn
mass of selected particle
void SetInitialGaussianSigma(Double_t sigma)
void SetInitialGaussianMean(Double_t mean)
Double_t minMass
Double_t fSigmaUncertainty
mass peak width from simultaneus fit
Double_t fRawYield
uncertainty on mass peak position from simultaneus fit
Int_t fSigmaFixed
number of parameters in vn bkg fit function
Double_t fMinRefl
refelction option
Int_t fVnBkgFuncType
type of mass bkg fit function
void Significance(Double_t nOfSigma, Double_t &significance, Double_t &errsignificance) const
Int_t fMassBkgFuncType
type of mass signal fit function
TF1 * fMassSgnFunc
mass bkg fit function (final, after simultaneus fit)
Double_t maxMass
bool Bool_t
Definition: External.C:53
Double_t fRawYieldHelp
switch for smoothing of reflection template
Bool_t fFixRflOverSig
reflection/signal
TString fRflOpt
mass bkg fit function plus reflections (final, after simultaneus fit)
Double_t fMassMin
mass fitter for mass prefit
Double_t MassBkg(Double_t *m, Double_t *pars)
Double_t fVnSecPeakUncertainty
flag to introduce second peak vn in the vn vs. mass fit
Bool_t fSmoothRfl
maximum for refelction histo
Double_t fVnSecPeak
flag to fix the width of the 2nd peak