1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FinalState.hh"
  4#include "Rivet/Projections/FastJets.hh"
  5#include "Rivet/Projections/HeavyHadrons.hh"
  6#include "Rivet/Projections/IdentifiedFinalState.hh"
  7
  8namespace Rivet {
  9
 10
 11  class D0_2000_I499943 : public Analysis {
 12  public:
 13
 14    /// Constructor
 15    RIVET_DEFAULT_ANALYSIS_CTOR(D0_2000_I499943);
 16
 17
 18    /// @name Analysis methods
 19    //@{
 20
 21    /// Book histograms and initialise projections before the run
 22    void init() {
 23      FinalState fs;
 24      IdentifiedFinalState muons(Cuts::abseta < 0.8 && Cuts::pT > 4.0*GeV);
 25      muons.acceptIdPair(PID::MUON);
 26      declare(muons, "Muons");
 27
 28      FastJets jetproj(fs, FastJets::D0ILCONE, 0.7);
 29      jetproj.useInvisibles();
 30      declare(jetproj, "Jets");
 31
 32      // Book histograms
 33      book(_h_pt_leading_mu ,1, 1, 1);
 34      book(_h_dphi_mumu ,3, 1, 1);
 35    }
 36
 37
 38    /// Perform the per-event analysis
 39    void analyze(const Event& event) {
 40      const Jets& jets = apply<FastJets>(event, "Jets").jetsByPt(12*GeV);
 41      if (jets.size() < 2) vetoEvent;
 42
 43      const Particles& muons = apply<IdentifiedFinalState>(event, "Muons").particlesByPt();
 44      if (muons.size() < 2) vetoEvent;
 45
 46      // Muon selection: require the muons to be *close* to jets, not the usual overlap vetoing!
 47      Particles cand_mu;
 48      for (const Particle& mu : muons) {
 49        // Ignore muons in "bad" region 80 < phi < 110 degrees
 50        /// @todo Is this really not corrected for?!
 51        if (inRange(mu.phi(), 1.4, 1.92)) continue;
 52
 53        // A muon is a good candidate if within R = 0.8 of a jet
 54        for (const Jet& jet : jets) {
 55          if (deltaR(mu, jet) < 0.8) {
 56            cand_mu.push_back(mu);
 57            break;
 58          }
 59        }
 60      }
 61
 62      // Must find at least two jet-matched muons in the event
 63      if (cand_mu.size() < 2) vetoEvent;
 64
 65      /// @todo Is this cut needed? Does space angle mean dR or 3D opening angle in lab frame?
 66      // Remove muon pairs closer than 165 deg in space angle (cosmic veto)
 67      // double dR_mumu = deltaR(cand_mu[0].momentum(), cand_mu[1].momentum());
 68      // if (dR_mumu < 165*degree) vetoEvent;
 69
 70      // Selecting muon pairs with 6 < mass < 35 GeV (we use the two with highest pT)
 71      double m_mumu = (cand_mu[0].momentum() + cand_mu[1].momentum()).mass();
 72      if (!inRange(m_mumu, 6*GeV, 35*GeV)) vetoEvent;
 73
 74      // Get phi angle between muons in degrees
 75      double dphi_mumu = deltaPhi(cand_mu[0], cand_mu[1]) * 180/M_PI;
 76
 77      // Fill histos
 78      _h_pt_leading_mu->fill(cand_mu[0].pt()/GeV);
 79      _h_dphi_mumu->fill(dphi_mumu);
 80    }
 81
 82
 83    // Normalise histograms to cross-section
 84    void finalize() {
 85      scale(_h_pt_leading_mu, crossSection()/sumOfWeights()/nanobarn);
 86      scale(_h_dphi_mumu, crossSection()/sumOfWeights()/nanobarn);
 87    }
 88
 89    //@}
 90
 91
 92  private:
 93
 94    /// @name Histograms
 95    //@{
 96    Histo1DPtr _h_pt_leading_mu, _h_dphi_mumu;
 97    //@}
 98
 99  };
100
101
102  // The hook for the plugin system
103  RIVET_DECLARE_PLUGIN(D0_2000_I499943);
104
105}