AliPhysics  master (3d17d9d)
ProjectCombinHFAndFit.C
Go to the documentation of this file.
1 #if !defined(__CINT__) || defined(__MAKECINT__)
2 #include <TString.h>
3 #include <TMath.h>
4 #include <TFile.h>
5 #include <TCanvas.h>
6 #include <TH1F.h>
7 #include <TH2F.h>
8 #include <TH3F.h>
9 #include <TH1D.h>
10 #include <TF1.h>
11 #include <TStyle.h>
12 #include <TLatex.h>
13 #include <TPaveText.h>
14 #include <TLegend.h>
15 #include <TLegendEntry.h>
16 #include <TDatabasePDG.h>
17 #include "AliHFInvMassFitter.h"
18 #include "AliVertexingHFUtils.h"
20 #endif
21 
23 
24 TString configFileName="configfile4lowptanalysis.txt";
25 
26 // input files and pt binning
31 
32 TString meson="Dzero";
33 const Int_t maxPtBins=30;
35 Double_t binLims[maxPtBins+1]={0.,1.,2.,3.,4.,5.,6.,8.,12.};
36 
37 // fit configuration
38 Int_t rebin[maxPtBins]={4,6,7,8,9,10,10,12};
39 Double_t minMass4Fit[maxPtBins]={1.72,1.72,1.72,1.72,1.72,1.72,1.72,1.72};
40 Double_t maxMass4Fit[maxPtBins]={2.04,2.04,2.04,2.04,2.04,2.04,2.04,2.04};
41 Int_t fixSigmaConf=1; // 0= all free, 1=all fixed, 2=use values per pt bin
42 Bool_t fixSigma[maxPtBins]={kFALSE,kFALSE,kFALSE,kFALSE,kFALSE,kFALSE,kFALSE,kFALSE};
43 Double_t tuneSigmaOnData=-1.; // scaling factor for the Gaussian sigma, if -1. use MC out of the box
44 Double_t sigmas[maxPtBins]={0.006,0.008,0.009,0.010,0.011,0.012,0.013,0.013};
45 Int_t fixMeanConf=0; // 0= all free, 1=all fixed, 2=use values per pt bin
46 Bool_t fixMean[maxPtBins]={kFALSE,kFALSE,kFALSE,kFALSE,kFALSE,kFALSE,kFALSE,kFALSE};
47 Double_t tuneMeanOnData=1.868; // if -1. use PDG value
48 // objects and options related to side-band fit method
50 Double_t fitSBrangelow[maxPtBins]={1.74,1.74,1.74,1.72,1.72,1.72,1.72,1.72};
51 Double_t fitSBrangeup[maxPtBins]={2.0,2.0,2.0,2.0,2.0,2.0,2.0,2.04};
52 Int_t nDegreeBackPolSB[maxPtBins]={4,4,4,2,2,2,2,2}; // degree of polynomial function describing the background
53 
60 
63 Int_t nDegreeBackPol[maxPtBins]={2,2,2,2,2,2,2,2}; // degree of polynomial function describing the background
64 // reflection option
66 TString reflopt="2gaus";
68 
70 Double_t nsigmaBinCounting=4.; // defines the mass interval over which the signal is bin counted
72 
73 // outputfiles
75 Int_t saveCanvasAsEps=1; //0=none, 1=main ones, 2=all
76 
77 
83 
86 TH1F* FitMCInvMassSpectra(TList* lMC, TString var);
87 Bool_t ReadConfig(TString configName);
88 void PrintConfig();
89 
90 AliHFInvMassFitter* ConfigureFitter(TH1D* histo, Int_t iPtBin, Int_t backcase, Double_t minFit, Double_t maxFit, Bool_t isDirect=kFALSE){
91  TH1F* histof=(TH1F*)histo->Clone(Form("%s_Fl",histo->GetName()));
92 
93 
94  AliHFInvMassFitter* fitter=new AliHFInvMassFitter(histof,minFit,maxFit,backcase,0);
95  if(backcase==6){
97  if(isDirect) fitter->SetPolDegreeForBackgroundFit(nDegreeBackPolSB[iPtBin]);
98  }
99  fitter->SetUseChi2Fit();
100  fitter->SetFitOption(fitoption.Data());
101  fitter->SetInitialGaussianMean(massD);
102  fitter->SetInitialGaussianSigma(sigmas[iPtBin]);
103  if(fixSigmaConf==1 || (fixSigmaConf==2 && fixSigma[iPtBin])) fitter->SetFixGaussianSigma(sigmas[iPtBin]);
104  if(fixMeanConf==1 || (fixMeanConf==2 && fixMean[iPtBin])) fitter->SetFixGaussianMean(massD);
105  if(correctForRefl){
106  TCanvas *cTest=new TCanvas("cTest","cTest",1600,800);
107  cTest->Divide(2,1);
108  cTest->cd(1);
109  TH1F *hmasstemp=fitter->GetHistoClone();
110  TH1F *hReflModif=(TH1F*)AliVertexingHFUtils::AdaptTemplateRangeAndBinning(hMCReflPtBin,hmasstemp,minFit,maxFit);
111  TH1F *hSigModif=(TH1F*)AliVertexingHFUtils::AdaptTemplateRangeAndBinning(hMCSigPtBin,hmasstemp,minFit,maxFit);
112  hReflModif->SetLineColor(kRed);
113  hSigModif->SetLineColor(kBlue);
114  hSigModif->Draw();
115  hMCSigPtBin->SetLineColor(kRed-9);
116  hMCSigPtBin->Draw("same");
117  hMCReflPtBin->SetLineColor(kGray+1);
118  hMCReflPtBin->Draw("same");
119  hReflModif->Draw("same");
120  delete hmasstemp;
121  Double_t fixSoverRefAt=rOverSmodif*(hReflModif->Integral(hReflModif->FindBin(minFit*1.0001),hReflModif->FindBin(maxFit*0.999))/hSigModif->Integral(hSigModif->FindBin(minFit*1.0001),hSigModif->FindBin(maxFit*0.999)));
122  TH1F* hrfl=fitter->SetTemplateReflections(hReflModif,reflopt,-1.,-1.);
123  cTest->cd(2);
124  hReflModif->Draw();
125  hrfl->Draw("same");
126  cTest->SaveAs(Form("figures/ReflectionConfig_PtBin%d.eps",iPtBin));
127  if(!hrfl){
128  Printf("SOMETHING WENT WRONG WHILE SETTINGS REFLECTIONS TEMPLATE");
129  delete hReflModif;
130  delete hSigModif;
131  delete cTest;
132  return 0x0;
133  }
134  fitter->SetFixReflOverS(fixSoverRefAt);
135  delete cTest;
136  delete hReflModif;
137  delete hSigModif;
138  }
139  return fitter;
140 }
141 
142 
143 
144 void PrintGausParams(TH1F* hPulls){
145  TF1* fgfit=(TF1*)hPulls->GetListOfFunctions()->FindObject("gaus");
146  TLatex* tg1=new TLatex(0.2,0.8,Form("Gauss mean = %.2f#pm%.2f",fgfit->GetParameter(1),fgfit->GetParError(1)));
147  tg1->SetNDC();
148  tg1->Draw();
149  TLatex* tg2=new TLatex(0.2,0.7,Form("Gauss #sigma = %.2f#pm%.2f",fgfit->GetParameter(2),fgfit->GetParError(2)));
150  tg2->SetNDC();
151  tg2->Draw();
152 
153 }
154 
155 
156 Bool_t QuadraticSmooth(TH1 *h,Int_t ntimes=1){// quadratic fit of 5 points
157  ntimes--;
158  TF1 *fp2=new TF1("fp2","pol2",h->GetXaxis()->GetBinLowEdge(1),h->GetXaxis()->GetBinUpEdge(h->GetNbinsX()));
159  TH1D *htmp=(TH1D*)h->Clone("htmp");
160  for(Int_t i=3;i<=h->GetNbinsX()-3;i++){ // first intermediate bins
161  htmp->Fit(fp2,"RLEMN0","",h->GetXaxis()->GetBinLowEdge(i-2),h->GetXaxis()->GetBinUpEdge(i+2));
162  h->SetBinContent(i,(fp2->Integral(h->GetXaxis()->GetBinLowEdge(i),h->GetXaxis()->GetBinUpEdge(i)))/h->GetBinWidth(i));// change only content, errors unchanged
163  }
164  // now initial and final bins
165  htmp->Fit(fp2,"RLEMN0","",h->GetXaxis()->GetBinLowEdge(1),h->GetXaxis()->GetBinUpEdge(5));// 5 pts, asymmetric fit region
166  h->SetBinContent(2,(fp2->Integral(h->GetXaxis()->GetBinLowEdge(2),h->GetXaxis()->GetBinUpEdge(2)))/h->GetBinWidth(2));// change only content, errors unchanged
167 
168  htmp->Fit(fp2,"RLEMN0","",h->GetXaxis()->GetBinLowEdge(1),h->GetXaxis()->GetBinUpEdge(4));// use only 4 pts
169  h->SetBinContent(1,(fp2->Integral(h->GetXaxis()->GetBinLowEdge(1),h->GetXaxis()->GetBinUpEdge(1)))/h->GetBinWidth(1));// change only content, errors unchanged
170 
171 
172  htmp->Fit(fp2,"RLEMN0","",h->GetXaxis()->GetBinLowEdge(h->GetNbinsX()-4),h->GetXaxis()->GetBinUpEdge(h->GetNbinsX()));// 5 pts, asymmetric fit region
173  h->SetBinContent(h->GetNbinsX()-1,(fp2->Integral(h->GetXaxis()->GetBinLowEdge(h->GetNbinsX()-1),h->GetXaxis()->GetBinUpEdge(h->GetNbinsX()-1)))/h->GetBinWidth(h->GetNbinsX()-1));// change only content, errors unchanged
174 
175  htmp->Fit(fp2,"RLEMN0","",h->GetXaxis()->GetBinLowEdge(h->GetNbinsX()-3),h->GetXaxis()->GetBinUpEdge(h->GetNbinsX()));
176  h->SetBinContent(h->GetNbinsX(),(fp2->Integral(h->GetXaxis()->GetBinLowEdge(h->GetNbinsX()),h->GetXaxis()->GetBinUpEdge(h->GetNbinsX())))/h->GetBinWidth(h->GetNbinsX()));// change only content, errors unchanged
177  if(ntimes>0)QuadraticSmooth(h,ntimes);
178  delete htmp;
179  return kTRUE;
180 }
181 
182 
183 void SetStyleHisto(TH1 *h,Int_t method,Int_t isXpt=-1){
184  if(isXpt==1)h->GetXaxis()->SetTitle("#it{p}_{T} (GeV/#it{c})");// just for convenience, most of the time it's pt :)
185 
186  if(method==kRot){
187  h->SetMarkerStyle(21);
188  h->GetYaxis()->SetTitleOffset(1.8);
189  h->SetMarkerColor(kBlack);
190  h->SetMarkerSize(1.0);
191  h->SetLineColor(kBlack);
192  return;
193  }
194 
195  if(method==kME){
196  h->SetMarkerStyle(25);
197  h->GetYaxis()->SetTitleOffset(1.8);
198  h->SetMarkerColor(4);
199  h->SetMarkerSize(1.0);
200  h->SetLineColor(4);
201  return;
202  }
203 
204  if(method==kLS){
205  h->SetMarkerStyle(22);
206  h->GetYaxis()->SetTitleOffset(1.8);
207  h->SetMarkerColor(kGreen+2);
208  h->SetMarkerSize(1.0);
209  h->SetLineColor(kGreen+2);
210  return;
211  }
212 
213  if(method==kSB){
214  h->SetMarkerStyle(27);
215  h->GetYaxis()->SetTitleOffset(1.8);
216  h->SetMarkerColor(6);
217  h->SetMarkerSize(1.0);
218  h->SetLineColor(6);
219  return;
220  }
221 
222  return;
223 }
224 
225 void DivideCanvas(TCanvas *c,Int_t ndivisions){
226  if(ndivisions>1){
227  if(ndivisions<3){
228  c->Divide(2,1);
229  }
230  else if(ndivisions<5){
231  c->Divide(2,2);
232  }
233  else if(ndivisions<7){
234  c->Divide(3,2);
235 
236  }
237  else if(ndivisions<=8){
238  c->Divide(4,2);
239 
240  }
241  else if(ndivisions<10){
242  c->Divide(3,3);
243 
244  }
245  else if(ndivisions<13){
246  c->Divide(4,3);
247 
248  }
249  else if(ndivisions<17){
250  c->Divide(4,4);
251 
252  }
253  else {
254  c->Divide(5,5);
255 
256  }
257  }else{
258  c->Divide(1,1);
259  }
260  return;
261 }
262 
263 
264 
265 TF1 *GausPlusLine(Double_t minRange=1.72,Double_t maxRange=2.05){
266  TF1 *f=new TF1("GausPlusLine","[0]*1./(TMath::Sqrt(2.*TMath::Pi())*[2])*TMath::Exp(-(x-[1])*(x-[1])/(2.*[2]*[2]))+[3]+[4]*x",minRange,maxRange);
267 
268  f->SetParameter(1,massD);
269  f->SetParLimits(0,0.,1.e+7);
270  f->SetParLimits(1,1.8,1.92);
271  f->SetParameter(2,0.010);
272  f->SetParLimits(2,0,0.1);
273  f->SetParameter(3,0.);
274  // f->SetParLimits(3,-10000,10000);
275  return f;
276 
277 }
278 
279 
280 
282  //
283  Double_t norm=hRatio->GetMaximum();
284  Double_t massForNorm=0;
285 
286  if(optForNorm==0){
287  Int_t binl=hRatio->GetXaxis()->FindBin(minMass);
288  Int_t binh=hRatio->GetXaxis()->FindBin(maxMass);
289  Double_t norml=hRatio->GetBinContent(1);
290  if(binl>0 && binl<=hRatio->GetNbinsX()) norml=hRatio->GetBinContent(binl);
291  Double_t normh=hRatio->GetBinContent(hRatio->GetNbinsX());
292  if(binh>0 && binh<=hRatio->GetNbinsX()) normh=hRatio->GetBinContent(binh);
293  if(norml>normh){
294  norm=norml;
295  massForNorm=hRatio->GetBinCenter(binl);
296  }else{
297  norm=normh;
298  massForNorm=hRatio->GetBinCenter(binh);
299  }
300  }else if(optForNorm==1){
301  norm=0.0001;
302  for(Int_t iMassBin=1; iMassBin<hRatio->GetNbinsX(); iMassBin++){
303  Double_t bce=hRatio->GetBinCenter(iMassBin);
304  if(bce>minMass && bce<maxMass){
305  Double_t bco=hRatio->GetBinContent(iMassBin);
306  if(bco>norm){
307  norm=bco;
308  massForNorm=bce;
309  }
310  }
311  }
312  }else if(optForNorm==2){
313  norm=0.0001;
314  for(Int_t iMassBin=1; iMassBin<hRatio->GetNbinsX(); iMassBin++){
315  Double_t bce=hRatio->GetBinCenter(iMassBin);
316  Double_t bco=hRatio->GetBinContent(iMassBin);
317  if(bco>norm){
318  norm=bco;
319  massForNorm=bce;
320  }
321  }
322  }else if(optForNorm==3){
323  Int_t binl=hRatio->GetXaxis()->FindBin(minMass);
324  Int_t binh=hRatio->GetXaxis()->FindBin(maxMass);
325  Double_t norml=0;
326  Double_t normReb=0;
327  Int_t iFirstBin=TMath::Max(binl-reb/2,1);
328  Int_t iLastBin=iFirstBin+reb;
329  if(iLastBin>hRatio->GetNbinsX()){
330  iLastBin=hRatio->GetNbinsX();
331  iFirstBin=iLastBin-reb;
332  }
333  for(Int_t jjj=iFirstBin; jjj<=iLastBin; jjj++){
334  norml+=hRatio->GetBinContent(jjj);
335  normReb+=1;
336  }
337  if(normReb>=1) norml/=normReb;
338  Double_t normh=0;
339  normReb=0;
340  iFirstBin=TMath::Max(binh-reb/2,1);
341  iLastBin=iFirstBin+reb;
342  if(iLastBin>hRatio->GetNbinsX()){
343  iLastBin=hRatio->GetNbinsX();
344  iFirstBin=iLastBin-reb;
345  }
346  for(Int_t jjj=iFirstBin; jjj<=iLastBin; jjj++){
347  normh+=hRatio->GetBinContent(jjj);
348  normReb+=1;
349  }
350  if(normReb>=1) normh/=normReb;
351  if(norml>normh){
352  norm=norml;
353  massForNorm=hRatio->GetBinCenter(binl);
354  }else{
355  norm=normh;
356  massForNorm=hRatio->GetBinCenter(binh);
357  }
358  }else if(optForNorm==4){
359  norm=0.0001;
360  for(Int_t iMassBin=1; iMassBin<hRatio->GetNbinsX(); iMassBin++){
361  Double_t bce=hRatio->GetBinCenter(iMassBin);
362  if(bce>minMass && bce<maxMass){
363  Double_t bco=0;
364  Double_t normReb=0;
365  Int_t iFirstBin=TMath::Max(iMassBin-reb/2,1);
366  Int_t iLastBin=iFirstBin+reb;
367  if(iLastBin>hRatio->GetNbinsX()){
368  iLastBin=hRatio->GetNbinsX();
369  iFirstBin=iLastBin-reb;
370  }
371  for(Int_t jjj=iFirstBin; jjj<=iLastBin; jjj++){
372  bco+=hRatio->GetBinContent(jjj);
373  normReb+=1;
374  }
375  if(normReb>=1){
376  bco/=normReb;
377  if(bco>norm){
378  norm=bco;
379  massForNorm=bce;
380  }
381  }
382  }
383  }
384  }else if(optForNorm==5){
385  hRatio->Fit("pol0","","",minMass,minMass+rangeForNorm);
386  TF1* func0=(TF1*)hRatio->GetListOfFunctions()->FindObject("pol0");
387  Double_t norml=func0->GetParameter(0);
388  hRatio->Fit("pol0","","",maxMass-rangeForNorm,maxMass);
389  func0=(TF1*)hRatio->GetListOfFunctions()->FindObject("pol0");
390  Double_t normh=func0->GetParameter(0);
391  if(norml>normh){
392  norm=norml;
393  massForNorm=maxMass-0.5*rangeForNorm;
394  }else{
395  norm=normh;
396  massForNorm=minMass+0.5*rangeForNorm;
397  }
398  }
399 
400  printf("Normalization factor = %f --> taken at mass=%f\n",norm,massForNorm);
401 
402  return norm;
403 }
404 
406 
407  if(configFileName.Length()>0){
408  if(gSystem->Exec(Form("ls -l %s > /dev/null 2>&1",configFileName.Data()))==0){
409  printf("Read configuration from file %s\n",configFileName.Data());
410  Bool_t readOK=ReadConfig(configFileName);
411  if(!readOK){
412  printf("Error in reading configuration file\n");
413  return;
414  }
415  }
416  }
417  printf("***************************************************\n");
418  printf("*** This is the configuration that will be used ***\n");
419  PrintConfig();
420  printf("*** ***\n");
421  printf("***************************************************\n");
422 
423  TString dirName=Form("PWG3_D2H_InvMass%sLowPt%s",meson.Data(),suffix.Data());
424  TString lstName=Form("coutput%s%s",meson.Data(),suffix.Data());
425  TString dirNameMC=Form("PWG3_D2H_InvMass%sLowPt%s",meson.Data(),suffixMC.Data());
426  TString lstNameMC=Form("coutput%s%s",meson.Data(),suffixMC.Data());
427 
428  if(correctForRefl) suffix.Prepend("Refl_");
429  if(fileName.Contains("FAST") && !fileName.Contains("wSDD")){
430  suffix.Prepend("FAST_");
431  }else if(!fileName.Contains("FAST") && fileName.Contains("wSDD")){
432  suffix.Prepend("wSDD_");
433  }
434  if(fileNameMC.Contains("_G3")) suffix.Append("_Geant3MC");
435  else if(fileNameMC.Contains("_G4")) suffix.Append("_Geant4MC");
436  if(smoothLS!=0) suffix.Append(Form("_smoothLS%d",smoothLS));
437  if(costhstcut<1.) suffix.Append(Form("_CosthSt%d",TMath::Nint(costhstcut*100)));
438  if(configFileName.Contains("coarse")) suffix.Append("_CoarsePt");
439  if(useEMwithLS) suffix.Append("_EMwithLS");
440 
441  if(tuneMeanOnData<0){
442  if(meson=="Dplus") massD=TDatabasePDG::Instance()->GetParticle(411)->Mass();
443  else if(meson=="Dzero") massD=TDatabasePDG::Instance()->GetParticle(421)->Mass();
444  }else{
446  }
447 
448  Int_t nBinsWithFixSig=0;
449  Int_t nBinsWithFixMean=0;
450  Int_t patSig=0;
451  Int_t patMean=0;
452  for(Int_t iPtBin=0; iPtBin<nPtBins; iPtBin++){
453  if(fixSigmaConf==1 || (fixSigmaConf==2 && fixSigma[iPtBin])){
454  nBinsWithFixSig++;
455  patSig+=1<<iPtBin;
456  }
457  if(fixMeanConf==1 || (fixMeanConf==2 && fixMean[iPtBin])){
458  nBinsWithFixMean++;
459  patMean+=1<<iPtBin;
460  }
461  }
462 
463  if(gSystem->Exec(Form("ls -l %s > /dev/null 2>&1",fileName.Data())) !=0){
464  printf("File %s with raw data results does not exist -> exiting\n",fileName.Data());
465  return;
466  }
467  TFile* fil=new TFile(fileName.Data());
468  if(!fil){
469  printf("File %s with raw data not opened -> exiting\n",fileName.Data());
470  return;
471  }
472  TDirectoryFile* df=(TDirectoryFile*)fil->Get(dirName.Data());
473  if(!df){
474  printf("TDirectoryFile %s not found in TFile\n",dirName.Data());
475  fil->ls();
476  return;
477  }
478 
479  AliNormalizationCounter *nC=(AliNormalizationCounter*)df->Get("NormalizationCounter");
480  if(!nC){
481  printf("AliNormalizationCounter object missing -> exiting\n");
482  return;
483  }
484  TH1F* hnEv=new TH1F("hEvForNorm","events for normalization",1,0,1);
485  printf("Number of Ev. for norm=%d\n",(Int_t)nC->GetNEventsForNorm());
486  hnEv->SetBinContent(1,nC->GetNEventsForNorm());
487 
488  TH1F* hRebin=new TH1F("hRebin","",nPtBins,binLims);
489  TH1F* hBkgFitFunc=new TH1F("hBkgFitFunc","",nPtBins,binLims);
490  TH1F* hBkgFitFuncSB=new TH1F("hBkgFitFuncSB","",nPtBins,binLims);
491 
492  TList* l=(TList*)df->Get(lstName.Data());
493 
494  TString var="Y";
495  if(costhstcut<1.) var="CosthSt";
496 
497  TH3F* h3d=(TH3F*)l->FindObject(Form("hMassVsPtVs%s",var.Data()));
498  if(!h3d){
499  printf("Histogram hMassVsPtVsCosthSt does not exist-> use Y\n");
500  var="Y";
501  h3d=(TH3F*)l->FindObject(Form("hMassVsPtVs%s",var.Data()));
502  if(!h3d){
503  printf("Histogram hMassVsPtVsY does not exist-> exit\n");
504  }
505  }
506  Int_t zbin1=0;
507  Int_t zbin2=h3d->GetZaxis()->GetNbins()+1;
508  if(var=="CosthSt") zbin2=h3d->GetZaxis()->FindBin(costhstcut-0.000001);
509 
510  printf("Binning in %s: %d %d\n",var.Data(),zbin1,zbin2);
511 
512  TH3F* h3dr=(TH3F*)l->FindObject(Form("hMassVsPtVs%sRot",var.Data()));
513  TH3F* h3dme=(TH3F*)l->FindObject(Form("hMassVsPtVs%sME",var.Data()));
514  TH3F* h3dmepp=(TH3F*)l->FindObject(Form("hMassVsPtVs%sMELSpp",var.Data()));
515  TH3F* h3dmemm=(TH3F*)l->FindObject(Form("hMassVsPtVs%sMELSmm",var.Data()));
516  TH2F* hpoolMix=(TH2F*)l->FindObject("hMixingsPerPool");
517  TH2F* hpoolEv=(TH2F*)l->FindObject("hEventsPerPool");
518  TH3F* h3dlsp=(TH3F*)l->FindObject(Form("hMassVsPtVs%sLSpp",var.Data()));
519  TH3F* h3dlsm=(TH3F*)l->FindObject(Form("hMassVsPtVs%sLSmm",var.Data()));
520 
521  TH3F* h3drefl=0x0;
522  TH3F* h3dmcsig=0x0;
523  TH1F* hSigmaMC=0x0;
524  if(gSystem->Exec(Form("ls -l %s > /dev/null 2>&1",fileNameMC.Data())) !=0){
525  printf("File %s with MC results does not exist -> exiting\n",fileNameMC.Data());
526  return;
527  }
528  TFile* filMC=new TFile(fileNameMC.Data());
529  if(filMC && filMC->IsOpen()){
530  TDirectoryFile* dfMC=(TDirectoryFile*)filMC->Get(dirNameMC.Data());
531  if(!dfMC){
532  printf("TDirectoryFile %s not found in TFile for MC\n",dirNameMC.Data());
533  filMC->ls();
534  return;
535  }
536  TList* lMC=(TList*)dfMC->Get(lstNameMC.Data());
537  hSigmaMC=FitMCInvMassSpectra(lMC,var);
538  if(nBinsWithFixSig>0 && !hSigmaMC){
539  printf("Fit to MC inv. mass spectra failed\n");
540  return;
541  }
542  if(correctForRefl){
543  h3drefl=(TH3F*)lMC->FindObject(Form("hMassVsPtVs%sRefl",var.Data()));
544  h3dmcsig=(TH3F*)lMC->FindObject(Form("hMassVsPtVs%sSig",var.Data()));
545  }
546  }
547 
548  TString sigConf="FixedSigma";
549  if(nBinsWithFixSig==0) sigConf="FreeSigma";
550  else if(nBinsWithFixSig==nPtBins) sigConf="FixedSigmaAll";
551  else sigConf=Form("FixedSigma%d",patSig);
552  if(nBinsWithFixSig>0 && tuneSigmaOnData>0.) sigConf+=Form("%d",TMath::Nint(tuneSigmaOnData*100.));
553  if(nBinsWithFixMean==nPtBins) sigConf+="FixedMeanAll";
554  else if(nBinsWithFixMean>0) sigConf+=Form("FixedMean%d",patMean);
555 
556  TCanvas* cem=new TCanvas("cem","Pools",1200,600);
557  cem->Divide(2,1);
558  cem->cd(1);
559  gPad->SetRightMargin(0.12);
560  gPad->SetLeftMargin(0.12);
561  hpoolMix->GetXaxis()->SetTitle("zVertex");
562  hpoolMix->GetYaxis()->SetTitle("Multiplicity");
563  hpoolMix->GetYaxis()->SetTitleOffset(1.4);
564  hpoolMix->Draw("colztext");
565  cem->cd(2);
566  gPad->SetRightMargin(0.12);
567  gPad->SetLeftMargin(0.12);
568  hpoolEv->Draw("colztext");
569  hpoolEv->GetXaxis()->SetTitle("zVertex");
570  hpoolEv->GetYaxis()->SetTitle("Multiplicity");
571  hpoolEv->GetYaxis()->SetTitleOffset(1.4);
572 
573 
574  TCanvas* c1=new TCanvas("c1","Mass",1200,800);
575  DivideCanvas(c1,nPtBins);
576 
577  TCanvas* c2=new TCanvas("c2","Mass-Bkg Rot",1200,800);
578  DivideCanvas(c2,nPtBins);
579  TCanvas* c2sub=new TCanvas("c2sub","Mass-Bkg-Fit Rot",1200,800);
580  DivideCanvas(c2sub,nPtBins);
581  TCanvas* c2pulls=new TCanvas("c2pulls","Mass-Bkg Rot pulls",1200,800);
582  DivideCanvas(c2pulls,nPtBins);
583  TCanvas* c2residuals=new TCanvas("c2residuals","Mass-Bkg Rot residuals",1200,800);
584  DivideCanvas(c2residuals,nPtBins);
585  TCanvas* c2residualTrend=new TCanvas("c2residualTrend","Mass-Bkg Rot residual trend vs. mass",1200,800);
586  DivideCanvas(c2residualTrend,nPtBins);
587  TCanvas* c2pullTrend=new TCanvas("c2pullTrend","Mass-Bkg Rot pull trend vs. mass",1200,800);
588  DivideCanvas(c2pullTrend,nPtBins);
589 
590  TCanvas* c3=new TCanvas("c3","Mass-Bkg LS",1200,800);
591  DivideCanvas(c3,nPtBins);
592  TCanvas* c3sub=new TCanvas("c3sub","Mass-Bkg-Fit LS",1200,800);
593  DivideCanvas(c3sub,nPtBins);
594  TCanvas* c3pulls=new TCanvas("c3pulls","Mass-Bkg LS pulls",1200,800);
595  DivideCanvas(c3pulls,nPtBins);
596  TCanvas* c3residuals=new TCanvas("c3residuals","Mass-Bkg LS residuals",1200,800);
597  DivideCanvas(c3residuals,nPtBins);
598  TCanvas* c3residualTrend=new TCanvas("c3residualTrend","Mass-Bkg LS residual trend vs. mass",1200,800);
599  DivideCanvas(c3residualTrend,nPtBins);
600  TCanvas* c3pullTrend=new TCanvas("c3pullTrend","Mass-Bkg LS pull trend vs. mass",1200,800);
601  DivideCanvas(c3pullTrend,nPtBins);
602 
603  TCanvas* c4=new TCanvas("c4","Mass-Bkg EM",1200,800);
604  DivideCanvas(c4,nPtBins);
605  TCanvas* c4sub=new TCanvas("c4sub","Mass-Bkg-Fit EM",1200,800);
606  DivideCanvas(c4sub,nPtBins);
607  TCanvas* c4pulls=new TCanvas("c4pulls","Mass-Bkg EM pulls",1200,800);
608  DivideCanvas(c4pulls,nPtBins);
609  TCanvas* c4residuals=new TCanvas("c4residuals","Mass-Bkg EM residuals",1200,800);
610  DivideCanvas(c4residuals,nPtBins);
611  TCanvas* c4residualTrend=new TCanvas("c4residualTrend","Mass-Bkg EM residual trend vs. mass",1200,800);
612  DivideCanvas(c4residualTrend,nPtBins);
613  TCanvas* c4pullTrend=new TCanvas("c4pullTrend","Mass-Bkg EM pull trend vs. mass",1200,800);
614  DivideCanvas(c4pullTrend,nPtBins);
615 
616  TCanvas *c5=0x0;
617  TCanvas *c5sub=0x0;
618  TCanvas *c5pulls=0x0;
619  TCanvas *c5residuals=0x0;
620  TCanvas *c5residualTrend=0x0;
621  TCanvas *c5pullTrend=0x0;
622 
623  if(tryDirectFit){
624  c5=new TCanvas("c5","Mass SB Fit",1200,800);
625  DivideCanvas(c5,nPtBins);
626  c5sub=new TCanvas("c5sub","Mass-Bkg SB",1200,800);
627  DivideCanvas(c5sub,nPtBins);
628  c5pulls=new TCanvas("c5pulls","Mass-Bkg SB pulls",1200,800);
629  DivideCanvas(c5pulls,nPtBins);
630  c5residuals=new TCanvas("c5residuals","Mass-Bkg SB residuals",1200,800);
631  DivideCanvas(c5residuals,nPtBins);
632  c5residualTrend=new TCanvas("c5residualTrend","Mass-Bkg SB residual trend vs. mass",1200,800);
633  DivideCanvas(c5residualTrend,nPtBins);
634  c5pullTrend=new TCanvas("c5pullTrend","Mass-Bkg SB pull trend vs. mass",1200,800);
635  DivideCanvas(c5pullTrend,nPtBins);
636  }
637 
638  TCanvas *cCompareResidualTrends=new TCanvas("cCompareResidualTrends","cCompareResidualTrends",1200,800);
639  DivideCanvas(cCompareResidualTrends,nPtBins);
640 
641 
642  AliHFInvMassFitter* fitterRot[nPtBins];
643  AliHFInvMassFitter* fitterLS[nPtBins];
644  AliHFInvMassFitter* fitterME[nPtBins];
645  AliHFInvMassFitter* fitterSB[nPtBins];
646 
647  TH1F* hRawYieldRot=new TH1F("hRawYieldRot"," ; p_{T} (GeV/c) ; Raw yield",nPtBins,binLims);
648  TH1F* hRawYieldLS=new TH1F("hRawYieldLS"," ; p_{T} (GeV/c) ; Raw yield",nPtBins,binLims);
649  TH1F* hRawYieldME=new TH1F("hRawYieldME"," ; p_{T} (GeV/c) ; Raw yield",nPtBins,binLims);
650  TH1F* hRawYieldSB=new TH1F("hRawYieldSB"," ; p_{T} (GeV/c) ; Raw yield",nPtBins,binLims);
651 
652  TH1F* hRelStatRot=new TH1F("hRelStatRot"," ; p_{T} (GeV/c) ; Relative stat. unc.",nPtBins,binLims);
653  TH1F* hRelStatLS=new TH1F("hRelStatLS"," ; p_{T} (GeV/c) ; Relative stat. unc.",nPtBins,binLims);
654  TH1F* hRelStatME=new TH1F("hRelStatME"," ; p_{T} (GeV/c) ; Relative stat. unc.",nPtBins,binLims);
655  TH1F* hRelStatSB=new TH1F("hRelStatSB"," ; p_{T} (GeV/c) ; Relative stat. unc.",nPtBins,binLims);
656 
657  TH1F* hSignifRot=new TH1F("hSignifRot","",nPtBins,binLims);
658  TH1F* hSignifLS=new TH1F("hSignifLS","",nPtBins,binLims);
659  TH1F* hSignifME=new TH1F("hSignifME","",nPtBins,binLims);
660  TH1F* hSignifSB=new TH1F("hSignifSB","",nPtBins,binLims);
661 
662  TH1F* hSoverBRot=new TH1F("hSoverBRot","",nPtBins,binLims);
663  TH1F* hSoverBLS=new TH1F("hSoverBLS","",nPtBins,binLims);
664  TH1F* hSoverBME=new TH1F("hSoverBME","",nPtBins,binLims);
665  TH1F* hSoverBSB=new TH1F("hSoverBSB","",nPtBins,binLims);
666 
667  TH1F* hGausMeanRot=new TH1F("hGausMeanRot","",nPtBins,binLims);
668  TH1F* hGausMeanLS=new TH1F("hGausMeanLS","",nPtBins,binLims);
669  TH1F* hGausMeanME=new TH1F("hGausMeanME","",nPtBins,binLims);
670  TH1F* hGausMeanSB=new TH1F("hGausMeanSB","",nPtBins,binLims);
671 
672  TH1F* hGausSigmaRot=new TH1F("hGausSigmaRot","",nPtBins,binLims);
673  TH1F* hGausSigmaLS=new TH1F("hGausSigmaLS","",nPtBins,binLims);
674  TH1F* hGausSigmaME=new TH1F("hGausSigmaME","",nPtBins,binLims);
675  TH1F* hGausSigmaSB=new TH1F("hGausSigmaSB","",nPtBins,binLims);
676 
677  TH1F* hChiSqRot=new TH1F("hChiSqRot","",nPtBins,binLims);
678  TH1F* hChiSqLS=new TH1F("hChiSqLS","",nPtBins,binLims);
679  TH1F* hChiSqME=new TH1F("hChiSqME","",nPtBins,binLims);
680  TH1F* hChiSqSB=new TH1F("hChiSqSB","",nPtBins,binLims);
681 
682  TH1F* hNdfRot=new TH1F("hNdfRot","",nPtBins,binLims);
683  TH1F* hNdfLS=new TH1F("hNdfLS","",nPtBins,binLims);
684  TH1F* hNdfME=new TH1F("hNdfME","",nPtBins,binLims);
685  TH1F* hNdfSB=new TH1F("hNdfSB","",nPtBins,binLims);
686 
687  TH1F* hRawYieldRotBC=new TH1F("hRawYieldRotBC","BC yield (rotational background)",nPtBins,binLims);
688  TH1F* hRawYieldLSBC=new TH1F("hRawYieldLSBC","BC yield (like-sign background)",nPtBins,binLims);
689  TH1F* hRawYieldMEBC=new TH1F("hRawYieldMEBC","BC yield (mixed-event background)",nPtBins,binLims);
690  TH1F* hRawYieldSBBC=new TH1F("hRawYieldSBBC","BC yield (side-band fit background)",nPtBins,binLims);
691 
692  TH1F* hInvMassHistoBinWidthRot=new TH1F("hInvMassHistoBinWidthRot","",nPtBins,binLims);
693  TH1F* hInvMassHistoBinWidthLS=new TH1F("hInvMassHistoBinWidthLS","",nPtBins,binLims);
694  TH1F* hInvMassHistoBinWidthME=new TH1F("hInvMassHistoBinWidthME","",nPtBins,binLims);
695  TH1F* hInvMassHistoBinWidthSB=new TH1F("hInvMassHistoBinWidthSB","",nPtBins,binLims);
696 
697  TLatex* tME=new TLatex(0.65,0.82,"MixEv +- pairs");
698  tME->SetNDC();
699  TLatex* tMEpp=new TLatex(0.65,0.75,"MixEv ++ pairs");
700  tMEpp->SetNDC();
701  TLatex* tMEmm=new TLatex(0.65,0.68,"MixEv -- pairs");
702  tMEmm->SetNDC();
703 
704  // TF1 *fpeak=new TF1("fpeak","[0]*1./(TMath::Sqrt(2.*TMath::Pi())*[2])*TMath::Exp(-(x-[1])*(x-[1])/(2.*[2]*[2]))",minMass,maxMass);
705 
706 
707  TDirectory *current = gDirectory;
708  TFile* fout=new TFile(Form("outputMassFits_%s_%s.root",sigConf.Data(),suffix.Data()),"recreate");
709  current->cd();
710 
711 
712  for(Int_t iPtBin=0; iPtBin<nPtBins; iPtBin++){
713 
714  printf("\n---------- pt interval %d (%.1f-%.1f)\n",iPtBin,binLims[iPtBin],binLims[iPtBin+1]);
715  minMass=minMass4Fit[iPtBin];
716  maxMass=maxMass4Fit[iPtBin];
717  Int_t bin1=h3d->GetYaxis()->FindBin(binLims[iPtBin]);
718  Int_t bin2=h3d->GetYaxis()->FindBin(binLims[iPtBin+1]-0.0001);
719  printf("Bin %d Pt range=%f %f\n",iPtBin,h3d->GetYaxis()->GetBinLowEdge(bin1),h3d->GetYaxis()->GetBinUpEdge(bin2));
720  TH1D* hMassPtBin=h3d->ProjectionX(Form("hMassPtBin%d",iPtBin),bin1,bin2,zbin1,zbin2);
721  TH1D* hMassPtBinr=h3dr->ProjectionX(Form("hMassPtBinr%d",iPtBin),bin1,bin2,zbin1,zbin2);
722  TH1D* hMassPtBinme=h3dme->ProjectionX(Form("hMassPtBinme%d",iPtBin),bin1,bin2,zbin1,zbin2);
723  TH1D* hMassPtBinmeLSpp=h3dmepp->ProjectionX(Form("hMassPtBinmeLSpp%d",iPtBin),bin1,bin2,zbin1,zbin2);
724  TH1D* hMassPtBinmeLSmm=h3dmemm->ProjectionX(Form("hMassPtBinmeLSmm%d",iPtBin),bin1,bin2,zbin1,zbin2);
725 
726  if(correctForRefl){
727  Int_t bin1MC=h3drefl->GetYaxis()->FindBin(binLims[iPtBin]);
728  Int_t bin2MC=h3drefl->GetYaxis()->FindBin(binLims[iPtBin+1]-0.0001);
729  hMCReflPtBin=h3drefl->ProjectionX(Form("hMCReflPtBin%d",iPtBin),bin1MC,bin2MC,zbin1,zbin2);
730  hMCSigPtBin=h3dmcsig->ProjectionX(Form("hMCSigPtBin%d",iPtBin),bin1MC,bin2MC,zbin1,zbin2);
731  }
732 
733  TH1D* hMassPtBinlsp=0x0;
734  TH1D* hMassPtBinlsm=0x0;
735  TH1D* hMassPtBinls=0x0;
736  if(h3dlsp){
737  hMassPtBinlsp=h3dlsp->ProjectionX(Form("hMassPtBinlsp%d",iPtBin),bin1,bin2,zbin1,zbin2);
738  hMassPtBinlsm=h3dlsm->ProjectionX(Form("hMassPtBinlsm%d",iPtBin),bin1,bin2,zbin1,zbin2);
739  hMassPtBinls=(TH1D*)hMassPtBinlsp->Clone(Form("hMassPtBinls%d",iPtBin));
740  hMassPtBinls->Reset();
741  for(Int_t iBin=1; iBin<=hMassPtBinlsp->GetNbinsX(); iBin++){
742  Double_t np=hMassPtBinlsp->GetBinContent(iBin);
743  Double_t nm=hMassPtBinlsm->GetBinContent(iBin);
744  Double_t enp=hMassPtBinlsp->GetBinError(iBin);
745  Double_t enm=hMassPtBinlsm->GetBinError(iBin);
746  Double_t tt=0;
747  Double_t ett=0;
748  if(useGeomMeanLS){
749  tt=2*TMath::Sqrt(np*nm);
750  if(tt>0) ett=2./tt*TMath::Sqrt(np*np*enm*enm+nm*nm*enp*enp);
751  }else{
752  tt=0.5*(np+nm);
753  ett=0.5*TMath::Sqrt(enm*enm+enp*enp);
754  }
755  hMassPtBinls->SetBinContent(iBin,tt);
756  hMassPtBinls->SetBinError(iBin,ett);
757  }
758  if(smoothLS<-0.5)hMassPtBinls->Smooth(-1*smoothLS);
759  else if(smoothLS>0.5)QuadraticSmooth(hMassPtBinls,smoothLS);
760  hMassPtBinls->SetLineColor(kGreen+1);
761  }
762  // hMassPtBin->Sumw2();
763  // hMassPtBinr->Sumw2();
764  // hMassPtBinme->Sumw2();
765  // hMassPtBinlsp->Sumw2();
766  // hMassPtBinls->Sumw2();
767  hMassPtBin->SetLineColor(1);
768  hMassPtBinr->SetLineColor(4);
769  hMassPtBinme->SetLineColor(kOrange+1);
770  hMassPtBinmeLSpp->SetLineColor(kGreen+2);
771  hMassPtBinmeLSmm->SetLineColor(kCyan);
772  hMassPtBin->SetTitle(Form("%.1f<p_{T}<%.1f GeV/c",binLims[iPtBin],binLims[iPtBin+1]));
773  hMassPtBinr->SetTitle(Form("%.1f<p_{T}<%.1f GeV/c",binLims[iPtBin],binLims[iPtBin+1]));
774  hMassPtBinme->SetTitle(Form("%.1f<p_{T}<%.1f GeV/c",binLims[iPtBin],binLims[iPtBin+1]));
775  TH1D* hRatioMEpp=(TH1D*)hMassPtBinmeLSpp->Clone(Form("hRatioPtBinmeLSpp%d",iPtBin));
776  hRatioMEpp->Divide(hMassPtBinme);
777  TH1D* hRatioMEmm=(TH1D*)hMassPtBinmeLSmm->Clone(Form("hRatioPtBinmeLSmm%d",iPtBin));
778  hRatioMEmm->Divide(hMassPtBinme);
779  TH1D* hMassPtBinmeAll=(TH1D*)hMassPtBinme->Clone(Form("hRatioPtBinmeAll%d",iPtBin));
780  hMassPtBinmeAll->Add(hMassPtBinmeLSpp);
781  hMassPtBinmeAll->Add(hMassPtBinmeLSmm);
782  hMassPtBinmeAll->SetLineColor(kRed);
783  TH1D* hRatioME=(TH1D*)hMassPtBinme->Clone(Form("hRatioME%d",iPtBin));
784  hRatioME->Divide(hMassPtBin);
785  TH1D* hRatioMEAll=(TH1D*)hMassPtBinmeAll->Clone(Form("hRatioMEAll%d",iPtBin));
786  hRatioMEAll->Divide(hMassPtBin);
787 
788 
789  TCanvas* c0=new TCanvas(Form("CBin%d",iPtBin),Form("Bin%d norm",iPtBin),1300,700);
790  c0->Divide(3,2);
791  c0->cd(1);
792  TH1D* hRatio=(TH1D*)hMassPtBinr->Clone(Form("hRatioFormNorm%d",iPtBin));
793  hRatio->Divide(hMassPtBin);
794  hRatio->Draw();
795  hRatio->GetYaxis()->SetTitle("Rotational/All");
796  hRatio->GetXaxis()->SetTitle("Invariant mass (GeV/c^{2})");
797  Double_t normRot=GetBackgroundNormalizationFactor(hRatio,rebin[iPtBin]);
798  hMassPtBinr->Scale(1./normRot);
799  TLatex* tnr=new TLatex(0.2,0.2,Form("Normaliz. factor = %f",normRot));
800  tnr->SetNDC();
801  tnr->Draw();
802  c0->cd(2);
803  TH1D* hRatioLS=(TH1D*)hMassPtBinls->Clone(Form("hRatioFormNormLS%d",iPtBin));
804  hRatioLS->Divide(hMassPtBin);
805  hRatioLS->Draw();
806  hRatioLS->GetYaxis()->SetTitle("LS/All");
807  hRatioLS->GetXaxis()->SetTitle("Invariant mass (GeV/c^{2})");
808  if(renormLS){
809  Double_t normLS=GetBackgroundNormalizationFactor(hRatioLS,rebin[iPtBin]);
810  hMassPtBinls->Scale(1./normLS);
811  TLatex* tnl=new TLatex(0.2,0.2,Form("Normaliz. factor = %f",normLS));
812  tnl->SetNDC();
813  tnl->Draw();
814  }
815  c0->cd(3);
816  hMassPtBinme->GetYaxis()->SetTitle("Entries (EvMix)");
817  hMassPtBinme->GetXaxis()->SetTitle("Invariant mass (GeV/c^{2})");
818  hMassPtBinme->SetMinimum(0);
819  hMassPtBinme->DrawCopy();
820  hMassPtBinmeLSpp->Draw("same");
821  hMassPtBinmeLSmm->Draw("same");
822  tME->SetTextColor(hMassPtBinme->GetLineColor());
823  tMEpp->SetTextColor(hMassPtBinmeLSpp->GetLineColor());
824  tMEmm->SetTextColor(hMassPtBinmeLSmm->GetLineColor());
825  tME->Draw();
826  tMEpp->Draw();
827  tMEmm->Draw();
828  c0->cd(4);
829  hRatioME->Draw();
830  hRatioME->GetYaxis()->SetTitle("EvMix (+-)/All");
831  hRatioME->GetXaxis()->SetTitle("Invariant mass (GeV/c^{2})");
832  c0->cd(5);
833  hRatioMEAll->Draw();
834  hRatioMEAll->GetYaxis()->SetTitle("EvMix (+-,++,--)/All");
835  hRatioMEAll->GetXaxis()->SetTitle("Invariant mass (GeV/c^{2})");
836  c0->cd(6);
837  hRatioMEpp->Draw();
838  hRatioMEpp->SetMinimum(0.49);
839  hRatioMEpp->SetMaximum(0.51);
840  hRatioMEpp->GetYaxis()->SetTitle("ME with LS / ME with OS");
841  hRatioMEpp->GetXaxis()->SetTitle("Invariant mass (GeV/c^{2})");
842  hRatioMEmm->Draw("same");
843 
844  Double_t normME=GetBackgroundNormalizationFactor(hRatioME,rebin[iPtBin]);
845  Double_t normMEAll=GetBackgroundNormalizationFactor(hRatioMEAll,rebin[iPtBin]);
846  hMassPtBinme->Scale(1./normME);
847  hMassPtBinmeAll->Scale(1./normMEAll);
848  printf("Background norm bin %d DONE\n",iPtBin);
849 
850  c1->cd(iPtBin+1);
851  hMassPtBin->GetXaxis()->SetRangeUser(minMass,maxMass);
852  hMassPtBin->SetMinimum(0);
853  hMassPtBin->Draw();
854  hMassPtBin->GetYaxis()->SetTitle("Counts");
855  hMassPtBin->GetYaxis()->SetTitleOffset(2.);
856  hMassPtBin->GetXaxis()->SetTitle("Invariant mass (GeV/c^{2})");
857  hMassPtBinr->Draw("same");
858  if(!useEMwithLS) hMassPtBinme->Draw("same");
859  else hMassPtBinmeAll->Draw("same");
860  if(hMassPtBinls) hMassPtBinls->Draw("same");
861  if(iPtBin==0){
862  TLegend* leg=new TLegend(0.5,0.6,0.89,0.89);
863  leg->SetFillStyle(0);
864  leg->AddEntry(hMassPtBin,"All candidates","L")->SetTextColor(hMassPtBin->GetLineColor());
865  leg->AddEntry(hMassPtBinr,"Background (rotations)","L")->SetTextColor(hMassPtBinr->GetLineColor());
866  if(!useEMwithLS) leg->AddEntry(hMassPtBinme,"Background (ME)","L")->SetTextColor(hMassPtBinme->GetLineColor());
867  else leg->AddEntry(hMassPtBinmeAll,"Background (ME)","L")->SetTextColor(hMassPtBinmeAll->GetLineColor());
868  if(hMassPtBinls) leg->AddEntry(hMassPtBinls,"Like-sign","L")->SetTextColor(hMassPtBinls->GetLineColor());
869  leg->Draw();
870  }
871  gPad->Update();
872 
873  TH1D* hMassSubRot=(TH1D*)hMassPtBin->Clone(Form("hMassSubRot_bin%d",iPtBin));
874  hMassSubRot->Add(hMassPtBinr,-1);
875  hMassSubRot->SetTitle(Form("%.1f<p_{T}<%.1f GeV/c -- Rotational",binLims[iPtBin],binLims[iPtBin+1]));
876  TH1D* hMassSubME=(TH1D*)hMassPtBin->Clone(Form("hMassSubME_bin%d",iPtBin));
877  if(!useEMwithLS) hMassSubME->Add(hMassPtBinme,-1);
878  else hMassSubME->Add(hMassPtBinmeAll,-1);
879  hMassSubME->SetTitle(Form("%.1f<p_{T}<%.1f GeV/c -- Mixed Ev",binLims[iPtBin],binLims[iPtBin+1]));
880  TH1D* hMassSubLS=0x0;
881  if(hMassPtBinls){
882  hMassSubLS=(TH1D*)hMassPtBin->Clone(Form("hMassSubLS_bin%d",iPtBin));
883  hMassSubLS->Add(hMassPtBinls,-1);
884  hMassSubLS->SetTitle(Form("%.1f<p_{T}<%.1f GeV/c -- Like Sign",binLims[iPtBin],binLims[iPtBin+1]));
885  }
886 
887  fout->cd();
888  hMassPtBin->Write();
889  hMassSubRot->Write();
890  hMassSubME->Write();
891  if(hMassPtBinls) hMassSubLS->Write();
892  if(hMCReflPtBin) hMCReflPtBin->Write();
893  if(hMCSigPtBin) hMCSigPtBin->Write();
894  current->cd();
895 
896  hMassSubRot=AliVertexingHFUtils::RebinHisto(hMassSubRot,rebin[iPtBin]);
897  hMassSubME=AliVertexingHFUtils::RebinHisto(hMassSubME,rebin[iPtBin]);
898  if(hMassPtBinls) hMassSubLS=AliVertexingHFUtils::RebinHisto(hMassSubLS,rebin[iPtBin]);
899 
900  hRebin->SetBinContent(iPtBin+1,rebin[iPtBin]);
901  Int_t bkgToFill=typeb;
902  if(typeb==6) bkgToFill=typeb+nDegreeBackPol[iPtBin];
903  Int_t bkgToFillSB=6+nDegreeBackPol[iPtBin];
904 
905  hBkgFitFunc->SetBinContent(iPtBin+1,bkgToFill);
906  hBkgFitFuncSB->SetBinContent(iPtBin+1,bkgToFillSB);
907 
908  fitterRot[iPtBin]=ConfigureFitter(hMassSubRot,iPtBin,typeb,minMass,maxMass);
909  if(hMassPtBinls) fitterLS[iPtBin]=ConfigureFitter(hMassSubLS,iPtBin,typeb,minMass,maxMass);
910  fitterME[iPtBin]=ConfigureFitter(hMassSubME,iPtBin,typeb,minMass,maxMass);
911 
912  Bool_t out1=fitterRot[iPtBin]->MassFitter(0);
913  Bool_t out2=kFALSE;
914  if(hMassPtBinls) out2=fitterLS[iPtBin]->MassFitter(0);
915  Bool_t out3=fitterME[iPtBin]->MassFitter(0);
916 
917  Double_t background=999999999.;
918  Bool_t out4=kFALSE;
919  if(tryDirectFit){
920  TH1D *hMassDirectFit=(TH1D*)hMassPtBin->Clone(Form("hMassDirectFit_bin%d",iPtBin));
921  hMassDirectFit=AliVertexingHFUtils::RebinHisto(hMassDirectFit,rebin[iPtBin]);
922  fitterSB[iPtBin]=ConfigureFitter(hMassDirectFit,iPtBin,6,fitSBrangelow[iPtBin],fitSBrangeup[iPtBin],kTRUE);
923  out4=fitterSB[iPtBin]->MassFitter(0);//DirectFit(hMassDirectFit,iPtBin,hRawYieldSB);
924 
925  Double_t background,ebkg;
926 
927  if(out4 && fitterSB[iPtBin]->GetMassFunc()){
928  c5->cd(iPtBin+1);
929  fitterSB[iPtBin]->DrawHere(gPad,3,0);
930  if(fitterSB[iPtBin]->GetRawYield()>0 && fitterSB[iPtBin]->GetReducedChiSquare()>0 && fitterSB[iPtBin]->GetReducedChiSquare()<5){
931  hRawYieldSB->SetBinContent(iPtBin+1,fitterSB[iPtBin]->GetRawYield());
932  hRawYieldSB->SetBinError(iPtBin+1,fitterSB[iPtBin]->GetRawYieldError());
933  hRelStatSB->SetBinContent(iPtBin+1,fitterSB[iPtBin]->GetRawYieldError()/fitterSB[iPtBin]->GetRawYield());
934  hRelStatSB->SetBinError(iPtBin+1,0.00000001);
935  fitterSB[iPtBin]->Background(3.,background,ebkg);
936  hSignifSB->SetBinContent(iPtBin+1,fitterSB[iPtBin]->GetRawYield()/TMath::Sqrt(background+fitterSB[iPtBin]->GetRawYield()));
937  hSignifSB->SetBinError(iPtBin+1,0.00000001);
938  hSoverBSB->SetBinContent(iPtBin+1,fitterSB[iPtBin]->GetRawYield()/background);
939  hSoverBSB->SetBinError(iPtBin+1,0.00000001);
940  hGausMeanSB->SetBinContent(iPtBin+1,fitterSB[iPtBin]->GetMean());
941  hGausMeanSB->SetBinError(iPtBin+1,fitterSB[iPtBin]->GetMeanUncertainty());
942  hGausSigmaSB->SetBinContent(iPtBin+1,fitterSB[iPtBin]->GetSigma());
943  hGausSigmaSB->SetBinError(iPtBin+1,fitterSB[iPtBin]->GetSigmaUncertainty());
944  hChiSqSB->SetBinContent(iPtBin+1,fitterSB[iPtBin]->GetReducedChiSquare());
945  hChiSqSB->SetBinError(iPtBin+1,0.00001); // very small number, for graphics
946  hNdfSB->SetBinContent(iPtBin+1,fitterSB[iPtBin]->GetMassFunc()->GetNDF());
947  hNdfSB->SetBinError(iPtBin+1,0.00001); // very small number, for graphics
948  }
949  hInvMassHistoBinWidthSB->SetBinContent(iPtBin+1,hMassDirectFit->GetBinWidth(1));
950  // if(!correctForRefl)
951  WriteFitInfo(fitterSB[iPtBin],hMassPtBin);
952  fout->cd();
953  WriteFitFunctionsToFile(fitterSB[iPtBin],"SB",iPtBin);
954  current->cd();
955 
956  c5sub->cd(iPtBin+1);
957  TH1F* hsubTemp=(TH1F*)hMassDirectFit->Clone(Form("%sSubBack%d",hMassDirectFit->GetName(),iPtBin));
958  TH1F* hsubTempAllRange=(TH1F*)hMassDirectFit->Clone(Form("%sSubBackAllRange%d",hMassDirectFit->GetName(),iPtBin));
959  TF1* funcAll=fitterSB[iPtBin]->GetMassFunc();
960  TF1* funcBkg=fitterSB[iPtBin]->GetBackgroundRecalcFunc();
961  // TF1* funcBkg2=fitterSB[iPtBin]->GetBackgroundFullRangeFunc();
962  for(Int_t jst=1;jst<=hsubTemp->GetNbinsX();jst++){
963  Double_t backg=funcBkg->Integral(hsubTemp->GetBinLowEdge(jst),hsubTemp->GetBinLowEdge(jst)+hsubTemp->GetBinWidth(jst))/hsubTemp->GetBinWidth(jst);
964  Double_t tot=funcAll->Integral(hsubTempAllRange->GetBinLowEdge(jst),hsubTempAllRange->GetBinLowEdge(jst)+hsubTempAllRange->GetBinWidth(jst))/hsubTempAllRange->GetBinWidth(jst);
965  hsubTemp->SetBinContent(jst,hsubTemp->GetBinContent(jst)-backg);
966  hsubTempAllRange->SetBinContent(jst,hsubTempAllRange->GetBinContent(jst)-tot);
967  }
968  hsubTemp->SetLineColor(kBlue);
969  hsubTempAllRange->SetLineColor(kGray+2);
970  Double_t ymin=0;
971  Double_t ymax=1;
972  for(Int_t ibs=1; ibs<hsubTemp->GetNbinsX(); ibs++){
973  Double_t binc=hsubTemp->GetBinCenter(ibs);
974  if(binc>fitSBrangelow[iPtBin] && binc<fitSBrangeup[iPtBin]){
975  Double_t yl=hsubTemp->GetBinContent(ibs)-hsubTemp->GetBinError(ibs);
976  Double_t yu=hsubTemp->GetBinContent(ibs)+hsubTemp->GetBinError(ibs);
977  if(yl<ymin) ymin=yl;
978  if(yu>ymax) ymax=yu;
979  }
980  }
981  if(ymax>0) ymax*=1.2;
982  else ymax*=0.8;
983  if(ymin<0) ymin*=1.2;
984  else ymin*=0.8;
985  hsubTemp->GetXaxis()->SetRangeUser(fitSBrangelow[iPtBin],fitSBrangeup[iPtBin]);
986  hsubTemp->SetMinimum(ymin);
987  hsubTemp->SetMaximum(ymax);
988  hsubTemp->SetMarkerStyle(20);
989  hsubTemp->SetMarkerColor(hsubTemp->GetLineColor());
990  hsubTemp->DrawCopy();
991  hsubTempAllRange->DrawCopy("same");
992  hsubTemp->DrawCopy("same");
993  TF1 *fpeak=fitterSB[iPtBin]->GetSignalFunc();
994  fpeak->DrawCopy("same");
995 
996  Double_t errbc;
998  hRawYieldSBBC->SetBinContent(iPtBin+1,bc);
999  hRawYieldSBBC->SetBinError(iPtBin+1,errbc);
1000 
1001  c5pulls->cd(iPtBin+1);
1002  TH1F *hPulls=new TH1F();//the name is changed afterwards, histo must not be deleted
1003  TH1F *hPullsTrend=new TH1F();// the name is changed afterwards, histo must not be deleted
1004  TH1F *hResidualTrend=new TH1F();// the name is changed afterwards, histo must not be deleted
1005  TH1F *hResiduals=fitterSB[iPtBin]->GetResidualsAndPulls(hPulls,hResidualTrend,hPullsTrend);
1006  hResiduals->SetName(Form("hResidualsSB_PtBin%d",iPtBin));
1007 
1008  hPulls->Draw();
1009  PrintGausParams(hPulls);
1010 
1011  c5residuals->cd(iPtBin+1);
1012  hResiduals->Draw();
1013 
1014  SetStyleHisto(hResidualTrend,kSB);
1015  cCompareResidualTrends->cd(iPtBin+1);
1016  hResidualTrend->Draw();
1017 
1018  c5residualTrend->cd(iPtBin+1);
1019  hResidualTrend->Draw();
1020 
1021  c5pullTrend->cd(iPtBin+1);
1022  hPullsTrend->Draw();
1023  //TPaveStats *tp=(TPaveStats*)hPulls->FindObject("stats");
1024  // tp->SetOptStat("remnpcev");
1025  delete hMassDirectFit;
1026  }
1027  }
1028 
1029  c2->cd(iPtBin+1);
1030  if(out1 && fitterRot[iPtBin]->GetMassFunc()){
1031  fitterRot[iPtBin]->DrawHere(gPad,3,0);
1032  gPad->Update();
1033  if(fitterRot[iPtBin]->GetRawYield()>0 && fitterRot[iPtBin]->GetReducedChiSquare()>0 && fitterRot[iPtBin]->GetReducedChiSquare()<5){
1034  hRawYieldRot->SetBinContent(iPtBin+1,fitterRot[iPtBin]->GetRawYield());
1035  hRawYieldRot->SetBinError(iPtBin+1,fitterRot[iPtBin]->GetRawYieldError());
1036  Double_t minBinBkg=hMassPtBin->FindBin(fitterRot[iPtBin]->GetMean()-3.*fitterRot[iPtBin]->GetSigma());
1037  Double_t maxBinBkg=hMassPtBin->FindBin(fitterRot[iPtBin]->GetMean()+3.*fitterRot[iPtBin]->GetSigma());
1038  background=hMassPtBin->Integral(minBinBkg,maxBinBkg);
1039  hSoverBRot->SetBinContent(iPtBin+1,fitterRot[iPtBin]->GetRawYield()/background);
1040  hSoverBRot->SetBinError(iPtBin+1,0.000001);
1041  hRelStatRot->SetBinContent(iPtBin+1,fitterRot[iPtBin]->GetRawYieldError()/fitterRot[iPtBin]->GetRawYield());
1042  hRelStatRot->SetBinError(iPtBin+1,0.000001);
1043  hSignifRot->SetBinContent(iPtBin+1,fitterRot[iPtBin]->GetRawYield()/TMath::Sqrt(background+fitterRot[iPtBin]->GetRawYield()));
1044  hSignifRot->SetBinError(iPtBin+1,0.00000001);
1045  hGausMeanRot->SetBinContent(iPtBin+1,fitterRot[iPtBin]->GetMean());
1046  hGausMeanRot->SetBinError(iPtBin+1,fitterRot[iPtBin]->GetMeanUncertainty());
1047  hGausSigmaRot->SetBinContent(iPtBin+1,fitterRot[iPtBin]->GetSigma());
1048  hGausSigmaRot->SetBinError(iPtBin+1,fitterRot[iPtBin]->GetSigmaUncertainty());
1049  hChiSqRot->SetBinContent(iPtBin+1,fitterRot[iPtBin]->GetReducedChiSquare());
1050  hChiSqRot->SetBinError(iPtBin+1,0.00001); // very small number, for graphics
1051  hNdfRot->SetBinContent(iPtBin+1,fitterRot[iPtBin]->GetMassFunc()->GetNDF());
1052  hNdfRot->SetBinError(iPtBin+1,0.00001); // very small number, for graphics
1053  }
1054  hInvMassHistoBinWidthRot->SetBinContent(iPtBin+1,hMassSubRot->GetBinWidth(1));
1055  // if(!correctForRefl)
1056  WriteFitInfo(fitterRot[iPtBin],hMassPtBin);
1057  fout->cd();
1058  WriteFitFunctionsToFile(fitterRot[iPtBin],"Rot",iPtBin);
1059  current->cd();
1060 
1061  Double_t errbc;
1062  Double_t bc=fitterRot[iPtBin]->GetRawYieldBinCounting(errbc,nsigmaBinCounting,optBkgBinCount);
1063  hRawYieldRotBC->SetBinContent(iPtBin+1,bc);
1064  hRawYieldRotBC->SetBinError(iPtBin+1,errbc);
1065 
1066  c2sub->cd(iPtBin+1);
1067  fitterRot[iPtBin]->DrawHistoMinusFit(gPad);
1068 
1069  c2pulls->cd(iPtBin+1);
1070  TH1F *hPullsTrend=new TH1F();// the name is changed afterwards, histo must not be deleted
1071  TH1F *hPulls=new TH1F();//the name is changed afterwards, histo must not be deleted
1072  TH1F *hResidualTrend=new TH1F();// the name is changed afterwards, histo must not be deleted
1073  TH1F *hResiduals=fitterRot[iPtBin]->GetResidualsAndPulls(hPulls,hResidualTrend,hPullsTrend);
1074  hResiduals->SetName(Form("hResidualsRot_PtBin%d",iPtBin));
1075 
1076  hPulls->Draw();
1077  PrintGausParams(hPulls);
1078 
1079  c2residuals->cd(iPtBin+1);
1080  hResiduals->Draw();
1081 
1082  SetStyleHisto(hResidualTrend,kRot);
1083  c2residualTrend->cd(iPtBin+1);
1084  hResidualTrend->Draw();
1085  cCompareResidualTrends->cd(iPtBin+1);
1086  if(out4)hResidualTrend->Draw("same");
1087  else hResidualTrend->Draw();
1088 
1089  c2pullTrend->cd(iPtBin+1);
1090  hPullsTrend->Draw();
1091 
1092  }
1093  else hMassSubRot->Draw("");
1094 
1095 
1096  c3->cd(iPtBin+1);
1097  if(out2 && fitterLS[iPtBin]->GetMassFunc()){
1098  fitterLS[iPtBin]->DrawHere(gPad,3,0);
1099  if(fitterLS[iPtBin]->GetRawYield()>0 && fitterLS[iPtBin]->GetReducedChiSquare()>0 && fitterLS[iPtBin]->GetReducedChiSquare()<5){
1100  hRawYieldLS->SetBinContent(iPtBin+1,fitterLS[iPtBin]->GetRawYield());
1101  hRawYieldLS->SetBinError(iPtBin+1,fitterLS[iPtBin]->GetRawYieldError());
1102  Double_t minBinBkg=hMassPtBin->FindBin(fitterLS[iPtBin]->GetMean()-3.*fitterLS[iPtBin]->GetSigma());
1103  Double_t maxBinBkg=hMassPtBin->FindBin(fitterLS[iPtBin]->GetMean()+3.*fitterLS[iPtBin]->GetSigma());
1104  background=hMassPtBin->Integral(minBinBkg,maxBinBkg);
1105  hSoverBLS->SetBinContent(iPtBin+1,fitterLS[iPtBin]->GetRawYield()/background);
1106  hSoverBLS->SetBinError(iPtBin+1,0.000001);
1107  hRelStatLS->SetBinContent(iPtBin+1,fitterLS[iPtBin]->GetRawYieldError()/fitterLS[iPtBin]->GetRawYield());
1108  hRelStatLS->SetBinError(iPtBin+1,0.0000001);
1109  hSignifLS->SetBinContent(iPtBin+1,fitterLS[iPtBin]->GetRawYield()/TMath::Sqrt(background+fitterLS[iPtBin]->GetRawYield()));
1110  hSignifLS->SetBinError(iPtBin+1,0.00000001);
1111  hGausMeanLS->SetBinContent(iPtBin+1,fitterLS[iPtBin]->GetMean());
1112  hGausMeanLS->SetBinError(iPtBin+1,fitterLS[iPtBin]->GetMeanUncertainty());
1113  hGausSigmaLS->SetBinContent(iPtBin+1,fitterLS[iPtBin]->GetSigma());
1114  hGausSigmaLS->SetBinError(iPtBin+1,fitterLS[iPtBin]->GetSigmaUncertainty());
1115  hChiSqLS->SetBinContent(iPtBin+1,fitterLS[iPtBin]->GetReducedChiSquare());
1116  hChiSqLS->SetBinError(iPtBin+1,0.00001); // very small number, for graphics
1117  hNdfLS->SetBinContent(iPtBin+1,fitterLS[iPtBin]->GetMassFunc()->GetNDF());
1118  hNdfLS->SetBinError(iPtBin+1,0.00001); // very small number, for graphics
1119  }
1120  hInvMassHistoBinWidthLS->SetBinContent(iPtBin+1,hMassSubLS->GetBinWidth(1));
1121  // if(!correctForRefl)
1122  WriteFitInfo(fitterLS[iPtBin],hMassPtBin);
1123  fout->cd();
1124  WriteFitFunctionsToFile(fitterLS[iPtBin],"LS",iPtBin);
1125  current->cd();
1126 
1127 
1128  Double_t errbc;
1130  hRawYieldLSBC->SetBinContent(iPtBin+1,bc);
1131  hRawYieldLSBC->SetBinError(iPtBin+1,errbc);
1132 
1133  c3sub->cd(iPtBin+1);
1134  fitterLS[iPtBin]->DrawHistoMinusFit(gPad);
1135 
1136  c3pulls->cd(iPtBin+1);
1137  TH1F *hPullsTrend=new TH1F();// the name is changed afterwards, histo must not be deleted
1138  TH1F *hPulls=new TH1F();//the name is changed afterwards, histo must not be deleted
1139  TH1F *hResidualTrend=new TH1F();// the name is changed afterwards, histo must not be deleted
1140  TH1F *hResiduals=fitterLS[iPtBin]->GetResidualsAndPulls(hPulls,hResidualTrend,hPullsTrend);
1141  hResiduals->SetName(Form("hResidualsLS_PtBin%d",iPtBin));
1142 
1143  hPulls->Draw();
1144  PrintGausParams(hPulls);
1145 
1146  c3residuals->cd(iPtBin+1);
1147  hResiduals->Draw();
1148 
1149  SetStyleHisto(hResidualTrend,kLS);
1150  c3residualTrend->cd(iPtBin+1);
1151  hResidualTrend->Draw();
1152 
1153  cCompareResidualTrends->cd(iPtBin+1);
1154  if(out4||out1)hResidualTrend->Draw("same");
1155  else hResidualTrend->Draw();
1156 
1157  c3pullTrend->cd(iPtBin+1);
1158  hPullsTrend->Draw();
1159  }
1160  else if(hMassPtBinls) hMassSubLS->Draw("");
1161 
1162  c4->cd(iPtBin+1);
1163  if(out3 && fitterME[iPtBin]->GetMassFunc()){
1164  fitterME[iPtBin]->DrawHere(gPad,3,0);
1165  if(fitterME[iPtBin]->GetRawYield()>0 && fitterME[iPtBin]->GetReducedChiSquare()>0 && fitterME[iPtBin]->GetReducedChiSquare()<5){
1166  hRawYieldME->SetBinContent(iPtBin+1,fitterME[iPtBin]->GetRawYield());
1167  hRawYieldME->SetBinError(iPtBin+1,fitterME[iPtBin]->GetRawYieldError());
1168  Double_t minBinBkg=hMassPtBin->FindBin(fitterME[iPtBin]->GetMean()-3.*fitterME[iPtBin]->GetSigma());
1169  Double_t maxBinBkg=hMassPtBin->FindBin(fitterME[iPtBin]->GetMean()+3.*fitterME[iPtBin]->GetSigma());
1170  background=hMassPtBin->Integral(minBinBkg,maxBinBkg);
1171  hSoverBME->SetBinContent(iPtBin+1,fitterME[iPtBin]->GetRawYield()/background);
1172  hSoverBME->SetBinError(iPtBin+1,0.000001);
1173  hRelStatME->SetBinContent(iPtBin+1,fitterME[iPtBin]->GetRawYieldError()/fitterME[iPtBin]->GetRawYield());
1174  hRelStatME->SetBinError(iPtBin+1,0.000001);
1175  hSignifME->SetBinContent(iPtBin+1,fitterME[iPtBin]->GetRawYield()/TMath::Sqrt(background+fitterME[iPtBin]->GetRawYield()));
1176  hSignifME->SetBinError(iPtBin+1,0.00000001);
1177  hGausMeanME->SetBinContent(iPtBin+1,fitterME[iPtBin]->GetMean());
1178  hGausMeanME->SetBinError(iPtBin+1,fitterME[iPtBin]->GetMeanUncertainty());
1179  hGausSigmaME->SetBinContent(iPtBin+1,fitterME[iPtBin]->GetSigma());
1180  hGausSigmaME->SetBinError(iPtBin+1,fitterME[iPtBin]->GetSigmaUncertainty());
1181  hChiSqME->SetBinContent(iPtBin+1,fitterME[iPtBin]->GetReducedChiSquare());
1182  hChiSqME->SetBinError(iPtBin+1,0.00001); // very small number, for graphics
1183  hNdfME->SetBinContent(iPtBin+1,fitterME[iPtBin]->GetMassFunc()->GetNDF());
1184  hNdfME->SetBinError(iPtBin+1,0.00001); // very small number, for graphics
1185  }
1186  hInvMassHistoBinWidthME->SetBinContent(iPtBin+1,hMassSubME->GetBinWidth(1));
1187  // if(!correctForRefl)
1188  WriteFitInfo(fitterME[iPtBin],hMassPtBin);
1189  fout->cd();
1190  WriteFitFunctionsToFile(fitterME[iPtBin],"ME",iPtBin);
1191  current->cd();
1192 
1193  Double_t errbc;
1195  hRawYieldMEBC->SetBinContent(iPtBin+1,bc);
1196  hRawYieldMEBC->SetBinError(iPtBin+1,errbc);
1197 
1198  c4sub->cd(iPtBin+1);
1199  fitterME[iPtBin]->DrawHistoMinusFit(gPad);
1200 
1201  c4pulls->cd(iPtBin+1);
1202  TH1F *hPullsTrend=new TH1F();// the name is changed afterwards, histo must not be deleted
1203  TH1F *hPulls=new TH1F();//the name is changed afterwards, histo must not be deleted
1204  TH1F *hResidualTrend=new TH1F();// the name is changed afterwards, histo must not be deleted
1205  TH1F *hResiduals=fitterME[iPtBin]->GetResidualsAndPulls(hPulls,hResidualTrend,hPullsTrend);
1206  hResiduals->SetName(Form("hResidualsME_PtBin%d",iPtBin));
1207  hPulls->Draw();
1208  PrintGausParams(hPulls);
1209 
1210  c4residuals->cd(iPtBin+1);
1211  hResiduals->Draw();
1212 
1213  SetStyleHisto(hResidualTrend,kME);
1214  c4residualTrend->cd(iPtBin+1);
1215  hResidualTrend->Draw();
1216 
1217  cCompareResidualTrends->cd(iPtBin+1);
1218  if(out4||out1||out2)hResidualTrend->Draw("same");
1219  else hResidualTrend->Draw();
1220 
1221  c4pullTrend->cd(iPtBin+1);
1222  hPullsTrend->Draw();
1223  }
1224  else hMassSubME->Draw("");
1225  if(correctForRefl){
1226  delete hMCReflPtBin;
1227  delete hMCSigPtBin;
1228  }
1229 
1230  }
1231  TString path(gSystem->pwd());
1232  path.Append("/figures");
1233  if(gSystem->AccessPathName(path.Data())){
1234  gROOT->ProcessLine(Form(".!mkdir -p %s",path.Data()));
1235  }
1236 
1237  if(saveCanvasAsEps>0){
1238  c1->SaveAs(Form("figures/InvMassSpectra_%s_%s_NoBkgSub.eps",sigConf.Data(),suffix.Data()));
1239  c2->SaveAs(Form("figures/InvMassSpectra_%s_%s_Rot.eps",sigConf.Data(),suffix.Data()));
1240  c3->SaveAs(Form("figures/InvMassSpectra_%s_%s_LS.eps",sigConf.Data(),suffix.Data()));
1241  c4->SaveAs(Form("figures/InvMassSpectra_%s_%s_EM.eps",sigConf.Data(),suffix.Data()));
1242  c2sub->SaveAs(Form("figures/InvMassSpectraFitSubtr_%s_%s_Rot.eps",sigConf.Data(),suffix.Data()));
1243  c3sub->SaveAs(Form("figures/InvMassSpectraFitSubtr_%s_%s_LS.eps",sigConf.Data(),suffix.Data()));
1244  c4sub->SaveAs(Form("figures/InvMassSpectraFitSubtr_%s_%s_EM.eps",sigConf.Data(),suffix.Data()));
1245  if(tryDirectFit){
1246  c5->SaveAs(Form("figures/InvMassSpectra_%s_%s_SB.eps",sigConf.Data(),suffix.Data()));
1247  c5sub->SaveAs(Form("figures/InvMassSpectraFitSubtr_%s_%s_SB.eps",sigConf.Data(),suffix.Data()));
1248  }
1249 
1250  if(saveCanvasAsEps>1){
1251  c2residuals->SaveAs(Form("figures/ResidualDistribution_%s_%s_Rot.eps",sigConf.Data(),suffix.Data()));
1252  c3residuals->SaveAs(Form("figures/ResidualDistribution_%s_%s_LS.eps",sigConf.Data(),suffix.Data()));
1253  c4residuals->SaveAs(Form("figures/ResidualDistribution_%s_%s_EM.eps",sigConf.Data(),suffix.Data()));
1254 
1255  c2residualTrend->SaveAs(Form("figures/residualTrendvsMass_%s_%s_Rot.eps",sigConf.Data(),suffix.Data()));
1256  c3residualTrend->SaveAs(Form("figures/residualTrendvsMass_%s_%s_LS.eps",sigConf.Data(),suffix.Data()));
1257  c4residualTrend->SaveAs(Form("figures/residualTrendvsMass_%s_%s_EM.eps",sigConf.Data(),suffix.Data()));
1258 
1259  c2pulls->SaveAs(Form("figures/PullDistribution_%s_%s_Rot.eps",sigConf.Data(),suffix.Data()));
1260  c3pulls->SaveAs(Form("figures/PullDistribution_%s_%s_LS.eps",sigConf.Data(),suffix.Data()));
1261  c4pulls->SaveAs(Form("figures/PullDistribution_%s_%s_EM.eps",sigConf.Data(),suffix.Data()));
1262 
1263  c2pullTrend->SaveAs(Form("figures/pullTrendvsMass_%s_%s_Rot.eps",sigConf.Data(),suffix.Data()));
1264  c3pullTrend->SaveAs(Form("figures/pullTrendvsMass_%s_%s_LS.eps",sigConf.Data(),suffix.Data()));
1265  c4pullTrend->SaveAs(Form("figures/pullTrendvsMass_%s_%s_EM.eps",sigConf.Data(),suffix.Data()));
1266 
1267  if(tryDirectFit){
1268  c5residuals->SaveAs(Form("figures/ResidualDistribution_%s_%s_SB.eps",sigConf.Data(),suffix.Data()));
1269  c5residualTrend->SaveAs(Form("figures/residualTrendvsMass_%s_%s_SB.eps",sigConf.Data(),suffix.Data()));
1270  c5pulls->SaveAs(Form("figures/PullDistribution_%s_%s_SB.eps",sigConf.Data(),suffix.Data()));
1271  c5pullTrend->SaveAs(Form("figures/pullTrendvsMass_%s_%s_SB.eps",sigConf.Data(),suffix.Data()));
1272  }
1273  }
1274  }
1275 
1276  // save also .root
1277  if(saveCanvasAsRoot){
1278  c2->SaveAs(Form("figures/InvMassSpectra_%s_%s_Rot.root",sigConf.Data(),suffix.Data()));
1279  c3->SaveAs(Form("figures/InvMassSpectra_%s_%s_LS.root",sigConf.Data(),suffix.Data()));
1280  c4->SaveAs(Form("figures/InvMassSpectra_%s_%s_EM.root",sigConf.Data(),suffix.Data()));
1281  c2sub->SaveAs(Form("figures/InvMassSpectraFitSubtr_%s_%s_Rot.root",sigConf.Data(),suffix.Data()));
1282  c3sub->SaveAs(Form("figures/InvMassSpectraFitSubtr_%s_%s_LS.root",sigConf.Data(),suffix.Data()));
1283  c4sub->SaveAs(Form("figures/InvMassSpectraFitSubtr_%s_%s_EM.root",sigConf.Data(),suffix.Data()));
1284 
1285  c2residuals->SaveAs(Form("figures/ResidualDistribution_%s_%s_Rot.root",sigConf.Data(),suffix.Data()));
1286  c3residuals->SaveAs(Form("figures/ResidualDistribution_%s_%s_LS.root",sigConf.Data(),suffix.Data()));
1287  c4residuals->SaveAs(Form("figures/ResidualDistribution_%s_%s_EM.root",sigConf.Data(),suffix.Data()));
1288 
1289  c2residualTrend->SaveAs(Form("figures/residualTrendvsMass_%s_%s_Rot.root",sigConf.Data(),suffix.Data()));
1290  c3residualTrend->SaveAs(Form("figures/residualTrendvsMass_%s_%s_LS.root",sigConf.Data(),suffix.Data()));
1291  c4residualTrend->SaveAs(Form("figures/residualTrendvsMass_%s_%s_EM.root",sigConf.Data(),suffix.Data()));
1292 
1293  c2pulls->SaveAs(Form("figures/PullDistribution_%s_%s_Rot.root",sigConf.Data(),suffix.Data()));
1294  c3pulls->SaveAs(Form("figures/PullDistribution_%s_%s_LS.root",sigConf.Data(),suffix.Data()));
1295  c4pulls->SaveAs(Form("figures/PullDistribution_%s_%s_EM.root",sigConf.Data(),suffix.Data()));
1296 
1297  c2pullTrend->SaveAs(Form("figures/pullTrendvsMass_%s_%s_Rot.root",sigConf.Data(),suffix.Data()));
1298  c3pullTrend->SaveAs(Form("figures/pullTrendvsMass_%s_%s_LS.root",sigConf.Data(),suffix.Data()));
1299  c4pullTrend->SaveAs(Form("figures/pullTrendvsMass_%s_%s_EM.root",sigConf.Data(),suffix.Data()));
1300 
1301  if(tryDirectFit){
1302  c5->SaveAs(Form("figures/InvMassSpectra_%s_%s_SB.root",sigConf.Data(),suffix.Data()));
1303  c5sub->SaveAs(Form("figures/InvMassSpectraFitSubtr_%s_%s_SB.root",sigConf.Data(),suffix.Data()));
1304  c5residuals->SaveAs(Form("figures/ResidualDistribution_%s_%s_SB.root",sigConf.Data(),suffix.Data()));
1305  c5residualTrend->SaveAs(Form("figures/residualTrendvsMass_%s_%s_SB.root",sigConf.Data(),suffix.Data()));
1306  c5pulls->SaveAs(Form("figures/PullDistribution_%s_%s_SB.root",sigConf.Data(),suffix.Data()));
1307  c5pullTrend->SaveAs(Form("figures/pullTrendvsMass_%s_%s_SB.root",sigConf.Data(),suffix.Data()));
1308  }
1309  }
1310 
1311 
1312  hRawYieldRot->SetMarkerStyle(21);
1313  hRawYieldLS->SetMarkerStyle(22);
1314  hRawYieldLS->SetMarkerColor(kGreen+2);
1315  hRawYieldLS->SetLineColor(kGreen+2);
1316  hRawYieldME->SetMarkerStyle(25);
1317  hRawYieldME->SetMarkerColor(4);
1318  hRawYieldME->SetLineColor(4);
1319  hRawYieldSB->SetMarkerStyle(27);
1320  hRawYieldSB->SetMarkerColor(6);
1321  hRawYieldSB->SetLineColor(6);
1322 
1323 
1324  TCanvas* cry=new TCanvas("cry","RawYield",800,700);
1325  cry->SetLeftMargin(0.15);
1326  hRawYieldRot->Draw("P");
1327  hRawYieldRot->SetMinimum(0);
1328  hRawYieldRot->GetYaxis()->SetTitleOffset(1.8);
1329  Double_t max=hRawYieldRot->GetMaximum();
1330  if(hRawYieldLS->GetMaximum()>max)max=hRawYieldLS->GetMaximum();
1331  if(hRawYieldME->GetMaximum()>max)max=hRawYieldME->GetMaximum();
1332  if(tryDirectFit){
1333  if(hRawYieldSB->GetMaximum()>max)max=hRawYieldSB->GetMaximum();
1334  }
1335  hRawYieldRot->SetMaximum(max*1.2);
1336  hRawYieldLS->Draw("PZSAME");
1337  hRawYieldME->Draw("PSAME");
1338  if(tryDirectFit) hRawYieldSB->Draw("PSAME");
1339  TLegend* legry=new TLegend(0.7,0.7,0.89,0.89);
1340  legry->SetFillStyle(0);
1341  legry->SetBorderSize(0);
1342  legry->AddEntry(hRawYieldRot,"Rotational","PL")->SetTextColor(1);
1343  legry->AddEntry(hRawYieldLS,"Like Sign","PL")->SetTextColor(kGreen+2);
1344  legry->AddEntry(hRawYieldME,"Ev Mix","PL")->SetTextColor(4);
1345  if(tryDirectFit){
1346  legry->AddEntry(hRawYieldSB,"Side-Band Fit","PL")->SetTextColor(6);
1347  }
1348  legry->Draw();
1349  if(saveCanvasAsEps>0) cry->SaveAs(Form("figures/RawYield_%s_%s.eps",sigConf.Data(),suffix.Data()));
1350  if(saveCanvasAsRoot) cry->SaveAs(Form("figures/RawYield_%s_%s.root",sigConf.Data(),suffix.Data()));
1351 
1352 
1353 
1354  TCanvas* cch2=new TCanvas("cch2","Chi2",800,700);
1355  cch2->SetLeftMargin(0.15);
1356  hChiSqRot->SetMarkerStyle(21);
1357  hChiSqRot->Draw("P");
1358  hChiSqRot->SetMinimum(0);
1359  hChiSqRot->GetXaxis()->SetTitle("p_{T} (GeV/c)");
1360  hChiSqRot->GetYaxis()->SetTitle("#chi^{2}/ndf");
1361  hChiSqRot->GetYaxis()->SetTitleOffset(1.8);
1362  Double_t maxc=hChiSqRot->GetMaximum();
1363  if(hChiSqLS->GetMaximum()>maxc)maxc=hChiSqLS->GetMaximum();
1364  if(hChiSqME->GetMaximum()>maxc)maxc=hChiSqME->GetMaximum();
1365  if(tryDirectFit){
1366  if(hChiSqSB->GetMaximum()>maxc)maxc=hChiSqSB->GetMaximum();
1367  }
1368  hChiSqRot->SetMaximum(maxc*1.2);
1369 
1370  hChiSqLS->SetMarkerStyle(22);
1371  hChiSqLS->SetMarkerColor(kGreen+2);
1372  hChiSqLS->SetLineColor(kGreen+2);
1373  hChiSqLS->Draw("PZSAME");
1374  hChiSqME->SetMarkerStyle(25);
1375  hChiSqME->SetMarkerColor(4);
1376  hChiSqME->SetLineColor(4);
1377  hChiSqME->Draw("PSAME");
1378  if(tryDirectFit){
1379  hChiSqSB->SetMarkerStyle(27);
1380  hChiSqSB->SetMarkerColor(6);
1381  hChiSqSB->SetLineColor(6);
1382  hChiSqSB->Draw("PSAME");
1383  }
1384  legry->Draw();
1385  if(saveCanvasAsEps>0) cch2->SaveAs(Form("figures/ChiSq_%s_%s.eps",sigConf.Data(),suffix.Data()));
1386  if(saveCanvasAsRoot) cch2->SaveAs(Form("figures/ChiSq_%s_%s.root",sigConf.Data(),suffix.Data()));
1387 
1388  TH1F* hRatioLSToME=(TH1F*)hRawYieldLS->Clone("hRatioLStoME");
1389  TH1F* hRatioRotToME=(TH1F*)hRawYieldRot->Clone("hRatioRottoME");
1390  TH1F* hRatioMEToME=(TH1F*)hRawYieldME->Clone("hRatioMEtoME");
1391  TH1F* hRatioSBToME=(TH1F*)hRawYieldSB->Clone("hRatioSBtoME");
1392  for(Int_t ib=1; ib<=hRawYieldME->GetNbinsX(); ib++){
1393  Double_t yme=hRawYieldME->GetBinContent(ib);
1394  if(yme>0.){
1395  hRatioLSToME->SetBinContent(ib,hRawYieldLS->GetBinContent(ib)/yme);
1396  hRatioLSToME->SetBinError(ib,hRawYieldLS->GetBinError(ib)/yme);
1397  hRatioMEToME->SetBinContent(ib,hRawYieldME->GetBinContent(ib)/yme);
1398  hRatioMEToME->SetBinError(ib,hRawYieldME->GetBinError(ib)/yme);
1399  hRatioRotToME->SetBinContent(ib,hRawYieldRot->GetBinContent(ib)/yme);
1400  hRatioRotToME->SetBinError(ib,hRawYieldRot->GetBinError(ib)/yme);
1401  hRatioSBToME->SetBinContent(ib,hRawYieldSB->GetBinContent(ib)/yme);
1402  hRatioSBToME->SetBinError(ib,hRawYieldSB->GetBinError(ib)/yme);
1403  }
1404  }
1405 
1406  TCanvas* cry2=new TCanvas("cry2","RawYield+Ratios",1400,700);
1407  cry2->Divide(2,1);
1408  cry2->cd(1);
1409  gPad->SetLeftMargin(0.15);
1410  gPad->SetRightMargin(0.05);
1411  hRawYieldRot->Draw("P");
1412  hRawYieldLS->Draw("PZSAME");
1413  hRawYieldME->Draw("PSAME");
1414  if(tryDirectFit){
1415  hRawYieldSB->Draw("PSAME");
1416  }
1417  legry->Draw();
1418  cry2->cd(2);
1419  hRatioLSToME->SetStats(0);
1420  hRatioLSToME->SetMinimum(0.3);
1421  hRatioLSToME->SetMaximum(1.7);
1422  hRatioLSToME->GetYaxis()->SetTitle("Ratio To EvMix");
1423  hRatioLSToME->Draw("same");
1424  hRatioRotToME->Draw("same");
1425  hRatioMEToME->Draw("same");
1426  hRatioSBToME->Draw("same");
1427  if(saveCanvasAsEps>0) cry2->SaveAs(Form("figures/RawYieldAndRatios_%s_%s.eps",sigConf.Data(),suffix.Data()));
1428  if(saveCanvasAsRoot) cry2->SaveAs(Form("figures/RawYieldAndRatios_%s_%s.root",sigConf.Data(),suffix.Data()));
1429 
1430  hRelStatRot->GetYaxis()->SetTitleOffset(1.8);
1431  hRelStatRot->SetMarkerStyle(21);
1432  hRelStatLS->SetMarkerStyle(22);
1433  hRelStatLS->SetMarkerColor(kGreen+2);
1434  hRelStatLS->SetLineColor(kGreen+2);
1435  hRelStatME->SetMarkerStyle(25);
1436  hRelStatME->SetMarkerColor(4);
1437  hRelStatME->SetLineColor(4);
1438  hRelStatSB->SetMarkerStyle(27);
1439  hRelStatSB->SetMarkerColor(6);
1440  hRelStatSB->SetLineColor(6);
1441 
1442  TCanvas* cry3=new TCanvas("cry3","RawYield+Ratios+Unc",1800,600);
1443  cry3->Divide(3,1);
1444  cry3->cd(1);
1445  gPad->SetLeftMargin(0.15);
1446  gPad->SetRightMargin(0.05);
1447  hRawYieldRot->Draw("P");
1448  hRawYieldLS->Draw("PZSAME");
1449  hRawYieldME->Draw("PSAME");
1450  if(tryDirectFit){
1451  hRawYieldSB->Draw("PSAME");
1452  }
1453  legry->Draw();
1454  cry3->cd(2);
1455  hRatioLSToME->SetStats(0);
1456  hRatioLSToME->SetMinimum(0.3);
1457  hRatioLSToME->SetMaximum(1.7);
1458  hRatioLSToME->GetYaxis()->SetTitle("Ratio To EvMix");
1459  hRatioLSToME->Draw("same");
1460  hRatioRotToME->Draw("same");
1461  hRatioMEToME->Draw("same");
1462  hRatioSBToME->Draw("same");
1463  cry3->cd(3);
1464  gPad->SetLeftMargin(0.15);
1465  gPad->SetRightMargin(0.05);
1466  hRelStatRot->SetStats(0);
1467  hRelStatRot->SetMinimum(0.04);
1468  hRelStatRot->SetMaximum(0.4);
1469  hRelStatRot->Draw();
1470  hRelStatLS->Draw("same");
1471  hRelStatME->Draw("same");
1472  if(tryDirectFit) hRelStatSB->Draw("same");
1473 
1474  if(saveCanvasAsEps>0) cry3->SaveAs(Form("figures/RawYieldRatiosAndUnc_%s_%s.eps",sigConf.Data(),suffix.Data()));
1475  if(saveCanvasAsRoot) cry3->SaveAs(Form("figures/RawYieldRatiosAndUnc_%s_%s.root",sigConf.Data(),suffix.Data()));
1476 
1477 
1478 
1479 
1480  cCompareResidualTrends->cd(1);
1481  TLegend *legRT=new TLegend(*legry);
1482  legRT->SetX1NDC(0.12);
1483  legRT->SetY1NDC(0.7);
1484  legRT->SetX2NDC(0.3);
1485  legRT->SetY2NDC(0.9);
1486  legRT->Draw();
1487 
1488  // NOW COMPARE YIELDS OBTAINED WITH BC FROM DIFFERENT APPROACHES
1489  TCanvas* cryBC=new TCanvas("cryBC","RawYield with BC",800,700);
1490  cryBC->SetLeftMargin(0.15);
1491  hRawYieldRot->Draw("P");
1492  hRawYieldME->Draw("PSAME");
1493  hRawYieldLS->Draw("PZSAME");
1494  if(tryDirectFit){
1495  hRawYieldSB->Draw("PSAME");
1496  }
1497 
1498  if(hRawYieldRotBC->GetMaximum()>max)max=hRawYieldRotBC->GetMaximum();
1499  if(hRawYieldLSBC->GetMaximum()>max)max=hRawYieldLSBC->GetMaximum();
1500  if(hRawYieldMEBC->GetMaximum()>max)max=hRawYieldMEBC->GetMaximum();
1501  if(tryDirectFit){
1502  if(hRawYieldSBBC->GetMaximum()>max)max=hRawYieldSBBC->GetMaximum();
1503  }
1504  hRawYieldRot->SetMaximum(max*1.2);
1505 
1506  hRawYieldRotBC->SetMarkerStyle(21);
1507  hRawYieldRotBC->SetLineStyle(2);
1508  hRawYieldRotBC->Draw("PSAME");
1509 
1510 
1511  hRawYieldLSBC->SetMarkerStyle(22);
1512  hRawYieldLSBC->SetMarkerColor(kGreen+2);
1513  hRawYieldLSBC->SetLineColor(kGreen+2);
1514  hRawYieldLSBC->SetLineStyle(2);
1515  hRawYieldLSBC->Draw("PZSAME");
1516 
1517  hRawYieldMEBC->SetMarkerStyle(25);
1518  hRawYieldMEBC->SetMarkerColor(4);
1519  hRawYieldMEBC->SetLineColor(4);
1520  hRawYieldMEBC->SetLineStyle(2);
1521  hRawYieldMEBC->Draw("PSAME");
1522  if(tryDirectFit){
1523  hRawYieldSBBC->SetMarkerStyle(27);
1524  hRawYieldSBBC->SetMarkerColor(6);
1525  hRawYieldSBBC->SetLineColor(6);
1526  hRawYieldSBBC->SetLineStyle(2);
1527  hRawYieldSBBC->Draw("PSAME");
1528  }
1529 
1530  TLegend* legryBC=new TLegend(0.7,0.7,0.89,0.89);
1531  legryBC->SetFillStyle(0);
1532  legryBC->SetBorderSize(0);
1533  legryBC->AddEntry(hRawYieldRot,"Rotational","PL")->SetTextColor(1);
1534  legryBC->AddEntry(hRawYieldRotBC,"Rotational BC","PL")->SetTextColor(1);
1535  legryBC->AddEntry(hRawYieldLS,"Like Sign ","PL")->SetTextColor(kGreen+2);
1536  legryBC->AddEntry(hRawYieldLSBC,"Like Sign BC","PL")->SetTextColor(kGreen+2);
1537  legryBC->AddEntry(hRawYieldME,"Ev Mix","PL")->SetTextColor(4);
1538  legryBC->AddEntry(hRawYieldMEBC,"Ev Mix BC","PL")->SetTextColor(4);
1539  if(tryDirectFit){
1540  legryBC->AddEntry(hRawYieldSB,"Side-Band Fit","PL")->SetTextColor(6);
1541  legryBC->AddEntry(hRawYieldSBBC,"Side-Band Fit BC","PL")->SetTextColor(6);
1542  }
1543  legryBC->Draw();
1544  if(saveCanvasAsEps>0) cryBC->SaveAs(Form("figures/RawYieldBC_%s_%s.eps",sigConf.Data(),suffix.Data()));
1545  if(saveCanvasAsRoot) cryBC->SaveAs(Form("figures/RawYieldBC_%s_%s.root",sigConf.Data(),suffix.Data()));
1546 
1547 
1548 
1549 
1550  fout->cd();
1551  hRawYieldRot->Write();
1552  hRawYieldLS->Write();
1553  hRawYieldME->Write();
1554  hRawYieldSB->Write();
1555  hRelStatRot->Write();
1556  hRelStatLS->Write();
1557  hRelStatME->Write();
1558  hRelStatSB->Write();
1559  hSignifRot->Write();
1560  hSignifLS->Write();
1561  hSignifME->Write();
1562  hSignifSB->Write();
1563  hSoverBRot->Write();
1564  hSoverBLS->Write();
1565  hSoverBME->Write();
1566  hSoverBSB->Write();
1567  hGausMeanRot->Write();
1568  hGausMeanLS->Write();
1569  hGausMeanME->Write();
1570  hGausMeanSB->Write();
1571  hGausSigmaRot->Write();
1572  hGausSigmaLS->Write();
1573  hGausSigmaME->Write();
1574  hGausSigmaSB->Write();
1575  hChiSqRot->Write();
1576  hChiSqLS->Write();
1577  hChiSqME->Write();
1578  hChiSqSB->Write();
1579  hNdfRot->Write();
1580  hNdfLS->Write();
1581  hNdfME->Write();
1582  hNdfSB->Write();
1583  hInvMassHistoBinWidthRot->Write();
1584  hInvMassHistoBinWidthLS->Write();
1585  hInvMassHistoBinWidthME->Write();
1586  hInvMassHistoBinWidthSB->Write();
1587  hRebin->Write();
1588  hBkgFitFunc->Write();
1589  hBkgFitFuncSB->Write();
1590  hnEv->Write();
1591  if(hSigmaMC) hSigmaMC->Write();
1592  fout->Close();
1593 
1594  return;
1595 }
1596 
1598  TString nameh;
1599  TF1* fTot=fitter->GetMassFunc();
1600  if(fTot){
1601  nameh.Form("FitFuncTot_%s_PtBin%d",meth.Data(),iPtBin);
1602  fTot->SetRange(1.6,2.2);
1603  fTot->SetNpx(500);
1604  fTot->Write(nameh.Data());
1605  }
1606  TF1* fSig=fitter->GetSignalFunc();
1607  nameh.Form("FitFuncSig_%s_PtBin%d",meth.Data(),iPtBin);
1608  if(fSig){
1609  fSig->SetRange(1.6,2.2);
1610  fSig->SetNpx(500);
1611  fSig->Write(nameh.Data());
1612  }
1613  TF1* fBkg=fitter->GetBackgroundRecalcFunc();
1614  nameh.Form("FitFuncBkg_%s_PtBin%d",meth.Data(),iPtBin);
1615  if(fBkg){
1616  fBkg->SetRange(1.6,2.2);
1617  fBkg->SetNpx(500);
1618  fBkg->Write(nameh.Data());
1619  }
1620  if(meson=="Dzero"){
1621  TF1* fBkgR=fitter->GetBkgPlusReflFunc();
1622  nameh.Form("FitFuncBkgRefl_%s_PtBin%d",meth.Data(),iPtBin);
1623  if(fBkgR){
1624  fBkgR->SetRange(1.6,2.2);
1625  fBkgR->SetNpx(500);
1626  fBkgR->Write(nameh.Data());
1627  }
1628  }
1629  return;
1630 }
1631 
1633  Double_t sig=fitter->GetRawYield();
1634  Double_t esig=fitter->GetRawYieldError();
1635  Double_t mean=fitter->GetMean();
1636  Double_t emean=fitter->GetMeanUncertainty();
1637  Double_t sigma=fitter->GetSigma();
1638  Double_t esigma=fitter->GetSigmaUncertainty();
1639  Double_t minBin=histo->FindBin(mean-3.*sigma);
1640  Double_t maxBin=histo->FindBin(mean+3.*sigma);
1641  Double_t back=histo->Integral(minBin,maxBin);
1642  TPaveText* tpar=new TPaveText(0.5,0.7,0.89,.87,"NDC");
1643  tpar->SetBorderSize(0);
1644  tpar->SetFillStyle(0);
1645  tpar->AddText(Form("Mean = %.3f #pm %.3f",mean,emean));
1646  tpar->AddText(Form("Sigma = %.3f #pm %.3f",sigma,esigma));
1647  tpar->SetTextColor(4);
1648  tpar->Draw();
1649 
1650  TPaveText* tss=new TPaveText(0.15,0.15,0.5,0.4,"NDC");
1651  tss->SetBorderSize(0);
1652  tss->SetFillStyle(0);
1653  tss->AddText(Form("S = %.0f #pm %.0f",sig,esig));
1654  tss->AddText(Form("B(3#sigma) = %.3g",back));
1655  tss->AddText(Form("S/B (3#sigma) = %.4g",sig/back));
1656  if(correctForRefl) tss->AddText(Form("Refl/Sig = %.3f #pm %.3f ",fitter->GetReflOverSig(),fitter->GetReflOverSigUncertainty()));
1657  tss->AddText(Form("Significance(3#sigma) = %.2f",sig/TMath::Sqrt(back+sig)));
1658  tss->SetTextColor(1);
1659  tss->Draw();
1660 
1661 }
1662 
1664  TH3F* h3d=(TH3F*)lMC->FindObject(Form("hMassVsPtVs%sSig",var.Data()));
1665  if(!h3d){
1666  printf("hMassVsPtVsYSig not found\n");
1667  return 0x0;
1668  }
1669  TH3F* h3dr=(TH3F*)lMC->FindObject(Form("hMassVsPtVs%sRefl",var.Data()));
1670  if(!h3dr){
1671  printf("hMassVsPtVsYRefl not found\n");
1672  }
1673  TH1F* hSigmaMC=new TH1F("hSigmaMC","",nPtBins,binLims);
1674  Int_t zbin1=0;
1675  Int_t zbin2=h3d->GetZaxis()->GetNbins()+1;
1676  if(var=="CosthSt" && costhstcut<1.){
1677  zbin2=h3d->GetZaxis()->FindBin(costhstcut-0.000001);
1678  }
1679 
1680  TCanvas* cmc1=new TCanvas("InvMassMC","InvMassMC",1200,800);
1681  DivideCanvas(cmc1,nPtBins);
1682 
1683  gStyle->SetOptFit(0);
1684  gStyle->SetOptStat(0);
1685  for(Int_t iPtBin=0; iPtBin<nPtBins; iPtBin++){
1686  Int_t bin1=h3d->GetYaxis()->FindBin(binLims[iPtBin]);
1687  Int_t bin2=h3d->GetYaxis()->FindBin(binLims[iPtBin+1]-0.0001);
1688  Double_t ptmed=(binLims[iPtBin]+binLims[iPtBin+1])*0.5;
1689  Double_t pthalfwid=(binLims[iPtBin+1]-binLims[iPtBin])*0.5;
1690  printf("Bin %d Pt range=%f %f\n",iPtBin,h3d->GetYaxis()->GetBinLowEdge(bin1),h3d->GetYaxis()->GetBinUpEdge(bin2));
1691  TH1D* hMassMCPtBin=h3d->ProjectionX(Form("hMassMCPtBin%d",iPtBin),bin1,bin2,zbin1,zbin2);
1692  hMassMCPtBin->SetTitle(Form("%.1f<p_{T}<%.1f GeV/c",binLims[iPtBin],binLims[iPtBin+1]));
1693  hMassMCPtBin->GetXaxis()->SetTitle("Invariant mass (GeV/c^{2})");
1694  hMassMCPtBin->GetXaxis()->SetRangeUser(1.65,2.06);
1695  cmc1->cd(iPtBin+1);
1696  gPad->SetLogy();
1697  hMassMCPtBin->Fit("gaus","");
1698  TF1* fg=(TF1*)hMassMCPtBin->GetListOfFunctions()->FindObject("gaus");
1699  TLatex* tsig=new TLatex(0.65,0.82,"Signal");
1700  tsig->SetNDC();
1701  tsig->SetTextColor(kBlue+1);
1702  tsig->Draw();
1703 
1704  TPaveText* t1=new TPaveText(0.15,0.55,0.4,0.88,"ndc");
1705  t1->SetFillStyle(0);
1706  t1->SetBorderSize(0);
1707  t1->AddText(Form("#mu = %.4f GeV/c^{2}",fg->GetParameter(1)));
1708  t1->AddText(Form("#sigma = %.4f GeV/c^{2}",fg->GetParameter(2)));
1709  t1->AddText(Form("Integral = %.0f",fg->Integral(1.7,2.1)/hMassMCPtBin->GetBinWidth(1)));
1710  t1->AddText(Form("Entries = %.0f",hMassMCPtBin->GetEntries()));
1711 
1712  if(h3dr){
1713  TH1D* hReflPtBin=h3dr->ProjectionX(Form("hReflPtBin%d",iPtBin),bin1,bin2,zbin1,zbin2);
1714  hReflPtBin->SetTitle(Form("%.1f<p_{T}<%.1f GeV/c",binLims[iPtBin],binLims[iPtBin+1]));
1715  hReflPtBin->GetXaxis()->SetTitle("Invariant mass (GeV/c^{2})");
1716  hReflPtBin->SetLineColor(2);
1717  Int_t bin1=hReflPtBin->FindBin(fg->GetParameter(1)-3.*fg->GetParameter(2));
1718  Int_t bin2=hReflPtBin->FindBin(fg->GetParameter(1)+3.*fg->GetParameter(2));
1719  TLatex* tref=new TLatex(0.65,0.75,"Reflections");
1720  tref->SetNDC();
1721  tref->SetTextColor(2);
1722  tref->Draw();
1723  t1->AddText(Form("Refl(3#sigma) = %.0f",hReflPtBin->Integral(bin1,bin2)));
1724  hReflPtBin->Draw("same");
1725  t1->Draw();
1726  }
1727  sigmas[iPtBin]=fg->GetParameter(2);
1728  Double_t errMCsigma=fg->GetParError(2);
1729  if(tuneSigmaOnData>0.){
1730  sigmas[iPtBin]*=tuneSigmaOnData;
1731  errMCsigma*=tuneSigmaOnData;
1732  }
1733  hSigmaMC->SetBinContent(iPtBin+1,sigmas[iPtBin]);
1734  hSigmaMC->SetBinError(iPtBin+1,errMCsigma);
1735  }
1736  return hSigmaMC;
1737 }
1739  printf("Meson: %s\n",meson.Data());
1740  printf("Data file: %s\n", fileName.Data());
1741  printf("Suffix Data: %s\n", suffix.Data());
1742  printf("MC file: %s\n", fileNameMC.Data());
1743  printf("Suffix MC: %s\n", suffixMC.Data());
1744  printf("Number of pt bins = %d\n",nPtBins);
1745  for(Int_t j=0; j<nPtBins+1; j++) printf(" %.1f",binLims[j]);
1746  printf("\n");
1747  for(Int_t j=0; j<nPtBins; j++) printf(" %d",rebin[j]);
1748  printf("\n");
1749  for(Int_t j=0; j<nPtBins; j++) printf(" %.2f",minMass4Fit[j]);
1750  printf("\n");
1751  for(Int_t j=0; j<nPtBins; j++) printf(" %.2f",maxMass4Fit[j]);
1752  printf("\n");
1753  printf("Gaussian sigma option: %d (tune on data = %f)\n",fixSigmaConf,tuneSigmaOnData);
1754  for(Int_t j=0; j<nPtBins; j++) printf(" %d",fixSigma[j]);
1755  printf("\n");
1756  for(Int_t j=0; j<nPtBins; j++) printf(" %.3f",sigmas[j]);
1757  printf("\n");
1758  printf("Gaussian mean option: %d (tune on data = %f)\n",fixMeanConf,tuneMeanOnData);
1759  for(Int_t j=0; j<nPtBins; j++) printf(" %d",fixMean[j]);
1760  printf("\n");
1761  printf("Background options: Norm=%d NormRange=%.2f UseLSinEvMex=%d SmoothLS=%d\n",optForNorm,rangeForNorm,useEMwithLS,smoothLS);
1762  printf("Direct fit: %d\n",tryDirectFit);
1763  for(Int_t j=0; j<nPtBins; j++) printf(" %d",nDegreeBackPol[j]);
1764  printf("\n");
1765  for(Int_t j=0; j<nPtBins; j++) printf(" %.2f",fitSBrangelow[j]);
1766  printf("\n");
1767  for(Int_t j=0; j<nPtBins; j++) printf(" %.2f",fitSBrangeup[j]);
1768  printf("\n");
1769  printf("Fit Configuration: Option=%s Background Func=%d\n",fitoption.Data(),typeb);
1770  printf("Reflections: Use=%d Func=%s rOverSmodif=%f\n",correctForRefl,reflopt.Data(),rOverSmodif);
1771  printf("BinCounting: Option=%d nsigma=%.2f\n",optBkgBinCount,nsigmaBinCounting);
1772 }
1773 
1775  FILE* confFil=fopen(configName.Data(),"r");
1776  char line[50];
1777  char name[200];
1778  int n;
1779  float x;
1780  bool readok;
1781  while(!feof(confFil)){
1782  readok=fscanf(confFil,"%s:",line);
1783  if(strstr(line,"DataFile")){
1784  readok=fscanf(confFil,"%s",name);
1785  fileName=name;
1786  }
1787  else if(strstr(line,"SuffixData")){
1788  readok=fscanf(confFil,"%s",name);
1789  suffix=name;
1790  }
1791  else if(strstr(line,"MCFile")){
1792  readok=fscanf(confFil,"%s",name);
1793  fileNameMC=name;
1794  }
1795  else if(strstr(line,"SuffixMC")){
1796  readok=fscanf(confFil,"%s",name);
1797  suffixMC=name;
1798  }
1799  else if(strstr(line,"Meson")){
1800  readok=fscanf(confFil,"%s",name);
1801  meson=name;
1802  }
1803  else if(strstr(line,"NumOfPtBins")){
1804  readok=fscanf(confFil,"%d",&n);
1805  nPtBins=n;
1806  }
1807  else if(strstr(line,"BinLimits")){
1808  readok=fscanf(confFil," [ ");
1809  for(int j=0; j<nPtBins; j++){
1810  readok=fscanf(confFil,"%f,",&x);
1811  binLims[j]=x;
1812  if(j>0 && binLims[j]<=binLims[j-1]){
1813  printf("ERROR in array of pt bin limits\n");
1814  return kFALSE;
1815  }
1816  }
1817  readok=fscanf(confFil,"%f",&x);
1818  binLims[nPtBins]=x;
1819  readok=fscanf(confFil," ] ");
1820  }
1821  else if(strstr(line,"Rebin")){
1822  readok=fscanf(confFil," [ ");
1823  for(int j=0; j<nPtBins; j++){
1824  if(j<nPtBins-1) readok=fscanf(confFil,"%d,",&n);
1825  else readok=fscanf(confFil,"%d",&n);
1826  rebin[j]=n;
1827  if(rebin[j]<=0 || rebin[j]>20){
1828  printf("ERROR in array of rebin values\n");
1829  return kFALSE;
1830  }
1831  }
1832  readok=fscanf(confFil," ] ");
1833  }
1834  else if(strstr(line,"MinFit")){
1835  readok=fscanf(confFil," [ ");
1836  for(int j=0; j<nPtBins; j++){
1837  if(j<nPtBins-1) readok=fscanf(confFil,"%f,",&x);
1838  else readok=fscanf(confFil,"%f",&x);
1839  minMass4Fit[j]=x;
1840  if(minMass4Fit[j]<1.5 || minMass4Fit[j]>1.86){
1841  printf("ERROR in array of min mass values\n");
1842  return kFALSE;
1843  }
1844  }
1845  readok=fscanf(confFil," ] ");
1846  }
1847  else if(strstr(line,"MaxFit")){
1848  readok=fscanf(confFil," [ ");
1849  for(int j=0; j<nPtBins; j++){
1850  if(j<nPtBins-1) readok=fscanf(confFil,"%f,",&x);
1851  else readok=fscanf(confFil,"%f",&x);
1852  maxMass4Fit[j]=x;
1853  if(maxMass4Fit[j]<1.86 || maxMass4Fit[j]>2.2){
1854  printf("ERROR in array of max mass values\n");
1855  return kFALSE;
1856  }
1857  }
1858  readok=fscanf(confFil," ] ");
1859  }
1860  else if(strstr(line,"FixSigmaConf")){
1861  readok=fscanf(confFil,"%d",&n);
1862  fixSigmaConf=n;
1863  }
1864  else if(strstr(line,"FixSigmaPerBin")){
1865  readok=fscanf(confFil," [ ");
1866  for(int j=0; j<nPtBins; j++){
1867  if(j<nPtBins-1) readok=fscanf(confFil,"%d,",&n);
1868  else readok=fscanf(confFil,"%d",&n);
1869  fixSigma[j]=n;
1870  if(n<0 || n>1){
1871  printf("ERROR in array of fix-sigma settings\n");
1872  return kFALSE;
1873  }
1874  }
1875  readok=fscanf(confFil," ] ");
1876  }
1877  else if(strstr(line,"TuneSigmaOnData")){
1878  readok=fscanf(confFil,"%f",&x);
1879  tuneSigmaOnData=x;
1880  }
1881  else if(strstr(line,"SigmaManual")){
1882  readok=fscanf(confFil," [ ");
1883  for(int j=0; j<nPtBins; j++){
1884  if(j<nPtBins-1) readok=fscanf(confFil,"%f,",&x);
1885  else readok=fscanf(confFil,"%f",&x);
1886  sigmas[j]=x;
1887  if(sigmas[j]<0.001 || sigmas[j]>0.04){
1888  printf("ERROR in array of sigma values\n");
1889  return kFALSE;
1890  }
1891  }
1892  readok=fscanf(confFil," ] ");
1893  }
1894  else if(strstr(line,"FixMeanConf")){
1895  readok=fscanf(confFil,"%d",&n);
1896  fixMeanConf=n;
1897  }
1898  else if(strstr(line,"FixMeanPerBin")){
1899  readok=fscanf(confFil," [ ");
1900  for(int j=0; j<nPtBins; j++){
1901  if(j<nPtBins-1) readok=fscanf(confFil,"%d,",&n);
1902  else readok=fscanf(confFil,"%d",&n);
1903  fixMean[j]=n;
1904  if(n<0 || n>1){
1905  printf("ERROR in array of fix-sigma settings\n");
1906  return kFALSE;
1907  }
1908  }
1909  readok=fscanf(confFil," ] ");
1910  }
1911  else if(strstr(line,"TuneMeanOnData")){
1912  readok=fscanf(confFil,"%f",&x);
1913  tuneMeanOnData=x;
1914  }
1915  else if(strstr(line,"NormalizationOption")){
1916  readok=fscanf(confFil,"%d",&n);
1917  optForNorm=n;
1918  }
1919  else if(strstr(line,"RangeForNormalization")){
1920  readok=fscanf(confFil,"%f",&x);
1921  rangeForNorm=x;
1922  }
1923  else if(strstr(line,"UseEvMixWithLS")){
1924  readok=fscanf(confFil,"%d",&n);
1925  if(n<0 || n>1){
1926  printf("ERROR in UseEvMixWithLS setting\n");
1927  return kFALSE;
1928  }
1929  useEMwithLS=n;
1930  }
1931  else if(strstr(line,"UseGeomMeanLS")){
1932  readok=fscanf(confFil,"%d",&n);
1933  if(n<0 || n>1){
1934  printf("ERROR in UseGeomMeanLS setting\n");
1935  return kFALSE;
1936  }
1937  useGeomMeanLS=n;
1938  }
1939  else if(strstr(line,"RenormalizeLS")){
1940  readok=fscanf(confFil,"%d",&n);
1941  if(n<0 || n>1){
1942  printf("ERROR in RenormalizeLS setting\n");
1943  return kFALSE;
1944  }
1945  renormLS=n;
1946  }
1947  else if(strstr(line,"SmoothLS")){
1948  readok=fscanf(confFil,"%d",&n);
1949  smoothLS=n;
1950  }
1951  else if(strstr(line,"TryDirectFit")){
1952  readok=fscanf(confFil,"%d",&n);
1953  if(n<0 || n>1){
1954  printf("ERROR in TryDirectFit setting\n");
1955  return kFALSE;
1956  }
1957  tryDirectFit=n;
1958  }
1959  else if(strstr(line,"MinSBFit")){
1960  readok=fscanf(confFil," [ ");
1961  for(int j=0; j<nPtBins; j++){
1962  if(j<nPtBins-1) readok=fscanf(confFil,"%f,",&x);
1963  else readok=fscanf(confFil,"%f",&x);
1964  fitSBrangelow[j]=x;
1965  if(fitSBrangelow[j]<1.5 || fitSBrangelow[j]>1.86){
1966  printf("ERROR in array of min mass SB values\n");
1967  return kFALSE;
1968  }
1969  }
1970  readok=fscanf(confFil," ] ");
1971  }
1972  else if(strstr(line,"MaxSBFit")){
1973  readok=fscanf(confFil," [ ");
1974  for(int j=0; j<nPtBins; j++){
1975  if(j<nPtBins-1) readok=fscanf(confFil,"%f,",&x);
1976  else readok=fscanf(confFil,"%f",&x);
1977  fitSBrangeup[j]=x;
1978  if(fitSBrangeup[j]<1.86 || fitSBrangeup[j]>2.2){
1979  printf("ERROR in array of max mass SB values\n");
1980  return kFALSE;
1981  }
1982  }
1983  readok=fscanf(confFil," ] ");
1984  }
1985  else if(strstr(line,"PolDegreeSB")){
1986  readok=fscanf(confFil," [ ");
1987  for(int j=0; j<nPtBins; j++){
1988  if(j<nPtBins-1) readok=fscanf(confFil,"%d,",&n);
1989  else readok=fscanf(confFil,"%d",&n);
1990  nDegreeBackPolSB[j]=n;
1991  if(nDegreeBackPolSB[j]<=0 || nDegreeBackPolSB[j]>10){
1992  printf("ERROR in array of polynomial degree values\n");
1993  return kFALSE;
1994  }
1995  }
1996  readok=fscanf(confFil," ] ");
1997  }
1998  else if(strstr(line,"FitOption")){
1999  readok=fscanf(confFil,"%s",name);
2000  fitoption=name;
2001  }
2002  else if(strstr(line,"BackgroundFunction")){
2003  readok=fscanf(confFil,"%d",&n);
2004  typeb=n;
2005  }
2006  else if(strstr(line,"PolDegreeFit")){
2007  readok=fscanf(confFil," [ ");
2008  for(int j=0; j<nPtBins; j++){
2009  if(j<nPtBins-1) readok=fscanf(confFil,"%d,",&n);
2010  else readok=fscanf(confFil,"%d",&n);
2011  nDegreeBackPol[j]=n;
2012  if(nDegreeBackPol[j]<=0 || nDegreeBackPol[j]>10){
2013  printf("ERROR in array of polynomial degree values\n");
2014  return kFALSE;
2015  }
2016  }
2017  readok=fscanf(confFil," ] ");
2018  }
2019  else if(strstr(line,"CorrectForReflections")){
2020  readok=fscanf(confFil,"%d",&n);
2021  correctForRefl=n;
2022  }
2023  else if(strstr(line,"ReflectionOption")){
2024  readok=fscanf(confFil,"%s",name);
2025  reflopt=name;
2026  }
2027  else if(strstr(line,"ReflOverSignalModif")){
2028  readok=fscanf(confFil,"%f",&x);
2029  rOverSmodif=x;
2030  }
2031  else if(strstr(line,"BackgroundFuncOptionForBinCount")){
2032  readok=fscanf(confFil,"%d",&n);
2033  optBkgBinCount=n;
2034  }
2035  else if(strstr(line,"NumOfSigmaForBinCount")){
2036  readok=fscanf(confFil,"%f",&x);
2038  }
2039  else if(strstr(line,"CosThetaStar")){
2040  readok=fscanf(confFil,"%f",&x);
2041  costhstcut=x;
2042  }
2043  }
2044  return kTRUE;
2045 }
Double_t fitSBrangeup[maxPtBins]
Double_t rOverSmodif
Double_t GetMeanUncertainty() const
Double_t GetReflOverSig() const
Double_t maxMass4Fit[maxPtBins]
TH1F * GetResidualsAndPulls(TH1 *hPulls=0x0, TH1 *hResidualTrend=0x0, TH1 *hPullsTrend=0x0, Double_t minrange=0, Double_t maxrange=-1, Int_t option=0)
static TH1D * RebinHisto(TH1 *hOrig, Int_t reb, Int_t firstUse=-1)
Rebinning of invariant mass histograms.
Bool_t ReadConfig(TString configName)
Double_t sigmas[maxPtBins]
const Double_t ymax
Definition: AddTaskCFDStar.C:7
double Double_t
Definition: External.C:58
Definition: External.C:260
Int_t nDegreeBackPolSB[maxPtBins]
void SetFixReflOverS(Double_t rovers)
Definition: External.C:236
Int_t saveCanvasAsEps
TString reflopt
TString configFileName
void DivideCanvas(TCanvas *c, Int_t ndivisions)
TString fileName
TSystem * gSystem
Int_t MassFitter(Bool_t draw=kTRUE)
AliHFInvMassFitter class for the fit of invariant mass distribution of charm hadrons.
void PrintGausParams(TH1F *hPulls)
void SetFitOption(TString opt)
TCanvas * c
Definition: TestFitELoss.C:172
Int_t optBkgBinCount
static TH1 * AdaptTemplateRangeAndBinning(const TH1 *hRef, TH1 *hData, Double_t minFit, Double_t maxFit)
Bool_t renormLS
Int_t smoothLS
Double_t GetSigma() const
Int_t rebin[maxPtBins]
Double_t GetBackgroundNormalizationFactor(TH1D *hRatio, Int_t reb=1)
TString fileNameMC
void SetPolDegreeForBackgroundFit(Int_t deg)
Bool_t fixMean[maxPtBins]
TFile * fil
Definition: InvMassFit.C:61
void ProjectCombinHFAndFit()
void SetFixGaussianSigma(Double_t sigma)
Double_t fitSBrangelow[maxPtBins]
TString suffix
Double_t * sigma
Double_t tuneMeanOnData
void SetStyleHisto(TH1 *h, Int_t method, Int_t isXpt=-1)
Double_t massD
Bool_t tryDirectFit
Bool_t correctForRefl
void SetFixGaussianMean(Double_t mean)
Int_t fixMeanConf
int Int_t
Definition: External.C:63
Double_t GetMean() const
Int_t method
TString suffixMC
void WriteFitFunctionsToFile(AliHFInvMassFitter *fitter, TString meth, Int_t iPtBin)
Definition: External.C:212
Bool_t QuadraticSmooth(TH1 *h, Int_t ntimes=1)
Double_t nsigmaBinCounting
TH1D * hMCReflPtBin
Double_t GetRawYieldBinCounting(Double_t &errRyBC, Double_t nSigma=3., Int_t option=0, Int_t pdgCode=0) const
Double_t binLims[maxPtBins+1]
TH1F * SetTemplateReflections(const TH1 *h, TString opt, Double_t minRange, Double_t maxRange)
TList * fitter
Definition: DrawAnaELoss.C:26
Bool_t useEMwithLS
AliHFInvMassFitter * ConfigureFitter(TH1D *histo, Int_t iPtBin, Int_t backcase, Double_t minFit, Double_t maxFit, Bool_t isDirect=kFALSE)
Double_t GetReflOverSigUncertainty() const
Bool_t saveCanvasAsRoot
void DrawHistoMinusFit(TVirtualPad *c, Int_t writeFitInfo=1)
Double_t costhstcut
TString fitoption
Double_t GetSigmaUncertainty() const
Double_t tuneSigmaOnData
Double_t GetRawYieldError() const
TH1D * hMCSigPtBin
void Background(Double_t nOfSigma, Double_t &background, Double_t &errbackground) const
void SetInitialGaussianSigma(Double_t sigma)
void SetInitialGaussianMean(Double_t mean)
Int_t fixSigmaConf
Double_t minMass
Int_t optForNorm
const Double_t ymin
Definition: AddTaskCFDStar.C:6
Bool_t useGeomMeanLS
Double_t maxMass
Double_t minMass4Fit[maxPtBins]
bool Bool_t
Definition: External.C:53
TFile * fout
input train file
TF1 * GausPlusLine(Double_t minRange=1.72, Double_t maxRange=2.05)
Int_t nPtBins
Int_t nDegreeBackPol[maxPtBins]
TH1F * FitMCInvMassSpectra(TList *lMC, TString var)
const Int_t maxPtBins
TString meson
Int_t typeb
Double_t rangeForNorm
void PrintConfig()
Bool_t fixSigma[maxPtBins]
TH1F * GetHistoClone() const
Definition: External.C:196
void DrawHere(TVirtualPad *c, Double_t nsigma=3, Int_t writeFitInfo=2)
Double_t GetRawYield() const
void WriteFitInfo(AliHFInvMassFitter *fitter, TH1D *histo)