11.2. Gaussian Fits

Table of Contents

Fitting histogram peaks

Gaussian fits are used when your histogram peak is shaped like a normal distribution, for example a detector energy peak or a time-of-flight signal. In ROOT, the typical workflow is to have a 1D histogram, then fit a Gaussian function to a region around the peak in order to extract its position and width.

Assume you already have a histogram, for example a TH1F *h filled with your data. The simplest way to fit a Gaussian peak is to call the Fit method on the histogram with the built-in function name "gaus":

cpp

h->Fit("gaus");

This creates, behind the scenes, a temporary TF1 with the Gaussian shape and fits it to the whole histogram range. After the fit, you can draw the histogram and the fitted function on a canvas:

cpp

h->Draw();

ROOT automatically overlays the Gaussian curve and prints the fit results in the canvas statistics box if the style allows it.

The underlying Gaussian function used by "gaus" has the form

$$
f(x) = p_0 \exp\left[-\frac{1}{2} \left(\frac{x - p_1}{p_2}\right)^2 \right],
$$

where the three parameters are:

Parameter index | Meaning
--------------- | -------
$p_0$ | Amplitude (peak height)
$p_1$ | Mean (peak position)
$p_2$ | Sigma (standard deviation, peak width)

To have more control, you can explicitly create a TF1:

cpp

TF1 *g = new TF1("g", "gaus", xMin, xMax);
h->Fit(g, "R");  // "R" tells ROOT to respect the range of g

This is very useful for fitting only one peak in a complex spectrum or for restricting the fit to a region that is approximately Gaussian. The function range and the initial guesses of parameters influence how well the fit converges, which you will see in the section on the fit range.

After the fit, all relevant information is stored both in the histogram and in the TF1 object. The histogram keeps a pointer to the last fitted function, accessible with

cpp

TF1 *fitFunc = h->GetFunction("gaus");  // or "g" if that was its name

You usually extract the peak position and width from this function, as described in the next sections.

A Gaussian fit with "gaus" always has three parameters:

$p_0$ amplitude.
$p_1$ mean.
$p_2$ sigma.
Always interpret the parameters in this order and use the corresponding index when reading them from the fitted function.

Mean

In the context of a Gaussian fit, the mean parameter $p_1$ gives the position of the center of the peak. This is often the main physical quantity you want, for example a calibrated energy, a mass or a time offset.

Once you have performed a Gaussian fit on a histogram h, you can retrieve the fitted mean using the GetParameter method of the fit function. If you used the automatic "gaus" fit,

cpp

h->Fit("gaus");
TF1 *fgaus = h->GetFunction("gaus");
double mean  = fgaus->GetParameter(1);
double emean = fgaus->GetParError(1);

Here mean is the fitted peak position and emean is the statistical uncertainty given by the fit. The index 1 is always the mean for the Gaussian built-in function.

In many analyses, you will compare the Gaussian mean to a known or expected value. For example, in a calibration run you might fit peaks from known reference lines and then use the fitted means as input to a calibration function. The uncertainty emean then propagates into the total calibration uncertainty.

If you create your own Gaussian TF1, the interpretation of parameter 1 as the mean remains the same as long as you use the "gaus" formula:

cpp

TF1 *g = new TF1("g", "gaus", xMin, xMax);
h->Fit(g, "R");
double mean  = g->GetParameter(1);
double emean = g->GetParError(1);

If you later define more complex models, for example a Gaussian plus a constant background, the mean will no longer be parameter 1 by default, so you should keep track of your parameter ordering. For the simple built-in "gaus" function, however, the mean is always parameter 1.

For the ROOT built-in Gaussian "gaus":

The mean of the peak is always parameter index 1.
Use GetParameter(1) for the mean value and GetParError(1) for its uncertainty.

Sigma

The sigma parameter $p_2$ is the standard deviation of the Gaussian and describes how wide the peak is. In physics applications, sigma is directly related to the detector resolution or to the spread of the measured quantity.

After fitting, you can read sigma and its uncertainty just like the mean:

cpp

h->Fit("gaus");
TF1 *fgaus = h->GetFunction("gaus");
double sigma  = fgaus->GetParameter(2);
double esigma = fgaus->GetParError(2);

The width of the peak is sometimes more conveniently expressed in terms of the full width at half maximum, FWHM. For a perfect Gaussian, FWHM is related to sigma by the formula

$$
\text{FWHM} = 2 \sqrt{2 \ln 2}\,\sigma \approx 2.355\,\sigma.
$$

You can compute this after the fit:

cpp

double fwhm = 2.35482 * sigma;

Sigma allows you to quantify how well your detector or measurement technique performs. A smaller sigma indicates better resolution. By comparing sigma for different settings or different detectors, you can evaluate performance improvements.

If your histogram bins are in physical units, then sigma and FWHM are also in these units. For example, if the x axis is calibrated in MeV, then sigma represents the energy resolution in MeV.

As with the mean, if you create your own multi-parameter model, the index 2 might no longer correspond to sigma, but for a pure "gaus" function it always does.

For a Gaussian distribution:

The sigma (standard deviation) of the peak is parameter index 2 for "gaus".
The full width at half maximum is related by
$$\text{FWHM} \approx 2.355\,\sigma.$$

Fit range

The choice of fit range is crucial for a reliable Gaussian fit. If you fit too wide a range, background structures or non-Gaussian tails can distort the result. If you fit too narrow a range, you may not capture enough of the peak shape for a stable fit.

By default,

cpp

h->Fit("gaus");

uses the full histogram x range. In many real data situations this is not optimal, especially when several peaks or a sloping background are present.

To specify a fit range explicitly, define a TF1 with a given x range and then fit with the "R" option so ROOT respects that range:

cpp

double xlow  = 100.0;
double xhigh = 140.0;
TF1 *g = new TF1("g", "gaus", xlow, xhigh);
h->Fit(g, "R");

Here, only bins whose centers lie between xlow and xhigh are used in the fit. This is the standard way to isolate a single peak.

A common practical strategy is:

Inspect the histogram visually to estimate the peak region.
Choose a range that covers about plus or minus two to three sigmas around the peak center.
Set that as the initial function range and perform the fit.
Optionally refine the range based on the fitted mean and sigma.

If the peak position is not known beforehand, you can still start with a rough range that clearly covers the whole peak and a bit of the nearby background. In some cases it helps to initialize the Gaussian parameters before fitting:

cpp

TF1 *g = new TF1("g", "gaus", xlow, xhigh);
g->SetParameters(h->GetMaximum(), peakGuess, widthGuess);
h->Fit(g, "R");

Here peakGuess and widthGuess are your approximate estimates of the mean and sigma. Good initial values and a sensible fit range together make the fit more stable and less likely to converge to an unphysical solution.

You can also change the fit range of an existing TF1 without recreating it:

cpp

g->SetRange(newLow, newHigh);
h->Fit(g, "R");

Finally, remember that the statistical quality of the fit depends on the number of entries in the fit range. If you choose a very narrow range around a small peak, the uncertainty on the mean and sigma will be large. In that case, collecting more data or widening the range slightly can help.

When fitting a Gaussian peak:

Do not blindly fit the full histogram range.
Define a reasonable fit range that covers the peak but avoids unrelated structures.
Use a TF1 with the desired range and the "R" fit option to constrain the fit.

11.2. Gaussian Fits

Fitting histogram peaks

Mean

Sigma

Fit range

Comments

Where to Move