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