Thrust.hh

// -*- C++ -*-
#ifndef RIVET_Thrust_HH
#define RIVET_Thrust_HH
 
#include "Rivet/Event.hh"
#include "Rivet/Projection.hh"
#include "Rivet/Projections/AxesDefinition.hh"
#include "Rivet/Projections/FinalState.hh"
 
namespace Rivet {
 

  class Thrust : public AxesDefinition {
  public:
 
    using AxesDefinition::operator=;

    Thrust() { }
 
    Thrust(const FinalState& fsp) {
      setName("Thrust");
      declare(fsp, "FS");
    }

    RIVET_DEFAULT_PROJ_CLONE(Thrust);

    using Projection::operator=;
 
 
  protected:

    void project(const Event& e) {
      const vector<Particle> ps = apply<FinalState>(e, "FS").particles();
      calc(ps);
    }

    CmpState compare(const Projection& p) const {
      return mkNamedPCmp(p, "FS");
    }
 
 
  public:

    double thrust() const {
      return _thrusts[0];
    }

    double thrustMajor() const {
      return _thrusts[1];
    }

    double thrustMinor() const {
      return _thrusts[2];
    }

    double oblateness() const {
      return _thrusts[1] - _thrusts[2];
    }


    const Vector3& thrustAxis() const {
      return _thrustAxes[0];
    }

    const Vector3& thrustMajorAxis() const {
      return _thrustAxes[1];
    }

    const Vector3& thrustMinorAxis() const {
      return _thrustAxes[2];
    }


    const Vector3& axis1() const {
      return thrustAxis();
    }
    const Vector3& axis2() const {
      return thrustMajorAxis();
    }
    const Vector3& axis3() const {
      return thrustMinorAxis();
    }

 
 
  public:


    void calc(const FinalState& fs);

    void calc(const vector<Particle>& fsparticles);

    void calc(const vector<FourMomentum>& fsmomenta);

    void calc(const vector<Vector3>& threeMomenta);

 
 
  protected:

    vector<double> _thrusts;

    vector<Vector3> _thrustAxes;
 
 
  protected:

    void _calcThrust(const vector<Vector3>& fsmomenta);
  };
 
}
 
#endif