Particle physics @utopian-io - exercise 1c solution

Repository

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

Introduction

This is a third article in my quest to ~~help Reptilians take over the World~~ better understand Madanalysis physics toolbox. This will hopefully end with a code to analyse results from Large Hadron Collider and test some particle physics hypothesis. Isn't that cool? Open source project helping to expand our knowledge about the world. The code I will describe below is a solution for exercise written by @lemouth.

[ credits: McCauley, Thomas; Taylor, Lucas; for the CMS Collaboration

Theory

However sophisticated the detectors of LHC are, it is sometimes hard to correctly identify what kind of particle was detected. To ensure that only good quality data is used for further analysis, a process called isolation is often implemented. This process makes sure that the particle of interest is well separated from other particles, so we have a better chance of not mistaking one for another.

In the exercise we use a method used also for CMS search for dark matter. This method consists of checking how much activity happens around examined photon. To measure the activity, we sum transverse momenta of all nearby particles. If they are small enough, we can use the photon for further analysis.

Implementation

A (rather lengthy) signature of a function that gets isolated photons looks like this:

  std::vector<RecPhotonFormat> GetIsolatedPhotons(const EventFormat& event,
                                                  MAfloat32 he_threshold,
                                                  MAfloat32 abseta_threshold,
                                                  const cms_constants::PhotonIsolationThresholds &isolation_thresholds
                                                  );

PhotonIsolationThresholds is a structure that contains the parameters needed to calculate momenta thresholds. Interesting parts of the function body are:

Rejection of the photons for which we don't have to calculate surrounding particles momenta, because they don't meet other criteria:

    if (photon.HEoverEE() >= he_threshold ||
        photon.abseta() >= abseta_threshold)
      continue;

Getting momenta sums:

    auto Ipi = PHYSICS->Isol->eflow->sumIsolation(photon, event.rec(), 0.3, 0.,
                                                    IsolationEFlow::TRACK_COMPONENT);
    auto In = PHYSICS->Isol->eflow->sumIsolation(photon, event.rec(), 0.3, 0.,
                                                   IsolationEFlow::NEUTRAL_COMPONENT);
    auto Igamma = PHYSICS->Isol->eflow->sumIsolation(photon, event.rec(), 0.3, 0.,
                                                     IsolationEFlow::PHOTON_COMPONENT);

Finally, comparing those sums to thresholds:

    if (Ipi < it.barrel_Ipi_const &&
        In < it.barrel_In_const + it.barrel_In_pt_multiplier * photon.pt() &&
        Igamma - photon.pt() < it.barrel_Igamma_const + it.barrel_Igamma_pt_multiplier * photon.pt()) {
      isolated_photons.push_back(photon);
    }

Histogramms

One of the exercise objectives was also to draw the histograms of photons and jets momenta and pseudo-rapidity. It was even a task request to implement a tool that does it in Python. I'm rather poor Python programmer and was running short on time, so I used a tool implemented by @crokkon - https://github.com/crokkon/saf2png.

Resulting plots:
histos.jets_eta.png
histos.jets_pt.png
histos.photons_eta.png
histos.photons_pt.png
@lemouth - I hope those are correct. Doing the exercise I got wondering: what is the meaning of y axis here? I was expecting it be the number of particles with a given momentum, but obviously it's not the case...