CaloCluster_v1.h

Go to the documentation of this file.

// -*- c++ -*-
/*
  Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
*/
#ifndef XAODCALOEVENT_VERSIONS_CALOCLUSTER_V1_H
#define XAODCALOEVENT_VERSIONS_CALOCLUSTER_V1_H
 
// System include(s):
extern "C" {
#   include <stdint.h>
}
#include <bitset>
#include <memory>
#include <vector>
#include <limits>
#include <bit>
 
#include <cmath>
 
// xAOD include(s):
#include "xAODBase/IParticle.h"
 
//Definition of CaloSamples (enum)
#include "CaloGeoHelpers/CaloSampling.h"
 
#include "xAODCaloEvent/CaloClusterBadChannelData.h"
#include "xAODCaloEvent/CaloClusterContainerFwd.h"
#include "AthLinks/ElementLink.h"
 
 
#if !(defined(SIMULATIONBASE) || defined(XAOD_ANALYSIS))
#include "CaloEvent/CaloClusterCellLinkContainer.h"
#include "CaloEvent/CaloRecoStatus.h"
#else
class CaloClusterCellLink {};
typedef unsigned CaloRecoStatus;
#endif // not (defined(SIMULATIONBASE) || defined(XAOD_ANALYSIS))
 
// ROOT include(s):
#include "Math/Vector4D.h"
 
//Already include the DataVector specialization for this type
#include "xAODCaloEvent/CaloClusterContainerFwd.h"
 
class CaloClusterChangeSignalState;
 
namespace xAOD {
   class CaloCluster_v1;
   typedef DataVector<CaloCluster_v1> CaloClusterContainer_v1;
 
   class CaloCluster_v1 : public IParticle {
     friend class ::CaloClusterChangeSignalState;
 
   public:
 
     typedef float flt_t;
 
     typedef CaloSampling::CaloSample CaloSample;
 
     typedef uint8_t ncells_t;
     typedef std::vector<uint16_t> ncells_store_t;
 
 
      enum ClusterSize {
         // electrons
         SW_55ele   = 1,
         SW_35ele   = 2,
         SW_37ele   = 3,
         // photons
         SW_55gam   = 4,
         SW_35gam   = 5,
         SW_37gam   = 6,
         // early converted photons
         SW_55Econv = 7,
         SW_35Econv = 8,
         SW_37Econv = 9,
         // soft electrons
         SW_softe   = 10,
         // topological clusters
         Topo_420   = 11,
         Topo_633   = 12,
         // transient cluster for AODCellContainer
         SW_7_11    = 13,
         //New (2016) egamma cluster
         SuperCluster=14,
         //New (2020) cluster representation of towers
         Tower_01_01 = 15,
         Tower_005_005 = 16,
         Tower_fixed_area = 17,
         CSize_Unknown = 99
      };
 
      enum MomentType {
         FIRST_PHI         = 101, 
         FIRST_ETA         = 102, 
         SECOND_R          = 201, 
         SECOND_LAMBDA     = 202, 
         DELTA_PHI         = 301,
         DELTA_THETA       = 302,
         DELTA_ALPHA       = 303, 
         CENTER_X          = 401, 
         CENTER_Y          = 402, 
         CENTER_Z          = 403, 
         CENTER_MAG        = 404,
         CENTER_LAMBDA     = 501, 
         LATERAL           = 601, 
         LONGITUDINAL      = 602, 
         ENG_FRAC_EM       = 701, 
         ENG_FRAC_MAX      = 702, 
         ENG_FRAC_CORE     = 703,
         FIRST_ENG_DENS    = 804, 
         SECOND_ENG_DENS   = 805, 
         ISOLATION         = 806,
         ENG_BAD_CELLS     = 807,
         N_BAD_CELLS       = 808, 
         N_BAD_CELLS_CORR  = 809,
         BAD_CELLS_CORR_E  = 813,
         BADLARQ_FRAC      = 821,
         ENG_POS           = 822, 
         SIGNIFICANCE      = 823, 
         CELL_SIGNIFICANCE = 824,
         CELL_SIG_SAMPLING = 825,
         AVG_LAR_Q         = 826,
         AVG_TILE_Q        = 827,
         ENG_BAD_HV_CELLS  = 828,
         N_BAD_HV_CELLS    = 829, 
         PTD               = 830,
         MASS              = 831,
         EM_PROBABILITY    = 900, 
         HAD_WEIGHT        = 901, 
         OOC_WEIGHT        = 902, 
         DM_WEIGHT         = 903, 
         TILE_CONFIDENCE_LEVEL = 904,
     SECOND_TIME      = 910,
     NCELL_SAMPLING   = 920,
 
         VERTEX_FRACTION = 1000, 
         NVERTEX_FRACTION = 1001, 
         ETACALOFRAME  = 1100, 
         PHICALOFRAME  = 1101, 
         ETA1CALOFRAME = 1102, 
         PHI1CALOFRAME = 1103, 
         ETA2CALOFRAME = 1104, 
         PHI2CALOFRAME = 1105, 
 
         ENG_CALIB_TOT     = 10001,
         ENG_CALIB_OUT_L   = 10010, 
         ENG_CALIB_OUT_M   = 10011, 
         ENG_CALIB_OUT_T   = 10012, 
         ENG_CALIB_DEAD_L  = 10020, 
         ENG_CALIB_DEAD_M  = 10021, 
         ENG_CALIB_DEAD_T  = 10022,  
         ENG_CALIB_EMB0         = 10030,
         ENG_CALIB_EME0         = 10031,
         ENG_CALIB_TILEG3       = 10032,
         ENG_CALIB_DEAD_TOT     = 10040,
         ENG_CALIB_DEAD_EMB0    = 10041,
         ENG_CALIB_DEAD_TILE0   = 10042,
         ENG_CALIB_DEAD_TILEG3  = 10043,
         ENG_CALIB_DEAD_EME0    = 10044,
         ENG_CALIB_DEAD_HEC0    = 10045,
         ENG_CALIB_DEAD_FCAL    = 10046,
         ENG_CALIB_DEAD_LEAKAGE = 10047,
         ENG_CALIB_DEAD_UNCLASS = 10048,
         ENG_CALIB_FRAC_EM      = 10051,
         ENG_CALIB_FRAC_HAD     = 10052,
         ENG_CALIB_FRAC_REST    = 10053,
 
 
         ENERGY_DigiHSTruth         = 40101, 
         ETA_DigiHSTruth         = 401024, 
         PHI_DigiHSTruth         = 401034, 
         TIME_DigiHSTruth         = 40104, 
         ENERGY_CALIB_DigiHSTruth         = 40105, 
         ETA_CALIB_DigiHSTruth         = 40106, 
         PHI_CALIB_DigiHSTruth         = 40107, 
         TIME_CALIB_DigiHSTruth         = 40108, 
         FIRST_PHI_DigiHSTruth         = 50101, 
         FIRST_ETA_DigiHSTruth         = 50102, 
         SECOND_R_DigiHSTruth          = 50201, 
         SECOND_LAMBDA_DigiHSTruth     = 50202, 
         DELTA_PHI_DigiHSTruth         = 50301,
         DELTA_THETA_DigiHSTruth       = 50302,
         DELTA_ALPHA_DigiHSTruth       = 50303,
         CENTER_X_DigiHSTruth          = 50401, 
         CENTER_Y_DigiHSTruth          = 50402, 
         CENTER_Z_DigiHSTruth          = 50403, 
         CENTER_MAG_DigiHSTruth        = 50404,
         CENTER_LAMBDA_DigiHSTruth     = 50501, 
         LATERAL_DigiHSTruth           = 50601, 
         LONGITUDINAL_DigiHSTruth      = 50602, 
         ENG_FRAC_EM_DigiHSTruth       = 50701, 
         ENG_FRAC_MAX_DigiHSTruth      = 50702, 
         ENG_FRAC_CORE_DigiHSTruth     = 75003,
         FIRST_ENG_DENS_DigiHSTruth    = 50804, 
         SECOND_ENG_DENS_DigiHSTruth   = 50805, 
         ISOLATION_DigiHSTruth         = 50806,
         ENG_BAD_CELLS_DigiHSTruth     = 50807,
         N_BAD_CELLS_DigiHSTruth       = 50808, 
         N_BAD_CELLS_CORR_DigiHSTruth  = 50809,
         BAD_CELLS_CORR_E_DigiHSTruth  = 50813,
         BADLARQ_FRAC_DigiHSTruth      = 50821,
         ENG_POS_DigiHSTruth           = 50822, 
         SIGNIFICANCE_DigiHSTruth      = 50823, 
         CELL_SIGNIFICANCE_DigiHSTruth = 50824,
         CELL_SIG_SAMPLING_DigiHSTruth = 50825,
         AVG_LAR_Q_DigiHSTruth         = 50826,
         AVG_TILE_Q_DigiHSTruth        = 50827,
         ENG_BAD_HV_CELLS_DigiHSTruth  = 50828,
         N_BAD_HV_CELLS_DigiHSTruth    = 50829, 
         EM_PROBABILITY_DigiHSTruth    = 50900, 
         HAD_WEIGHT_DigiHSTruth        = 50901, 
         OOC_WEIGHT_DigiHSTruth        = 50902, 
         DM_WEIGHT_DigiHSTruth         = 50903  
      };
 
     enum State {
       UNKNOWN             = -1,
       UNCALIBRATED        =  0,
       CALIBRATED          =  1,
       ALTCALIBRATED       =  2,
       NSTATES             =  3
     };
 
 
     CaloCluster_v1();
 
     CaloCluster_v1(const CaloCluster_v1& other);
 
     virtual ~CaloCluster_v1();
 
     CaloCluster_v1& operator=(const xAOD::CaloCluster_v1& other);
 
 
     virtual double           pt() const;
     virtual double           eta() const;
     virtual double           phi() const;
     virtual double           m() const;
     virtual double           e() const;
     virtual double           rapidity() const;
 
     typedef IParticle::FourMom_t FourMom_t;
 
     typedef ROOT::Math::LorentzVector<ROOT::Math::PtEtaPhiM4D<double> > GenVecFourMom_t;
 
     GenVecFourMom_t genvecP4() const;
 
     GenVecFourMom_t genvecP4(const State s) const;
 
     virtual FourMom_t p4() const;
 
     FourMom_t p4(const State s) const;
 
     virtual Type::ObjectType type() const;
 
     double et() const;
 
 
   public:
 
     float eSample(const CaloSample sampling) const;
     float etaSample(const CaloSample sampling) const;
     float phiSample(const CaloSample sampling) const;
 
     float energy_max(const CaloSample sampling) const;
     float etamax(const CaloSample sampling) const;
     float phimax(const CaloSample sampling) const;
 
     float etasize(const CaloSample sampling) const;
     float phisize(const CaloSample sampling) const;
 
     int numberCellsInSampling(const CaloSample samp,bool isInnerWheel=false) const; 
 
     void setNumberCellsInSampling(CaloSampling::CaloSample samp, int ncells,bool isInnerWheel=false); 
 
     int numberCells() const; 
 
     template<class CDATA>
     bool getNumberCellsInSampling(CDATA& cdata) const {
       ncells_store_t clist;
       if ( !retrieveMoment(NCELL_SAMPLING,clist) ) { return false; }
       cdata.clear(); cdata.reserve(clist.size());
       for ( auto cn : clist ) { cdata.push_back(extractLowerCount<typename CDATA::value_type>(cn)); }
       return true;
     }
 
 
     float energyBE(const unsigned layer) const;
 
     float etaBE(const unsigned layer) const;
     float phiBE(const unsigned layer) const;
 
 
     bool setEnergy(const CaloSample sampling, const float e);
     bool setEta(const CaloSample sampling, const float eta);
     bool setPhi(const CaloSample sampling, const float phi );
 
     bool setEmax(const CaloSample sampling, const float eMax );
     bool setEtamax(const CaloSample sampling, const float etaMax );
     bool setPhimax(const CaloSample sampling, const float phiMax );
     bool setEtasize(const CaloSample sampling, const float etaSize );
     bool setPhisize(const CaloSample sampling, const float phiSize );
 
 
 
 
 
 
     void insertMoment( MomentType type, double value );
 
     bool retrieveMoment( MomentType type, double& value ) const;
 
     double getMomentValue( MomentType type) const;
 
     void insertMoment( MomentType type, const ncells_store_t& values);
 
     bool retrieveMoment( MomentType type, ncells_store_t& values) const;
 
 
 
 
     void setEta0(flt_t);
     flt_t eta0() const;
 
     void setPhi0(flt_t);
     flt_t phi0() const;
 
     void setTime(flt_t);
     flt_t time() const;
     void setSecondTime(flt_t stime); 
     flt_t secondTime() const; 
 
     unsigned samplingPattern() const;
     void  setSamplingPattern(const unsigned sp, const bool clearSamplingVars=false);
 
     void clearSamplingData();
 
     unsigned nSamples() const;
 
     bool hasSampling(const CaloSample s) const;
 
 
     ClusterSize clusterSize() const;
     void  setClusterSize(const ClusterSize);
 
     bool inBarrel() const;
     bool inEndcap() const;
 
 
     void setE(flt_t);
     void setEta(flt_t);
     void setPhi(flt_t);
     void setM(flt_t);
 
 
     flt_t rawE() const;
     void  setRawE(flt_t);
     flt_t rawEta() const;
     void  setRawEta(flt_t);
     flt_t rawPhi() const;
     void  setRawPhi(flt_t);
     flt_t rawM() const;
     void  setRawM(flt_t);
 
     flt_t altE() const;
     void  setAltE(flt_t);
     flt_t altEta() const;
     void  setAltEta(flt_t);
     flt_t altPhi() const;
     void  setAltPhi(flt_t);
     flt_t altM() const;
     void  setAltM(flt_t);
 
     flt_t calE() const;
     void  setCalE(flt_t);
     flt_t calEta() const;
     void  setCalEta(flt_t);
     flt_t calPhi() const;
     void  setCalPhi(flt_t);
     flt_t calM() const;
     void  setCalM(flt_t);
#if !(defined(SIMULATIONBASE) || defined(XAOD_ANALYSIS))
   private:
#endif //not defined(SIMULATIONBASE) || defined(XAOD_ANALYSIS)
     bool setSignalState(const State s) ;
   public:
     State signalState() const {return m_signalState;}
 
     double pt(const State s) const;
 
     double e(const State s) const;
 
     double eta(const State s) const;
 
     double phi(const State s) const;
 
     double m(const State s) const;
 
 
 
 
     unsigned int getClusterEtaSize() const;
 
     unsigned int getClusterPhiSize() const;
 
     void setBadChannelList(const CaloClusterBadChannelList& bcl);
     const CaloClusterBadChannelList& badChannelList() const;
 
     const CaloCluster_v1* getSisterCluster() const;
 
     const ElementLink<xAOD::CaloClusterContainer_v1>& getSisterClusterLink() const;
 
     bool setSisterClusterLink(const ElementLink<CaloClusterContainer_v1>& sister);
 
     //For debugging only...
     //std::vector<std::pair<std::string,float> > getAllMoments();
 
   private:
     unsigned m_samplingPattern;
 
     State m_signalState;
     std::unique_ptr<CaloClusterCellLink> m_cellLinks;
 
     CaloRecoStatus m_recoStatus;
 
     double m_secondTime = { -1. };
 
     template<class UNSIGNED> UNSIGNED extractLowerCount(ncells_store_t::value_type cdata)  const { return static_cast<UNSIGNED>(cdata & 0x00ff); }        
     template<class UNSIGNED> UNSIGNED extractUpperCount(ncells_store_t::value_type cdata)  const { return static_cast<UNSIGNED>( (cdata & 0xff00)>>8 ); } 
     template<class UNSIGNED> ncells_store_t::value_type setLowerCount(ncells_store_t::value_type cdata,UNSIGNED clower)                            
     { return static_cast<ncells_store_t::value_type>((cdata & 0xff00)|(clower & 0x00ff)); }
     template<class UNSIGNED> ncells_store_t::value_type setUpperCount(ncells_store_t::value_type cdata,UNSIGNED cupper)                            
     { return static_cast<ncells_store_t::value_type>((cdata & 0x00ff)|((cupper & 0x00ff)<<8)); }
     template<class UNSIGNED> UNSIGNED lowerCellCountBound() const { return static_cast<UNSIGNED>(std::numeric_limits<ncells_t>::lowest()); }  
     template<class UNSIGNED> UNSIGNED upperCellCountBound() const { return static_cast<UNSIGNED>(std::numeric_limits<ncells_t>::max()); }     
     template<class UNSIGNED,class CCTYPE> CCTYPE adjustToRange(UNSIGNED count) const                                                          
     { return static_cast<CCTYPE>(std::min(std::max(count,lowerCellCountBound< UNSIGNED >()),upperCellCountBound< UNSIGNED >())); }
 
     unsigned sampVarIdx(const CaloSample) const;
 
     float getSamplVarFromAcc(const Accessor<std::vector<float > >& acc,
                  const CaloSample sampling, const float errorvalue=-999) const; //FIXME find a better errorcode
 
     bool setSamplVarFromAcc(const Accessor<std::vector<float> >& acc,
                 const CaloSample sampling, const float value);
 
   public:
#if !(defined(SIMULATIONBASE) || defined(XAOD_ANALYSIS))
 
     void addCellLink(CaloClusterCellLink* CCCL) {
       m_cellLinks.reset(CCCL);
     }
 
     void addCellLink(std::unique_ptr<CaloClusterCellLink> CCCL) {
       m_cellLinks=std::move(CCCL);
     }
     //bool createCellElemLink(const std::string& CCCL_key, const size_t index);
 
     bool setLink(CaloClusterCellLinkContainer* CCCL,
                  IProxyDict* sg = nullptr);
 
 
     bool setLink(CaloClusterCellLinkContainer* CCCL,
                  const EventContext& ctx);
 
 
     const CaloClusterCellLink* getCellLinks() const;
 
     CaloClusterCellLink* getOwnCellLinks() {
       return m_cellLinks.get();
     }
 
     bool addCell(const unsigned index, const double weight) {
       if (!m_cellLinks) {
         return false;
       }
       return m_cellLinks->addCell(index, weight);
     }
 
      bool removeCell(const CaloCell* ptr);
 
 
      size_t size() const;
 
      typedef CaloClusterCellLink::const_iterator const_cell_iterator;
      // Fixme: Check ret-val of getCellLinks (might be NULL);
      const_cell_iterator cell_cbegin() const
      {
        const CaloClusterCellLink* links = getCellLinks();
        if (!links) {
          return CaloClusterCellLink::dummyIt;
        }
        return links->begin();
      }
     const_cell_iterator cell_cend() const {
       const CaloClusterCellLink* links=getCellLinks();
       if (!links) {
         return CaloClusterCellLink::dummyIt;
       }
         return getCellLinks()->end();
     }
 
     const_cell_iterator cell_begin() const { return cell_cbegin(); }
     const_cell_iterator cell_end() const { return cell_cend(); }
 
     typedef CaloClusterCellLink::iterator cell_iterator;
     //Fixme: Check ret-val of getCellLinks (might be NULL);
     cell_iterator cell_begin() { return getOwnCellLinks()->begin();}
     cell_iterator cell_end() { return getOwnCellLinks()->end();}
 
     typedef const_cell_iterator const_iterator;
     typedef cell_iterator iterator;
     const_iterator begin() const { return cell_cbegin(); }
     const_iterator end() const { return cell_cend(); }
     const_iterator cbegin() const { return cell_cbegin(); }
     const_iterator cend() const { return cell_cend(); }
     iterator begin() { return cell_begin(); }
     iterator end() { return cell_end(); }
 
     void reweightCell(cell_iterator it, const double weight) {it.reweight(weight);}
 
     CaloRecoStatus& recoStatus() {return m_recoStatus;}
 
     const CaloRecoStatus& recoStatus() const {return m_recoStatus;}
 
#endif // not defined(SIMULATIONBASE) || defined(XAOD_ANALYSIS)
 
      void toPersistent();
 
   }; // class CaloCluster_v1
 
 
  inline  double CaloCluster_v1::et() const {
    if (this->m() == 0) {
      return this->pt();
    }
 
      return this->p4().Et();
  }
 
  inline unsigned CaloCluster_v1::samplingPattern() const {
    return m_samplingPattern;
  }
 
  inline unsigned CaloCluster_v1::sampVarIdx(const CaloCluster_v1::CaloSample s) const {
    const uint32_t& pattern= m_samplingPattern;
    //std::cout << "Pattern=" << std::hex << pattern << std::dec << ", Sampling=" << s << std::endl;
    if ((pattern & (0x1U << s)) == 0) {
      return CaloSampling::Unknown;
    }
    if (s == 0) {
      return 0;
    } // shifting a 32-bit int by 32 bits is undefined behavior!
    return std::popcount(pattern << (32 - s));
    // Explanation: Need to get the number of bit (=samples) before the sampling in question
    // Shift to the left, so bits after the sampling in question fall off the 32bit integer
    // Then use  popcount to count the numbers of 1 in the rest
  }
 
 
  inline unsigned CaloCluster_v1::nSamples() const {
    const uint32_t pattern=samplingPattern();
    return  std::popcount(pattern);
  }
 
 
  inline bool CaloCluster_v1::hasSampling(const CaloSample s) const {
   const unsigned pattern=samplingPattern();
   return (pattern & (0x1U<<(uint32_t)s));
 }
 
 
  inline bool CaloCluster_v1::inBarrel() const {
       return (m_samplingPattern & CaloSampling::barrelPattern());
  }
 
 
  inline bool CaloCluster_v1::inEndcap() const {
    return (m_samplingPattern & CaloSampling::endcapPattern());
  }
 
 
  inline double CaloCluster_v1::getMomentValue( CaloCluster_v1::MomentType type) const {
    double value=-999;
    (void)this->retrieveMoment(type,value);
    return value;
  }
 
} // namespace xAOD
 
// Finish declaration of IParticle as a base class of CaloCluster_v1:
DATAVECTOR_BASE_FIN( xAOD::CaloCluster_v1, xAOD::IParticle );
 
 
#endif // XAODCALOEVENT_VERSIONS_CALOCLUSTER_V1_H