16.3. Simple Monte Carlo Simulation

Table of Contents

Generating events

A simple Monte Carlo simulation in ROOT follows a clear pattern. You define a model for how events should look, you generate many random events according to that model, and then you study the resulting distributions with histograms.

In ROOT, random numbers are produced by classes such as TRandom and TRandom3. For simple exercises it is common to create one global random generator and use it for the whole simulation. A minimal setup in a macro might look like:

cpp

void simple_mc()
{
   TRandom3 rng(12345);  // fixed seed for reproducibility
   const Int_t nEvents = 100000;
   // ... create histograms here ...
   for (Int_t i = 0; i < nEvents; ++i) {
      // generate one event here
   }
}

The important idea is that an "event" is one independent trial of your process. In a physics context one event could be a collision, a decay, or a single measurement. In a more abstract context an event could be one draw from a probability distribution.

Within the event loop, you use methods of rng to draw random numbers from the distributions that describe your model. For instance:

cpp

double x_uniform  = rng.Uniform(0.0, 1.0);     // uniform in [0, 1)
double x_gauss    = rng.Gaus(0.0, 1.0);        // Gaussian with mean 0, sigma 1
int    n_poisson  = rng.Poisson(3.0);          // Poisson with mean 3
double t_exp      = rng.Exp(2.0);              // exponential with mean 2

You can combine several random draws to represent more complex situations. For example, you might model a measurement as a "true" physical quantity plus a random Gaussian fluctuation that represents detector resolution:

cpp

double true_energy   = rng.Exp(1.0);      // true energy, exponential distribution
double measured_energy = rng.Gaus(true_energy, 0.1);  // smear with sigma = 0.1

In a simple Monte Carlo, each pass through the event loop corresponds to generating such quantities once and then recording them. The number of events is under your control. In most cases you choose it so that statistical fluctuations become small enough for your purpose. If the number of events is $N$, the typical relative statistical uncertainty on simple counts scales as $1/\sqrt{N}$.

In a Monte Carlo simulation each event must be statistically independent and generated from the correct underlying distribution. Reusing the same random number or introducing correlations that are not part of the model will distort the simulated distributions.

Filling histograms

Once you generate values event by event, you usually fill histograms to build the corresponding distributions. This lets you visualize the model and compare it with data later.

To prepare for this, create histograms before the event loop:

cpp

TH1D *h_true   = new TH1D("h_true",   "True energy;E_{true};Events",     100, 0.0, 5.0);
TH1D *h_meas   = new TH1D("h_meas",   "Measured energy;E_{meas};Events", 100, 0.0, 5.0);
TH1D *h_resid  = new TH1D("h_resid",  "Residuals;E_{meas} - E_{true};Events", 100, -1.0, 1.0);

Inside the event loop you fill them:

cpp

for (Int_t i = 0; i < nEvents; ++i) {
   double true_energy    = rng.Exp(1.0);
   double measured_energy = rng.Gaus(true_energy, 0.1);
   double residual        = measured_energy - true_energy;
   h_true->Fill(true_energy);
   h_meas->Fill(measured_energy);
   h_resid->Fill(residual);
}

Each call to Fill() corresponds to one event contributing to the appropriate bin. The bin content is a count of how many events fall in that range. If you use weights, for example to represent variable event importance or cross sections, you can pass a second argument:

cpp

double weight = 0.5;   // example weight
h_meas->Fill(measured_energy, weight);

Weights are common when simulated events represent different amounts of real data. If you use weights, you should enable storage of proper bin errors before filling:

cpp

h_meas->Sumw2();

After your event loop, you can draw the histograms on a canvas:

cpp

TCanvas *c1 = new TCanvas("c1", "Simple Monte Carlo", 800, 600);
h_meas->SetLineColor(kRed);
h_true->SetLineColor(kBlue);
h_true->Draw();
h_meas->Draw("SAME");
c1->Update();

The whole structure of a simple Monte Carlo simulation is summarized as:

Step	Example in ROOT
Create random generator	`TRandom3 rng(seed);`
Define histograms	`new TH1D("name", "title", nbins, min, max);`
Event loop	`for (int i = 0; i < nEvents; ++i)`
Generate event quantities	`x = rng.Gaus(mean, sigma);`
Fill histograms	`hist->Fill(x);`
Visualize results	`hist->Draw();`

Always create histograms before the event loop and fill them consistently inside the loop. Changing histogram binning or range after filling is not equivalent to regenerating the events and can misrepresent your simulated distribution.

Comparing simulated distributions

A key purpose of Monte Carlo is to compare your simulated distributions with expectations or with real experimental data. In simple exercises, this often means overlaying histograms or checking whether the shape matches an analytical probability density function.

To compare two histograms that represent the same physical quantity, such as a simulated distribution and a reference, you usually draw them on the same canvas and make sure they are normalized in a meaningful way. For example:

cpp

TCanvas *c2 = new TCanvas("c2", "Comparison", 800, 600);
h_true->SetLineColor(kBlue);
h_meas->SetLineColor(kRed);
// Option 1: compare raw counts
h_true->Draw("HIST");
h_meas->Draw("HIST SAME");

If the total number of events is different in the two histograms, it is often better to normalize them to unit area before comparing shapes:

cpp

h_true->Scale(1.0 / h_true->Integral());
h_meas->Scale(1.0 / h_meas->Integral());

After scaling, each histogram represents an estimate of the underlying probability density, up to bin width. If you want to make this explicit, you can divide by bin width as well, but for many beginner uses comparing unit-area histograms is sufficient.

You can also compare simulated distributions to known theoretical forms. For example, if you expect the measured energy residuals to follow a Gaussian, you can fit a Gaussian to the h_resid histogram and inspect the parameters:

cpp

h_resid->Fit("gaus");

This gives you estimates of the mean and sigma from the simulation and lets you check whether they match the values you used as input.

More formal comparisons can rely on statistical tests, such as chi-square tests, that compare binned distributions or test a histogram against a known function. Even without formal tests, you can gain intuition by visually inspecting whether the simulated histogram follows the expected shape, where deviations occur, and how fluctuations decrease as you increase the number of events.

A simple workflow might be:

Generate events according to a model.
Fill histograms for observables of interest.
Normalize histograms if you want to compare shapes.
Overlay several histograms on the same axes.
Optionally fit analytical functions to the histograms.
Adjust your model if the simulated distributions do not match expectations.

When comparing simulated distributions, always ensure that you are comparing quantities with the same definition, units, and normalization. Misaligned definitions or inconsistent normalization can lead to incorrect conclusions about agreement or disagreement between models and data.

16.3. Simple Monte Carlo Simulation

Generating events

Filling histograms

Comparing simulated distributions

Comments

Where to Move