AliPhysics  251aa1e (251aa1e)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
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  fNumOfFirstBinSteps(1),
40  fNumOfLowLimFitSteps(6),
41  fLowLimFitSteps(0x0),
42  fNumOfUpLimFitSteps(6),
43  fUpLimFitSteps(0x0),
44  fNumOfnSigmaBinCSteps(11),
45  fnSigmaBinCSteps(0x0),
46  fSigmaGausMC(0.010),
47  fSigmaMCVariation(0.15),
48  fMassD(1.86484),
49  fSuffix(""),
50  fFitOption(0),
51  fUseExpoBkg(kTRUE),
52  fUseLinBkg(kTRUE),
53  fUsePol2Bkg(kTRUE),
54  fUsePol3Bkg(kTRUE),
55  fUsePol4Bkg(kTRUE),
56  fUsePol5Bkg(kFALSE),
57  fUsePowLawBkg(kFALSE),
58  fUsePowLawTimesExpoBkg(kFALSE),
59  fUseFixSigUpFreeMean(kTRUE),
60  fUseFixSigDownFreeMean(kTRUE),
61  fUseFreeS(kTRUE),
62  fUseFixedMeanFreeS(kTRUE),
63  fUseFixSigFreeMean(kTRUE),
64  fUseFixSigFixMean(kTRUE),
65  fDrawIndividualFits(kFALSE),
66  fHistoRawYieldDistAll(0x0),
67  fHistoRawYieldTrialAll(0x0),
68  fHistoSigmaTrialAll(0x0),
69  fHistoMeanTrialAll(0x0),
70  fHistoChi2TrialAll(0x0),
71  fHistoSignifTrialAll(0x0),
72  fHistoRawYieldDistBinCAll(0x0),
73  fHistoRawYieldTrialBinCAll(0x0),
74  fHistoRawYieldDist(0x0),
75  fHistoRawYieldTrial(0x0),
76  fHistoSigmaTrial(0x0),
77  fHistoMeanTrial(0x0),
78  fHistoChi2Trial(0x0),
79  fHistoSignifTrial(0x0),
80  fHistoRawYieldDistBinC(0x0),
81  fHistoRawYieldTrialBinC(0x0),
82  fhTemplRefl(0x0),
83  fFixRefloS(0),
84  fNtupleMultiTrials(0x0),
85  fMinYieldGlob(0),
86  fMaxYieldGlob(0)
87 {
88  // constructor
89  Int_t rebinStep[4]={3,4,5,6};
90  Double_t minMassStep[6]={1.68,1.70,1.72,1.74,1.76,1.78};
91  Double_t maxMassStep[6]={2.06,2.04,2.02,2.00,1.98,1.96};
92  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};
93  ConfigureRebinSteps(4,rebinStep);
94  ConfigureLowLimFitSteps(6,minMassStep);
95  ConfigureUpLimFitSteps(6,maxMassStep);
96  ConfigurenSigmaBinCSteps(11,nSigmasBC);
97 }
98 
99 //________________________________________________________________________
101  // destructor
102  delete [] fRebinSteps;
103  delete [] fLowLimFitSteps;
104  delete [] fUpLimFitSteps;
105  if(fhTemplRefl) delete fhTemplRefl;
106 }
107 
108 //________________________________________________________________________
110  // creates output histograms
111 
112  const Int_t nCases=kNBkgFuncCases*kNFitConfCases;
113 
114  TString funcBkg[kNBkgFuncCases]={"Expo","Lin","Pol2","Pol3","Pol4","Pol5","PowLaw","PowLawExpo"};
115  TString gausSig[kNFitConfCases]={"FixedS","FixedSp20","FixedSm20","FreeS","FixedMeanFixedS","FixedMeanFreeS"};
116 
118 
119  fHistoRawYieldDistAll = new TH1F(Form("hRawYieldDistAll%s",fSuffix.Data())," ; Raw Yield",5000,0.,50000.);
120  fHistoRawYieldTrialAll = new TH1F(Form("hRawYieldTrialAll%s",fSuffix.Data())," ; Trial # ; Raw Yield",nCases*totTrials,-0.5,nCases*totTrials-0.5);
121  fHistoSigmaTrialAll = new TH1F(Form("hSigmaTrialAll%s",fSuffix.Data())," ; Trial # ; Sigma (GeV/c^{2})",nCases*totTrials,-0.5,nCases*totTrials-0.5);
122  fHistoMeanTrialAll = new TH1F(Form("hMeanTrialAll%s",fSuffix.Data())," ; Trial # ; Mean (GeV/c^{2})",nCases*totTrials,-0.5,nCases*totTrials-0.5);
123  fHistoChi2TrialAll = new TH1F(Form("hChi2TrialAll%s",fSuffix.Data())," ; Trial # ; #chi^{2}",nCases*totTrials,-0.5,nCases*totTrials-0.5);
124  fHistoSignifTrialAll = new TH1F(Form("hSignifTrialAll%s",fSuffix.Data())," ; Trial # ; Significance",nCases*totTrials,-0.5,nCases*totTrials-0.5);
125 
126  fHistoRawYieldDistBinCAll = new TH1F(Form("hRawYieldDistBinCAll%s",fSuffix.Data())," ; Raw Yield (bin count)",5000,0.,50000.);
127  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);
128 
129  fHistoRawYieldDist = new TH1F*[nCases];
130  fHistoRawYieldTrial = new TH1F*[nCases];
131  fHistoSigmaTrial = new TH1F*[nCases];
132  fHistoMeanTrial = new TH1F*[nCases];
133  fHistoChi2Trial = new TH1F*[nCases];
134  fHistoSignifTrial = new TH1F*[nCases];
135 
136  fHistoRawYieldDistBinC = new TH1F*[nCases];
137  fHistoRawYieldTrialBinC = new TH2F*[nCases];
138 
139  for(Int_t ib=0; ib<kNBkgFuncCases; ib++){
140  for(Int_t igs=0; igs<kNFitConfCases; igs++){
141  Int_t theCase=igs*kNBkgFuncCases+ib;
142  fHistoRawYieldDist[theCase]=new TH1F(Form("hRawYieldDist%s%s%s",funcBkg[ib].Data(),gausSig[igs].Data(),fSuffix.Data())," ; Raw Yield",5000,0.,50000.);
143  fHistoRawYieldDistBinC[theCase]=new TH1F(Form("hRawYieldDistBinC%s%s%s",funcBkg[ib].Data(),gausSig[igs].Data(),fSuffix.Data())," ; Raw Yield (bin count)",5000,0.,50000.);
144  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);
145  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);
146  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);
147  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);
148  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);
149  fHistoSignifTrial[theCase]=new TH1F(Form("hSignifTrial%s%s%s",funcBkg[ib].Data(),gausSig[igs].Data(),fSuffix.Data())," ; Trial # ; Significance",totTrials,-0.5,totTrials-0.5);
150 
151  fHistoChi2Trial[theCase]->SetMarkerStyle(7);
152  fHistoSignifTrial[theCase]->SetMarkerStyle(7);
153  fHistoSigmaTrial[theCase]->SetMarkerStyle(7);
154  fHistoMeanTrial[theCase]->SetMarkerStyle(7);
155 
156  }
157  }
158  fNtupleMultiTrials = new TNtuple(Form("ntuMultiTrial%s",fSuffix.Data()),Form("ntuMultiTrial%s",fSuffix.Data()),"rebin:firstb:minfit:maxfit:bkgfunc:confsig:confmean:chi2:signif:mean:emean:sigma:esigma:rawy:erawy",128000);
159  return kTRUE;
160 
161 }
162 
163 //________________________________________________________________________
164 Bool_t AliHFMultiTrials::DoMultiTrials(TH1D* hInvMassHisto, TPad* thePad){
165  // perform the multiple fits
166 
167  Bool_t hOK=CreateHistos();
168  if(!hOK) return kFALSE;
169 
170  Int_t itrial=0;
171  Int_t types=0;
172  Int_t itrialBC=0;
174 
175  fMinYieldGlob=999999.;
176  fMaxYieldGlob=0.;
177  Float_t xnt[15];
178 
179  for(Int_t ir=0; ir<fNumOfRebinSteps; ir++){
180  Int_t rebin=fRebinSteps[ir];
181  for(Int_t iFirstBin=1; iFirstBin<=fNumOfFirstBinSteps; iFirstBin++) {
182  TH1F* hRebinned=0x0;
183  if(fNumOfFirstBinSteps==1) hRebinned=RebinHisto(hInvMassHisto,rebin,-1);
184  else hRebinned=RebinHisto(hInvMassHisto,rebin,iFirstBin);
185  for(Int_t iMinMass=0; iMinMass<fNumOfLowLimFitSteps; iMinMass++){
187  Double_t hmin=TMath::Max(minMassForFit,hRebinned->GetBinLowEdge(2));
188  for(Int_t iMaxMass=0; iMaxMass<fNumOfUpLimFitSteps; iMaxMass++){
190  Double_t hmax=TMath::Min(maxMassForFit,hRebinned->GetBinLowEdge(hRebinned->GetNbinsX()));
191  ++itrial;
192  for(Int_t typeb=0; typeb<kNBkgFuncCases; typeb++){
193  if(typeb==kExpoBkg && !fUseExpoBkg) continue;
194  if(typeb==kLinBkg && !fUseLinBkg) continue;
195  if(typeb==kPol2Bkg && !fUsePol2Bkg) continue;
196  if(typeb==kPol3Bkg && !fUsePol3Bkg) continue;
197  if(typeb==kPol4Bkg && !fUsePol4Bkg) continue;
198  if(typeb==kPol5Bkg && !fUsePol5Bkg) continue;
199  if(typeb==kPowBkg && !fUsePowLawBkg) continue;
201  for(Int_t igs=0; igs<kNFitConfCases; igs++){
202  if (igs==kFixSigUpFreeMean && !fUseFixSigUpFreeMean) continue;
203  if (igs==kFixSigDownFreeMean && !fUseFixSigDownFreeMean) continue;
204  if (igs==kFreeSigFixMean && !fUseFixedMeanFreeS) continue;
205  if (igs==kFreeSigFreeMean && !fUseFreeS) continue;
206  if (igs==kFixSigFreeMean && !fUseFixSigFreeMean) continue;
207  if (igs==kFixSigFixMean && !fUseFixSigFixMean) continue;
208  Int_t theCase=igs*kNBkgFuncCases+typeb;
209  Int_t globBin=itrial+theCase*totTrials;
210  for(Int_t j=0; j<15; j++) xnt[j]=0.;
211 
213  if(typeb<=kPol2Bkg){
214  fitter=new AliHFMassFitterVAR(hRebinned,hmin, hmax,1,typeb,types);
215  }else if(typeb==kPowBkg){
216  fitter=new AliHFMassFitterVAR(hRebinned,hmin, hmax,1,4,types);
217  }else if(typeb==kPowTimesExpoBkg){
218  fitter=new AliHFMassFitterVAR(hRebinned,hmin, hmax,1,5,types);
219  }else{
220  fitter=new AliHFMassFitterVAR(hRebinned,hmin, hmax,1,6,types);
221  if(typeb==kPol3Bkg) fitter->SetBackHighPolDegree(3);
222  if(typeb==kPol4Bkg) fitter->SetBackHighPolDegree(4);
223  if(typeb==kPol5Bkg) fitter->SetBackHighPolDegree(5);
224  }
225  fitter->SetReflectionSigmaFactor(0);
226  //if D0 Reflection
227  if(fhTemplRefl){
228  fitter=new AliHFMassFitterVAR(hRebinned,hmin,hmax,1,typeb,2);
230  fitter->SetFixReflOverS(fFixRefloS,kTRUE);
231  }
232  if(fFitOption==1) fitter->SetUseChi2Fit();
235  xnt[0]=rebin;
236  xnt[1]=iFirstBin;
237  xnt[2]=minMassForFit;
238  xnt[3]=maxMassForFit;
239  xnt[4]=typeb;
240  xnt[6]=0;
241  if(igs==kFixSigFreeMean){
242  fitter->SetFixGaussianSigma(fSigmaGausMC,kTRUE);
243  xnt[5]=1;
244  }else if(igs==kFixSigUpFreeMean){
246  xnt[5]=2;
247  }else if(igs==kFixSigDownFreeMean){
249  xnt[5]=3;
250  }else if(igs==kFreeSigFreeMean){
251  xnt[5]=0;
252  }else if(igs==kFixSigFixMean){
253  fitter->SetFixGaussianSigma(fSigmaGausMC,kTRUE);
254  fitter->SetFixGaussianMean(fMassD,kTRUE);
255  xnt[5]=1;
256  xnt[6]=1;
257  }else if(igs==kFreeSigFixMean){
258  fitter->SetFixGaussianMean(fMassD,kTRUE);
259  xnt[5]=0;
260  xnt[6]=1;
261  }
262  Bool_t out=kFALSE;
263  Double_t chisq=-1.;
264  Double_t sigma=0.;
265  Double_t esigma=0.;
266  Double_t pos=.0;
267  Double_t epos=.0;
268  Double_t ry=.0;
269  Double_t ery=.0;
270  Double_t significance=0.;
271  Double_t erSignif=0.;
272  TF1* fB1=0x0;
273  if(typeb<kNBkgFuncCases){
274  printf("****** START FIT OF HISTO %s WITH REBIN %d FIRST BIN %d MASS RANGE %f-%f BACKGROUND FIT FUNCTION=%d CONFIG SIGMA/MEAN=%d\n",hInvMassHisto->GetName(),rebin,iFirstBin,minMassForFit,maxMassForFit,typeb,igs);
275  out=fitter->MassFitter(0);
276  chisq=fitter->GetReducedChiSquare();
277  fitter->Significance(3,significance,erSignif);
278  sigma=fitter->GetSigma();
279  pos=fitter->GetMean();
280  esigma=fitter->GetSigmaUncertainty();
281  if(esigma<0.00001) esigma=0.0001;
282  epos=fitter->GetMeanUncertainty();
283  if(epos<0.00001) epos=0.0001;
284  ry=fitter->GetRawYield();
285  ery=fitter->GetRawYieldError();
286  fB1=fitter->GetBackgroundFullRangeFunc();
287  if(out && fDrawIndividualFits && thePad){
288  thePad->Clear();
289  fitter->DrawHere(thePad);
290  }
291  }
292  // else{
293  // out=DoFitWithPol3Bkg(hRebinned,hmin,hmax,igs);
294  // if(out && thePad){
295  // thePad->Clear();
296  // hRebinned->Draw();
297  // TF1* fSB=(TF1*)hRebinned->GetListOfFunctions()->FindObject("fSB");
298  // fB1=new TF1("fB1","[0]+[1]*x+[2]*x*x+[3]*x*x*x",hmin,hmax);
299  // for(Int_t j=0; j<4; j++) fB1->SetParameter(j,fSB->GetParameter(3+j));
300  // fB1->SetLineColor(2);
301  // fB1->Draw("same");
302  // fSB->SetLineColor(4);
303  // fSB->Draw("same");
304  // thePad->Update();
305  // chisq=fSB->GetChisquare()/fSB->GetNDF();;
306  // sigma=fSB->GetParameter(2);
307  // esigma=fSB->GetParError(2);
308  // if(esigma<0.00001) esigma=0.0001;
309  // pos=fSB->GetParameter(1);
310  // epos=fSB->GetParError(1);
311  // if(epos<0.00001) epos=0.0001;
312  // ry=fSB->GetParameter(0)/hRebinned->GetBinWidth(1);
313  // ery=fSB->GetParError(0)/hRebinned->GetBinWidth(1);
314  // }
315  // }
316  xnt[7]=chisq;
317  if(out && chisq>0. && sigma>0.5*fSigmaGausMC && sigma<2.0*fSigmaGausMC){
318  xnt[8]=significance;
319  xnt[9]=pos;
320  xnt[10]=epos;
321  xnt[11]=sigma;
322  xnt[12]=esigma;
323  xnt[13]=ry;
324  xnt[14]=ery;
325  fHistoRawYieldDistAll->Fill(ry);
326  fHistoRawYieldTrialAll->SetBinContent(globBin,ry);
327  fHistoRawYieldTrialAll->SetBinError(globBin,ery);
328  fHistoSigmaTrialAll->SetBinContent(globBin,sigma);
329  fHistoSigmaTrialAll->SetBinError(globBin,esigma);
330  fHistoMeanTrialAll->SetBinContent(globBin,pos);
331  fHistoMeanTrialAll->SetBinError(globBin,epos);
332  fHistoChi2TrialAll->SetBinContent(globBin,chisq);
333  fHistoChi2TrialAll->SetBinError(globBin,0.00001);
334  fHistoSignifTrialAll->SetBinContent(globBin,significance);
335  fHistoSignifTrialAll->SetBinError(globBin,erSignif);
336  if(ry<fMinYieldGlob) fMinYieldGlob=ry;
337  if(ry>fMaxYieldGlob) fMaxYieldGlob=ry;
338  fHistoRawYieldDist[theCase]->Fill(ry);
339  fHistoRawYieldTrial[theCase]->SetBinContent(itrial,ry);
340  fHistoRawYieldTrial[theCase]->SetBinError(itrial,ery);
341  fHistoSigmaTrial[theCase]->SetBinContent(itrial,sigma);
342  fHistoSigmaTrial[theCase]->SetBinError(itrial,esigma);
343  fHistoMeanTrial[theCase]->SetBinContent(itrial,pos);
344  fHistoMeanTrial[theCase]->SetBinError(itrial,epos);
345  fHistoChi2Trial[theCase]->SetBinContent(itrial,chisq);
346  fHistoChi2Trial[theCase]->SetBinError(itrial,0.00001);
347  fHistoSignifTrial[theCase]->SetBinContent(itrial,significance);
348  fHistoSignifTrial[theCase]->SetBinError(itrial,erSignif);
349  for(Int_t iStepBC=0; iStepBC<fNumOfnSigmaBinCSteps; iStepBC++){
350  Double_t minMassBC=fMassD-fnSigmaBinCSteps[iStepBC]*sigma;
351  Double_t maxMassBC=fMassD+fnSigmaBinCSteps[iStepBC]*sigma;
352  if(minMassBC>minMassForFit &&
353  maxMassBC<maxMassForFit &&
354  minMassBC>(hRebinned->GetXaxis()->GetXmin()) &&
355  maxMassBC<(hRebinned->GetXaxis()->GetXmax())){
356  Double_t cnts,ecnts;
357  BinCount(hRebinned,fB1,1,minMassBC,maxMassBC,cnts,ecnts);
358  ++itrialBC;
359  fHistoRawYieldDistBinCAll->Fill(cnts);
360  fHistoRawYieldTrialBinCAll->SetBinContent(globBin,iStepBC+1,cnts);
361  fHistoRawYieldTrialBinCAll->SetBinError(globBin,iStepBC+1,ecnts);
362  fHistoRawYieldTrialBinC[theCase]->SetBinContent(itrial,iStepBC+1,cnts);
363  fHistoRawYieldTrialBinC[theCase]->SetBinError(itrial,iStepBC+1,ecnts);
364  fHistoRawYieldDistBinC[theCase]->Fill(cnts);
365  }
366  }
367  }
368  delete fitter;
369  // if(typeb>4) delete fB1;
370  fNtupleMultiTrials->Fill(xnt);
371  }
372  }
373  }
374  }
375  delete hRebinned;
376  }
377  }
378  return kTRUE;
379 }
380 
381 //________________________________________________________________________
383  // save histos in a root file for further analysis
384  const Int_t nCases=kNBkgFuncCases*kNFitConfCases;
385  TFile* outHistos=new TFile(fileName.Data(),option.Data());
386  fHistoRawYieldTrialAll->Write();
387  fHistoSigmaTrialAll->Write();
388  fHistoMeanTrialAll->Write();
389  fHistoChi2TrialAll->Write();
390  fHistoSignifTrialAll->Write();
391  fHistoRawYieldDistBinCAll->Write();
392  fHistoRawYieldTrialBinCAll->Write();
393  for(Int_t ic=0; ic<nCases; ic++){
394  fHistoRawYieldTrial[ic]->Write();
395  fHistoSigmaTrial[ic]->Write();
396  fHistoMeanTrial[ic]->Write();
397  fHistoChi2Trial[ic]->Write();
398  fHistoSignifTrial[ic]->Write();
399  fHistoRawYieldTrialBinC[ic]->Write();
400  fHistoRawYieldDistBinC[ic]->Write();
401  }
402  fNtupleMultiTrials->Write();
403  outHistos->Close();
404  delete outHistos;
405 }
406 
407 //________________________________________________________________________
408 void AliHFMultiTrials::DrawHistos(TCanvas* cry) const{
409  // draw histos
410  cry->Divide(2,2);
411  cry->cd(1);
412  fHistoSigmaTrialAll->Draw();
413  cry->cd(3);
414  fHistoChi2TrialAll->Draw();
415  cry->cd(2);
416  fHistoRawYieldTrialAll->Draw();
417  cry->cd(4);
418  fHistoRawYieldDistAll->Draw();
419  TLatex* tmean=new TLatex(0.15,0.8,Form("mean=%.1f",fHistoRawYieldDistAll->GetMean()));
420  tmean->SetNDC();
421  tmean->Draw();
422  TLatex* thrms=new TLatex(0.15,0.72,Form("rms=%.1f",fHistoRawYieldDistAll->GetRMS()));
423  thrms->SetNDC();
424  thrms->Draw();
425  TLatex* tmax=new TLatex(0.6,0.8,Form("max=%.1f",fMaxYieldGlob));
426  tmax->SetNDC();
427  tmax->Draw();
428  TLatex* tmin=new TLatex(0.6,0.72,Form("min=%.1f",fMinYieldGlob));
429  tmin->SetNDC();
430  tmin->Draw();
431  TLatex* trms=new TLatex(0.6,0.64,Form("(max-min)/sqrt(12)=%.1f",(fMaxYieldGlob-fMinYieldGlob)/sqrt(12)));
432  trms->SetNDC();
433  trms->Draw();
434 
435 }
436 //________________________________________________________________________
437 TH1F* AliHFMultiTrials::RebinHisto(TH1D* hOrig, Int_t reb, Int_t firstUse) const{
438  // Rebin histogram, from bin firstUse to lastUse
439  // Use all bins if firstUse=-1
440 
441  Int_t nBinOrig=hOrig->GetNbinsX();
442  Int_t firstBinOrig=1;
443  Int_t lastBinOrig=nBinOrig;
444  Int_t nBinOrigUsed=nBinOrig;
445  Int_t nBinFinal=nBinOrig/reb;
446  if(firstUse>=1){
447  firstBinOrig=firstUse;
448  nBinFinal=(nBinOrig-firstUse+1)/reb;
449  nBinOrigUsed=nBinFinal*reb;
450  lastBinOrig=firstBinOrig+nBinOrigUsed-1;
451  }else{
452  Int_t exc=nBinOrigUsed%reb;
453  if(exc!=0){
454  nBinOrigUsed-=exc;
455  firstBinOrig+=exc/2;
456  lastBinOrig=firstBinOrig+nBinOrigUsed-1;
457  }
458  }
459 
460  printf("Rebin from %d bins to %d bins -- Used bins=%d in range %d-%d\n",nBinOrig,nBinFinal,nBinOrigUsed,firstBinOrig,lastBinOrig);
461  Float_t lowLim=hOrig->GetXaxis()->GetBinLowEdge(firstBinOrig);
462  Float_t hiLim=hOrig->GetXaxis()->GetBinUpEdge(lastBinOrig);
463  TH1F* hRebin=new TH1F(Form("%s-rebin%d_%d",hOrig->GetName(),reb,firstUse),hOrig->GetTitle(),nBinFinal,lowLim,hiLim);
464  Int_t lastSummed=firstBinOrig-1;
465  for(Int_t iBin=1;iBin<=nBinFinal; iBin++){
466  Float_t sum=0.;
467  Float_t sum2=0.;
468  for(Int_t iOrigBin=0;iOrigBin<reb;iOrigBin++){
469  sum+=hOrig->GetBinContent(lastSummed+1);
470  sum2+=hOrig->GetBinError(lastSummed+1)*hOrig->GetBinError(lastSummed+1);
471  lastSummed++;
472  }
473  hRebin->SetBinContent(iBin,sum);
474  hRebin->SetBinError(iBin,TMath::Sqrt(sum2));
475  }
476  return hRebin;
477 }
478 
479 //________________________________________________________________________
480 void AliHFMultiTrials::BinCount(TH1F* h, TF1* fB, Int_t rebin, Double_t minMass, Double_t maxMass, Double_t& count, Double_t& ecount) const{
481  // compute yield with bin couting
482  Int_t minBinSum=h->FindBin(minMass);
483  Int_t maxBinSum=h->FindBin(maxMass);
484  Double_t cntSig=0.;
485  Double_t cntErr=0.;
486  for(Int_t iMB=minBinSum; iMB<=maxBinSum; iMB++){
487  Double_t bkg=fB ? fB->Eval(h->GetBinCenter(iMB))/(Double_t)rebin : 0;
488  cntSig+=(h->GetBinContent(iMB)-bkg);
489  cntErr+=(h->GetBinError(iMB)*h->GetBinError(iMB));
490  }
491  count=cntSig;
492  ecount=TMath::Sqrt(cntErr);
493 }
494 //________________________________________________________________________
496  Int_t iCase){
497  //
498 
499  TH1F *hCutTmp=(TH1F*)histoToFit->Clone("hCutTmp");
500  for(Int_t ib=1; ib<=hCutTmp->GetNbinsX(); ib++){
501  Double_t xc=hCutTmp->GetBinCenter(ib);
502  if(xc>(fMassD-5.*fSigmaGausMC) && xc<(fMassD+5.*fSigmaGausMC)){
503  hCutTmp->SetBinContent(ib,0.);
504  hCutTmp->SetBinError(ib,0.);
505  }
506  }
507 
508  hCutTmp->Fit("pol2","E0","",hmin,hmax);
509  TF1* f2=(TF1*)hCutTmp->GetListOfFunctions()->FindObject("pol2");
510  TF1* f3=new TF1("myPol3","pol3");
511  for(Int_t i=0; i<3;i++) f3->SetParameter(i,f2->GetParameter(i));
512  hCutTmp->Fit(f3,"E0","",hmin,hmax);
513  Double_t quickCount=0.;
514  for(Int_t ib=1; ib<=histoToFit->GetNbinsX(); ib++){
515  Double_t xc=hCutTmp->GetBinCenter(ib);
516  if(xc>(fMassD-3.*fSigmaGausMC) && xc<(fMassD+3.*fSigmaGausMC)){
517  quickCount+=(histoToFit->GetBinContent(ib)-f3->Eval(xc));
518  }
519  }
520  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");
521  fSB->SetParameter(0,quickCount);
522  fSB->SetParameter(1,fMassD);
523  fSB->SetParameter(2,fSigmaGausMC);
524  for(Int_t j=0; j<4; j++) fSB->SetParameter(j+3,f3->GetParameter(j));
525  if(iCase==0) fSB->FixParameter(2,fSigmaGausMC);
526  else if(iCase==1) fSB->FixParameter(2,fSigmaGausMC*(1.+fSigmaMCVariation));
527  else if(iCase==2) fSB->FixParameter(2,fSigmaGausMC*(1.-fSigmaMCVariation));
528  else if(iCase==4){
529  fSB->FixParameter(1,fMassD);
530  fSB->FixParameter(2,fSigmaGausMC);
531  } else if(iCase==5){
532  fSB->FixParameter(1,fMassD);
533  }
534  histoToFit->Fit(fSB,"ME0","",hmin,hmax);
535  // quality cuts
536  if(fSB->GetParError(0)<0.01*fSB->GetParameter(0)) return kFALSE;
537  if(fSB->GetParError(0)>0.6*fSB->GetParameter(0)) return kFALSE;
538 
539  delete hCutTmp;
540  return kTRUE;
541 }
542 //__________________________________________________________________________________
543 TH1F* AliHFMultiTrials::SetTemplateRefl(const TH1F *h) {
544  fhTemplRefl=(TH1F*)h->Clone("hTemplRefl");
545  return fhTemplRefl;
546 }
547 
void Significance(Double_t nOfSigma, Double_t &significance, Double_t &errsignificance) const
backgournd in (min, max) with error
void SetFixReflOverS(Double_t val, Bool_t fixpar=kTRUE)
double Double_t
Definition: External.C:58
Double_t fMaxYieldGlob
minimum yield
Bool_t fUsePowLawBkg
switch for pol5 background
Definition: External.C:236
TH1F * fhTemplRefl
histo with bin counts from subsamples of trials
TH1F * fHistoSigmaTrialAll
histo with yield from all trials
TNtuple * fNtupleMultiTrials
Double_t fMinYieldGlob
tree
TH1F * SetTemplateRefl(const TH1F *hTemplRefl)
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)
Bool_t fUseFixSigFixMean
switch for FixSigFreeMean
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
Int_t fNumOfLowLimFitSteps
number of steps in the first bin for rebin
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 ** fHistoSignifTrial
histo with chi2 from subsamples of trials
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
int Int_t
Definition: External.C:63
TH1F * fHistoSignifTrialAll
histo with chi2 from all trials
void ConfigureLowLimFitSteps(Int_t nSteps, Double_t *values)
Bool_t fDrawIndividualFits
switch for FixSigFixMean
Double_t GetMeanUncertainty() const
Bool_t fUsePol2Bkg
switch for linear background
Bool_t fUsePol5Bkg
switch for pol4 background
float Float_t
Definition: External.C:68
Double_t GetMean() const
Bool_t fUseFreeS
switch for FixSigDownFreeMean
TH1F ** fHistoRawYieldDistBinC
histo with chi2 from subsamples of trials
Double_t GetReducedChiSquare() const
Bool_t fUseFixedMeanFreeS
switch for FreeSigma
void SetReflectionSigmaFactor(Int_t constant)
Definition: External.C:212
Bool_t fUseFixSigUpFreeMean
switch for power law background
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
TList * fitter
Definition: DrawAnaELoss.C:26
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
TH1F * SetTemplateReflections(const TH1 *h, TString option="templ", Double_t minRange=1.72, Double_t maxRange=2.05)
Bool_t fUsePol3Bkg
switch for pol2 background
void SetFixGaussianSigma(Double_t sigma=0.012, Bool_t fixpar=kTRUE)
Bool_t fUseFixSigDownFreeMean
switch for FixSigUpFreeMean
Bool_t fUsePol4Bkg
switch for pol3 background
ClassImp(AliAnalysisTaskCRC) AliAnalysisTaskCRC
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
Bool_t fUseFixSigFreeMean
switch for FixedMeanFreeS
Double_t maxMass
Double_t GetSigma() const
bool Bool_t
Definition: External.C:53
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
Definition: External.C:196
Float_t fFixRefloS
template of reflection contribution
Double_t minMassForFit
TH1F * fHistoChi2TrialAll
histo with gauss mean from all trials
Double_t maxMassForFit
Double_t GetSigmaUncertainty() const
Bool_t fUsePowLawTimesExpoBkg
switch for power law background