ZmmYCommon.h

Go to the documentation of this file.

#pragma once
 
#include <vector>
 
#include "Core/Objects/interface/Event.h"
#include "Core/Objects/interface/Photon.h"
#include "Core/Objects/interface/Di.h"
 
namespace DAS::Unfolding::ZmmY {
 
#if !defined(__CLING__) || defined(__ROOTCLING__)
enum class EventCategory {
    NotSelected, 
    MuMu,        
    MuMuGamma,   
};
 
template<class Muon, class Photon>
struct EventInfo
{
    // Phase space cuts.
    constexpr static double minMuon0Pt = 20,    
                            minMuon1Pt = 15,    
                            maxMuonEta = 2.4,   
                            minPhotonPt = 20,   
                            maxPhotonEta = 2.4, 
                            minZMass = 76,      
                            maxZMass = 106,     
                            minMuMuMass = 40;   
    static_assert(minMuon0Pt >= minMuon1Pt, "leading pT always >= subleading");
    static_assert(minZMass < maxZMass, "min mass always < max mass");
 
    EventCategory category;
 
    Di<const Muon, const Muon> dimuon;
 
    float mMuMu;
 
    std::vector<Photon> selectedPhotons;
 
    float mMuMuGamma;
 
    float weight;
 
    EventInfo ()
        : category(EventCategory::NotSelected)
        , selectedPhotons{}
        , mMuMu(-1)
        , mMuMuGamma(-1)
        , weight(0)
    {
    }
 
    EventInfo (const GenEvent* gEv, 
               const RecEvent* rEv, 
               const std::vector<Muon>& muons,
               const std::vector<Photon>& photons,
               const Uncertainties::Variation& v);
 
private:
    bool IsMassOnZPeak (float mass) const
    {
        return mass > minZMass && mass < maxZMass;
    }
 
    DAS::Di<const Muon, const Muon> SelectMuons (
        const std::vector<Muon>& muons, const Uncertainties::Variation& v) const
    {
        if (muons.size() < 2) return {};
 
        const auto& m0 = muons[0],
                  & m1 = muons[1];
 
        const auto p0 = m0.CorrP4(v),
                   p1 = m1.CorrP4(v);
 
        if (m0.Q * m1.Q != -1 // opposite charges
            || p0.Pt() < minMuon0Pt
            || p1.Pt() < minMuon1Pt
            || abs(p0.Eta()) >= maxMuonEta
            || abs(p1.Eta()) >= maxMuonEta)
            return {};
 
        return m0 + m1;
    }
 
    bool IsGoodPhoton (const RecPhoton& photon,
                       const Uncertainties::Variation& v) const
    {
        auto p4 = photon.CorrP4(v);
        return photon.Weight(v) != 0
            && p4.Pt() > minPhotonPt && std::abs(p4.Eta()) < maxPhotonEta;
    }
 
    bool IsGoodPhoton (const GenPhoton& photon,
                       const Uncertainties::Variation& v) const
    {
        auto p4 = photon.CorrP4(v);
        return photon.Weight(v) != 0
            && photon.zAncestor
            && p4.Pt() > minPhotonPt && std::abs(p4.Eta()) < maxPhotonEta;
    }
};
 
template<class Muon, class Photon>
EventInfo<Muon, Photon>::EventInfo (const GenEvent* gEv,
                                    const RecEvent* rEv,
                                    const std::vector<Muon>& muons,
                                    const std::vector<Photon>& photons,
                                    const Uncertainties::Variation& v)
    : EventInfo()
{
    weight = 1;
    if (rEv) weight *= rEv->Weight(v);
    if (gEv) weight *= gEv->Weight(v);
 
    // Start by looking at muons.
    dimuon = SelectMuons(muons, v);
    if (!dimuon) return;
 
    mMuMu = dimuon.CorrP4(v).M();
 
    // The dimuon system is on the Z peak. Don't consider photons.
    if (IsMassOnZPeak(mMuMu)) {
        category = EventCategory::MuMu;
        weight *= dimuon.Weight(v);
        return;
    }
 
    // Do not consider muon pairs below minMuMuMass. This reduces
    // backgrounds (e.g. vector meson + gamma).
    // We also remove here events with mMuMu above the Z peak, as there is
    // no way adding photons can bring them back to the peak.
    if (mMuMu < minMuMuMass || mMuMu > minZMass)
        return;
 
    // Now get photons. We must copy them because the photons collection could
    // be reused with another variation.
    std::copy_if(photons.begin(), photons.end(),
                    std::back_inserter(selectedPhotons),
                    [&] (auto photon) { return IsGoodPhoton(photon, v); });
 
    if (selectedPhotons.empty()) return;
 
    mMuMuGamma = (dimuon.CorrP4(v) + selectedPhotons[0].CorrP4(v)).M();
    if (IsMassOnZPeak(mMuMuGamma)) {
        category = EventCategory::MuMuGamma;
        weight *= dimuon.Weight(v) * selectedPhotons[0].Weight(v);
    }
}
#endif
 
} // end of DAS::Unfolding::ZmmY namespace