AliPhysics  c923f52 (c923f52)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
DrawFpromptVsRaaElossHypoCombined.C
Go to the documentation of this file.
1 #if !defined(__CINT__) || defined(__MAKECINT__)
2 #include <Riostream.h>
3 #include "TH1D.h"
4 #include "TH1.h"
5 #include "TH2F.h"
6 #include "TNtuple.h"
7 #include "TFile.h"
8 #include "TGraphAsymmErrors.h"
9 #include "TCanvas.h"
10 #include "TROOT.h"
11 #include "TStyle.h"
12 #include "TLegend.h"
13 #include "TMath.h"
14 #endif
15 
16 /******************************************************
17  *
18  * Macro to compute fprompt for pPb and PbPb data
19  * using the outputs of the HFPtSpectrumRaa maco
20  * Questions and complains to Z.Conesa del Valle
21  *
22  *******************************************************/
23 
25 
26 TGraphAsymmErrors * DrawFpromptVsRaaElossHypo(const char *infile="", //const char *outfile="",
27  Double_t MinHypo=1./3., Double_t MaxHypo=3.0, Double_t CentralHypo = 1.0, Bool_t isYieldUnc=true);
28 
29 
30 //_______________________________________________________________________
31 void DrawFpromptVsRaaElossHypoCombined(const char *infileNb="", const char *infilefc="", const char *outfile="", Double_t MinHypo=1./3., Double_t MaxHypo=3.0, Double_t CentralHypo = 1.0, Bool_t isYieldUnc=true)
32 {
33 
34  TH2F *hdraw = new TH2F("hdraw","fprompt vs pt; p_{T} (GeV/c); f_{prompt}",36,0.,36.,10.,0.,1);
35  hdraw->Draw();
36 
37  cout<<endl<<endl<<" ********* Checking Nb method ********* "<<endl<<endl;
38  TGraphAsymmErrors *gfPromptNb = DrawFpromptVsRaaElossHypo(infileNb,MinHypo,MaxHypo,CentralHypo,isYieldUnc);
39  gfPromptNb->SetNameTitle("gfPromptNb","fprompt Nb method, Eloss range considered");
40  gfPromptNb->SetLineColor(kBlack);
41  gfPromptNb->SetMarkerColor(kGreen+2);
42  gfPromptNb->SetMarkerStyle(25);
43  gfPromptNb->Draw("P");
44 
45  cout<<endl<<endl<<" ********* Checking fc method ********* "<<endl<<endl;
46  TGraphAsymmErrors *gfPromptfc = DrawFpromptVsRaaElossHypo(infilefc,MinHypo,MaxHypo,CentralHypo,isYieldUnc);
47  gfPromptfc->SetNameTitle("gfPromptfc","fprompt fc method, Eloss range considered");
48  gfPromptfc->SetLineColor(kBlack);
49  gfPromptfc->SetMarkerColor(kRed);
50  gfPromptfc->SetMarkerStyle(21);
51  gfPromptfc->Draw("P");
52 
53  cout<<endl<<endl<<" ********* Computing the envelope ******** "<<endl<<endl;
54  TGraphAsymmErrors *gFPromptCombined = new TGraphAsymmErrors(0);
55  gFPromptCombined->SetNameTitle("gFPromptCombined","fprompt Nb + fc, Eloss range considered");
56 
57  Int_t nbins=gfPromptNb->GetN();
58  for(Int_t i=0; i<nbins; i++) {
59  Double_t pt=0.,value=0.,ptfc=0.,valuefc=0.,exlow=0.,eylow=0.,eyNblow=0.,eyfclow=0.,eyhigh=0.,eyNbhigh=0.,eyfchigh=0.;
60  gfPromptNb->GetPoint(i,pt,value);
61  gfPromptfc->GetPoint(i,ptfc,valuefc);
62  exlow = gfPromptNb->GetErrorXlow(i);
63  eyNblow = gfPromptNb->GetErrorYlow(i);
64  eyNbhigh = gfPromptNb->GetErrorYhigh(i);
65  eyfclow = gfPromptfc->GetErrorYlow(i);
66  eyfchigh = gfPromptfc->GetErrorYhigh(i);
67  Double_t uncertainties[5]={ value, value-eyNblow, value+eyNbhigh, valuefc-eyfclow, valuefc+eyfchigh };
68  eylow = value - TMath::MinElement(5,uncertainties);
69  eyhigh = TMath::MaxElement(5,uncertainties) - value;
70  gFPromptCombined->SetPoint(i,pt,value);
71  gFPromptCombined->SetPointError(i,exlow,exlow,eylow,eyhigh);
72  }
73 
74  gFPromptCombined->SetFillStyle(0);
75  gFPromptCombined->SetLineColor(kBlack);
76  gFPromptCombined->SetMarkerColor(kBlack);
77  gFPromptCombined->SetMarkerStyle(20);
78  gFPromptCombined->Draw("2P");
79 
80  cout<<endl<<endl;
81 
82  TFile* fout= new TFile(outfile,"recreate");
83  gfPromptNb->Write();
84  gfPromptfc->Write();
85  gFPromptCombined->Write();
86  fout->Write();
87 }
88 
89 //_______________________________________________________________________
91  //const char *outfile="",
92  Double_t MinHypo, Double_t MaxHypo, Double_t CentralHypo, Bool_t isYieldUnc)
93 {
94 
95  Float_t pt=0.,TAB=0.,sigmaPP=0.,invyieldAB=0.,RaaCharm=0.,RaaBeauty=0., fc=0.;
96  Float_t yieldFdHigh=0., yieldFdLow=0.,RaaCharmFdHigh=0.,RaaCharmFdLow=0.;
97 
98  TFile *fRaa = new TFile(infile,"read");
99  TH1D *hRABvsPt = (TH1D*)fRaa->Get("hRABvsPt");
100  TNtuple *ntRaa = (TNtuple*)fRaa->Get("ntupleRAB");
101  ntRaa->SetBranchAddress("pt",&pt);
102  ntRaa->SetBranchAddress("TAB",&TAB);
103  ntRaa->SetBranchAddress("sigmaPP",&sigmaPP);
104  ntRaa->SetBranchAddress("invyieldAB",&invyieldAB);
105  ntRaa->SetBranchAddress("RABCharm",&RaaCharm);
106  ntRaa->SetBranchAddress("RABBeauty",&RaaBeauty);
107  ntRaa->SetBranchAddress("RABCharmFDHigh",&RaaCharmFdHigh);
108  ntRaa->SetBranchAddress("RABCharmFDLow",&RaaCharmFdLow);
109  ntRaa->SetBranchAddress("invyieldABFDHigh",&yieldFdHigh);
110  ntRaa->SetBranchAddress("invyieldABFDLow",&yieldFdLow);
111  ntRaa->SetBranchAddress("fc",&fc);
112 
113  // define the binning
114  Int_t nbins = hRABvsPt->GetNbinsX();
115 
116  //
117  //
118  // Search the central value of the energy loss hypothesis Rb = Rc (bin)
119  //
120  Double_t ElossCentral[nbins+1];
121  Double_t fcV[nbins+1], fcVmax[nbins+1], fcVmin[nbins+1];
122  Double_t fdMax[nbins+1], fdMin[nbins+1];
123  for(Int_t i=0; i<=nbins; i++) {
124  ElossCentral[i]=0.;
125  fcV[i]=0.; fcVmax[i]=0.; fcVmin[i]=6.;
126  fdMax[i]=0.; fdMin[i]=6.;
127  }
128  //
129  Int_t netries = ntRaa->GetEntries();
130  //
131  for(Int_t ientry=0; ientry<=netries; ientry++){
132  ntRaa->GetEntry(ientry);
133  Double_t ElossHypo = 1. / (RaaCharm / RaaBeauty) ;
134  //
135  // Find the bin for the central Eloss hypo
136  //
137  if( TMath::Abs( ElossHypo - CentralHypo ) < 0.075 ){
138  Int_t hABbin = hRABvsPt->FindBin( pt );
139  Double_t DeltaIni = TMath::Abs( ElossCentral[ hABbin ] - CentralHypo );
140  Double_t DeltaV = TMath::Abs( ElossHypo - CentralHypo );
141  // cout << " pt " << ptTot << " ECentral Tot " << ElossCentralTot[ hABbin ] << " Ehypo "<< ElossHypo ;
142  if ( DeltaV < DeltaIni ) ElossCentral[ hABbin ] = ElossHypo;
143  // cout << " final ECentral " << ElossCentralTot[ hABbin ] << endl;
144  }
145  }
146 
147  //
148  for(Int_t ientry=0; ientry<=netries; ientry++){
149 
150  ntRaa->GetEntry(ientry);
151  Double_t ElossHypo = 1. / (RaaCharm / RaaBeauty) ;
152  Int_t hABbin = hRABvsPt->FindBin( pt );
153 
154  // Get the central value of Raa
155  if ( ElossHypo == ElossCentral[ hABbin ] ) {
156 
157  fcV[hABbin] = fc;
158 
159  // Now look at FONLL band for the central value
160  Double_t elems[3]= {0., 0., 0.};
161  // fdMax[hABbin]= (RaaCharmFdHigh-RaaCharm)/RaaCharm;
162  // fdMin[hABbin]= (RaaCharm-RaaCharmFdLow)/RaaCharm;
163  if(isYieldUnc) {
164  elems[0]=invyieldAB;
165  elems[1]=yieldFdHigh;
166  elems[2]=yieldFdLow;
167  } else {
168  elems[0]=RaaCharm;
169  elems[1]=RaaCharmFdHigh;
170  elems[2]=RaaCharmFdLow;
171  }
172  fdMax[hABbin]= (TMath::MaxElement(3,elems)-elems[0])/elems[0];
173  fdMin[hABbin]= (elems[0]-TMath::MinElement(3,elems))/elems[0];
174  // cout<<" Check FONLL band giving for pt="<< pt <<" +"<< fdMax[hABbin]<< " -"<< fdMin[hABbin]<<" (relative variation)"<<endl;
175  }
176 
177  // Here check the Eloss band
178  if( RaaCharm>0 && ElossHypo >= MinHypo && ElossHypo <=MaxHypo ) {
179  if( fc < fcVmin[hABbin] ) fcVmin[hABbin] = fc;
180  if( fc > fcVmax[hABbin] ) fcVmax[hABbin] = fc;
181  }
182  }
183 
184  TGraphAsymmErrors *gFPrompt = new TGraphAsymmErrors(0);
185  Int_t j=0;
186  for(Int_t i=1; i<=nbins; i++) {
187  Double_t xpt = hRABvsPt->GetBinCenter(i);
188  Double_t width = hRABvsPt->GetBinWidth(i)/2;
189  Double_t elossLow = (fcV[i]-fcVmin[i])/fcV[i];
190  Double_t elossHigh = (fcVmax[i]-fcV[i])/fcV[i];
191  Double_t fcLow = fcV[i]-fcVmin[i];
192  Double_t fcHigh = fcVmax[i]-fcV[i];
193  // cout<<" pt="<<xpt<<" , fprompt="<<fcV[i]<<" ["<< fcVmin[i] <<","<<fcVmax[i] <<"] ( Eloss) "<<endl;
194 
195  if(elossFDQuadSum){ // add quadratically
196  fcLow = TMath::Sqrt( fdMin[i]*fdMin[i] + elossLow*elossLow )*fcV[i];
197  fcHigh = TMath::Sqrt( fdMax[i]*fdMax[i] + elossHigh*elossHigh )*fcV[i];
198  } else { // add linearly
199  fcLow = ( fdMin[i] + elossLow )*fcV[i];
200  fcHigh = ( fdMax[i] + elossHigh )*fcV[i];
201  }
202 
203  gFPrompt->SetPoint(j,xpt,fcV[i]);
204  // gFPrompt->SetPointError(j,width,width,fcV[i]-fcVmin[i],fcVmax[i]-fcV[i]);
205  gFPrompt->SetPointError(j,width,width,fcLow,fcHigh);
206  j++;
207 
208  cout<<" pt="<<xpt<<" , fprompt="<<fcV[i]<<" +"<<fcHigh<<" -"<<fcLow<<" ( splitted it is +"<<elossLow*fcV[i]<<" -"<<elossHigh*fcV[i]<<" Eloss, +"<<fdMax[i]*fcV[i]<<" -"<<fdMin[i]*fcV[i]<<" FONLL-band )"<<endl;
209  }
210  // gFPrompt->Draw("AP");
211 
212  // if(strcmp(outfile,"")!=0) {
213  // TFile *out = new TFile(outfile,"recreate");
214  // gFPrompt->Write();
215  // }
216 
217  fRaa->Close();
218 
219  return gFPrompt;
220 
221 }
double Double_t
Definition: External.C:58
Definition: External.C:236
TGraphAsymmErrors * DrawFpromptVsRaaElossHypo(const char *infile="", Double_t MinHypo=1./3., Double_t MaxHypo=3.0, Double_t CentralHypo=1.0, Bool_t isYieldUnc=true)
void DrawFpromptVsRaaElossHypoCombined(const char *infileNb="", const char *infilefc="", const char *outfile="", Double_t MinHypo=1./3., Double_t MaxHypo=3.0, Double_t CentralHypo=1.0, Bool_t isYieldUnc=true)
int Int_t
Definition: External.C:63
float Float_t
Definition: External.C:68
Definition: External.C:212
const Int_t nbins
bool Bool_t
Definition: External.C:53