ECLShower.h

/**************************************************************************
 * basf2 (Belle II Analysis Software Framework)                           *
 * Author: The Belle II Collaboration                                     *
 *                                                                        *
 * See git log for contributors and copyright holders.                    *
 * This file is licensed under LGPL-3.0, see LICENSE.md.                  *
 **************************************************************************/
 
#pragma once
 
/* Basf2 headers. */
#include <framework/datastore/RelationsObject.h>
 
/* ROOT headers. */
#include <Math/Vector3D.h>
#include <TMatrixDSym.h>
 
/* C++ headers. */
#include <cmath>
 
namespace Belle2 {
  class ECLShower : public RelationsObject {
  public:
 
    enum Hypothesis : unsigned int {
      c_muonNPhotons = 1,
      c_chargedHadron = 2,
      c_electronNPhotons = 3,
      c_nPhotons = 5,
      c_neutralHadron = 6,
      c_mergedPi0 = 7
    };
 
    enum StatusBit : unsigned int  {
      c_hasDeadCrystal = 1 << 0,
 
      c_hasHotCrystal = 1 << 1,
 
      c_hasProblematicCrystal = c_hasDeadCrystal | c_hasHotCrystal,
 
      c_hasPulseShapeDiscrimination = 1 << 2,
 
    };
 
    ECLShower() :
      m_isTrk(false),                          
      m_status(0),                             
      m_showerId(0),                           
      m_connectedRegionId(0),                  
      m_hypothesisId(0),                       
      m_centralCellId(0),                      
      m_energy(0.0),                           
      m_energyRaw(0.0),                        
      m_theta(0.0),                            
      m_phi(0.0),                              
      m_r(0.0),                                
      m_Error{0.0, 0.0, 0.0, 0.0, 0.0, 0.0},   
      m_time(0),                               
      m_deltaTime99(0),                        
      m_energyHighestCrystal(0.0),             
      m_lateralEnergy(0.0),                    
      m_minTrkDistance(0.0),                   
      m_trkDepth(0.0),                         
      m_showerDepth(0.0),                      
      m_numberOfCrystals(0.0),                 
      m_absZernikeMoments{ -999.0},            
      m_zernikeMVA(0.0),                       
      m_secondMoment(0.0),                     
      m_E1oE9(0.0),                            
      m_E9oE21(0.0),                           
      m_ShowerHadronIntensity(0.0),            
      m_PulseShapeDiscriminationMVA(0.5),      
      m_NumberOfHadronDigits(0.0),             
      m_numberOfCrystalsForEnergy(0.0),        
      m_nominalNumberOfCrystalsForEnergy(0.0), 
      m_nOptimalGroupIndex(-1),                
      m_nOptimalEnergyBin(-1),                      
      m_nOptimalEnergy(0.0),                    
      m_listOfCrystalsForEnergy{},             
      m_listOfCrystalEnergyRankAndQuality{}    
    { }
 
    void setIsTrack(bool val) { m_isTrk = val; }
 
    void setStatus(int Status) { m_status = Status; }
 
    void setShowerId(int ShowerId) { m_showerId = ShowerId; }
 
    void setConnectedRegionId(int connectedRegionId) { m_connectedRegionId = connectedRegionId; }
 
    void setHypothesisId(int hypothesisId) { m_hypothesisId = hypothesisId; }
 
    void setCentralCellId(int centralCellId) { m_centralCellId = centralCellId; }
 
    void setEnergy(double Energy) { m_energy = Energy; }
 
    void setEnergyRaw(double EnergySum) { m_energyRaw = EnergySum; }
 
    void setTheta(double Theta) { m_theta = Theta; }
 
    void setPhi(double Phi) { m_phi = Phi; }
 
    void setR(double R) { m_r = R; }
 
    void setCovarianceMatrix(double covArray[6])
    {
      for (unsigned int i = 0; i < 6; i++) {
        m_Error[i] = covArray[i];
      }
    }
 
    void setTime(double Time) { m_time = Time; }
 
    void setDeltaTime99(double TimeReso) { m_deltaTime99 = TimeReso; }
 
    void setEnergyHighestCrystal(double HighestEnergy) { m_energyHighestCrystal = HighestEnergy; }
 
    void setLateralEnergy(double lateralEnergy) { m_lateralEnergy = lateralEnergy; }
 
    void setMinTrkDistance(double dist) { m_minTrkDistance = dist; }
 
    void setTrkDepth(double trkDepth) { m_trkDepth = trkDepth; }
 
    void setShowerDepth(double showerDepth) { m_showerDepth = showerDepth; }
 
    void setNumberOfCrystals(double nofCrystals) { m_numberOfCrystals = nofCrystals; }
 
    void setAbsZernikeMoment(unsigned int n, unsigned int m, double absZernikeMoment) { m_absZernikeMoments[(n * (n + 1)) / 2 + m - 1] = absZernikeMoment; }
 
    void setZernikeMVA(double zernikeMVA) {m_zernikeMVA = zernikeMVA; }
 
    void setSecondMoment(double secondMoment) { m_secondMoment = secondMoment; }
 
    void setE1oE9(double E1oE9) { m_E1oE9 = E1oE9; }
 
    void setE9oE21(double E9oE21) { m_E9oE21 = E9oE21; }
 
    void setShowerHadronIntensity(double hadronIntensity) { m_ShowerHadronIntensity = hadronIntensity; }
 
    void setPulseShapeDiscriminationMVA(double mvaVal) { m_PulseShapeDiscriminationMVA = mvaVal; }
 
    void setNumberOfHadronDigits(double NumberOfHadronDigits) { m_NumberOfHadronDigits = NumberOfHadronDigits; }
 
    void setNumberOfCrystalsForEnergy(double numberOfCrystalsForEnergy) { m_numberOfCrystalsForEnergy = numberOfCrystalsForEnergy; }
 
    void setNominalNumberOfCrystalsForEnergy(double nominalNumberOfCrystalsForEnergy) { m_nominalNumberOfCrystalsForEnergy = nominalNumberOfCrystalsForEnergy; }
 
    void setNOptimalGroupIndex(int nOptimalGroupIndex) { m_nOptimalGroupIndex = nOptimalGroupIndex; }
 
    void setNOptimalEnergyBin(int nOptimalEnergyBin) { m_nOptimalEnergyBin = nOptimalEnergyBin; }
 
    void setNOptimalEnergy(double nOptimalEnergy) { m_nOptimalEnergy = nOptimalEnergy; }
 
    void setListOfCrystalsForEnergy(const std::vector<unsigned int>& listofcrystals) { m_listOfCrystalsForEnergy = listofcrystals;}
 
    void setListOfCrystalEnergyRankAndQuality(std::vector<std::pair<unsigned int, bool>>
                                              listOfCrystalEnergyRankAndQuality) {m_listOfCrystalEnergyRankAndQuality = listOfCrystalEnergyRankAndQuality;}
 
    bool getIsTrack() const { return m_isTrk; }
 
    int getStatus() const { return m_status; }
 
    int getShowerId() const { return m_showerId; }
 
    int getConnectedRegionId() const { return m_connectedRegionId; }
 
    int getHypothesisId() const { return m_hypothesisId; }
 
    int getCentralCellId() const { return m_centralCellId; }
 
    double getEnergy() const { return m_energy; }
 
    double getEnergyRaw() const { return m_energyRaw; }
 
    double getTheta() const { return m_theta; }
 
    double getPhi() const { return m_phi; }
 
    double getR() const { return m_r; }
 
    void getCovarianceMatrixAsArray(double covArray[6]) const
    {
      for (unsigned int i = 0; i < 6; i++) {
        covArray[i] = m_Error[i];
      }
    }
 
    double getUncertaintyEnergy() const { return sqrt(m_Error[0]);}
 
    double getUncertaintyTheta() const { return  sqrt(m_Error[5]);}
 
    double getUncertaintyPhi() const {return sqrt(m_Error[2]);}
 
    double getTime() const { return m_time; }
 
    double getDeltaTime99() const { return m_deltaTime99; }
 
    double getEnergyHighestCrystal() const { return m_energyHighestCrystal; }
 
    double getLateralEnergy() const { return m_lateralEnergy; }
 
    double getMinTrkDistance() const { return m_minTrkDistance; }
 
    double getTrkDepth() const { return m_trkDepth; }
 
    double getShowerDepth() const { return m_showerDepth; }
 
    double getNumberOfCrystals() const { return m_numberOfCrystals; }
 
    double getAbsZernikeMoment(unsigned int n, unsigned int m) const { return m_absZernikeMoments[(n * (n + 1)) / 2 + m - 1]; }
 
    double getZernikeMVA() const {return m_zernikeMVA; }
 
    double getSecondMoment() const { return m_secondMoment; }
 
    double getE1oE9() const { return m_E1oE9; }
 
    double getE9oE21() const { return m_E9oE21; }
 
    double getShowerHadronIntensity() const { return m_ShowerHadronIntensity; }
 
    double getPulseShapeDiscriminationMVA() const { return m_PulseShapeDiscriminationMVA; }
 
    double getNumberOfHadronDigits() const { return m_NumberOfHadronDigits; }
 
    double getNumberOfCrystalsForEnergy() const { return m_numberOfCrystalsForEnergy; }
 
    double getNominalNumberOfCrystalsForEnergy() const { return m_nominalNumberOfCrystalsForEnergy; }
 
    int getNOptimalGroupIndex() const { return m_nOptimalGroupIndex; }
 
    int getNOptimalEnergyBin() const { return m_nOptimalEnergyBin; }
 
    double getNOptimalEnergy() const { return m_nOptimalEnergy; }
 
    std::vector<unsigned int>& getListOfCrystalsForEnergy()  { return m_listOfCrystalsForEnergy; }
 
    std::vector<std::pair<unsigned int, bool>> getListOfCrystalEnergyRankAndQuality() {return m_listOfCrystalEnergyRankAndQuality; }
 
 
 
    ROOT::Math::XYZVector getMomentum() const
    {
      return ROOT::Math::XYZVector(
               m_energy * std::sin(m_theta) * std::cos(m_phi),
               m_energy * std::sin(m_theta) * std::sin(m_phi),
               m_energy * std::cos(m_theta)
             );
    }
 
    TMatrixDSym getCovarianceMatrix3x3() const
    {
      TMatrixDSym covmatecl(3);
      covmatecl(0, 0) = m_Error[0];
      covmatecl(1, 0) = m_Error[1];
      covmatecl(1, 1) = m_Error[2];
      covmatecl(2, 0) = m_Error[3];
      covmatecl(2, 1) = m_Error[4];
      covmatecl(2, 2) = m_Error[5];
 
      //make symmetric
      for (int i = 0; i < 3; i++)
        for (int j = 0; j < i ; j++)
          covmatecl(j, i) = covmatecl(i, j);
      return covmatecl;
    }
 
    int getDetectorRegion() const
    {
      const double theta = getTheta();
 
      if (theta < 0.2164208) return 0;   // < 12.4deg
      if (theta < 0.5480334) return 1;   // < 31.4deg
      if (theta < 0.561996) return 11;   // < 32.2deg
      if (theta < 2.2462387) return 2;   // < 128.7deg
      if (theta < 2.2811453) return 13;   // < 130.7deg
      if (theta < 2.7070057) return 3;   // < 155.1deg
      else return 0;
    }
 
    bool hasStatus(unsigned short int bitmask) const { return (m_status & bitmask) == bitmask; }
 
    void addStatus(unsigned short int bitmask) { m_status |= bitmask; }
 
    bool hasHotCrystal() const;
 
    bool hasDeadCrystal() const;
 
    bool hasProblematicCrystal() const;
 
    bool hasPulseShapeDiscrimination() const;
 
    int getUniqueId() const
    {
      return 100000 * m_connectedRegionId + 1000 * m_hypothesisId + m_showerId;
    }
 
 
  private:
    bool m_isTrk;                   
    int m_status;                   
    int m_showerId;                 
    int m_connectedRegionId;        
    int m_hypothesisId;             
    int m_centralCellId;            
    Double32_t m_energy;            
    Double32_t m_energyRaw;         
    Double32_t m_theta;             
    Double32_t m_phi;               
    Double32_t m_r;                 
    Double32_t m_Error[6];          
    Double32_t m_time;              
    Double32_t m_deltaTime99;       
    Double32_t m_energyHighestCrystal;     
    Double32_t m_lateralEnergy;     
    Double32_t m_minTrkDistance;    
    Double32_t m_trkDepth;          
    Double32_t m_showerDepth;       
    Double32_t m_numberOfCrystals;       
    Double32_t m_absZernikeMoments[20];  
    Double32_t m_zernikeMVA;        
    Double32_t m_secondMoment;      
    Double32_t m_E1oE9;             
    Double32_t m_E9oE21;            
    Double32_t
    m_ShowerHadronIntensity;        
    Double32_t
    m_PulseShapeDiscriminationMVA;        
    Double32_t
    m_NumberOfHadronDigits;         
    Double32_t m_numberOfCrystalsForEnergy; 
    Double32_t m_nominalNumberOfCrystalsForEnergy; 
    int m_nOptimalGroupIndex; 
    int m_nOptimalEnergyBin; 
    Double32_t m_nOptimalEnergy; 
    std::vector<unsigned int> m_listOfCrystalsForEnergy; 
    std::vector<std::pair<unsigned int, bool>>
                                            m_listOfCrystalEnergyRankAndQuality; 
    // 2: added uniqueID and highestE (TF)
    // 3: added LAT and distance to closest track and trk match flag (GDN)
    // 4: added time resolution (TF)
    // 5: clean up, float to Double32_t, and new variables (TF)
    // 6: changed names of Zernike moment variables/getters/setters to indicate that they are the absolute values of the moments (TF and AH)
    // 7: added centralCellId (TF)
    // 8: added zernikeMVA, removed absZernike20, 42, 53 (AH)
    // 9: renamed variables according to the new mdst scheme (TF)
    // 10: added getUniqueId()
    // 11: added m_ShowerHadronIntensity and m_NumberOfHadronDigits variables (SL)
    // 12: added m_PulseShapeDiscriminationMVA.  Noted m_ShowerHadronIntensity will be removed in release-04 (SL)
    // 13: made enums strongly typed
    // 14: added m_numberOfCrystalsForEnergy of crystals for energy determination
    // 15: added m_listOfCrystalsForEnergy, m_nominalNumberOfCrystalsForEnergy
    // 16: removed m_absZernike40 and 51, added m_absZernikeMoments, m_listOfCrystalEnergyRankAndQuality (MH)
    // 17: added m_nOptimalGroupIndex, m_nOptimalEnergyBin, and m_nOptimalEnergy (CH)
    ClassDef(ECLShower, 17);
  };
 
  inline bool ECLShower::hasHotCrystal() const
  {
    return hasStatus(c_hasHotCrystal);
  }
 
  inline bool ECLShower::hasDeadCrystal() const
  {
    return hasStatus(c_hasDeadCrystal);
  }
 
  inline bool ECLShower::hasProblematicCrystal() const
  {
    return hasStatus(c_hasProblematicCrystal);
  }
 
  inline bool ECLShower::hasPulseShapeDiscrimination() const
  {
    return hasStatus(c_hasPulseShapeDiscrimination);
  }
 
} // end namespace Belle2