Rivet Analyses Reference

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137

// -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/FinalState.hh"
#include "Rivet/Projections/IdentifiedFinalState.hh"
#include "Rivet/Projections/DressedLeptons.hh"

namespace Rivet {


  class ATLAS_2015_I1351916 : public Analysis {
  public:

    /// Constructor
    ATLAS_2015_I1351916(const string name="ATLAS_2015_I1351916", size_t mode=0,
                        const string ref_data="ATLAS_2015_I1351916") : Analysis(name) {
      _mode = mode; // pick electron channel by default
      setRefDataName(ref_data);
    }


    /// @name Analysis methods
    //@{

    /// Book histograms and initialise projections before the run
    void init() {

      const FinalState fs;

      IdentifiedFinalState bareleptons(fs);
      bareleptons.acceptIdPair(_mode? PID::MUON : PID::ELECTRON);

      const Cut cuts = (_mode == 0) ? (Cuts::pT > 25*GeV && Cuts::abseta < 4.9) : (Cuts::pT > 20*GeV && Cuts::abseta < 2.47);
      DressedLeptons leptons(fs, bareleptons, 0.1, cuts, true);
      declare(leptons, "leptons");


      // Book dummy histograms for heterogeneous merging
      const Scatter2D& ref = refData(_mode? 4 : 2, 1, 2);
      book(_h["NCC_pos"], "_ncc_pos", ref);
      book(_h["NCC_neg"], "_ncc_neg", ref);
      book(_s["CC"], _mode ? 4 : 2, 1, 2, true);

      if (_mode == 0) { // electron-channel only
        const Scatter2D& ref_cf = refData(3, 1, 2);
        book(_h["NCF_pos"], "_ncf_pos", ref_cf);
        book(_h["NCF_neg"], "_ncf_neg", ref_cf);
        book(_s["CF"], 3, 1, 2, true);
      }
    }


    /// Perform the per-event analysis
    void analyze(const Event& e) {

      // Get and cut on dressed leptons
      const vector<DressedLepton>& leptons = apply<DressedLeptons>(e, "leptons").dressedLeptons();
      if (leptons.size() != 2) vetoEvent; // require exactly two leptons
      if (leptons[0].charge3() * leptons[1].charge3() > 0) vetoEvent; // require opposite charge

      // Identify lepton vs antilepton
      const Particle& lpos = leptons[(leptons[0].charge3() > 0) ? 0 : 1];
      const Particle& lneg = leptons[(leptons[0].charge3() < 0) ? 0 : 1];

      string label = "N";
      if (_mode == 1) {// electron channel
        label += "CC"; // only central-central for muons
      } else { // electron channel
        const double eta1 = lpos.abseta();
        const double eta2 = lneg.abseta();
        if ( (eta1 < 2.47 && inRange(eta2, 2.5, 4.9)) || (eta2 < 2.47 && inRange(eta1, 2.5, 4.9)) )
          label += "CF"; // central-forward
        else if (eta1 < 2.47 && eta2 < 2.47)
          label += "CC"; // central-central
        else vetoEvent; // ain't no forward-forward
      }

      const double cosThetaStar = cosCollinsSoper(lneg, lpos);
      const double mll = (lpos.mom() + lneg.mom()).mass();
      label += cosThetaStar < 0.0?  "_neg" : "_pos";
      _h[label]->fill(mll/GeV);
    }


    /// Normalise histograms etc., after the run
    void finalize() {
      const double sf = crossSectionPerEvent() / picobarn;
      for (const auto& key_hist : _h) scale(key_hist.second, sf);
      divide(*_h["NCC_pos"] - *_h["NCC_neg"], *_h["NCC_pos"] + *_h["NCC_neg"], _s["CC"]);
      if (!_mode)  divide(*_h["NCF_pos"] - *_h["NCF_neg"], *_h["NCF_pos"] + *_h["NCF_neg"], _s["CF"]);
    }


    // Cosine of the decay angle in the Collins-Soper frame
    double cosCollinsSoper(const FourMomentum& l1, const FourMomentum& l2) {
      const FourMomentum ll = l1 + l2;
      const double nom  = (l1.E() + l1.pz()) * (l2.E() - l2.pz()) - (l1.E() - l1.pz()) * (l2.E() + l2.pz());
      const double denom = ll.mass() * sqrt( sqr(ll.mass()) + sqr(ll.pt()) );
      return sign(ll.pz()) * safediv(nom, denom); // protect against division by zero, you never know...
    }

    //@}


  protected:

    /// Electron or muon mode = 0 or 1, for use by derived _EL, _MU analysis classes
    size_t _mode;


  private:

    /// Histograms
    map<string, Histo1DPtr> _h;
    /// Asymmetries
    map<string, Scatter2DPtr> _s;

  };



  class ATLAS_2015_I1351916_EL : public ATLAS_2015_I1351916 {
  public:
    ATLAS_2015_I1351916_EL() : ATLAS_2015_I1351916("ATLAS_2015_I1351916_EL", 0) { }
  };


  class ATLAS_2015_I1351916_MU : public ATLAS_2015_I1351916 {
  public:
    ATLAS_2015_I1351916_MU() : ATLAS_2015_I1351916("ATLAS_2015_I1351916_MU", 1) { }
  };


  RIVET_DECLARE_PLUGIN(ATLAS_2015_I1351916);
  RIVET_DECLARE_PLUGIN(ATLAS_2015_I1351916_EL);
  RIVET_DECLARE_PLUGIN(ATLAS_2015_I1351916_MU);

}