EventShapeContainer.h

/**************************************************************************
 * basf2 (Belle II Analysis Software Framework)                           *
 * Author: The Belle II Collaboration                                     *
 *                                                                        *
 * See git log for contributors and copyright holders.                    *
 * This file is licensed under LGPL-3.0, see LICENSE.md.                  *
 **************************************************************************/
 
#pragma once
 
#include <framework/datastore/RelationsObject.h>
#include <framework/logging/Logger.h>
 
#include <Math/Vector3D.h>
#include <Math/Vector4D.h>
#include <TLorentzVector.h>  // needed for schema evolution
#include <vector>
 
namespace Belle2 {
  class EventShapeContainer : public RelationsObject {
 
  public:
 
    EventShapeContainer()
    {
      m_thrustAxis.SetXYZ(0., 0., 0.);
      m_spherocityAxis.SetXYZ(0., 0., 0.);
      m_forwardHemisphere4Momentum.SetXYZT(0., 0., 0., 0.);
      m_backwardHemisphere4Momentum.SetXYZT(0., 0., 0., 0.);
 
      m_thrustAxis.SetXYZ(0., 0., 0.);
      for (int i = 0; i < 3; i++)
        m_sphericityEigenvector[i].SetXYZ(0., 0., 0.);
    };
 
 
    // --------------------------------------
    // Setters
    // --------------------------------------
 
 
    void setSphericityEigenvector(short i, ROOT::Math::XYZVector eigenvector)
    {
      if (i < 0 || i > 2)
        B2ERROR("Invalid  EigenVector number  (n = " << i << "). You can set only 3 eigenvectors.");
      else
        m_sphericityEigenvector[i] = eigenvector;
    };
 
    void setSphericityEigenvalue(short i, float eigenvalue)
    {
      if (i < 0 || i > 2)
        B2ERROR("Invalid  EigenVector number  (n = " << i << "). You can set only 3 eigenvectors.");
      else
        m_sphericityEigenvalue[i] = eigenvalue;
    };
 
    void setFWMoment(short order, float moment)
    {
      if (order < 0 || order > 9)
        B2ERROR("Invalid Fox-Wolfram moment order (n = " << order << "). The order must be in the [0,9] range.");
      else
        m_foxWolframMoments[order] = moment;
    };
 
    void setThrust(float thrust)
    {
      m_thrust = thrust;
    };
 
    void setThrustAxis(ROOT::Math::XYZVector thrustAxis)
    {
      if (thrustAxis.R() < 1.E-10)
        B2WARNING("The thrust axis you are trying to set has magnitude numerically compatible with 0.");
      m_thrustAxis = thrustAxis.Unit();
    };
 
    void setSpherocityAxis(ROOT::Math::XYZVector spherocityAxis)
    {
      if (spherocityAxis.R() < 1.E-10)
        B2WARNING("The spherocity axis you are trying to set has magnitude numerically compatible with 0.");
      m_spherocityAxis = spherocityAxis.Unit();
    };
 
    void setHarmonicMomentThrust(short order, float moment)
    {
      if (order < 0 || order > 9)
        B2ERROR("Invalid harmonic moment order. It  must be in the [0,9] range.");
      else
        m_harmonicMomentsThrust[order] = moment;
    };
 
    void setHarmonicMomentCollision(short order, float moment)
    {
      if (order < 0 || order > 9)
        B2ERROR("Invalid harmonic moment order. It  must be in the [0,9] range.");
      else
        m_harmonicMomentsCollision[order] = moment;
    };
 
    void setCleoConeThrust(short order, float moment)
    {
      if (order < 0 || order > 9)
        B2ERROR("Invalid Cleo cone order. It  must be in the [0,9] range.");
      else
        m_cleoConesThrust[order] = moment;
    };
 
    void setCleoConeCollision(short order, float moment)
    {
      if (order < 0 || order > 9)
        B2ERROR("Invalid Cleo cone order. It must be in the [0,9] range.");
      else
        m_cleoConesCollision[order] = moment;
    };
 
    void setForwardHemisphere4Momentum(ROOT::Math::PxPyPzEVector mom)
    {
      m_forwardHemisphere4Momentum = mom;
    };
 
    void setBackwardHemisphere4Momentum(ROOT::Math::PxPyPzEVector mom)
    {
      m_backwardHemisphere4Momentum = mom;
    };
 
 
    // --------------------------------------
    // Getters
    // --------------------------------------
 
 
    float getThrust() const
    {
      return m_thrust;
    };
 
    ROOT::Math::XYZVector getThrustAxis() const
    {
      return m_thrustAxis;
    };
 
    ROOT::Math::XYZVector getSpherocityAxis() const
    {
      return m_spherocityAxis;
    };
 
    float getSphericityEigenvalue(short i) const
    {
      if (i < 0 || i > 2) {
        B2ERROR("Invalid  Eigenvalue number  (n = " << i << "). There are only 3 eigenvalues...");
        return std::numeric_limits<float>::quiet_NaN();
      } else
        return m_sphericityEigenvalue[i];
    };
 
 
    ROOT::Math::XYZVector getSphericityEigenvector(short i) const
    {
      if (i < 0 || i > 2) {
        B2ERROR("Invalid  Eigenvalue number  (n = " << i << "). There are only 3 eigenvalues...");
        return ROOT::Math::XYZVector(std::numeric_limits<float>::quiet_NaN(), std::numeric_limits<float>::quiet_NaN(),
                                     std::numeric_limits<float>::quiet_NaN());
      } else
        return m_sphericityEigenvector[i];
    };
 
 
 
    float getFWMoment(short order) const
    {
      if (order < 0 || order > 9) {
        B2ERROR("Invalid Fox-Wolfram moment order. It must be in the [0,9] range.");
        return std::numeric_limits<float>::quiet_NaN();
      } else
        return m_foxWolframMoments[order];
    };
 
    float getHarmonicMomentCollision(short order) const
    {
      if (order < 0 || order > 9) {
        B2ERROR("Invalid harmonic moment order. It  must be in the [0,9] range.");
        return std::numeric_limits<float>::quiet_NaN();
      } else
        return m_harmonicMomentsCollision[order];
    };
 
    float getHarmonicMomentThrust(short order) const
    {
      if (order < 0 || order > 9) {
        B2ERROR("Invalid harmonic moment order. The order must be in the [0,9] range.");
        return std::numeric_limits<float>::quiet_NaN();
      } else
        return m_harmonicMomentsThrust[order];
    };
 
 
    float getCleoConeCollision(short order) const
    {
      if (order < 0 || order > 9) {
        B2ERROR("Invalid CLEO cone order. The order must be in the [0,9] range.");
        return std::numeric_limits<float>::quiet_NaN();
      } else
        return m_cleoConesCollision[order];
    };
 
    float getCleoConeThrust(short order) const
    {
      if (order < 0 || order > 9) {
        B2ERROR("Invalid CLEO cone order. The order must be in the [0,9] range.");
        return std::numeric_limits<float>::quiet_NaN();
      } else
        return m_cleoConesThrust[order];
    };
 
 
    ROOT::Math::PxPyPzEVector getForwardHemisphere4Momentum() const
    {
      return m_forwardHemisphere4Momentum;
    };
 
    ROOT::Math::PxPyPzEVector getBackwardHemisphere4Momentum() const
    {
      return m_backwardHemisphere4Momentum;
    };
 
 
  private:
 
    // Axes
    ROOT::Math::XYZVector m_thrustAxis; 
    ROOT::Math::XYZVector m_spherocityAxis; 
    // Axis-independent quantities
    ROOT::Math::XYZVector m_sphericityEigenvector[3]; 
    float m_sphericityEigenvalue[3] = {0.}; 
    float m_foxWolframMoments[10] = {0.}; 
    // Axis-dependent quantities
    float m_thrust = 0; 
    float m_harmonicMomentsThrust[10] = {0.}; 
    float m_cleoConesThrust[10] = {0.}; 
    float m_harmonicMomentsCollision[10] = {0.}; 
    float m_cleoConesCollision[10] = {0.}; 
    // Hemisphere related quantities
    ROOT::Math::PxPyPzEVector m_forwardHemisphere4Momentum; 
    ROOT::Math::PxPyPzEVector m_backwardHemisphere4Momentum; 
    ClassDef(EventShapeContainer, 2) 
  };
 
 
} // end namespace Belle2