AliPhysics  97344c9 (97344c9)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliHFInvMassFitter.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-2019, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Author: The ALICE Off-line Project. *
5  * Contributors are mentioned in the code where appropriate. *
6  * *
7  * Permission to use, copy, modify and distribute this software and its *
8  * documentation strictly for non-commercial purposes is hereby granted *
9  * without fee, provided that the above copyright notice appears in all *
10  * copies and that both the copyright notice and this permission notice *
11  * appear in the supporting documentation. The authors make no claims *
12  * about the suitability of this software for any purpose. It is *
13  * provided "as is" without express or implied warranty. *
14  **************************************************************************/
15 
16 #include <TH1F.h>
17 #include <TF1.h>
18 #include <TMath.h>
19 #include <TVirtualFitter.h>
20 #include <TDatabasePDG.h>
21 #include <TCanvas.h>
22 #include <TStyle.h>
23 #include <TPaveText.h>
24 
25 #include "AliVertexingHFUtils.h"
26 #include "AliHFInvMassFitter.h"
27 
31 
47 
48 //__________________________________________________________________________
50  TNamed(),
51  fHistoInvMass(0x0),
52  fMinMass(0),
53  fMaxMass(5),
54  fTypeOfFit4Bkg(kExpo),
55  fPolDegreeBkg(4),
56  fCurPolDegreeBkg(-1),
57  fMassParticle(1.864),
58  fTypeOfFit4Sgn(kGaus),
59  fMass(1.865),
60  fMassErr(0.),
61  fSigmaSgn(0.012),
62  fSigmaSgnErr(0.),
63  fFixedMean(kFALSE),
64  fFixedSigma(kFALSE),
65  fFixedRawYield(-1.),
66  fNParsSig(3),
67  fNParsBkg(2),
68  fOnlySideBands(kFALSE),
69  fNSigma4SideBands(4.),
70  fFitOption("L,E"),
71  fRawYield(0.),
72  fRawYieldErr(0.),
73  fSigFunc(0x0),
74  fBkgFuncSb(0x0),
75  fBkgFunc(0x0),
76  fBkgFuncRefit(0x0),
77  fReflections(kFALSE),
78  fNParsRfl(0),
79  fRflOverSig(0),
80  fFixRflOverSig(kFALSE),
81  fHistoTemplRfl(0x0),
82  fSmoothRfl(kFALSE),
83  fRawYieldHelp(0),
84  fRflFunc(0x0),
85  fBkRFunc(0x0),
86  fSecondPeak(kFALSE),
87  fNParsSec(0),
88  fSecMass(-999.),
89  fSecWidth(9999.),
90  fFixSecMass(kFALSE),
91  fFixSecWidth(kFALSE),
92  fSecFunc(0x0),
93  fTotFunc(0x0)
94 {
96 }
97 
98 //__________________________________________________________________________
99 AliHFInvMassFitter::AliHFInvMassFitter(const TH1F *histoToFit, Double_t minvalue, Double_t maxvalue, Int_t fittypeb, Int_t fittypes):
100  TNamed(),
101  fHistoInvMass(0x0),
102  fMinMass(minvalue),
103  fMaxMass(maxvalue),
104  fTypeOfFit4Bkg(fittypeb),
105  fPolDegreeBkg(4),
106  fCurPolDegreeBkg(-1),
107  fMassParticle(1.864),
108  fTypeOfFit4Sgn(fittypes),
109  fMass(1.865),
110  fMassErr(0.),
111  fSigmaSgn(0.012),
112  fSigmaSgnErr(0.),
113  fFixedMean(kFALSE),
114  fFixedSigma(kFALSE),
115  fFixedRawYield(-1.),
116  fNParsSig(3),
117  fNParsBkg(2),
118  fOnlySideBands(kFALSE),
119  fNSigma4SideBands(4.),
120  fFitOption("L,E"),
121  fRawYield(0.),
122  fRawYieldErr(0.),
123  fSigFunc(0x0),
124  fBkgFuncSb(0x0),
125  fBkgFunc(0x0),
126  fBkgFuncRefit(0x0),
127  fReflections(kFALSE),
128  fNParsRfl(0),
129  fRflOverSig(0),
130  fFixRflOverSig(kFALSE),
131  fHistoTemplRfl(0x0),
132  fSmoothRfl(kFALSE),
133  fRawYieldHelp(0),
134  fRflFunc(0x0),
135  fBkRFunc(0x0),
136  fSecondPeak(kFALSE),
137  fNParsSec(0),
138  fSecMass(-999.),
139  fSecWidth(9999.),
140  fFixSecMass(kFALSE),
141  fFixSecWidth(kFALSE),
142  fSecFunc(0x0),
143  fTotFunc(0x0)
144 {
146  fHistoInvMass=(TH1F*)histoToFit->Clone("fHistoInvMass");
147  fHistoInvMass->SetDirectory(0);
149 }
150 //_________________________________________________________________________
152 
154 
155  delete fHistoInvMass;
156  delete fSigFunc;
157  delete fBkgFunc;
158  delete fBkgFuncSb;
159  delete fBkgFuncRefit;
160  delete fRflFunc;
161  delete fBkRFunc;
162  delete fSecFunc;
163  delete fTotFunc;
164  delete fHistoTemplRfl;
165 
166 }
167 //__________________________________________________________________________
171 
172  switch (fTypeOfFit4Bkg) {
173  case 0:
174  fNParsBkg=2;
175  break;
176  case 1:
177  fNParsBkg=2;
178  break;
179  case 2:
180  fNParsBkg=3;
181  break;
182  case 3:
183  fNParsBkg=1;
184  break;
185  case 4:
186  fNParsBkg=2;
187  break;
188  case 5:
189  fNParsBkg=3;
190  break;
191  case 6:
193  break;
194  default:
195  AliError("Error in computing fNParsBkg: check fTypeOfFit4Bkg");
196  break;
197  }
198 
199  switch (fTypeOfFit4Sgn) {
200  case 0:
201  fNParsSig=3;
202  break;
203  case 1:
204  fNParsSig=5;
205  break;
206  default:
207  AliError("Error in computing fNParsSig: check fTypeOfFit4Sgn");
208  break;
209  }
210 
211  if(fReflections) fNParsRfl=1;
212  else fNParsSec=0;
213 
214  if(fSecondPeak) fNParsSec=3;
215  else fNParsSec=0;
216 
217 }
218 //__________________________________________________________________________
222 
223  TVirtualFitter::SetDefaultFitter("Minuit");
224 
225  Double_t integralHisto=fHistoInvMass->Integral(fHistoInvMass->FindBin(fMinMass),fHistoInvMass->FindBin(fMaxMass),"width");
226 
227  fOnlySideBands = kTRUE;
228  fBkgFuncSb = CreateBackgroundFitFunction("funcbkgsb",integralHisto);
229  Int_t status=-1;
230  printf("\n--- First fit with only background on the side bands - Exclusion region = %.2f sigma ---\n",fNSigma4SideBands);
231  if(fTypeOfFit4Bkg==6){
233  fHistoInvMass->GetListOfFunctions()->Add(fBkgFuncSb);
234  fHistoInvMass->GetFunction(fBkgFuncSb->GetName())->SetBit(1<<9,kTRUE);
235  status=0;
236  }
237  }
238  else status=fHistoInvMass->Fit("funcbkgsb",Form("R,%s,+,0",fFitOption.Data()));
239  fBkgFuncSb->SetLineColor(kGray+1);
240  if (status != 0){
241  printf(" ---> Failed first fit with only background, minuit status = %d\n",status);
242  return kFALSE;
243  }
244 
245  fOnlySideBands = kFALSE;
246  if(!fBkgFunc){
247  fBkgFunc = CreateBackgroundFitFunction("funcbkg",integralHisto);
248  for(Int_t ipar=0; ipar<fNParsBkg; ipar++) fBkgFunc->SetParameter(ipar,fBkgFuncSb->GetParameter(ipar));
249  }
250  fBkgFunc->SetLineColor(kGray+1);
251 
252  printf("\n--- Estimate signal counts in the peak region ---\n");
253  Double_t estimSignal=CheckForSignal(fMass,fSigmaSgn);
254  if(estimSignal<0.){
255  if(draw) DrawFit();
256  return kTRUE;
257  }
258 
259  printf("\n--- Final fit with signal+background on the full range ---\n");
260  fRawYieldHelp=estimSignal; // needed for reflection normalization
261  if(!fBkgFuncRefit){
262  fBkgFuncRefit = CreateBackgroundFitFunction("funcbkgrefit",integralHisto);
263  for(Int_t ipar=0; ipar<fNParsBkg; ipar++) fBkgFuncRefit->SetParameter(ipar,fBkgFunc->GetParameter(ipar));
264  }
265  fBkgFuncRefit->SetLineColor(2);
266  fSigFunc = CreateSignalFitFunction("fsigfit",estimSignal);
267  if(fSecondPeak){
268  printf(" ---> Final fit includes a second inv. mass peak\n");
270  fSecFunc = CreateSecondPeakFunction("fsecpeak",estimSec);
271  }
272  if(fReflections){
273  printf(" ---> Final fit includes reflections\n");
274  fRflFunc = CreateReflectionFunction("freflect");
276  }
277  fTotFunc = CreateTotalFitFunction("funcmass");
278  status=fHistoInvMass->Fit("funcmass",Form("R,%s,+,0",fFitOption.Data()));
279  if (status != 0){
280  printf(" ---> Failed fit with signal+background, minuit status = %d\n",status);
281  return kFALSE;
282  }
283  for(Int_t ipar=0; ipar<fNParsBkg; ipar++) fBkgFuncRefit->SetParameter(ipar,fTotFunc->GetParameter(ipar));
284  for(Int_t ipar=0; ipar<fNParsSig; ipar++) fSigFunc->SetParameter(ipar,fTotFunc->GetParameter(ipar+fNParsBkg));
285  if(fSecondPeak){
286  for(Int_t ipar=0; ipar<fNParsSec; ipar++)fSecFunc->SetParameter(ipar,fTotFunc->GetParameter(ipar+fNParsBkg+fNParsSig));
287  fSecFunc->SetLineColor(kMagenta+1);
288  fSecFunc->SetLineStyle(3);
289  }
290  if(fReflections){
291  for(Int_t ipar=0; ipar<fNParsRfl; ipar++){
292  fRflFunc->SetParameter(ipar,fTotFunc->GetParameter(ipar+fNParsBkg+fNParsSig+fNParsSec));
293  fRflFunc->SetParError(ipar,fTotFunc->GetParError(ipar+fNParsBkg+fNParsSig+fNParsSec));
294  fBkRFunc->SetParameter(ipar+fNParsBkg,fTotFunc->GetParameter(ipar+fNParsBkg+fNParsSig+fNParsSec));
295  }
296  for(Int_t ipar=0; ipar<fNParsBkg; ipar++) fBkRFunc->SetParameter(ipar,fTotFunc->GetParameter(ipar));
297  fRflFunc->SetLineColor(kGreen+1);
298  fBkRFunc->SetLineColor(kRed+1);
299  fBkRFunc->SetLineStyle(7);
300  }
301  fMass=fSigFunc->GetParameter(1);
302  fMassErr=fSigFunc->GetParError(1);
303  fSigmaSgn=fSigFunc->GetParameter(2);
304  fSigmaSgnErr=fSigFunc->GetParError(2);
305  fTotFunc->SetLineColor(4);
306  fRawYield=fTotFunc->GetParameter(fNParsBkg)/fHistoInvMass->GetBinWidth(1);
307  fRawYieldErr=fTotFunc->GetParError(fNParsBkg)/fHistoInvMass->GetBinWidth(1);
309  if(draw) DrawFit();
310  return kTRUE;
311 }
312 //______________________________________________________________________________
316 
317  TCanvas* c0=new TCanvas("c0");
318  DrawHere(c0);
319 }
320 //______________________________________________________________________________
321 void AliHFInvMassFitter::DrawHere(TVirtualPad* c){
324 
325  gStyle->SetOptStat(0);
326  gStyle->SetCanvasColor(0);
327  gStyle->SetFrameFillColor(0);
328  c->cd();
329  fHistoInvMass->GetXaxis()->SetRangeUser(fMinMass,fMaxMass);
330  fHistoInvMass->SetMarkerStyle(20);
331  fHistoInvMass->SetMinimum(0.);
332  fHistoInvMass->Draw("PE");
333  if(fBkgFunc) fBkgFunc->Draw("same");
334  if(fBkgFuncRefit) fBkgFuncRefit->Draw("same");
335  if(fBkRFunc) fBkRFunc->Draw("same");
336  if(fRflFunc) fRflFunc->Draw("same");
337  if(fSecFunc) fSecFunc->Draw("same");
338  if(fTotFunc) fTotFunc->Draw("same");
339  TPaveText *pinfos=new TPaveText(0.12,0.65,0.47,0.89,"NDC");
340  TPaveText *pinfom=new TPaveText(0.6,0.7,1.,.87,"NDC");
341  pinfos->SetBorderSize(0);
342  pinfos->SetFillStyle(0);
343  pinfom->SetBorderSize(0);
344  pinfom->SetFillStyle(0);
345  if(fTotFunc){
346  pinfom->SetTextColor(kBlue);
347  for(Int_t ipar=1; ipar<fNParsSig; ipar++){
348  pinfom->AddText(Form("%s = %.3f #pm %.3f",fTotFunc->GetParName(ipar+fNParsBkg),fTotFunc->GetParameter(ipar+fNParsBkg),fTotFunc->GetParError(ipar+fNParsBkg)));
349  }
350  pinfom->Draw();
351 
352  Double_t nsigma=3;
353  Double_t bkg,errbkg;
354  Background(nsigma,bkg,errbkg);
355  Double_t signif,errsignif;
356  Significance(nsigma,signif,errsignif);
357 
358  pinfos->AddText(Form("S = %.0f #pm %.0f ",fRawYield,fRawYieldErr));
359  pinfos->AddText(Form("B (%.0f#sigma) = %.0f #pm %.0f",nsigma,bkg,errbkg));
360  pinfos->AddText(Form("S/B = (%.0f#sigma) %.4f ",nsigma,fRawYield/bkg));
361  if(fRflFunc) pinfos->AddText(Form("Refl/Sig = %.3f #pm %.3f ",fRflFunc->GetParameter(0),fRflFunc->GetParError(0)));
362  pinfos->AddText(Form("Signif (%.0f#sigma) = %.1f #pm %.1f ",nsigma,signif,errsignif));
363  pinfos->Draw();
364  }
365  c->Update();
366  return;
367 }
368 //______________________________________________________________________________
372 
373  Double_t minForSig=mean-4.*sigma;
374  Double_t maxForSig=mean+4.*sigma;
375  Int_t binForMinSig=fHistoInvMass->FindBin(minForSig);
376  Int_t binForMaxSig=fHistoInvMass->FindBin(maxForSig);
377  Double_t sum=0.;
378  Double_t sumback=0.;
379  fBkgFunc->Print();
380  for(Int_t ibin=binForMinSig; ibin<=binForMaxSig; ibin++){
381  sum+=fHistoInvMass->GetBinContent(ibin);
382  sumback+=fBkgFunc->Eval(fHistoInvMass->GetBinCenter(ibin));
383  }
384  Double_t diffUnderPeak=(sum-sumback);
385  printf(" ---> IntegralUnderFitFunc=%f IntegralUnderHisto=%f EstimatedSignal=%f\n",sum,sumback,diffUnderPeak);
386  if(diffUnderPeak/TMath::Sqrt(sum)<1.){
387  printf(" ---> (Tot-Bkg)/sqrt(Tot)=%f ---> Likely no signal/\n",diffUnderPeak/TMath::Sqrt(sum));
388  return -1;
389  }
390  return diffUnderPeak*fHistoInvMass->GetBinWidth(1);
391 }
392 
393 //______________________________________________________________________________
397 
399  TF1* funcbkg = new TF1(fname.Data(),this,&AliHFInvMassFitter::FitFunction4Bkg,fMinMass,fMaxMass,fNParsBkg,"AliHFInvMassFitter","FitFunction4Bkg");
400  switch (fTypeOfFit4Bkg) {
401  case 0: //gaus+expo
402  funcbkg->SetParNames("BkgInt","Slope");
403  funcbkg->SetParameters(integral,-2.);
404  break;
405  case 1:
406  funcbkg->SetParNames("BkgInt","Slope");
407  funcbkg->SetParameters(integral,-100.);
408  break;
409  case 2:
410  funcbkg->SetParNames("BkgInt","Coef1","Coef2");
411  funcbkg->SetParameters(integral,-10.,5);
412  break;
413  case 3:
414  funcbkg->SetParNames("Const");
415  funcbkg->SetParameter(0,0.);
416  funcbkg->FixParameter(0,0.);
417  break;
418  case 4:
419  funcbkg->SetParNames("BkgInt","Coef1");
420  funcbkg->SetParameters(integral,0.5);
421  break;
422  case 5:
423  funcbkg->SetParNames("BkgInt","Coef1","Coef2");
424  funcbkg->SetParameters(integral,-10.,5.);
425  break;
426  case 6:
427  for(Int_t j=0;j<fNParsBkg; j++){
428  funcbkg->SetParName(j,Form("Coef%d",j));
429  funcbkg->SetParameter(j,0);
430  }
431  funcbkg->SetParameter(0,integral);
432  break;
433  default:
434  AliError(Form("Wrong choice of fTypeOfFit4Bkg (%d)",fTypeOfFit4Bkg));
435  delete funcbkg;
436  return 0x0;
437  break;
438  }
439  funcbkg->SetLineColor(kBlue+3);
440  return funcbkg;
441 }
442 //______________________________________________________________________________
446 
447  TF1* funcsec = new TF1(fname.Data(),this,&AliHFInvMassFitter::FitFunction4SecPeak,fMinMass,fMaxMass,3,"AliHFInvMassFitter","FitFunction4SecPeak");
448  funcsec->SetParameter(0,integsig);
449  funcsec->SetParameter(1,fSecMass);
450  if(fFixSecMass) funcsec->FixParameter(1,fSecMass);
451  funcsec->SetParameter(2,fSecWidth);
452  if(fFixSecWidth) funcsec->FixParameter(2,fSecWidth);
453  funcsec->SetParNames("SecPeakInt","SecMean","SecSigma");
454  return funcsec;
455 }
456 //______________________________________________________________________________
460  TF1* funcrfl = new TF1(fname.Data(),this,&AliHFInvMassFitter::FitFunction4Refl,fMinMass,fMaxMass,1,"AliHFInvMassFitter","FitFunction4Refl");
461  funcrfl->SetParameter(0,fRflOverSig);
462  funcrfl->SetParLimits(0,0.,1.);
463  if(fFixRflOverSig) funcrfl->FixParameter(0,fRflOverSig);
464  funcrfl->SetParNames("ReflOverS");
465  return funcrfl;
466 }
467 //______________________________________________________________________________
472  Int_t totParams=fNParsBkg+fNParsRfl;
473  TF1* fbr=new TF1(fname.Data(),this,&AliHFInvMassFitter::FitFunction4BkgAndRefl,fMinMass,fMaxMass,totParams,"AliHFInvMassFitter","FitFunction4BkgAndRefl");
474  for(Int_t ipar=0; ipar<fNParsBkg; ipar++){
475  fbr->SetParameter(ipar,fBkgFunc->GetParameter(ipar));
476  fbr->SetParName(ipar,fBkgFunc->GetParName(ipar));
477  }
478  for(Int_t ipar=0; ipar<fNParsRfl; ipar++){
479  fbr->SetParameter(ipar+fNParsBkg,fRflFunc->GetParameter(ipar));
480  fbr->SetParName(ipar+fNParsBkg,fRflFunc->GetParName(ipar));
481  // par limits not set because this function is not used for fitting but only for drawing
482  }
483  return fbr;
484 }
485 //______________________________________________________________________________
489 
491  TF1* funcsig = new TF1(fname.Data(),this,&AliHFInvMassFitter::FitFunction4Sgn,fMinMass,fMaxMass,fNParsSig,"AliHFInvMassFitter","FitFunction4Sgn");
492  if(fTypeOfFit4Sgn==kGaus){
493  funcsig->SetParameter(0,integsig);
494  if(fFixedRawYield>-0.1) funcsig->FixParameter(0,fFixedRawYield);
495  funcsig->SetParameter(1,fMass);
496  if(fFixedMean) funcsig->FixParameter(1,fMass);
497  funcsig->SetParameter(2,fSigmaSgn);
498  if(fFixedSigma) funcsig->FixParameter(2,fSigmaSgn);
499  funcsig->SetParNames("SgnInt","Mean","Sigma");
500  }
501  if(fTypeOfFit4Sgn==k2Gaus){
502  funcsig->SetParameter(0,integsig);
503  if(fFixedRawYield>-0.1) funcsig->FixParameter(0,fFixedRawYield);
504  funcsig->SetParameter(1,fMass);
505  if(fFixedMean) funcsig->FixParameter(1,fMass);
506  funcsig->SetParameter(2,fSigmaSgn);
507  funcsig->SetParLimits(2,0.004,0.05);
508  if(fFixedSigma) funcsig->FixParameter(2,fSigmaSgn);
509  funcsig->SetParameter(3,0.2);
510  funcsig->SetParLimits(3,0.,1.);
511  funcsig->SetParameter(4,fSigmaSgn);
512  funcsig->SetParLimits(4,0.004,0.05);
513  funcsig->SetParNames("SgnInt","Mean","Sigma1","Frac","Sigma2");
514  }
515  return funcsig;
516 }
517 
518 //______________________________________________________________________________
522 
525  TF1* ftot=ftot=new TF1(fname.Data(),this,&AliHFInvMassFitter::FitFunction4Mass,fMinMass,fMaxMass,totParams,"AliHFInvMassFitter","FitFunction4Mass");
526  for(Int_t ipar=0; ipar<fNParsBkg; ipar++){
527  ftot->SetParameter(ipar,fBkgFunc->GetParameter(ipar));
528  ftot->SetParName(ipar,fBkgFunc->GetParName(ipar));
529  }
530  for(Int_t ipar=0; ipar<fNParsSig; ipar++){
531  ftot->SetParameter(ipar+fNParsBkg,fSigFunc->GetParameter(ipar));
532  ftot->SetParName(ipar+fNParsBkg,fSigFunc->GetParName(ipar));
533  Double_t parmin,parmax;
534  fSigFunc->GetParLimits(ipar,parmin,parmax);
535  ftot->SetParLimits(ipar+fNParsBkg,parmin,parmax);
536  }
537  if(fSecondPeak && fSecFunc){
538  for(Int_t ipar=0; ipar<fNParsSec; ipar++){
539  ftot->SetParameter(ipar+fNParsBkg+fNParsSig,fSecFunc->GetParameter(ipar));
540  ftot->SetParName(ipar+fNParsBkg+fNParsSig,fSecFunc->GetParName(ipar));
541  Double_t parmin,parmax;
542  fSecFunc->GetParLimits(ipar,parmin,parmax);
543  ftot->SetParLimits(ipar+fNParsBkg+fNParsSig,parmin,parmax);
544  }
545  }
546  if(fReflections && fRflFunc){
547  for(Int_t ipar=0; ipar<fNParsRfl; ipar++){
548  ftot->SetParameter(ipar+fNParsBkg+fNParsSig+fNParsSec,fRflFunc->GetParameter(ipar));
549  ftot->SetParName(ipar+fNParsBkg+fNParsSig+fNParsSec,fRflFunc->GetParName(ipar));
550  Double_t parmin,parmax;
551  fRflFunc->GetParLimits(ipar,parmin,parmax);
552  ftot->SetParLimits(ipar+fNParsBkg+fNParsSig+fNParsSec,parmin,parmax);
553  }
554  }
555  return ftot;
556 }
557 //__________________________________________________________________________
561 
563  if(fOnlySideBands && TMath::Abs(x[0]-fMass) < maxDeltaM) {
564  TF1::RejectPoint();
565  return 0;
566  }
567  if(fOnlySideBands && fSecondPeak && TMath::Abs(x[0]-fSecMass) < (fNSigma4SideBands*fSecWidth)){
568  TF1::RejectPoint();
569  return 0;
570  }
571  Double_t total=0;
572 
573  switch (fTypeOfFit4Bkg){
574  case 0:
575  //exponential
576  //exponential = A*exp(B*x) -> integral(exponential)=A/B*exp(B*x)](min,max)
577  //-> A = B*integral/(exp(B*max)-exp(B*min)) where integral can be written
578  //as integralTot- integralGaus (=par [2])
579  //Par:
580  // * [0] = integralBkg;
581  // * [1] = B;
582  //exponential = [1]*[0]/(exp([1]*max)-exp([1]*min))*exp([1]*x)
583  total = par[0]*par[1]/(TMath::Exp(par[1]*fMaxMass)-TMath::Exp(par[1]*fMinMass))*TMath::Exp(par[1]*x[0]);
584  // AliInfo("Background function set to: exponential");
585  break;
586  case 1:
587  //linear
588  //y=a+b*x -> integral = a(max-min)+1/2*b*(max^2-min^2) -> a = (integral-1/2*b*(max^2-min^2))/(max-min)=integral/(max-min)-1/2*b*(max+min)
589  // * [0] = integralBkg;
590  // * [1] = b;
591  total= par[0]/(fMaxMass-fMinMass)+par[1]*(x[0]-0.5*(fMaxMass+fMinMass));
592  // AliInfo("Background function set to: linear");
593  break;
594  case 2:
595  //parabola
596  //y=a+b*x+c*x**2 -> integral = a(max-min) + 1/2*b*(max^2-min^2) +
597  //+ 1/3*c*(max^3-min^3) ->
598  //a = (integral-1/2*b*(max^2-min^2)-1/3*c*(max^3-min^3))/(max-min)
599  // * [0] = integralBkg;
600  // * [1] = b;
601  // * [2] = c;
602  total = par[0]/(fMaxMass-fMinMass)+par[1]*(x[0]-0.5*(fMaxMass+fMinMass))+par[2]*(x[0]*x[0]-1/3.*(fMaxMass*fMaxMass*fMaxMass-fMinMass*fMinMass*fMinMass)/(fMaxMass-fMinMass));
603  // AliInfo("Background function set to: polynomial");
604  break;
605  case 3:
606  total=par[0];
607  break;
608  case 4:
609  //power function
610  //y=a(x-m_pi)^b -> integral = a/(b+1)*((max-m_pi)^(b+1)-(min-m_pi)^(b+1))
611  //
612  //a = integral*(b+1)/((max-m_pi)^(b+1)-(min-m_pi)^(b+1))
613  // * [0] = integralBkg;
614  // * [1] = b;
615  // a(power function) = [0]*([1]+1)/((max-m_pi)^([1]+1)-(min-m_pi)^([1]+1))*(x-m_pi)^[1]
616  {
617  Double_t mpi = TDatabasePDG::Instance()->GetParticle(211)->Mass();
618 
619  total = par[0]*(par[1]+1.)/(TMath::Power(fMaxMass-mpi,par[1]+1.)-TMath::Power(fMinMass-mpi,par[1]+1.))*TMath::Power(x[0]-mpi,par[1]);
620  // AliInfo("Background function set to: powerlaw");
621  }
622  break;
623  case 5:
624  //power function wit exponential
625  //y=a*Sqrt(x-m_pi)*exp(-b*(x-m_pi))
626  {
627  Double_t mpi = TDatabasePDG::Instance()->GetParticle(211)->Mass();
628 
629  total = par[1]*TMath::Sqrt(x[0] - mpi)*TMath::Exp(-1.*par[2]*(x[0]-mpi));
630  // AliInfo("Background function set to: wit exponential");
631  }
632  break;
633  case 6:
634  // the following comment must be removed
635  // // pol 3, following convention for pol 2
636  // //y=a+b*x+c*x**2+d*x**3 -> integral = a(max-min) + 1/2*b*(max^2-min^2) +
637  // //+ 1/3*c*(max^3-min^3) + 1/4 d * (max^4-min^4) ->
638  // //a = (integral-1/2*b*(max^2-min^2)-1/3*c*(max^3-min^3) - 1/4 d * (max^4-min^4) )/(max-min)
639  // // * [0] = integralBkg;
640  // // * [1] = b;
641  // // * [2] = c;
642  // // * [3] = d;
643  {
644  total=par[0];
645  for(Int_t it=1;it<=fPolDegreeBkg;it++){
646  total+=par[it]*TMath::Power(x[0]-fMassParticle,it)/TMath::Factorial(it);
647  }
648  }
649  break;
650  }
651  return total;
652 }
653 //_________________________________________________________________________
657 
658  // AliInfo("Signal function set to: Gaussian");
659  Double_t sigval=0;
660  switch (fTypeOfFit4Sgn){
661  case 0:
662  //gaussian = A/(sigma*sqrt(2*pi))*exp(-(x-mean)^2/2/sigma^2)
663  //Par:
664  // * [0] = integralSgn
665  // * [1] = mean
666  // * [2] = sigma
667  //gaussian = [0]/TMath::Sqrt(2.*TMath::Pi())/[2]*exp[-(x-[1])*(x-[1])/(2*[2]*[2])]
668  sigval=par[0]/TMath::Sqrt(2.*TMath::Pi())/par[2]*TMath::Exp(-(x[0]-par[1])*(x[0]-par[1])/2./par[2]/par[2]);
669  break;
670  case 1:
671  //double gaussian = A/(sigma*sqrt(2*pi))*exp(-(x-mean)^2/2/sigma^2)
672  //Par:
673  // * [0] = integralSgn
674  // * [1] = mean
675  // * [2] = sigma1
676  // * [3] = 2nd gaussian ratio
677  // * [4] = deltaSigma
678  //gaussian = [0]/TMath::Sqrt(2.*TMath::Pi())/[2]*exp[-(x-[1])*(x-[1])/(2*[2]*[2])]
679  Double_t g1=(1.-par[3])/TMath::Sqrt(2.*TMath::Pi())/par[2]*TMath::Exp(-(x[0]-par[1])*(x[0]-par[1])/2./par[2]/par[2]);
680  Double_t g2=par[3]/TMath::Sqrt(2.*TMath::Pi())/par[4]*TMath::Exp(-(x[0]-par[1])*(x[0]-par[1])/2./par[4]/par[4]);
681  sigval=par[0]*(g1+g2);
682  break;
683  }
684  fRawYieldHelp=par[0]/fHistoInvMass->GetBinWidth(1);
685  return sigval;
686 }
687 //_________________________________________________________________________
691  if(!fHistoTemplRfl) return 0;
692 
693  Int_t bin =fHistoTemplRfl->FindBin(x[0]);
694  Double_t value=fHistoTemplRfl->GetBinContent(bin);
695  Int_t binmin=fHistoTemplRfl->FindBin(fMinMass*1.00001);
696  Int_t binmax=fHistoTemplRfl->FindBin(fMaxMass*0.99999);
697  Double_t norm=fHistoTemplRfl->Integral(binmin,binmax)*fHistoTemplRfl->GetBinWidth(bin);
698  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
699  value+=fHistoTemplRfl->GetBinContent(bin-1)+fHistoTemplRfl->GetBinContent(bin+1);
700  value/=3.;
701  }
702  return par[0]*value/norm*fRawYieldHelp*fHistoInvMass->GetBinWidth(1);
703 }
704 //_________________________________________________________________________
708  Double_t bkg=FitFunction4Bkg(x,par);
709  Double_t refl=0;
710  if(fReflections) refl=FitFunction4Refl(x,&par[fNParsBkg]);
711  return bkg+refl;
712 }
713 //_________________________________________________________________________
717 
718  //gaussian = A/(sigma*sqrt(2*pi))*exp(-(x-mean)^2/2/sigma^2)
719  //Par:
720  // * [0] = integralSgn
721  // * [1] = mean
722  // * [2] = sigma
723  Double_t secgaval=par[0]/TMath::Sqrt(2.*TMath::Pi())/par[2]*TMath::Exp(-(x[0]-par[1])*(x[0]-par[1])/2./par[2]/par[2]);
724  return secgaval;
725 }
726 //_________________________________________________________________________
730 
731  Double_t bkg=FitFunction4Bkg(x,par);
732  Double_t sig=FitFunction4Sgn(x,&par[fNParsBkg]);
733  Double_t sec=0.;
734  if(fSecondPeak) sec=FitFunction4SecPeak(x,&par[fNParsBkg+fNParsSig]);
735  Double_t refl=0;
736  if(fReflections) refl=FitFunction4Refl(x,&par[fNParsBkg+fNParsSig+fNParsSec]);
737  return bkg+sig+sec+refl;
738 }
739 
740 //_________________________________________________________________________
741 void AliHFInvMassFitter::Signal(Double_t nOfSigma,Double_t &signal,Double_t &errsignal) const {
744 
745  Double_t min=fMass-nOfSigma*fSigmaSgn;
746  Double_t max=fMass+nOfSigma*fSigmaSgn;
747  Signal(min,max,signal,errsignal);
748  return;
749 }
750 
751 //_________________________________________________________________________
752 void AliHFInvMassFitter::Signal(Double_t min, Double_t max, Double_t &signal,Double_t &errsignal) const {
755 
756  signal=fSigFunc->Integral(min, max)/(Double_t)fHistoInvMass->GetBinWidth(1);
757  errsignal=(fRawYieldErr/fRawYield)*signal;/*assume relative error is the same as for total integral*/
758  return;
759 }
760 
761 //___________________________________________________________________________
762 void AliHFInvMassFitter::Background(Double_t nOfSigma,Double_t &background,Double_t &errbackground) const {
765 
766  Double_t min=fMass-nOfSigma*fSigmaSgn;
767  Double_t max=fMass+nOfSigma*fSigmaSgn;
768  Background(min,max,background,errbackground);
769  return;
770 
771 }
772 //___________________________________________________________________________
773 void AliHFInvMassFitter::Background(Double_t min, Double_t max, Double_t &background,Double_t &errbackground) const {
776 
777  TF1 *funcbkg=0x0;
778  if(fBkgFuncRefit) funcbkg=fBkgFuncRefit;
779  else if(fBkgFunc) funcbkg=fBkgFunc;
780  if(!funcbkg){
781  AliError("Bkg function not found!");
782  return;
783  }
784 
785  Double_t intB=funcbkg->GetParameter(0);
786  Double_t intBerr=funcbkg->GetParError(0);
787 
788  //relative error evaluation: from histo
789 
791  Int_t rightBand=fHistoInvMass->FindBin(fMass+fNSigma4SideBands*fSigmaSgn);
792  intB=fHistoInvMass->Integral(1,leftBand)+fHistoInvMass->Integral(rightBand,fHistoInvMass->GetNbinsX());
793  Double_t sum2=0;
794  for(Int_t i=1;i<=leftBand;i++){
795  sum2+=fHistoInvMass->GetBinError(i)*fHistoInvMass->GetBinError(i);
796  }
797  for(Int_t i=rightBand; i<=fHistoInvMass->GetNbinsX();i++){
798  sum2+=fHistoInvMass->GetBinError(i)*fHistoInvMass->GetBinError(i);
799  }
800 
801  intBerr=TMath::Sqrt(sum2);
802 
803  background=funcbkg->Integral(min,max)/(Double_t)fHistoInvMass->GetBinWidth(1);
804  errbackground=intBerr/intB*background;
805 
806  return;
807 
808 }
809 //__________________________________________________________________________
810 
811 void AliHFInvMassFitter::Significance(Double_t nOfSigma,Double_t &significance,Double_t &errsignificance) const {
814 
815  Double_t min=fMass-nOfSigma*fSigmaSgn;
816  Double_t max=fMass+nOfSigma*fSigmaSgn;
817  Significance(min, max, significance, errsignificance);
818 
819  return;
820 }
821 
822 //__________________________________________________________________________
823 
824 void AliHFInvMassFitter::Significance(Double_t min, Double_t max, Double_t &significance,Double_t &errsignificance) const {
827 
828  Double_t background,errbackground;
829  Background(min,max,background,errbackground);
830 
831  if (fRawYield+background <= 0.){
832  significance=-1;
833  errsignificance=0;
834  return;
835  }
836 
837  AliVertexingHFUtils::ComputeSignificance(fRawYield,fRawYieldErr,background,errbackground,significance,errsignificance);
838 
839  return;
840 }
841 //________________________________________________________________________
845 
846  TH1F *hCp=(TH1F*)fHistoInvMass->Clone("htemp");
847  Double_t estimatecent=0.5*(hCp->GetBinContent(hCp->FindBin(fMass-3.5*fSigmaSgn))+hCp->GetBinContent(hCp->FindBin(fMass+3.5*fSigmaSgn)));// just a first rough estimate
848  Double_t estimateslope=(hCp->GetBinContent(hCp->FindBin(fMass+3.5*fSigmaSgn))-hCp->GetBinContent(hCp->FindBin(fMass-3.5*fSigmaSgn)))/(7*fSigmaSgn);// first rough estimate
849 
851  TF1 *funcbkg,*funcPrev=0x0;
853  funcbkg = new TF1(Form("temp%d",fCurPolDegreeBkg),this,&AliHFInvMassFitter::BackFitFuncPolHelper,fMinMass,fMaxMass,fCurPolDegreeBkg+1,"AliHFInvMassFitter","BackFitFuncPolHelper");
854  if(funcPrev){
855  for(Int_t j=0;j<fCurPolDegreeBkg;j++){// now is +1 degree w.r.t. previous fit funct
856  funcbkg->SetParameter(j,funcPrev->GetParameter(j));
857  }
858  delete funcPrev;
859  }
860  else{
861  funcbkg->SetParameter(0,estimatecent);
862  funcbkg->SetParameter(1,estimateslope);
863  }
864  printf(" ---> Pre-fit of background with pol degree %d ---\n",fCurPolDegreeBkg);
865  hCp->Fit(funcbkg,"REMN","");
866  funcPrev=(TF1*)funcbkg->Clone("ftemp");
867  delete funcbkg;
869  }
870 
871  for(Int_t j=0;j<=fPolDegreeBkg;j++){
872  fback->SetParameter(j,funcPrev->GetParameter(j));
873  fback->SetParError(j,funcPrev->GetParError(j));
874  }
875  printf(" ---> Final background fit with pol degree %d ---\n",fPolDegreeBkg);
876  hCp->Fit(fback,"REMN","");// THIS IS JUST TO SET NOT ONLY THE PARAMETERS BUT ALSO chi2, etc...
877 
878 
879  // The following lines might be useful for debugging
880  // TCanvas *cDebug=new TCanvas();
881  // cDebug->cd();
882  // hCp->Draw();
883  // TString strout=Form("Test%d.root",(Int_t)fhistoInvMass->GetBinContent(fhistoInvMass->FindBin(fMass)));
884  // cDebug->Print(strout.Data());
885  // delete cDebug;
886 
887  delete funcPrev;
888  delete hCp;
889  return kTRUE;
890 
891 }
892 // _______________________________________________________________________
896 
898  if(fOnlySideBands && TMath::Abs(x[0]-fMass) < maxDeltaM) {
899  TF1::RejectPoint();
900  return 0;
901  }
902  Double_t back=par[0];
903  for(Int_t it=1;it<=fCurPolDegreeBkg;it++){
904  back+=par[it]*TMath::Power(x[0]-fMassParticle,it)/TMath::Factorial(it);
905  }
906  return back;
907 }
908 // _______________________________________________________________________
917 
918  if(!h){
919  fReflections=kFALSE;
920  return 0x0;
921  }
922  fHistoTemplRfl=(TH1F*)h->Clone("hTemplRfl");
923  opt.ToLower();
924  fReflections=kTRUE;
925  printf("\n--- Reflection templates from simulation ---\n");
926  if(opt.Contains("templ")){
927  printf(" ---> Reflection contribution using directly the histogram from simulation\n");
928  return fHistoTemplRfl;
929  }
930 
931  TF1 *f=0x0;
932  Bool_t isPoissErr=kTRUE;
933  Double_t xMinForFit=h->GetBinLowEdge(1);
934  Double_t xMaxForFit=h->GetXaxis()->GetBinUpEdge(h->GetNbinsX());
935  if(minRange>=0 && maxRange>=0){
936  xMinForFit=TMath::Max(minRange,h->GetBinLowEdge(1));
937  xMaxForFit=TMath::Min(maxRange,h->GetXaxis()->GetBinUpEdge(h->GetNbinsX()));
938  }
939  if(opt.EqualTo("1gaus") || opt.EqualTo("singlegaus")){
940  printf(" ---> Reflection contribution from single-Gaussian fit to histogram from simulation\n");
941  f=new TF1("mygaus","gaus",xMinForFit,xMaxForFit);
942  f->SetParameter(0,h->GetMaximum());
943  // f->SetParLimits(0,0,100.*h->Integral());
944  f->SetParameter(1,1.865);
945  f->SetParameter(2,0.050);
946  fHistoTemplRfl->Fit(f,"REM","");//,h->GetBinLowEdge(1),h->GetXaxis()->GetBinUpEdge(h->GetNbinsX()));
947  }
948  else if(opt.EqualTo("2gaus") || opt.EqualTo("doublegaus")){
949  printf(" ---> Reflection contribution from double-Gaussian fit to histogram from simulation\n");
950  f=new TF1("my2gaus","[0]*([3]/( TMath::Sqrt(2.*TMath::Pi())*[2])*TMath::Exp(-(x-[1])*(x-[1])/(2.*[2]*[2]))+(1.-[3])/( TMath::Sqrt(2.*TMath::Pi())*[5])*TMath::Exp(-(x-[4])*(x-[4])/(2.*[5]*[5])))",xMinForFit,xMaxForFit);
951  f->SetParameter(0,h->GetMaximum());
952  // f->SetParLimits(0,0,100.*h->Integral());
953  f->SetParLimits(3,0.,1.);
954  f->SetParameter(3,0.5);
955  f->SetParameter(1,1.84);
956  f->SetParameter(2,0.050);
957  f->SetParameter(4,1.88);
958  f->SetParameter(5,0.050);
959  fHistoTemplRfl->Fit(f,"REM","");//,h->GetBinLowEdge(1),h->GetXaxis()->GetBinUpEdge(h->GetNbinsX()));
960  }
961  else if(opt.EqualTo("pol3")){
962  printf(" ---> Reflection contribution from pol3 fit to histogram from simulation\n");
963  f=new TF1("mypol3","pol3",xMinForFit,xMaxForFit);
964  f->SetParameter(0,h->GetMaximum());
965  // f->SetParLimits(0,0,100.*h->Integral());
966  // Hard to initialize the other parameters...
967  fHistoTemplRfl->Fit(f,"REM","");
968  // Printf("We USED %d POINTS in the Fit",f->GetNumberFitPoints());
969  }
970  else if(opt.EqualTo("pol6")){
971  printf(" ---> Reflection contribution from pol6 fit to histogram from simulation\n");
972  f=new TF1("mypol6","pol6",xMinForFit,xMaxForFit);
973  f->SetParameter(0,h->GetMaximum());
974  // f->SetParLimits(0,0,100.*h->Integral());
975  // Hard to initialize the other parameters...
976  fHistoTemplRfl->Fit(f,"RLEMI","");//,h->GetBinLowEdge(1),h->GetXaxis()->GetBinUpEdge(h->GetNbinsX()));
977  }
978  else{
979  // no good option passed
980  printf(" ---> Bad option for reflection configuration -> reflections will not be included in the fit\n");
981  fReflections=kFALSE;
982  delete fHistoTemplRfl;
983  fHistoTemplRfl=0x0;
984  return 0x0;
985  }
986 
987  // Fill fHistoTemplRfl with values of fit function
988  if(f){
989  for(Int_t j=1;j<=fHistoTemplRfl->GetNbinsX();j++){
990  fHistoTemplRfl->SetBinContent(j,f->Integral(fHistoTemplRfl->GetBinLowEdge(j),fHistoTemplRfl->GetXaxis()->GetBinUpEdge(j))/fHistoTemplRfl->GetBinWidth(j));
991  if(fHistoTemplRfl->GetBinContent(j)>=0.&&TMath::Abs(h->GetBinError(j)*h->GetBinError(j)-h->GetBinContent(j))>0.1*h->GetBinContent(j))isPoissErr=kFALSE;
992  }
993  for(Int_t j=1;j<=fHistoTemplRfl->GetNbinsX();j++){
994  if(isPoissErr){
995  fHistoTemplRfl->SetBinError(j,TMath::Sqrt(fHistoTemplRfl->GetBinContent(j)));
996  }
997  else fHistoTemplRfl->SetBinError(j,0.001*fHistoTemplRfl->GetBinContent(j));
998  }
999  fReflections=kTRUE;
1000  return fHistoTemplRfl;
1001  }else{
1002  printf(" ---> Fit to MC template for reflection failed -> reflections will not be included in the fit\n");
1003  fReflections=kFALSE;
1004  delete fHistoTemplRfl;
1005  fHistoTemplRfl=0x0;
1006  return 0x0;
1007  }
1008  return 0x0;
1009 }
1010 // _______________________________________________________________________
1014  if(fBkgFunc){
1015  printf("--- Background function in 1st step fit to side bands ---\n");
1016  fBkgFunc->Print();
1017  }
1018  if(fTotFunc){
1019  printf("--- Total fit function from 2nd step fit ---\n");
1020  fTotFunc->Print();
1021  }
1022  if(fBkgFuncRefit){
1023  printf("--- Background function in 2nd step fit ---\n");
1024  fBkgFuncRefit->Print();
1025  }
1026  if(fRflFunc){
1027  printf("--- Reflections ---\n");
1028  fRflFunc->Print();
1029  }
1030  if(fBkRFunc){
1031  printf("--- Background + reflections ---\n");
1032  fBkRFunc->Print();
1033  }
1034  if(fSecFunc){
1035  printf("--- Additional Gaussian peak ---\n");
1036  fSecFunc->Print();
1037  }
1038 }
Int_t fTypeOfFit4Sgn
pdg value of particle mass
Double_t fSigmaSgnErr
signal gaussian sigma
Bool_t fFixSecWidth
flag to fix the position of the 2nd peak
Double_t fRawYieldHelp
switch for smoothing of reflection template
double Double_t
Definition: External.C:58
Double_t fRawYield
L, LW or Chi2.
TF1 * CreateReflectionFunction(TString fname)
Int_t fCurPolDegreeBkg
degree of polynomial expansion for back fit (option 6 for back)
Bool_t PrepareHighPolFit(TF1 *fback)
Bool_t fOnlySideBands
fit parameters in background fit function
virtual void Signal(Double_t nOfSigma, Double_t &signal, Double_t &errsignal) const
Double_t FitFunction4SecPeak(Double_t *x, Double_t *par)
Double_t fMass
signal fit func
TF1 * CreateTotalFitFunction(TString fname)
Double_t fRawYieldErr
signal gaussian integral
Int_t fNParsBkg
fit parameters in signal fit function
AliHFInvMassFitter class for the fit of invariant mass distribution of charm hadrons.
Int_t fNParsRfl
flag use/not use reflections
TCanvas * c
Definition: TestFitELoss.C:172
TF1 * CreateSecondPeakFunction(TString fname, Double_t integral)
TF1 * fTotFunc
fit function for second peak
Double_t fFixedRawYield
switch for fix Sigma of gaussian
Double_t FitFunction4Bkg(Double_t *x, Double_t *par)
Double_t fSigmaSgn
unc on signal gaussian mean value
Double_t fSecMass
fit parameters in 2nd peak fit function
Bool_t fFixSecMass
width of the 2nd peak
TF1 * fSecFunc
flag to fix the width of the 2nd peak
Double_t * sigma
int Int_t
Definition: External.C:63
Double_t fSecWidth
position of the 2nd peak
TString fFitOption
number of sigmas to veto the signal peak
TF1 * fBkRFunc
fit function for reflections
Int_t fNParsSec
switch off/on second peak (for D+->KKpi in Ds)
TF1 * CreateSignalFitFunction(TString fname, Double_t integral)
Bool_t fFixedMean
unc on signal gaussian sigma
Double_t fMassParticle
help variable
Double_t FitFunction4Refl(Double_t *x, Double_t *par)
TH1F * SetTemplateReflections(const TH1 *h, TString opt, Double_t minRange, Double_t maxRange)
Double_t nsigma
Int_t fTypeOfFit4Bkg
upper mass limit
TH1F * fHistoTemplRfl
switch for fix refl/signal
TF1 * fBkgFunc
background fit function (1st step, side bands only)
Double_t fMassErr
signal gaussian mean value
Double_t fRflOverSig
fit parameters in reflection fit function
Double_t fMinMass
histogram to fit
Bool_t fFixRflOverSig
reflection/signal
static void ComputeSignificance(Double_t signal, Double_t errsignal, Double_t background, Double_t errbackground, Double_t &significance, Double_t &errsignificance)
Significance calculator.
Double_t FitFunction4Mass(Double_t *x, Double_t *par)
void Significance(Double_t nOfSigma, Double_t &significance, Double_t &errsignificance) const
Double_t BackFitFuncPolHelper(Double_t *x, Double_t *par)
Bool_t MassFitter(Bool_t draw=kTRUE)
Bool_t draw[nPtBins]
Double_t fMaxMass
lower mass limit
TF1 * CreateBackgroundFitFunction(TString fname, Double_t integral)
TF1 * fRflFunc
internal variable for fit with reflections
ClassImp(AliAnalysisTaskCRC) AliAnalysisTaskCRC
Double_t fNSigma4SideBands
kTRUE = only side bands considered
void Background(Double_t nOfSigma, Double_t &background, Double_t &errbackground) const
Double_t FitFunction4Sgn(Double_t *x, Double_t *par)
Int_t fPolDegreeBkg
background fit func
bool Bool_t
Definition: External.C:53
TF1 * fBkgFuncRefit
background fit function (1st step, extended in peak region)
Double_t CheckForSignal(Double_t mean, Double_t sigma)
TF1 * fSigFunc
err on signal gaussian integral
TF1 * fBkgFuncSb
Signal fit function.
Bool_t fSecondPeak
fit function for reflections
void DrawHere(TVirtualPad *c)
Bool_t fSmoothRfl
histogram with reflection template
Definition: External.C:196
TF1 * CreateBackgroundPlusReflectionFunction(TString fname)
Double_t FitFunction4BkgAndRefl(Double_t *x, Double_t *par)
Int_t fNParsSig
initialization for wa yield
Bool_t fReflections
background fit function (2nd step)
Bool_t fFixedSigma
switch for fix mean of gaussian