NeutrinoFlux.h

Go to the documentation of this file.

#ifndef NEUTRINOFLUX_H
#define NEUTRINOFLUX_H

// IceTray specific stuff
#ifdef ICETRAY
#include "dataclasses/Utility.h"
#include "icetray/I3DefaultName.h"
#include <icetray/I3Logging.h>
#include <dataclasses/physics/I3Particle.h>
#else  //ICETRAY
namespace I3Particle{ 
  enum ParticleType {
    unknown = 0,
    NuE = 66,
    NuEBar = 67,
    NuMu = 68,
    NuMuBar = 69,
    TauPlus = 131,
    TauMinus = 132,
    NuTau = 133,
    NuTauBar = 134,
    };
}
#endif //ICETRAY

#include <iostream>
#include <fstream>
#include <math.h>

using namespace std;

class NeutrinoFlux
{
 public:
  
  virtual ~NeutrinoFlux() {};

  /* Methods */    
  virtual double getFlux(I3Particle::ParticleType neutType, double eNeut, double costheta)const=0;

  virtual double getFluxIntegral(I3Particle::ParticleType neutType, double eNeut)const=0;

  virtual string getModelName()const=0;


#ifdef ICETRAY

  double getFlux(const I3Particle& particle) const;

 SET_LOGGER("NeutrinoFlux");

#endif
  
};


class ConventionalNeutrinoFlux : public NeutrinoFlux
{
 public:

  ConventionalNeutrinoFlux(string model);

  virtual double getFlux(I3Particle::ParticleType neutType, double eNeut, double costheta)const;
  virtual double getFluxIntegral(I3Particle::ParticleType neutType, double eNeut)const;
  virtual string getModelName()const  {return model_;};

 private:

  string model_;

  static const int POLYT = 31;
  static const int nbinCT_ = 10;
  static const double CMIN=-1.;
  static const double CMAX=1.;
  
  // tables for nu, anu, respectively
  double par_conv_nu_[POLYT];
  double par_conv_anu_[POLYT];

};

class PromptNeutrinoFlux : public NeutrinoFlux
{
 public:
  PromptNeutrinoFlux(string model);

  double getFlux(I3Particle::ParticleType neutType, double eNeut, double costheta) const;
  double getFluxIntegral(I3Particle::ParticleType neutType, double eNeut) const;
  string getModelName() const {return model_;};

 private:

  string model_;

  static const int nbinCT_prompt_ = 10;
  static const double CMIN_PROMPT=-1.;
  static const double CMAX_PROMPT=1.; 

  static const int POLYT_PROMPT=15;
  double e_joint_nu_;
  double e_joint_anu_;
  //parameters for nu
  double par_prompt_nu_[POLYT_PROMPT];
  //nu parameter errors
  double pare_prompt_nu_[POLYT_PROMPT];
  //parameters for anu
  double par_prompt_anu_[POLYT_PROMPT];
  //anu parameter errors
  double pare_prompt_anu_[POLYT_PROMPT];

};


class AtmosphericNeutrinoFlux: public NeutrinoFlux
{
 public:

  AtmosphericNeutrinoFlux(string conv_model,string prompt_model);

  double getFlux(I3Particle::ParticleType neutType, double eNeut, double costheta) const;
  double getFluxIntegral(I3Particle::ParticleType neutType, double eNeut) const;
  string getModelName() const;

 private:
  ConventionalNeutrinoFlux conv_flux_;
  PromptNeutrinoFlux prompt_flux_;

};

class AstroPointNeutrinoFlux: public NeutrinoFlux
{
 public:
  AstroPointNeutrinoFlux(string model);

  double getFlux(I3Particle::ParticleType neutType, double eNeut, double costheta) const;
  double getFlux(I3Particle::ParticleType neutType, double eNeut) const;
  double getFluxIntegral(I3Particle::ParticleType neutType, double eNeut) const;
  double getDecli() const {return decli_;} ;
  string getModelName() const {return model_;} ;

 private:

  string model_;

  int piece_;
  double decli_;
  double Emin_[100];
  double Emax_[100];
  int nparam_[100];
  double param_[100][100];
  //1- energy limits
  static const double Em_=1.;
  static const double Ema_=1.0e+12;

};

class AstroNeutrinoFlux: public NeutrinoFlux
{
 public:
  AstroNeutrinoFlux(string model);

  double getFlux(I3Particle::ParticleType neutType, double eNeut, double costheta) const;
  double getFlux(I3Particle::ParticleType neutType, double eNeut) const;
  double getFluxIntegral(I3Particle::ParticleType neutType, double eNeut) const;
  string getModelName() const {return model_;} ;

 private:

  string model_;

  int piece_;
  double Emin_[100];
  double Emax_[100];
  int nparam_[100];
  double param_[100][100];

  //1- energy limits
  static const double Em_=1.;
  static const double Ema_=1.0e+12;

};

                                                                                                    
#ifdef ICETRAY
I3_POINTER_TYPEDEFS(NeutrinoFlux);
I3_DEFAULT_NAME(NeutrinoFlux);
#endif //ICETRAY
    
#endif // NEUTRINOFLUX_H

---