JShowerNPE.hh

Go to the documentation of this file.

#ifndef __JFIT__JSHOWERNPE__
#define __JFIT__JSHOWERNPE__

#include <cmath>

#include "JTools/JConstants.hh"
#include "JPhysics/JPDFToolkit.hh"
#include "JFit/JK40.hh"
#include "JFit/JFitToolkit.hh"


/**
 * \author mdejong, adomi
 */

namespace JFIT {}
namespace JPP { using namespace JFIT; }

namespace JFIT {
  
  /**
   * Auxiliary class for handling EM shower light yield.
   * Note that the effective light yield due to bremsstrahlung is proportional to the EM shower energy.
   */
  struct JShowerNPE :
    public JK40
  {
    
    using JK40::getP;

    
    /**
     * Default constructor.
     */
    JShowerNPE() :
      JK40(),
      __y(0.0)
    {}


    /**
     * Constructor.
     *
     * \param  k40                K40 light yields [npe]
     * \param  y                  light yield due to bremsstrahlung            [npe/GeV]
     */
    JShowerNPE(const JK40&  k40,
               const double y) :
      JK40(k40),
      __y(y)
    {}


    /**
     * Get light yield due to bremsstrahlung.
     *
     * \return                    light yield [npe/GeV]
     */
    double getY() const
    {
      return __y;
    }

    /**
     * Expected number of photo-electrons for random background hypothesis.
     *
     * \return                    light yield [npe]
     */
    inline double getH0() const
    {
      return this->getY0();
    }


    /**
     * Expected number of photo-electrons for random background hypothesis.
     *
     * \param  M                  multiplicity
     * \return                    light yield [npe]
     */
    inline double getH0(const size_t M) const
    {
      return this->getY1(M);
    }

    
    /**
     * Expected number of photo-electrons for EM shower hypothesis as a function of shower energy.
     *
     * \param  E_GeV              energy [GeV]
     * \return                    light yield [npe]
     */
    inline double getH1(const double E_GeV) const
    {
      using namespace JTOOLS;
      using namespace JPHYSICS;

      return E_GeV * this->getY();
    }


    /**
     * Get probability for observing a hit or not as a function of shower energy.
     *
     * \param  E_GeV              energy [GeV]
     * \param  hit                hit
     * \return                    probability
     */
    double getP(const double E_GeV, const bool hit) const
    {
      return JFIT::getP(this->getH1(E_GeV) + this->getH0(), hit);
    }


    /**
     * Get chi2 for observing a hit or not as a function of shower energy.
     *
     * \param  E_GeV              energy [GeV]
     * \param  hit                hit
     * \return                    probability
     */
    double getChi2(const double E_GeV, const bool hit) const
    {
      return JFIT::getChi2(this->getH1(E_GeV) + this->getH0(), hit);
    }

  protected:
    double __y;      //!< light yield       [npe/GeV]
  };
}

#endif