Histograms for Systematics

Last update: 05 Jul 2023 [History] [Edit]

The common CP algorithms provide systematic variations of every analysis object. In a typical analysis workflow, it is necessary to save information about the events with each systematic variation. However, this can be a significant bookkeeping challenge and the most common approach of creating separate trees for each systematic variation is very space inefficient.

tip There is an optimal approach to handling systematic variations in ntuples that results in files that are 1-2 orders of magnitude smaller than the common approach. This makes use of the TreeMaker, AsgxAODNTupleMaker and TreeFiller algorithms. However, this approach is beyond the scope of this tutorial. An example implementation can be found here.

Instead of dealing with ntuple bookkeeping, this tutorial will show you how to produce systematic variations of histograms.

Plot Sum Of Muon Pt

As mentioned earlier, the CP algorithms put their outputs in the event store, and those outputs can then be accessed as if they come directly from the input file. To begin with, we won’t deal with systematics just yet.

Let’s make a histogram of the sum pT of all muons in an event. First, we will create the histogram itself. This will be done in a slightly different way than was introduced in Making Histograms. This is because we will modify it in the next section when creating systematic variations.

First declare it in your class definition inside the header file (MyxAODAnalysis.h):

TH1 *m_sumPtHist {nullptr};

Next we will create the histogram. This is done in execute(), not initialize() so that we can pick up the systematics in the next section (systematics are not available during initialize()). So add these lines at the beginning of execute():

  if (m_sumPtHist == nullptr)
  {
    std::string name = "sumPtHist";
    ANA_CHECK (book (TH1F (name.c_str(), "pt", 20, 0, 200e3)));
    m_sumPtHist = hist (name);
  }

Now let’s retrieve the list of muons, calculate the sum of pT, and fill it into the histogram:

  const xAOD::MuonContainer *systMuons = nullptr;
  ANA_CHECK (evtStore()->retrieve (systMuons, "AnaMuons_NOSYS"));
  float sumPt = 0;
  for (const xAOD::Muon *mu : *systMuons)
    sumPt += mu->pt();
  m_sumPtHist->Fill (sumPt);

Now compile and run it and see if you get the newly defined histogram and if it makes sense. As a bonus exercise, create a histogram containing only the pT of the first muon and compare it to what we got from the histogram the muon chain produces intrinsically.

Some things to note here:

  • The name of the muon container is based on what you specified when you configured your muon analysis sequence. You can change this if you want, as long as you do it everywhere.

  • There is no specific ordering to the muons in the container, or more specifically they have the same ordering as the input container, i.e., they are ordered in pT before momentum corrections are applied. You can change this by setting the sortPt option on the algorithm creating the final deep copy. Still, in many ways it’s better to get used to the containers not being ordered (as that allows certain optimizations).

  • We don’t use any weights here. This can be excused by saying we are running on data. If you want to read out the scale factor, each muon has a scale factor attached as muon_eff_tight_NOSYS, etc. As an exercise you can collect them for an overall histogram weight.

  • We are using const containers and objects, as we do not intend to modify them at all. This will become more important in the next section.