FLSHit.h

Go to the documentation of this file.

////////////////////////////////////////////////////////////////////////
/// \brief   A single unit of energy deposition in the liquid scintillator
///
/// \author  messier@indiana.edu
/// \date
////////////////////////////////////////////////////////////////////////
#ifndef SIM_FLSHIT_H
#define SIM_FLSHIT_H
#include "GeometryObjects/CellUniqueId.h"

#include <iosfwd>
#include <iostream>
#include <iomanip>
#include <numeric>


namespace sim {
  /// A single unit of energy deposition in the liquid scintillator
  class FLSHit {
  public:
    FLSHit();
    ~FLSHit();
    
    bool operator==(const FLSHit& rhs) const;
    bool operator<(const FLSHit& rhs) const;

    /// Clear the FLS hit
    void Clear();
    
    /// Get total Energy deposited into the cell for the whole FLSHit
    float             GetEdep()         const {return getSum(fVectorEdep);}
    /// Get total Energy with Birks suppression deposited into the cell for the whole FLSHit
    float             GetEdepBirks()    const {return getSum(fVectorEdepBirks);}
    /// Get total N Cerenkov photons
    float             GetNCerenkov()    const {return getSum(fVectorNCerenkov);}
    /// Plane ID
    int               GetPlaneID()      const {return fPlaneId;}
    /// Cell ID
    int               GetCellID ()      const {return fCellId;}
    /// Unique Cell ID
    geo::CellUniqueId GetCellUniqueId() const {return fId;}
    /// PDG
    int               GetPDG    ()      const {return fPDG;}
    /// Track ID
    int               GetTrackID()      const {return fTrackId;}

    /// Entry point of the particle (position, time and energy)
    float GetEntryX()      const {return getFirst(fX);}
    float GetEntryY()      const {return getFirst(fY);}
    float GetEntryZ()      const {return getFirst(fZ);}
    float GetEntryT()      const {return getFirst(fT);}
    float GetEntryEnergy() const {return getFirst(fEnergy);}
    /// Exit point of the particle (position, time and energy)
    float GetExitX()       const {return getLast(fX);}
    float GetExitY()       const {return getLast(fY);}
    float GetExitZ()       const {return getLast(fZ);}
    float GetExitT()       const {return getLast(fT);}
    float GetExitEnergy()  const {return getLast(fEnergy);}

    /// Number of end points in FLSHit (could be Geant end points)
    int   GetNPoints() const { return fX.size();}
    /// Number of steps in FLSHit. Should be equal to GetNPoints() - 1
    int   GetNSteps () const { return fVectorEdep.size();}

    /// Get path length of a particular step
    float GetPathLength(const int step) const;
    /// Get path length of all steps in FLSHit
    float GetTotalPathLength() const;

    /// Get energy loss of a particular step
    float GetEnergyLoss(const int step) const;
    /// Get total energy loss in the entire FLSHit
    float GetTotalEnergyLoss() const;

    /// Get X-average for the step. This is in local coordinates.
    float GetXAverage(const int step) const {return GetXYZAverage(fX, step);}
    /// Get Y-average for the step. This is in local coordinates.
    float GetYAverage(const int step) const {return GetXYZAverage(fY, step);}
    /// Get Z-average for the step. This is in local coordinates.
    float GetZAverage(const int step) const {return GetXYZAverage(fZ, step);}
    /// Get X-average for the FLSHit. This is in local coordinates.
    float GetXAverage() const { return GetXYZAverage(fX);}
    /// Get Y-average for the FLSHit. This is in local coordinates.
    float GetYAverage() const { return GetXYZAverage(fY);}
    /// Get Z-average for the FLSHit. This is in local coordinates.
    float GetZAverage() const { return GetXYZAverage(fZ);}


    /// Get a point of the particle in the FLSHit (position, time and energy)
    float GetX        (const int ipoint) const {return getElement(fX              , ipoint);}
    float GetY        (const int ipoint) const {return getElement(fY              , ipoint);}
    float GetZ        (const int ipoint) const {return getElement(fZ              , ipoint);}
    float GetT        (const int ipoint) const {return getElement(fT              , ipoint);}
    float GetEnergy   (const int ipoint) const {return getElement(fEnergy         , ipoint);}
    float GetEdep     (const int ipoint) const {return getElement(fVectorEdep     , ipoint);}
    float GetEdepBirks(const int ipoint) const {return getElement(fVectorEdepBirks, ipoint);}
    float GetNCerenkov(const int ipoint) const {return getElement(fVectorNCerenkov, ipoint);}

    /// Set methods
    void  SetPDG         (const int pdg)               { fPDG = pdg; }
    void  SetTrackId     (const int trackid)           { fTrackId = trackid; }
    void  SetPlane       (const int ip)                { fPlaneId = ip; }
    void  SetCell        (const int ic)                { fCellId = ic; }
    void  SetCellUniqueId(const geo::CellUniqueId& id) { fId = id;}
    void  AddEdep        (const float edep)            { fVectorEdep.push_back(edep); }
    void  AddEdepBirks   (const float edepbirks)       { fVectorEdepBirks.push_back(edepbirks);}
    void  AddNCerenkov   (const float ncerenkov)       { fVectorNCerenkov.push_back(ncerenkov); }
    void  AddPos (const float x, const float y, const float z, const double t, const float energy) {fX.push_back(x);fY.push_back(y);fZ.push_back(z);fT.push_back(t);fEnergy.push_back(energy);}

    /// Performs various checks on the size of the data and its content.
    /// Has an option to print errors out, which is turned off by default.
    bool IsFLSHitReasonable(const bool    print_errors          = false
                           ,std::ostream& error_printint_stream = std::cout) const;

  private:

    /// Sums over the array of the vector
    template< class T>
    T getSum(const std::vector<T>& input_vector) const{
        const T sum = std::accumulate(input_vector.begin(), input_vector.end(), 0.0);
        return sum;
    }// end of getSum

    /// Gets i-th element of the vector
    template< class T>
    T getElement(const std::vector<T>& input_vector, const unsigned int ipoint) const{
        return input_vector.at(ipoint);
    }// end of getElement

    /// Gets the last element of the vector
    template< class T>
    T getFirst(const std::vector<T>& input_vector) const{
        return getElement<T>(input_vector, 0);
    }// end of getFirst
    
    /// Gets the last element of the vector
    template< class T>
    T getLast(const std::vector<T>& input_vector) const{
        const int n_size = input_vector.size();
        if(n_size == 0) return (T)(0.0);
        return input_vector[n_size - 1];
    }// end of getLast

    float GetXYZAverage(const std::vector<float>& xyz, const int step) const;
    float GetXYZAverage(const std::vector<float>& xyz) const;

    /// Checks whether the input vector increases monotonically
    template< class T>
    bool DoesIncreaseMonotonically(const std::vector<T>& input_vector) const {

        const int n_size = input_vector.size();
        
        /// if size is less than 2, then it's monotonic, obviously
        if(n_size < 2) return true;

        for(int i=1; i<n_size; ++i){
            if(input_vector[i-1] > input_vector[i])
                return false;
        }

        return true;
    }// end of DoesIncreaseMonotonically

    /// Checks whether the input vector decreases monotonically
    template< class T>
    bool DoesDecreaseMonotonically(const std::vector<T>& input_vector) const {

        const int n_size = input_vector.size();

        /// if size is less than 2, then it's monotonic, obviously
        if(n_size < 2) return true;

        for(int i=1; i<n_size; ++i){
            if(input_vector[i-1] < input_vector[i])
                return false;
        }

        return true;
    }// end of DoesDecreaseMonotonically


  private:

    geo::CellUniqueId   fId;              ///< Unique cell ID

    unsigned short int  fPlaneId;         ///< Plane number
    unsigned short int  fCellId;          ///< Cell number
    int                 fPDG;             ///< PDG code of parent particle
    int                 fTrackId;         ///< Index number of the Global track producing hit

    std::vector<float>  fX;               ///< Start or end of the step inside scintillator. X position, cm
    std::vector<float>  fY;               ///< Start or end of the step inside scintillator. Y position, cm
    std::vector<float>  fZ;               ///< Start or end of the step inside scintillator. Z position, cm
    std::vector<float>  fEnergy;          ///< Start or end of the step inside scintillator. Energy, GeV
    std::vector<double> fT;               ///< Start or end of the step inside scintillator. Time, ns

    std::vector<float>  fVectorEdep;      ///< Amount of energy deposited during a Geant4 step (GeV)
    std::vector<float>  fVectorEdepBirks; ///< Amount of energy deposited during a Geant4 step ; Birks weighted (GeV)
    std::vector<float>  fVectorNCerenkov; ///< Amount of Cerenkov light emitteded in scintillator
  };


  std::ostream& operator<< (std::ostream & o, const FLSHit & a);
  
}
#endif
////////////////////////////////////////////////////////////////////////