Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
fitConvolution.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_fit
3## \notebook
4## Tutorial for convolution of two functions
5##
6## \macro_image
7## \macro_output
8## \macro_code
9##
10## \author Jonas Rembser, Aurelie Flandi (C++ version)
11
12import ROOT
13import numpy as np
14
15# Construction of histogram to fit.
16h_ExpGauss = ROOT.TH1F("h_ExpGauss", "Exponential convoluted by Gaussian", 100, 0.0, 5.0)
17h_ExpGauss.Fill(
18 np.array(
19 [
20 # Gives a alpha of -0.3 in the exp.
21 # Probability density function of the addition of two variables is the
22 # convolution of two density functions.
23 ROOT.gRandom.Exp(1.0 / 0.3) + ROOT.gRandom.Gaus(0.0, 3.0)
24 for _ in range(1000000)
25 ]
26 )
27)
28f_conv = ROOT.TF1Convolution("expo", "gaus", -1, 6, True)
29f_conv.SetRange(-1.0, 6.0)
30f_conv.SetNofPointsFFT(1000)
31f = ROOT.TF1("f", f_conv, 0.0, 5.0, f_conv.GetNpar())
32f.SetParameters(1.0, -0.3, 0.0, 1.0)
33
34c1 = ROOT.TCanvas("c1", "c1", 800, 1000)
35
36# Fit and draw result of the fit
37h_ExpGauss.Fit("f")
38
39c1.SaveAs("fitConvolution.png")