Particle physics - exercise 1c solution (CORRECTED)

Repository

https://github.com/BFuks/mad5-utopian-exercises

The files described in this post are saved as:

ex1c-irelandscape/test_cms.cpp
ex1c-irelandscape/saf_reader.py
ex1c-irelandscape/pldoc.py
ex1c-irelandscape/ex1c.py
ex1c-irelandscape/xml_as_dict/

Introduction

This is my second post for providing a solution to @lemouth third exercise in his series about open-source implementation of LHC particle analysis.

It contains corrections to were kindly provided by @lemouth and I hope that the new results are correct.

My original answer to the exercise can be found at here

[Image credits - Pixabay]

Corrections

In my previous post I wrote:

This new isolation procedure is applied to photons in the exercise and consists in calculating the sum of all transverse momenta of charged hadrons (Iπ), neutral hadrons (Iη) and photons (Iγ) within an angular distance dR of the subject photon.

These sums are compared to the transverse momentum of the subject photon and if these are small enough the photon is tagged as a signal photon.

Correction: The photon that matches the above criteria is not yet selected as a signal photon.
Instead it means that the photon has been isolated.
That is, it is a photon with little activity from other particles within a certain angular distance.
This photon can then be further assessed to see if it should be considered as a signal photon for the analysis.

In my original C++ code, I attempted to clean the jets from the electrons as follows:

// Remove jets electrons overlap auto cleaned_jets = PHYSICS->Isol->JetCleaning(jets, electrons, 0.2);

Correction:: This is only partly correct and @lemouth will explain the full procedure in a future post.
As a result, that code is now removed.

In my original C++ code, I tried to isolate photons using the sum of transverse momenta using the following code:

    if (Igam >= (0.7 + (0.005 * photon_pt)) &&
        In >= (1.0  + (0.04 * photon_pt)) &&
        iPi >= 1.5 &&
        photon.HEoverEE() < 0.05 &&
        photon.abseta() < 1.44)
    {
        ...

Correction: The above check is incorrect. Firstly, the transverse momenta should be lower than the expression on the right of each condition.
Secondly the photon itself contributes to the gamma transverse momentum so its transverse momentum should be substracted from the Igamma element.

The corrected code is thus:

    if ((Igam - photon_pt) < (0.7 + (0.005 * photon_pt)) &&
        In < (1.0  + (0.04 * photon_pt)) &&
        iPi < 1.5 &&
        photon.HEoverEE() < 0.05 &&
        photon.abseta() < 1.44)
    {
        ...

In my original code I added every signal jet and photon to the histograms.
@lemouth asked to only add the first of each element since Madanalysis5 already provide jets and photons sorted by descending transverse momentum.

The corrected code is thus:

  // Extract jets (pt > 30GeV, |n| < 5)
  bool found = false;
  std::vector<RecJetFormat> signal_jets;
  for (auto &jet : event.rec()->jets())
  {
    if (jet.pt() > 30 && jet.abseta() < 5)
    {
      if (!found)
      {
        Manager()->FillHisto("pt_jets", jet.pt());
        Manager()->FillHisto("n_jets", jet.eta());
        found = true;
      }
      signal_jets.push_back(jet);
      cout << "Adding jet - Pt = " << jet.pt() << ", |n| = " << jet.eta() << endl;
    }
  }
...
    auto photon_pt = photon.pt();
    if ((Igam - photon_pt) < (0.7 + (0.005 * photon_pt)) &&
        In < (1.0  + (0.04 * photon_pt)) &&
        iPi < 1.5 &&
        photon.HEoverEE() < 0.05 &&
        photon.abseta() < 1.44)
    {
      cout << "Adding  photon: Pt = " << photon_pt << ", |n| = " << photon.eta() << endl;
      if (!found)
      {
        Manager()->FillHisto("pt_photons", photon_pt);
        Manager()->FillHisto("n_photons", photon.eta());
        found = true;
      }
    }
  }

New Result

After applying the above corrections and executing the python script, the following diagrams are generated:

Conclusions

I hope I have understood my mistakes and that the result is now correct.