AliPhysics  vAN-20150827 (3e81cbb)
 All Classes Namespaces Files Functions Variables Enumerations Enumerator Macros
AliHFMultiTrials.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 2008-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 <TMath.h>
17 #include <TPad.h>
18 #include <TCanvas.h>
19 #include <TH1F.h>
20 #include <TH2F.h>
21 #include <TNtuple.h>
22 #include <TF1.h>
23 #include <TLatex.h>
24 #include <TFile.h>
25 #include "AliHFMassFitter.h"
26 #include "AliHFMassFitterVAR.h"
27 #include "AliHFMultiTrials.h"
28 
32 
33 
34 //_________________________________________________________________________
36  TNamed(),
37  fNumOfRebinSteps(4),
38  fRebinSteps(0x0),
39  fNumOfLowLimFitSteps(6),
40  fLowLimFitSteps(0x0),
41  fNumOfUpLimFitSteps(6),
42  fUpLimFitSteps(0x0),
43  fNumOfnSigmaBinCSteps(11),
44  fnSigmaBinCSteps(0x0),
45  fSigmaGausMC(0.010),
46  fSigmaMCVariation(0.15),
47  fMassD(1.86484),
48  fSuffix(""),
49  fFitOption(0),
50  fUseExpoBkg(kTRUE),
51  fUseLinBkg(kTRUE),
52  fUsePol2Bkg(kTRUE),
53  fUsePol3Bkg(kTRUE),
54  fUsePol4Bkg(kTRUE),
55  fUsePol5Bkg(kFALSE),
56  fDrawIndividualFits(kFALSE),
57  fHistoRawYieldDistAll(0x0),
58  fHistoRawYieldTrialAll(0x0),
59  fHistoSigmaTrialAll(0x0),
60  fHistoMeanTrialAll(0x0),
61  fHistoChi2TrialAll(0x0),
62  fHistoRawYieldDistBinCAll(0x0),
63  fHistoRawYieldTrialBinCAll(0x0),
64  fHistoRawYieldDist(0x0),
65  fHistoRawYieldTrial(0x0),
66  fHistoSigmaTrial(0x0),
67  fHistoMeanTrial(0x0),
68  fHistoChi2Trial(0x0),
69  fHistoRawYieldDistBinC(0x0),
70  fHistoRawYieldTrialBinC(0x0),
71  fNtupleMultiTrials(0x0),
72  fMinYieldGlob(0),
73  fMaxYieldGlob(0)
74 {
75  // constructor
76  Int_t rebinStep[4]={3,4,5,6};
77  Double_t minMassStep[6]={1.68,1.70,1.72,1.74,1.76,1.78};
78  Double_t maxMassStep[6]={2.06,2.04,2.02,2.00,1.98,1.96};
79  Double_t nSigmasBC[11]={2.0,2.5,3.0,3.5,4.0,4.5,5.0,5.5,6.0,6.5,7.0};
80  ConfigureRebinSteps(4,rebinStep);
81  ConfigureLowLimFitSteps(6,minMassStep);
82  ConfigureUpLimFitSteps(6,maxMassStep);
83  ConfigurenSigmaBinCSteps(11,nSigmasBC);
84 }
85 
86 //________________________________________________________________________
88  // destructor
89  delete [] fRebinSteps;
90  delete [] fLowLimFitSteps;
91  delete [] fUpLimFitSteps;
92 }
93 
94 //________________________________________________________________________
96  // creates output histograms
97 
98  const Int_t nCases=kNBkgFuncCases*kNFitConfCases;
99 
100  TString funcBkg[kNBkgFuncCases]={"Expo","Lin","Pol2","Pol3","Pol4","Pol5"};
101  TString gausSig[kNFitConfCases]={"FixedS","FixedSp20","FixedSm20","FreeS","FixedMeanFixedS","FixedMeanFreeS"};
102 
104 
105  fHistoRawYieldDistAll = new TH1F(Form("hRawYieldDistAll%s",fSuffix.Data())," ; Raw Yield",5000,0.,50000.);
106  fHistoRawYieldTrialAll = new TH1F(Form("hRawYieldTrialAll%s",fSuffix.Data())," ; Trial # ; Raw Yield",nCases*totTrials,-0.5,nCases*totTrials-0.5);
107  fHistoSigmaTrialAll = new TH1F(Form("hSigmaTrialAll%s",fSuffix.Data())," ; Trial # ; Sigma (GeV/c^{2})",nCases*totTrials,-0.5,nCases*totTrials-0.5);
108  fHistoMeanTrialAll = new TH1F(Form("hMeanTrialAll%s",fSuffix.Data())," ; Trial # ; Mean (GeV/c^{2})",nCases*totTrials,-0.5,nCases*totTrials-0.5);
109  fHistoChi2TrialAll = new TH1F(Form("hChi2TrialAll%s",fSuffix.Data())," ; Trial # ; #chi^{2}",nCases*totTrials,-0.5,nCases*totTrials-0.5);
110 
111  fHistoRawYieldDistBinCAll = new TH1F(Form("hRawYieldDistBinCAll%s",fSuffix.Data())," ; Raw Yield (bin count)",5000,0.,50000.);
112  fHistoRawYieldTrialBinCAll = new TH2F(Form("hRawYieldTrialBinCAll%s",fSuffix.Data())," ; Trial # ; Range for count ; Raw Yield (bin count)",totTrials,-0.5,totTrials-0.5,fNumOfnSigmaBinCSteps,-0.5,fNumOfnSigmaBinCSteps-0.5);
113 
114  fHistoRawYieldDist = new TH1F*[nCases];
115  fHistoRawYieldTrial = new TH1F*[nCases];
116  fHistoSigmaTrial = new TH1F*[nCases];
117  fHistoMeanTrial = new TH1F*[nCases];
118  fHistoChi2Trial = new TH1F*[nCases];
119 
120  fHistoRawYieldDistBinC = new TH1F*[nCases];
121  fHistoRawYieldTrialBinC = new TH2F*[nCases];
122 
123  for(Int_t ib=0; ib<kNBkgFuncCases; ib++){
124  for(Int_t igs=0; igs<kNFitConfCases; igs++){
125  Int_t theCase=igs*kNBkgFuncCases+ib;
126  fHistoRawYieldDist[theCase]=new TH1F(Form("hRawYieldDist%s%s%s",funcBkg[ib].Data(),gausSig[igs].Data(),fSuffix.Data())," ; Raw Yield",5000,0.,50000.);
127  fHistoRawYieldDistBinC[theCase]=new TH1F(Form("hRawYieldDistBinC%s%s%s",funcBkg[ib].Data(),gausSig[igs].Data(),fSuffix.Data())," ; Raw Yield (bin count)",5000,0.,50000.);
128  fHistoRawYieldTrial[theCase]=new TH1F(Form("hRawYieldTrial%s%s%s",funcBkg[ib].Data(),gausSig[igs].Data(),fSuffix.Data())," ; Trial # ; Raw Yield",totTrials,-0.5,totTrials-0.5);
129  fHistoRawYieldTrialBinC[theCase]=new TH2F(Form("hRawYieldTrialBinC%s%s%s",funcBkg[ib].Data(),gausSig[igs].Data(),fSuffix.Data())," ; Trial # ; Range for count ; Raw Yield (bin count)",totTrials,-0.5,totTrials-0.5,fNumOfnSigmaBinCSteps,-0.5,fNumOfnSigmaBinCSteps-0.5);
130  fHistoSigmaTrial[theCase]=new TH1F(Form("hSigmaTrial%s%s%s",funcBkg[ib].Data(),gausSig[igs].Data(),fSuffix.Data())," ; Trial # ; Sigma (GeV/c^{2})",totTrials,-0.5,totTrials-0.5);
131  fHistoMeanTrial[theCase]=new TH1F(Form("hMeanTrial%s%s%s",funcBkg[ib].Data(),gausSig[igs].Data(),fSuffix.Data())," ; Trial # ; Mean (GeV/c^{2})",totTrials,-0.5,totTrials-0.5);
132  fHistoChi2Trial[theCase]=new TH1F(Form("hChi2Trial%s%s%s",funcBkg[ib].Data(),gausSig[igs].Data(),fSuffix.Data())," ; Trial # ; #chi^{2}",totTrials,-0.5,totTrials-0.5);
133 
134  fHistoChi2Trial[theCase]->SetMarkerStyle(7);
135  fHistoSigmaTrial[theCase]->SetMarkerStyle(7);
136  fHistoMeanTrial[theCase]->SetMarkerStyle(7);
137 
138  }
139  }
140  fNtupleMultiTrials = new TNtuple(Form("ntuMultiTrial%s",fSuffix.Data()),Form("ntuMultiTrial%s",fSuffix.Data()),"rebin:minfit:maxfit:bkgfunc:confsig:confmean:chi2:mean:emean:sigma:esigma:rawy:erawy");
141  return kTRUE;
142 
143 }
144 
145 //________________________________________________________________________
146 Bool_t AliHFMultiTrials::DoMultiTrials(TH1D* hInvMassHisto, TPad* thePad){
147  // perform the multiple fits
148 
149  Bool_t hOK=CreateHistos();
150  if(!hOK) return kFALSE;
151 
152  Int_t itrial=0;
153  Int_t types=0;
154  Int_t itrialBC=0;
156 
157  fMinYieldGlob=999999.;
158  fMaxYieldGlob=0.;
159  Float_t xnt[13];
160 
161  for(Int_t ir=0; ir<fNumOfRebinSteps; ir++){
162  Int_t rebin=fRebinSteps[ir];
163  TH1F* hRebinned=RebinHisto(hInvMassHisto,rebin,-1);
164  for(Int_t iMinMass=0; iMinMass<fNumOfLowLimFitSteps; iMinMass++){
165  Double_t minMassForFit=fLowLimFitSteps[iMinMass];
166  Double_t hmin=TMath::Max(minMassForFit,hRebinned->GetBinLowEdge(2));
167  for(Int_t iMaxMass=0; iMaxMass<fNumOfUpLimFitSteps; iMaxMass++){
168  Double_t maxMassForFit=fUpLimFitSteps[iMaxMass];
169  Double_t hmax=TMath::Min(maxMassForFit,hRebinned->GetBinLowEdge(hRebinned->GetNbinsX()));
170  ++itrial;
171  for(Int_t typeb=0; typeb<kNBkgFuncCases; typeb++){
172  if(typeb==kExpoBkg && !fUseExpoBkg) continue;
173  if(typeb==kLinBkg && !fUseLinBkg) continue;
174  if(typeb==kPol2Bkg && !fUsePol2Bkg) continue;
175  if(typeb==kPol3Bkg && !fUsePol3Bkg) continue;
176  if(typeb==kPol4Bkg && !fUsePol4Bkg) continue;
177  if(typeb==kPol5Bkg && !fUsePol5Bkg) continue;
178  for(Int_t igs=0; igs<kNFitConfCases; igs++){
179  Int_t theCase=igs*kNBkgFuncCases+typeb;
180  Int_t globBin=itrial+theCase*totTrials;
181  for(Int_t j=0; j<13; j++) xnt[j]=0.;
182 
183  AliHFMassFitterVAR* fitter=0x0;
184  if(typeb<=kPol2Bkg){
185  fitter=new AliHFMassFitterVAR(hRebinned,hmin, hmax,1,typeb,types);
186  }else{
187  fitter=new AliHFMassFitterVAR(hRebinned,hmin, hmax,1,6,types);
188  if(typeb==kPol3Bkg) fitter->SetBackHighPolDegree(3);
189  if(typeb==kPol4Bkg) fitter->SetBackHighPolDegree(4);
190  if(typeb==kPol5Bkg) fitter->SetBackHighPolDegree(5);
191  }
192  fitter->SetReflectionSigmaFactor(0);
193  if(fFitOption==1) fitter->SetUseChi2Fit();
196  xnt[0]=rebin;
197  xnt[1]=minMassForFit;
198  xnt[2]=maxMassForFit;
199  xnt[3]=typeb;
200  xnt[5]=0;
201  if(igs==kFixSigFreeMean){
202  fitter->SetFixGaussianSigma(fSigmaGausMC,kTRUE);
203  xnt[4]=1;
204  }else if(igs==kFixSigUpFreeMean){
206  xnt[4]=2;
207  }else if(igs==kFixSigDownFreeMean){
209  xnt[4]=3;
210  }else if(igs==kFreeSigFreeMean){
211  xnt[4]=0;
212  }else if(igs==kFixSigFixMean){
213  fitter->SetFixGaussianSigma(fSigmaGausMC,kTRUE);
214  fitter->SetFixGaussianMean(fMassD,kTRUE);
215  xnt[4]=1;
216  xnt[5]=1;
217  }else if(igs==kFreeSigFixMean){
218  fitter->SetFixGaussianMean(fMassD,kTRUE);
219  xnt[4]=0;
220  xnt[5]=1;
221  }
222  Bool_t out=kFALSE;
223  Double_t chisq=-1.;
224  Double_t sigma=0.;
225  Double_t esigma=0.;
226  Double_t pos=.0;
227  Double_t epos=.0;
228  Double_t ry=.0;
229  Double_t ery=.0;
230  TF1* fB1=0x0;
231  if(typeb<kNBkgFuncCases){
232  printf("****** START FIT OF HISTO %s WITH REBIN %d MASS RANGE %f-%f BACKGROUND FIT FUNCTION=%d CONFIG SIGMA/MEAN=%d\n",hInvMassHisto->GetName(),rebin,minMassForFit,maxMassForFit,typeb,igs);
233  out=fitter->MassFitter(0);
234  chisq=fitter->GetReducedChiSquare();
235  sigma=fitter->GetSigma();
236  pos=fitter->GetMean();
237  esigma=fitter->GetSigmaUncertainty();
238  if(esigma<0.00001) esigma=0.0001;
239  epos=fitter->GetMeanUncertainty();
240  if(epos<0.00001) epos=0.0001;
241  ry=fitter->GetRawYield();
242  ery=fitter->GetRawYieldError();
243  fB1=fitter->GetBackgroundFullRangeFunc();
244  if(fDrawIndividualFits && thePad){
245  thePad->Clear();
246  fitter->DrawHere(thePad);
247  }
248  }
249  // else{
250  // out=DoFitWithPol3Bkg(hRebinned,hmin,hmax,igs);
251  // if(out && thePad){
252  // thePad->Clear();
253  // hRebinned->Draw();
254  // TF1* fSB=(TF1*)hRebinned->GetListOfFunctions()->FindObject("fSB");
255  // fB1=new TF1("fB1","[0]+[1]*x+[2]*x*x+[3]*x*x*x",hmin,hmax);
256  // for(Int_t j=0; j<4; j++) fB1->SetParameter(j,fSB->GetParameter(3+j));
257  // fB1->SetLineColor(2);
258  // fB1->Draw("same");
259  // fSB->SetLineColor(4);
260  // fSB->Draw("same");
261  // thePad->Update();
262  // chisq=fSB->GetChisquare()/fSB->GetNDF();;
263  // sigma=fSB->GetParameter(2);
264  // esigma=fSB->GetParError(2);
265  // if(esigma<0.00001) esigma=0.0001;
266  // pos=fSB->GetParameter(1);
267  // epos=fSB->GetParError(1);
268  // if(epos<0.00001) epos=0.0001;
269  // ry=fSB->GetParameter(0)/hRebinned->GetBinWidth(1);
270  // ery=fSB->GetParError(0)/hRebinned->GetBinWidth(1);
271  // }
272  // }
273  xnt[6]=chisq;
274  if(out && chisq>0. && sigma>0.5*fSigmaGausMC && sigma<2.0*fSigmaGausMC){
275  xnt[7]=pos;
276  xnt[8]=epos;
277  xnt[9]=sigma;
278  xnt[10]=esigma;
279  xnt[11]=ry;
280  xnt[12]=ery;
281  fHistoRawYieldDistAll->Fill(ry);
282  fHistoRawYieldTrialAll->SetBinContent(globBin,ry);
283  fHistoRawYieldTrialAll->SetBinError(globBin,ery);
284  fHistoSigmaTrialAll->SetBinContent(globBin,sigma);
285  fHistoSigmaTrialAll->SetBinError(globBin,esigma);
286  fHistoMeanTrialAll->SetBinContent(globBin,pos);
287  fHistoMeanTrialAll->SetBinError(globBin,epos);
288  fHistoChi2TrialAll->SetBinContent(globBin,chisq);
289  fHistoChi2TrialAll->SetBinError(globBin,0.00001);
290  if(ry<fMinYieldGlob) fMinYieldGlob=ry;
291  if(ry>fMaxYieldGlob) fMaxYieldGlob=ry;
292  fHistoRawYieldDist[theCase]->Fill(ry);
293  fHistoRawYieldTrial[theCase]->SetBinContent(itrial,ry);
294  fHistoRawYieldTrial[theCase]->SetBinError(itrial,ery);
295  fHistoSigmaTrial[theCase]->SetBinContent(itrial,sigma);
296  fHistoSigmaTrial[theCase]->SetBinError(itrial,esigma);
297  fHistoMeanTrial[theCase]->SetBinContent(itrial,pos);
298  fHistoMeanTrial[theCase]->SetBinError(itrial,epos);
299  fHistoChi2Trial[theCase]->SetBinContent(itrial,chisq);
300  fHistoChi2Trial[theCase]->SetBinError(itrial,0.00001);
301  for(Int_t iStepBC=0; iStepBC<fNumOfnSigmaBinCSteps; iStepBC++){
302  Double_t minMassBC=fMassD-fnSigmaBinCSteps[iStepBC]*sigma;
303  Double_t maxMassBC=fMassD+fnSigmaBinCSteps[iStepBC]*sigma;
304  if(minMassBC>minMassForFit &&
305  maxMassBC<maxMassForFit &&
306  minMassBC>(hRebinned->GetXaxis()->GetXmin()) &&
307  maxMassBC<(hRebinned->GetXaxis()->GetXmax())){
308  Double_t cnts,ecnts;
309  BinCount(hRebinned,fB1,1,minMassBC,maxMassBC,cnts,ecnts);
310  ++itrialBC;
311  fHistoRawYieldDistBinCAll->Fill(cnts);
312  fHistoRawYieldTrialBinCAll->SetBinContent(globBin,iStepBC+1,cnts);
313  fHistoRawYieldTrialBinCAll->SetBinError(globBin,iStepBC+1,ecnts);
314  fHistoRawYieldTrialBinC[theCase]->SetBinContent(itrial,iStepBC+1,cnts);
315  fHistoRawYieldTrialBinC[theCase]->SetBinError(itrial,iStepBC+1,ecnts);
316  fHistoRawYieldDistBinC[theCase]->Fill(cnts);
317  }
318  }
319  }
320  delete fitter;
321  // if(typeb>4) delete fB1;
322  fNtupleMultiTrials->Fill(xnt);
323  }
324  }
325  }
326  }
327  delete hRebinned;
328  }
329  return kTRUE;
330 }
331 
332 //________________________________________________________________________
333 void AliHFMultiTrials::SaveToRoot(TString fileName, TString option) const{
334  // save histos in a root file for further analysis
335  const Int_t nCases=kNBkgFuncCases*kNFitConfCases;
336  TFile* outHistos=new TFile(fileName.Data(),option.Data());
337  fHistoRawYieldTrialAll->Write();
338  fHistoSigmaTrialAll->Write();
339  fHistoMeanTrialAll->Write();
340  fHistoChi2TrialAll->Write();
341  fHistoRawYieldDistBinCAll->Write();
342  fHistoRawYieldTrialBinCAll->Write();
343  for(Int_t ic=0; ic<nCases; ic++){
344  fHistoRawYieldTrial[ic]->Write();
345  fHistoSigmaTrial[ic]->Write();
346  fHistoMeanTrial[ic]->Write();
347  fHistoChi2Trial[ic]->Write();
348  fHistoRawYieldTrialBinC[ic]->Write();
349  fHistoRawYieldDistBinC[ic]->Write();
350  }
351  fNtupleMultiTrials->Write();
352  outHistos->Close();
353  delete outHistos;
354 }
355 
356 //________________________________________________________________________
357 void AliHFMultiTrials::DrawHistos(TCanvas* cry) const{
358  // draw histos
359  cry->Divide(2,2);
360  cry->cd(1);
361  fHistoSigmaTrialAll->Draw();
362  cry->cd(3);
363  fHistoChi2TrialAll->Draw();
364  cry->cd(2);
365  fHistoRawYieldTrialAll->Draw();
366  cry->cd(4);
367  fHistoRawYieldDistAll->Draw();
368  TLatex* tmean=new TLatex(0.15,0.8,Form("mean=%.1f",fHistoRawYieldDistAll->GetMean()));
369  tmean->SetNDC();
370  tmean->Draw();
371  TLatex* thrms=new TLatex(0.15,0.72,Form("rms=%.1f",fHistoRawYieldDistAll->GetRMS()));
372  thrms->SetNDC();
373  thrms->Draw();
374  TLatex* tmax=new TLatex(0.6,0.8,Form("max=%.1f",fMaxYieldGlob));
375  tmax->SetNDC();
376  tmax->Draw();
377  TLatex* tmin=new TLatex(0.6,0.72,Form("min=%.1f",fMinYieldGlob));
378  tmin->SetNDC();
379  tmin->Draw();
380  TLatex* trms=new TLatex(0.6,0.64,Form("(max-min)/sqrt(12)=%.1f",(fMaxYieldGlob-fMinYieldGlob)/sqrt(12)));
381  trms->SetNDC();
382  trms->Draw();
383 
384 }
385 //________________________________________________________________________
386 TH1F* AliHFMultiTrials::RebinHisto(TH1D* hOrig, Int_t reb, Int_t firstUse) const{
387  // Rebin histogram, from bin firstUse to lastUse
388  // Use all bins if firstUse=-1
389 
390  Int_t nBinOrig=hOrig->GetNbinsX();
391  Int_t firstBinOrig=1;
392  Int_t lastBinOrig=nBinOrig;
393  Int_t nBinOrigUsed=nBinOrig;
394  Int_t nBinFinal=nBinOrig/reb;
395  if(firstUse>=1){
396  firstBinOrig=firstUse;
397  nBinFinal=(nBinOrig-firstUse+1)/reb;
398  nBinOrigUsed=nBinFinal*reb;
399  lastBinOrig=firstBinOrig+nBinOrigUsed-1;
400  }else{
401  Int_t exc=nBinOrigUsed%reb;
402  if(exc!=0){
403  nBinOrigUsed-=exc;
404  firstBinOrig+=exc/2;
405  lastBinOrig=firstBinOrig+nBinOrigUsed-1;
406  }
407  }
408 
409  printf("Rebin from %d bins to %d bins -- Used bins=%d in range %d-%d\n",nBinOrig,nBinFinal,nBinOrigUsed,firstBinOrig,lastBinOrig);
410  Float_t lowLim=hOrig->GetXaxis()->GetBinLowEdge(firstBinOrig);
411  Float_t hiLim=hOrig->GetXaxis()->GetBinUpEdge(lastBinOrig);
412  TH1F* hRebin=new TH1F(Form("%s-rebin%d_%d",hOrig->GetName(),reb,firstUse),hOrig->GetTitle(),nBinFinal,lowLim,hiLim);
413  Int_t lastSummed=firstBinOrig-1;
414  for(Int_t iBin=1;iBin<=nBinFinal; iBin++){
415  Float_t sum=0.;
416  Float_t sum2=0.;
417  for(Int_t iOrigBin=0;iOrigBin<reb;iOrigBin++){
418  sum+=hOrig->GetBinContent(lastSummed+1);
419  sum2+=hOrig->GetBinError(lastSummed+1)*hOrig->GetBinError(lastSummed+1);
420  lastSummed++;
421  }
422  hRebin->SetBinContent(iBin,sum);
423  hRebin->SetBinError(iBin,TMath::Sqrt(sum2));
424  }
425  return hRebin;
426 }
427 
428 //________________________________________________________________________
429 void AliHFMultiTrials::BinCount(TH1F* h, TF1* fB, Int_t rebin, Double_t minMass, Double_t maxMass, Double_t& count, Double_t& ecount) const{
430  // compute yield with bin couting
431  Int_t minBinSum=h->FindBin(minMass);
432  Int_t maxBinSum=h->FindBin(maxMass);
433  Double_t cntSig=0.;
434  Double_t cntErr=0.;
435  for(Int_t iMB=minBinSum; iMB<=maxBinSum; iMB++){
436  Double_t bkg=fB ? fB->Eval(h->GetBinCenter(iMB))/(Double_t)rebin : 0;
437  cntSig+=(h->GetBinContent(iMB)-bkg);
438  cntErr+=(h->GetBinError(iMB)*h->GetBinError(iMB));
439  }
440  count=cntSig;
441  ecount=TMath::Sqrt(cntErr);
442 }
443 //________________________________________________________________________
444 Bool_t AliHFMultiTrials::DoFitWithPol3Bkg(TH1F* histoToFit, Double_t hmin, Double_t hmax,
445  Int_t iCase){
446  //
447 
448  TH1F *hCutTmp=(TH1F*)histoToFit->Clone("hCutTmp");
449  for(Int_t ib=1; ib<=hCutTmp->GetNbinsX(); ib++){
450  Double_t xc=hCutTmp->GetBinCenter(ib);
451  if(xc>(fMassD-5.*fSigmaGausMC) && xc<(fMassD+5.*fSigmaGausMC)){
452  hCutTmp->SetBinContent(ib,0.);
453  hCutTmp->SetBinError(ib,0.);
454  }
455  }
456 
457  hCutTmp->Fit("pol2","E0","",hmin,hmax);
458  TF1* f2=(TF1*)hCutTmp->GetListOfFunctions()->FindObject("pol2");
459  TF1* f3=new TF1("myPol3","pol3");
460  for(Int_t i=0; i<3;i++) f3->SetParameter(i,f2->GetParameter(i));
461  hCutTmp->Fit(f3,"E0","",hmin,hmax);
462  Double_t quickCount=0.;
463  for(Int_t ib=1; ib<=histoToFit->GetNbinsX(); ib++){
464  Double_t xc=hCutTmp->GetBinCenter(ib);
465  if(xc>(fMassD-3.*fSigmaGausMC) && xc<(fMassD+3.*fSigmaGausMC)){
466  quickCount+=(histoToFit->GetBinContent(ib)-f3->Eval(xc));
467  }
468  }
469  TF1* fSB=new TF1("fSB","[0]*1./(TMath::Sqrt(2.*TMath::Pi())*[2])*TMath::Exp(-(x-[1])*(x-[1])/(2.*[2]*[2]))+[3]+[4]*x+[5]*x*x+[6]*x*x*x");
470  fSB->SetParameter(0,quickCount);
471  fSB->SetParameter(1,fMassD);
472  fSB->SetParameter(2,fSigmaGausMC);
473  for(Int_t j=0; j<4; j++) fSB->SetParameter(j+3,f3->GetParameter(j));
474  if(iCase==0) fSB->FixParameter(2,fSigmaGausMC);
475  else if(iCase==1) fSB->FixParameter(2,fSigmaGausMC*(1.+fSigmaMCVariation));
476  else if(iCase==2) fSB->FixParameter(2,fSigmaGausMC*(1.-fSigmaMCVariation));
477  else if(iCase==4){
478  fSB->FixParameter(1,fMassD);
479  fSB->FixParameter(2,fSigmaGausMC);
480  } else if(iCase==5){
481  fSB->FixParameter(1,fMassD);
482  }
483  histoToFit->Fit(fSB,"ME0","",hmin,hmax);
484  // quality cuts
485  if(fSB->GetParError(0)<0.01*fSB->GetParameter(0)) return kFALSE;
486  if(fSB->GetParError(0)>0.6*fSB->GetParameter(0)) return kFALSE;
487 
488  delete hCutTmp;
489  return kTRUE;
490 }
ClassImp(AliAnalysisTaskTriggerRates) AliAnalysisTaskTriggerRates
Double_t fMaxYieldGlob
minimum yield
TH1F * fHistoSigmaTrialAll
histo with yield from all trials
TNtuple * fNtupleMultiTrials
histo with bin counts from subsamples of trials
Double_t fMinYieldGlob
tree
TH1F * fHistoRawYieldDistAll
flag for drawing fits
Bool_t fUseLinBkg
switch for exponential background
TString fileName
void SetInitialGaussianMean(Double_t mean)
Double_t * fLowLimFitSteps
number of steps on the min. mass for fit
void DrawHere(TVirtualPad *pd, Double_t nsigma=3, Int_t writeFitInfo=1)
void SetBackHighPolDegree(Int_t deg)
void ConfigurenSigmaBinCSteps(Int_t nSteps, Double_t *values)
TString fSuffix
mass of D meson
TH1F * fHistoMeanTrialAll
histo with gauss sigma from all trials
Int_t types
Bool_t MassFitter(Bool_t draw=kTRUE)
TH1F ** fHistoRawYieldTrial
histo with yield from subsamples of trials
Double_t fMassD
relative variation of the sigma
TF1 * GetBackgroundFullRangeFunc()
Int_t * fRebinSteps
number of rebin steps
TH2F ** fHistoRawYieldTrialBinC
histo with bin counts from subsamples of trials
TH1F ** fHistoSigmaTrial
histo with yield from subsamples of trials
void DrawHistos(TCanvas *cry) const
Double_t GetRawYieldError() const
Int_t typeb
Double_t * sigma
AliHFMassFitterVAR for the fit of invariant mass distribution of charmed mesons.
TH1F ** fHistoRawYieldDist
histo with bin counts from all trials
Double_t * fnSigmaBinCSteps
number of steps on the bin counting
Int_t fFitOption
name to characterize analysis case
void ConfigureLowLimFitSteps(Int_t nSteps, Double_t *values)
Bool_t fDrawIndividualFits
switch for pol5 background
Double_t GetMeanUncertainty() const
Bool_t fUsePol2Bkg
switch for linear background
Bool_t fUsePol5Bkg
switch for pol4 background
Double_t GetMean() const
TH1F ** fHistoRawYieldDistBinC
histo with chi2 from subsamples of trials
Double_t GetReducedChiSquare() const
void SetReflectionSigmaFactor(Int_t constant)
Int_t totTrials
TH1F * fHistoRawYieldDistBinCAll
histo with chi2 from all trials
Bool_t DoFitWithPol3Bkg(TH1F *histoToFit, Double_t hmin, Double_t hmax, Int_t theCase)
Double_t GetRawYield() const
void SetInitialGaussianSigma(Double_t sigma)
change the default value of the mean
TH1F ** fHistoMeanTrial
histo with gauss sigma from subsamples of trials
Bool_t Data(TH1F *h, Double_t *rangefit, Bool_t writefit, Double_t &sgn, Double_t &errsgn, Double_t &bkg, Double_t &errbkg, Double_t &sgnf, Double_t &errsgnf, Double_t &sigmafit, Int_t &status)
TH1F * RebinHisto(TH1D *hOrig, Int_t reb, Int_t firstUse) const
void BinCount(TH1F *h, TF1 *fB, Int_t rebin, Double_t minMass, Double_t maxMass, Double_t &count, Double_t &ecount) const
TH2F * fHistoRawYieldTrialBinCAll
histo with bin counts from all trials
Double_t fSigmaMCVariation
sigma of D meson peak from MC
Bool_t fUsePol3Bkg
switch for pol2 background
void SetFixGaussianSigma(Double_t sigma=0.012, Bool_t fixpar=kTRUE)
Bool_t fUsePol4Bkg
switch for pol3 background
Double_t minMass
Int_t rebin
TH1F ** fHistoChi2Trial
histo with gauss mean from subsamples of trials
TH1F * fHistoRawYieldTrialAll
histo with yield from all trials
void ConfigureRebinSteps(Int_t nSteps, Int_t *values)
void ConfigureUpLimFitSteps(Int_t nSteps, Double_t *values)
void SaveToRoot(TString fileName, TString option="recreate") const
Double_t maxMass
Double_t GetSigma() const
virtual ~AliHFMultiTrials()
Bool_t fUseExpoBkg
LL or chi2 fit.
Bool_t DoMultiTrials(TH1D *hInvMassHisto, TPad *thePad=0x0)
void SetFixGaussianMean(Double_t mean=1.865, Bool_t fixpar=kTRUE)
Double_t * fUpLimFitSteps
number of steps on the max. mass for fit
Double_t minMassForFit
TH1F * fHistoChi2TrialAll
histo with gauss mean from all trials
Double_t maxMassForFit
Double_t GetSigmaUncertainty() const