JOscProbToolkit.hh

Go to the documentation of this file.

#ifndef __JOSCPROB__JOSCPROBTOOLKIT__
#define __JOSCPROB__JOSCPROBTOOLKIT__

#include "km3net-dataformat/offline/Evt.hh"
#include "km3net-dataformat/offline/Trk.hh"

#include "Jeep/JMessage.hh"

#include "JLang/JException.hh"

#include "JTools/JRange.hh"


/**
 * \author bjung
 * Auxiliary methods for oscillation probabilities.
 */

namespace JOSCPROB {}
namespace JPP { using namespace JOSCPROB; }

namespace JOSCPROB {

  using JLANG::JValueOutOfRange;
  
  using JTOOLS::JRange;
  

  /**
   * Allowed ranges of neutrino oscillation parameters.
   */
  static const JRange<double>   SINSQTH12Range(0.250,     0.354);
  static const JRange<double>   DM21Range     (6.93e-5,   7.96e-5);   // [eV^2]
  static const JRange<double>   SINSQTH13Range(0.019,     0.025);
  static const JRange<double>   DM31Range     (+2.45e-3, +2.69e-3);   // [eV^2]
  static const JRange<double>   SINSQTH23Range(0.381,     0.636);
  static const JRange<double>   DCPRange      (0.000,     2.000);     // [pi]

  
  /**
   * Earth radii.
   */
  static const double R_EARTH_BELOW = 6368; //!< Earth radius, up until detector layer        [km]
  static const double R_EARTH_ABOVE = 6386; //!< Earth radius, including layer above detector [km]


  /**
   * Get baseline for a given cosine zenith angle.
   *
   * \param  costh            cosine zenith angle
   * \return                  baseline [km]
   */
  inline double getL(const double costh)
  {
    using namespace std;
    using namespace JPP;

    static const double epsilon = 1e-6;
    
    if (fabs(costh) < 1.0 || fabs(costh) - 1.0 < epsilon) {

      return (-R_EARTH_BELOW * costh + sqrt(R_EARTH_ABOVE*R_EARTH_ABOVE -
                                            R_EARTH_BELOW*R_EARTH_BELOW * (1 - costh) * (1 + costh)));
      
    } else {

      THROW(JValueOutOfRange, "getL(): Invalid cosine zenith angle " << costh << endl);
    }
  }
}

#endif