AliPhysics  master (3d17d9d)
PlotUtils.C
Go to the documentation of this file.
1 
12 #if !defined(__CINT__) || defined(__MAKECINT__)
13 
14 #include <TH1D.h>
15 #include <TH2F.h>
16 #include <TF1.h>
17 #include <TMath.h>
18 #include <TGraphErrors.h>
19 
20 #endif
21 
25 //-----------------------------------------------------------------------------
27 {
28  Double_t N = par[0];
29  Double_t width = par[1];
30 
31  Double_t mean = par[2];
32  Double_t sigma = par[3];
33  Double_t alpha = par[4];
34  Double_t n = par[5];
35 
36  Double_t A = pow(n/fabs(alpha),n)*exp(-pow(alpha,2)/2);
37  Double_t B = n/fabs(alpha) - fabs(alpha);
38 
39  if ((x[0]-mean)/sigma>-alpha)
40  return N*width*TMath::Gaus(x[0],mean,sigma,1);
41  else
42  return N/(sqrt(TMath::TwoPi())*sigma)*width*A*pow(B-(x[0]-mean)/sigma,-n);
43 }
44 
48 //-----------------------------------------------------------------------------
50 {
51  Double_t gaus = par[0] * TMath::Exp( -(x[0]-par[1])*(x[0]-par[1]) /
52  (2*par[2]*par[2]) );
53  Double_t back = par[3] + par[4]*x[0] + par[5]*x[0]*x[0] + par[6]*x[0]*x[0]*x[0];
54  return gaus+back;
55 }
56 
60 //-----------------------------------------------------------------------------
62 {
63  Double_t gaus = par[0] * TMath::Exp( -(x[0]-par[1])*(x[0]-par[1]) /
64  (2*par[2]*par[2]) );
65  Double_t back = par[3] + par[4]*x[0] + par[5]*x[0]*x[0];
66  return gaus+back;
67 }
68 
72 //-----------------------------------------------------------------------------
74 {
75  Double_t gaus = par[0] * TMath::Exp( -(x[0]-par[1])*(x[0]-par[1]) /
76  (2*par[2]*par[2]) );
77  Double_t back = par[3] + par[4]*x[0];
78  return gaus+back;
79 }
80 
84 //-----------------------------------------------------------------------------
86 {
87  Double_t gaus = par[0] * TMath::Exp( -(x[0]-par[1])*(x[0]-par[1]) /
88  (2*par[2]*par[2]) );
89  Double_t back = 0;//par[3];
90  return gaus+back;
91 }
92 
96 //-----------------------------------------------------------------------------
98 {
99  if ( x[0] > 0.425 && x[0] < 0.65 )
100  {
101  TF1::RejectPoint();
102  return 0;
103  }
104 
105  return par[0] + par[1]*x[0];
106  // + x[0]*x[0]*par[2] + x[0]*x[0]*x[0]*par[3];//+ x[0]*x[0]*x[0]*x[0]*par[4];
107 }
108 
112 //-----------------------------------------------------------------------------
114 {
115  //if ( x[0] > 0.07 && x[0] < 0.2 )
116  if ( x[0] < 0.25 )
117  {
118  TF1::RejectPoint();
119  return 0;
120  }
121 
122  return par[0] + par[1]*x[0];
123  // + x[0]*x[0]*par[2] + x[0]*x[0]*x[0]*par[3] ;//+ x[0]*x[0]*x[0]*x[0]*par[4];
124 }
125 
126 //-----------------------------------------------------------------------------
131 //-----------------------------------------------------------------------------
133 (Double_t num , Double_t den,
134  Double_t numErr, Double_t denErr)
135 {
136  if ( num == 0 || den == 0 ) return 0.;
137 
138  //printf("\t num %e den %e numErr %e denErr %e\n",num,den,numErr,denErr);
139 
140  return num/den * TMath::Sqrt( ( numErr * numErr ) / ( num * num ) +
141  ( denErr * denErr ) / ( den * den ) );
142 }
143 
144 //-----------------------------------------------------------------------------
150 //-----------------------------------------------------------------------------
151 static TGraphErrors * DivideGraphs
153 {
154  if ( !gDen || !gNum )
155  {
156  printf("Graph num %p or graph den %p not available\n",gNum,gDen);
157  return 0x0;
158  }
159  const Int_t nBins = gNum->GetN();
160  if ( nBins != gDen->GetN() )
161  {
162  printf("Cannot divide %s with %d bins and %s with %d bins!\n",
163  gNum->GetName(),nBins,gDen->GetName(),gDen->GetN());
164  return 0x0;
165  }
166 
167  Double_t ratio [nBins];
168  Double_t ratioErr[nBins];
169  Double_t x [nBins];
170  Double_t xErr [nBins];
171  for (Int_t ibin = 0; ibin < nBins; ibin++)
172  {
173  Double_t num = gNum->GetY ()[ibin];
174  Double_t den = gDen->GetY ()[ibin];
175  Double_t numErr = gNum->GetEY()[ibin];
176  Double_t denErr = gDen->GetEY()[ibin];
177 
178  x [ibin] = gNum->GetX ()[ibin];
179  xErr[ibin] = gNum->GetEX()[ibin];
180 
181  if ( num == 0 || den == 0 )
182  {
183  ratio [ibin] = 0;
184  ratioErr[ibin] = 0;
185  continue;
186  }
187 
188  ratio [ibin] = num / den ;
189  ratioErr[ibin] = GetFractionError(num,den,numErr,denErr);
190 
191  // printf("bin %d, x %f (%f) num %f (%f), den %f (%f), ratio %f (%f) \n",
192  // ibin,x[ibin],xErr[ibin],num,numErr,den,denErr,ratio[ibin],ratioErr[ibin]);
193  } // do the ratio to sum
194 
195  return new TGraphErrors(nBins,x,ratio,xErr,ratioErr);
196 }
197 
198 
199 //-----------------------------------------------------------------------------
207 //-----------------------------------------------------------------------------
208 static void GetRangeIntegralAndError
209 (TH1D* h, Int_t binMin, Int_t binMax,
210  Double_t & integral, Double_t & integralErr )
211 {
212  integral = 0;
213  integralErr = 0;
214  for(Int_t ibin = binMin; ibin <= binMax; ibin++)
215  {
216  //if ( h->GetBinContent(ibin) == 0 ) continue ;
217 
218  integral += h->GetBinContent(ibin);
219  integralErr += h->GetBinError (ibin) * h->GetBinError(ibin);
220  }
221 // printf("\t bin range [%d,%d], n %d, sum %2.2e err %2.2e\n",
222 // binMin,binMax,n,integral,integralErr);
223  if ( integralErr > 0 ) integralErr = TMath::Sqrt(integralErr);
224 }
225 
226 //-----------------------------------------------------------------------------
234 //-----------------------------------------------------------------------------
235 static void GetRangeIntegralAndError
236 (TH1D* h, Float_t minX, Float_t maxX,
237  Double_t & integral, Double_t & integralErr )
238 {
239  Int_t binMin = h->FindBin(minX);
240  Int_t binMax = h->FindBin(maxX)-1;
241  return GetRangeIntegralAndError(h, binMin, binMax,
242  integral, integralErr );
243 }
244 
245 //-----------------------------------------------------------------------------
256 //-----------------------------------------------------------------------------
258 (TGraphErrors* graph,
259  Double_t min, Double_t max,
260  Float_t errFraction = 0.5,
261  Float_t value = -1, Float_t valueErr = 0 )
262 {
263  for(Int_t ibin = 0; ibin < graph->GetN(); ibin++)
264  {
265  Float_t x = graph->GetX()[ibin];
266  Float_t y = graph->GetY()[ibin];
267  Float_t yE = graph->GetEY()[ibin];
268  if(x < min || x > max )
269  {
270  graph->SetPoint (ibin,x, value );
271  graph->SetPointError(ibin,x, valueErr);
272  }
273  else
274  if(yE > errFraction*y) // Remove points with more than 50% error
275  {
276  //printf("ibin %d, y %f, yE %f\n",ibin,y,yE);
277  graph->SetPoint (ibin,x,value );
278  graph->SetPointError(ibin,x,valueErr);
279  }
280  }
281 }
282 
283 //-----------------------------------------------------------------------------
285 //-----------------------------------------------------------------------------
286 static void ScaleBinBySize(TH1D* h)
287 {
288  for(Int_t ibin = 1; ibin < h->GetNbinsX();ibin++)
289  {
290  Double_t width = h->GetBinWidth(ibin);
291  Double_t content = h->GetBinContent(ibin);
292  Double_t error = h->GetBinError(ibin);
293 
294  //printf("bin %d, width %f, content %e\n",ibin,width,content);
295  h->SetBinContent(ibin,content/width);
296  h->SetBinError (ibin, error/width);
297  }
298 }
299 
300 //-----------------------------------------------------------------------------
302 //-----------------------------------------------------------------------------
303 void ScaleXaxis2D(TH2F* h2D)
304 {
305  Int_t nbinsy = h2D->GetNbinsY();
306 
307  for(Int_t j = 1; j <= h2D->GetNbinsX(); j++)
308  {
309  TH1D* temp1 = h2D->ProjectionY(Form("Bin%d",j),j,j);
310 
311  Float_t scale1 = temp1 -> Integral(-1,-1);
312 
313  for(Int_t i = 1; i <= nbinsy; i++)
314  {
315  //printf("i %d, j %d; content %f / scale %f = %f\n",
316  //i,j,h2D->GetBinContent(j,i),scale2,h2D->GetBinContent(j,i)/scale2);
317 
318  if ( scale1 > 0 )
319  {
320  h2D->SetBinContent(j,i, h2D->GetBinContent(j,i)/scale1);
321  h2D->SetBinError (j,i, h2D->GetBinError (j,i)/scale1);
322  }
323  else
324  {
325  h2D->SetBinContent(j,i, 0);
326  h2D->SetBinError (j,i, 0);
327  }
328 
329  } // y bin loop
330  } // x bin loop
331 
332  h2D->SetZTitle("x bin norm. to integral");
333 }
334 
335 
336 //-----------------------------------------------------------------------------
343 //-----------------------------------------------------------------------------
344 void GetCanvasColRowNumber(Int_t npad, Int_t & ncol, Int_t & nrow,
345  Bool_t square = kFALSE)
346 {
347  if(npad <= 0)
348  {
349  ncol = 0;
350  nrow = 0;
351  return;
352  }
353 
354  ncol = TMath::Sqrt(npad);
355  nrow = ncol;
356 
357  if ( square )
358  return;
359 
360  if ( npad == 1 )
361  {
362  ncol = 1;
363  nrow = 1;
364  }
365  else if ( npad == 2 )
366  {
367  ncol = 2;
368  nrow = 1;
369  }
370  else if ( npad == 3 )
371  {
372  ncol = 3;
373  nrow = 1;
374  }
375  else if ( npad == 4 )
376  {
377  ncol = 2;
378  nrow = 2;
379  }
380  else if ( npad < 7 )
381  {
382  ncol = 3;
383  nrow = 2;
384  }
385  else if ( npad < 10 )
386  {
387  ncol = 3;
388  nrow = 3;
389  }
390  else if ( npad < 13 )
391  {
392  ncol = 4;
393  nrow = 3;
394  }
395  else if ( npad < 17 )
396  {
397  ncol = 4;
398  nrow = 4;
399  }
400  else if ( npad < 21 )
401  {
402  ncol = 5;
403  nrow = 4;
404  }
405  else if ( npad < 26 )
406  {
407  ncol = 5;
408  nrow = 5;
409  }
410  else if ( npad < 31 )
411  {
412  ncol = 6;
413  nrow = 5;
414  }
415  else if ( npad < 37 )
416  {
417  ncol = 6;
418  nrow = 6;
419  }
420 
421 }
422 
static Double_t FunctionCrystalBall(Double_t *x, Double_t *par)
Definition: PlotUtils.C:26
static Double_t FunctionTruncatedPolPi0(Double_t *x, Double_t *par)
Definition: PlotUtils.C:113
double Double_t
Definition: External.C:58
Definition: External.C:236
static void ScaleBinBySize(TH1D *h)
Scale histogram bins by its size.
Definition: PlotUtils.C:286
Int_t nbinsy
static Double_t FunctionGaussPol0(Double_t *x, Double_t *par)
Definition: PlotUtils.C:85
void ScaleXaxis2D(TH2F *h2D)
Scale 2D histogram X bins by its integral.
Definition: PlotUtils.C:303
static Double_t FunctionGaussPol1(Double_t *x, Double_t *par)
Definition: PlotUtils.C:73
static TGraphErrors * DivideGraphs(TGraphErrors *gNum, TGraphErrors *gDen)
Definition: PlotUtils.C:152
static Double_t GetFractionError(Double_t num, Double_t den, Double_t numErr, Double_t denErr)
Definition: PlotUtils.C:133
void RemovePointsOutOfRangeOrLargeErrorFromGraph(TGraphErrors *graph, Double_t min, Double_t max, Float_t errFraction=0.5, Float_t value=-1, Float_t valueErr=0)
Definition: PlotUtils.C:258
Double_t * sigma
int Int_t
Definition: External.C:63
float Float_t
Definition: External.C:68
Definition: External.C:212
static Double_t FunctionGaussPol3(Double_t *x, Double_t *par)
Definition: PlotUtils.C:49
static void GetRangeIntegralAndError(TH1D *h, Int_t binMin, Int_t binMax, Double_t &integral, Double_t &integralErr)
Definition: PlotUtils.C:209
static Double_t FunctionGaussPol2(Double_t *x, Double_t *par)
Definition: PlotUtils.C:61
static Double_t FunctionTruncatedPolEta(Double_t *x, Double_t *par)
Definition: PlotUtils.C:97
void GetCanvasColRowNumber(Int_t npad, Int_t &ncol, Int_t &nrow, Bool_t square=kFALSE)
Definition: PlotUtils.C:344
bool Bool_t
Definition: External.C:53