rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

MC_WPOL

Monte Carlo validation observables for $W$ polarisation
Experiment: ()
Status: VALIDATED
Authors:
  • Frank Siegert
No references listed
Beams: * *
Beam energies: ANY
Run details:
  • $W \to e \, \nu$ + jets.

Observables sensitive to the polarisation of the W boson: A0, ... A7, fR, fL, f0, separately for W+ and W-.

Source code: MC_WPOL.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FinalState.hh"
  4#include "Rivet/Projections/WFinder.hh"
  5#include "Rivet/Projections/Beam.hh"
  6
  7namespace Rivet {
  8
  9
 10  /// @brief MC validation analysis for W polarisation
 11  class MC_WPOL : public Analysis {
 12  public:
 13
 14    /// @name Constructors etc.
 15    //@{
 16
 17    /// Constructor
 18    MC_WPOL()
 19      : Analysis("MC_WPOL")
 20    {    }
 21
 22    //@}
 23
 24
 25  public:
 26
 27    /// @name Analysis methods
 28    //@{
 29
 30    /// Book histograms and initialise projections before the run
 31    void init() {
 32
 33      FinalState fs;
 34      WFinder wfinder(fs, Cuts::open(), PID::ELECTRON,
 35                      60.0*GeV, 100.0*GeV, 0.0*GeV, 0.0);
 36      declare(wfinder, "WFinder");
 37      Beam beams;
 38      declare(beams, "Beams");
 39
 40      vector<string> tags{"_wplus", "_wminus"};
 41      _h_dists.resize(tags.size());
 42      _h_histos.resize(tags.size());
 43      for (size_t i=0; i<tags.size(); ++i) {
 44        _h_dists[i].resize(11,Profile1DPtr());
 45        double sqrts = sqrtS()>0. ? sqrtS() : 14000.;
 46        book(_h_dists[i][0] ,"A0"+tags[i],logspace(100, 1.0, 0.5*sqrts));
 47        book(_h_dists[i][1] ,"A1"+tags[i],logspace(100, 1.0, 0.5*sqrts));
 48        book(_h_dists[i][2] ,"A2"+tags[i],logspace(100, 1.0, 0.5*sqrts));
 49        book(_h_dists[i][3] ,"A3"+tags[i],logspace(100, 1.0, 0.5*sqrts));
 50        book(_h_dists[i][4] ,"A4"+tags[i],logspace(100, 1.0, 0.5*sqrts));
 51        book(_h_dists[i][5] ,"A5"+tags[i],logspace(100, 1.0, 0.5*sqrts));
 52        book(_h_dists[i][6] ,"A6"+tags[i],logspace(100, 1.0, 0.5*sqrts));
 53        book(_h_dists[i][7] ,"A7"+tags[i],logspace(100, 1.0, 0.5*sqrts));
 54        book(_h_dists[i][8] ,"fL"+tags[i],logspace(100, 1.0, 0.5*sqrts));
 55        book(_h_dists[i][9] ,"fR"+tags[i],logspace(100, 1.0, 0.5*sqrts));
 56        book(_h_dists[i][10] ,"f0"+tags[i],logspace(100, 1.0, 0.5*sqrts));
 57        _h_histos[i].resize(4,Histo1DPtr());
 58        book(_h_histos[i][0] ,"thetastar"+tags[i],100,-1.0,1.0);
 59        book(_h_histos[i][1] ,"phistar"+tags[i],90,0.0,360.0);
 60        book(_h_histos[i][2] ,"thetastar_ptw20"+tags[i],100,-1.0,1.0);
 61        book(_h_histos[i][3] ,"phistar_ptw20"+tags[i],90,0.0,360.0);
 62      }
 63    }
 64
 65
 66    /// Perform the per-event analysis
 67    void analyze(const Event& event) {
 68      const double weight = 1.0;
 69
 70      const WFinder& wfinder = apply<WFinder>(event, "WFinder");
 71      if (wfinder.bosons().size() != 1) {
 72        vetoEvent;
 73      }
 74      const ParticlePair& beams = apply<Beam>(event, "Beams").beams();
 75
 76      FourMomentum pb1(beams.second.momentum()), pb2(beams.first.momentum());
 77      Particle lepton = wfinder.constituentLeptons()[0];
 78      FourMomentum pl(lepton.momentum());
 79      size_t idx = (PID::charge3(lepton.pid())>0 ? 0 : 1);
 80      FourMomentum plnu(wfinder.bosons()[0].momentum());
 81
 82      const LorentzTransform cms = LorentzTransform::mkFrameTransformFromBeta(plnu.betaVec());
 83      Matrix3 zrot(plnu.p3(), Vector3(0.0, 0.0, 1.0));
 84      pl=cms.transform(pl);
 85      pb1=cms.transform(pb1);
 86      pb2=cms.transform(pb2);
 87      Vector3 pl3=pl.p3();
 88      Vector3 pb13=pb1.p3();
 89      Vector3 pb23=pb2.p3();
 90      pl3=zrot*pl3;
 91      pb13=zrot*pb13;
 92      pb23=zrot*pb23;
 93      Vector3 xref(cos(pb13.theta())>cos(pb23.theta())?pb13:pb23);
 94      Matrix3 xrot(Vector3(xref.x(), xref.y(), 0.0), Vector3(1.0, 0.0, 0.0));
 95      pl3=xrot*pl3;
 96
 97      double ptw(wfinder.bosons()[0].pT()/GeV);
 98      double thetas(pl3.theta()), phis(pl3.phi());
 99      double costhetas(cos(thetas)), sinthetas(sin(thetas));
100      double cosphis(cos(phis)), sinphis(sin(phis));
101      if (phis<0.0) phis+=2.0*M_PI;
102
103      _h_histos[idx][0]->fill(costhetas,weight);
104      _h_histos[idx][1]->fill(phis*180.0/M_PI,weight);
105      if (ptw>20.0) {
106        _h_histos[idx][2]->fill(costhetas,weight);
107        _h_histos[idx][3]->fill(phis*180.0/M_PI,weight);
108      }
109      _h_dists[idx][0]->fill(ptw,10.0/3.0*(1.0-3.0*sqr(costhetas))+2.0/3.0,weight);
110      _h_dists[idx][1]->fill(ptw,10.0*sinthetas*costhetas*cosphis,weight);
111      _h_dists[idx][2]->fill(ptw,10.0*sqr(sinthetas)*(sqr(cosphis)-sqr(sinphis)),weight);
112      _h_dists[idx][3]->fill(ptw,4.0*sinthetas*cosphis,weight);
113      _h_dists[idx][4]->fill(ptw,4.0*costhetas,weight);
114      _h_dists[idx][5]->fill(ptw,4.0*sinthetas*sinphis,weight);
115      _h_dists[idx][6]->fill(ptw,10.0*costhetas*sinthetas*sinphis,weight);
116      _h_dists[idx][7]->fill(ptw,10.0*sqr(sinthetas)*cosphis*sinphis,weight);
117      _h_dists[idx][8]->fill(ptw,0.5*sqr(1.0-costhetas)-(1.0-2.0*sqr(costhetas)),weight);
118      _h_dists[idx][9]->fill(ptw,0.5*sqr(1.0+costhetas)-(1.0-2.0*sqr(costhetas)),weight);
119      _h_dists[idx][10]->fill(ptw,5.0*sqr(sinthetas)-3.0,weight);
120
121    }
122
123
124    /// Normalise histograms etc., after the run
125    void finalize() {
126
127      for (size_t i=0; i<_h_histos.size(); ++i) {
128        for (Histo1DPtr histo : _h_histos[i]) {
129          scale(histo, crossSection()/picobarn/sumOfWeights());
130        }
131      }
132
133    }
134
135    //@}
136
137
138  private:
139
140    /// @name Histograms
141    //@{
142
143    vector<vector<Profile1DPtr> > _h_dists;
144    vector<vector<Histo1DPtr> > _h_histos;
145    //@}
146
147
148  };
149
150
151
152  // The hook for the plugin system
153  RIVET_DECLARE_PLUGIN(MC_WPOL);
154
155}