1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FinalState.hh"
  4#include "Rivet/Projections/WFinder.hh"
  5#include "Rivet/Projections/LeadingParticlesFinalState.hh"
  6#include "Rivet/Projections/IdentifiedFinalState.hh"
  7
  8namespace Rivet {
  9
 10
 11  /// @brief D0 Run II measurement of W charge asymmetry
 12  ///
 13  /// @author Andy Buckley
 14  /// @author Gavin Hesketh
 15  class D0_2008_S7837160 : public Analysis {
 16  public:
 17
 18    RIVET_DEFAULT_ANALYSIS_CTOR(D0_2008_S7837160);
 19
 20
 21    /// @name Analysis methods
 22    /// @{
 23
 24    // Book histograms and set up projections
 25    void init() {
 26      // Projections
 27      FinalState fs;
 28      /// @todo Use separate pT and ETmiss cuts in WFinder
 29      const WFinder wfe(fs, Cuts::abseta < 5 && Cuts::pT > 25*GeV, PID::ELECTRON, 60.0*GeV, 100.0*GeV, 25.0*GeV, 0.2);
 30      declare(wfe, "WFe");
 31
 32      // Histograms (temporary +- charge histos and scatters to store the calculated asymmetries)
 33      for (size_t pmindex = 0; pmindex < 2; ++pmindex) {
 34        const string suffix = (pmindex == 0) ? "plus" : "minus";
 35        book(_hs_dsigpm_deta_25_35[pmindex] ,"TMP/dsigpm_deta_25_35_" + suffix, refData(1, 1, 1));
 36        book(_hs_dsigpm_deta_35[pmindex] ,"TMP/dsigpm_deta_35_" + suffix, refData(1, 1, 2));
 37        book(_hs_dsigpm_deta_25[pmindex] ,"TMP/dsigpm_deta_25_" + suffix, refData(1, 1, 3));
 38      }
 39      book(_h_asym1, 1, 1, 1);
 40      book(_h_asym2, 1, 1, 2);
 41      book(_h_asym3, 1, 1, 3);
 42    }
 43
 44
 45    /// Do the analysis
 46    void analyze(const Event & event) {
 47      const WFinder& wf = apply<WFinder>(event, "WFe");
 48      if (wf.bosons().size() == 0) {
 49        MSG_DEBUG("No W candidates found: vetoing");
 50        vetoEvent;
 51      }
 52
 53      // Get the e+- momentum, and an effective charge including the eta sign
 54      /// @todo Is it correct to multiply the eta sign into the charge to "fold" the plot?
 55      const FourMomentum p_e = wf.constituentLeptons()[0].momentum();
 56      const int chg_e = sign(p_e.eta()) * sign(charge(wf.constituentLeptons()[0]));
 57      assert(chg_e == 1 || chg_e == -1);
 58      MSG_TRACE("Charged lepton sign = " << chg_e);
 59
 60      // Fill histos with appropriate +- indexing
 61      const size_t pmindex = (chg_e > 0) ? 0 : 1;
 62      if (p_e.Et() < 35*GeV) _hs_dsigpm_deta_25_35[pmindex]->fill(fabs(p_e.eta()));
 63      else _hs_dsigpm_deta_35[pmindex]->fill(fabs(p_e.eta()));
 64      _hs_dsigpm_deta_25[pmindex]->fill(fabs(p_e.eta()));
 65    }
 66
 67
 68    /// @brief Finalize
 69    ///
 70    /// Construct asymmetry: (dsig+/deta - dsig-/deta) / (dsig+/deta + dsig-/deta) for each ET region
 71    void finalize() {
 72      calc_asymm(_hs_dsigpm_deta_25_35, _h_asym1);
 73      calc_asymm(_hs_dsigpm_deta_35, _h_asym2);
 74      calc_asymm(_hs_dsigpm_deta_25, _h_asym3);
 75      _h_asym1->scaleY(100.);
 76      _h_asym2->scaleY(100.);
 77      _h_asym3->scaleY(100.);
 78    }
 79
 80    /// @}
 81
 82
 83  private:
 84
 85    /// @name Helper functions for constructing asymmetry histograms in finalize()
 86    /// @{
 87    void calc_asymm(const Histo1DPtr plus, const Histo1DPtr minus, Scatter2DPtr target) {
 88      divide(*plus - *minus, *plus + *minus, target);
 89    }
 90    void calc_asymm(const Histo1DPtr histos[2], Scatter2DPtr target) {
 91      calc_asymm(histos[0], histos[1], target);
 92    }
 93    /// @}
 94
 95
 96    /// @name Histograms
 97    /// @{
 98    Histo1DPtr _hs_dsigpm_deta_25_35[2], _hs_dsigpm_deta_35[2], _hs_dsigpm_deta_25[2];
 99    Scatter2DPtr _h_asym1, _h_asym2, _h_asym3;
100    /// @}
101
102  };
103
104
105
106
107  RIVET_DECLARE_ALIASED_PLUGIN(D0_2008_S7837160, D0_2008_I791230);
108
109}