JNPETable1D.hh

Go to the documentation of this file.

#ifndef __JNPETABLE1D__
#define __JNPETABLE1D__

#include <cmath>

#include "JTools/JSpline.hh"
#include "JTools/JElement.hh"
#include "JTools/JGridCollection.hh"
#include "JTools/JTransformableMultiFunction.hh"
#include "JGeometry3D/JOmega3D.hh"
#include "JPhysics/JNPETable.hh"


/**
 * \author mdejong
 */

namespace JPHYSICS {}
namespace JPP { using namespace JPHYSICS; }

namespace JPHYSICS {

  using JTOOLS::JSplineFunction1D;
  using JTOOLS::JSplineElement2D;
  using JTOOLS::JGridCollection;
  using JTOOLS::JTransformableMultiFunction;


  /**
   * Custom class for the integrated values of the PDF of the arrival time of Cherenkov light.
   *
   * This class provides for the average number of photo-electrons as a function of
   * the leading parameter value (<tt>R</tt>).
   * The average is taken over the orientations of the PMT in the full solid angle. 
   */
  template<class JFunction1D_t>
  class JNPETable1D :
    public JFunction1D_t
  {
  public:

    typedef typename JFunction1D_t::argument_type                  argument_type;
    typedef typename JFunction1D_t::result_type                    result_type;


    /**
     * Default constructor.
     */
    JNPETable1D()
    {}


    /**
     * Constructor.
     *
     * \param  input               multi-dimensional PDF
     */
    template<class JPDF_t, class JPDFMaplist_t, class JPDFDistance_t>
    JNPETable1D(const JTransformableMultiFunction<JPDF_t, JPDFMaplist_t, JPDFDistance_t>& input)
    {
      using namespace JTOOLS;
      using namespace JGEOMETRY3D;

      typedef JTransformableMultiFunction<JPDF_t, JPDFMaplist_t, JPDFDistance_t>   JTransformableMultiFunction_t;

      const JNPETable<typename JPDF_t::abscissa_type, typename JPDF_t::ordinate_type, JPDFMaplist_t, JPDFDistance_t> intensity(input);

      const JOmega3D omega(0.03*PI);

      for (typename JTransformableMultiFunction_t::const_iterator i = input.begin(); i != input.end(); ++i) {

        const double R  = i->getX(); 

        // average intensity

        double y = 0.0;

        for (JOmega3D::const_iterator dir = omega.begin(); dir != omega.end(); ++dir) {
          y += intensity(R, dir->getTheta(), fabs(dir->getPhi()));
        }

        this->put(R,y/omega.size());
      }

      this->compile();
      this->setExceptionHandler(new typename JFunction1D_t::JDefaultResult(JMATH::zero));
    }
  };
}

#endif