CaloClusterDetails.h

Go to the documentation of this file.

/*
  Copyright (C) 2002-2025 CERN for the benefit of the ATLAS collaboration
*/
 
#ifndef XAODCALOEVENT_CALOCLUSTERDETAILS_H
#define XAODCALOEVENT_CALOCLUSTERDETAILS_H
 
#include <CaloGeoHelpers/CaloPhiRange.h>
#include <CaloGeoHelpers/CaloSampling.h>
#include <bit>
#include <cstdint>
#include <span>
 
namespace xAOD
{
  namespace CaloClusterDetails
  {
 
 
    using CaloSample = CaloSampling::CaloSample;
 
 
    constexpr float defaultErrorValue = -999;
 
 
    [[nodiscard]] inline bool hasSampling(const CaloSample s, const std::uint32_t samplingPattern) {
      return (samplingPattern & (0x1U<<(std::uint32_t)s));
    }
 
 
    [[nodiscard]] inline unsigned sampVarIdx(const CaloSample s, const std::uint32_t samplingPattern) {
      //std::cout << "Pattern=" << std::hex << pattern << std::dec << ", Sampling=" << s << std::endl;
      if ((samplingPattern & (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(samplingPattern << (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
    }
 
 
    [[nodiscard]] inline float getSamplVar(const CaloSample sampling, const std::uint32_t samplingPattern, const std::span<const float> vec, const float errorvalue=defaultErrorValue) {
      const unsigned idx=sampVarIdx(sampling, samplingPattern);
      if (idx<vec.size() ) {
        return vec[idx];
      }
  
        //std::cout <<Sampling " << sampling << ", Pattern=" << std::hex <<m_samplingPattern << std::dec << ", index=" << idx << " size=" << vec.size() << std::endl;
      return errorvalue;
    }
 
 
    [[nodiscard]] inline float energyBE(const unsigned sample, const std::uint32_t samplingPattern, const std::span<const float> e_sampl) {
      if (sample>3) return defaultErrorValue;
      const CaloSample barrelSample=(CaloSample)(CaloSampling::PreSamplerB+sample);
      const CaloSample endcapSample=(CaloSample)(CaloSampling::PreSamplerE+sample);
      double energy=0;
      if (hasSampling(barrelSample, samplingPattern)) {
        energy+=getSamplVar(barrelSample,samplingPattern,e_sampl); //Check for errorcode? Should not happen...
      }
      if (hasSampling(endcapSample, samplingPattern)) {
        energy+=getSamplVar(endcapSample,samplingPattern,e_sampl);
      }
      return energy;
    }
 
 
 
    [[nodiscard]] inline float etaBE(const unsigned sample, const std::uint32_t samplingPattern, const std::span<const float> e_sampl, const std::span<const float> eta_sampl) {
      if (sample>3) {return defaultErrorValue;}
      const CaloSample barrelSample=(CaloSample)(CaloSampling::PreSamplerB+sample);
      const CaloSample endcapSample=(CaloSample)(CaloSampling::PreSamplerE+sample);
      const bool haveBarrel=hasSampling(barrelSample, samplingPattern);
      const bool haveEndcap=hasSampling(endcapSample, samplingPattern);
      if (haveBarrel && haveEndcap) {
        //cluster spans barren and endcap
        float eBarrel=getSamplVar(barrelSample,samplingPattern,e_sampl);  //Check for errorcode? Should not happen...
        float eEndcap=getSamplVar(endcapSample,samplingPattern,e_sampl);
 
        float etaBarrel=getSamplVar(barrelSample,samplingPattern,eta_sampl);
        float etaEndcap=getSamplVar(endcapSample,samplingPattern,eta_sampl);
        float eSum=eBarrel + eEndcap;
        if (eSum > 100 /*MeV*/) {
          //E-weighted average ...
          if ((eBarrel > 0 && eEndcap > 0) || (eBarrel < 0 && eEndcap < 0))
            return (eBarrel * etaBarrel + eEndcap * etaEndcap) / eSum;
          else if (eBarrel > 0)
            return etaBarrel;
          else
            return etaEndcap;
         }//else eSum==0 case, should never happen
         return (0.5 * (etaBarrel + etaEndcap));
      }
      if (haveBarrel) {
        return getSamplVar(barrelSample,samplingPattern,eta_sampl);
      }
      if (haveEndcap) {
        return getSamplVar(endcapSample,samplingPattern,eta_sampl);
      }
  
      //Should never reach this point ...
      return defaultErrorValue;
    }
 
 
 
    [[nodiscard]] inline float phiBE(const unsigned sample, const std::uint32_t samplingPattern, const std::span<const float> e_sampl, const std::span<const float> phi_sampl) {
      if (sample>3) {return defaultErrorValue;}
      const CaloSample barrelSample=(CaloSample)(CaloSampling::PreSamplerB+sample);
      const CaloSample endcapSample=(CaloSample)(CaloSampling::PreSamplerE+sample);
      const bool haveBarrel=hasSampling(barrelSample, samplingPattern);
      const bool haveEndcap=hasSampling(endcapSample, samplingPattern);
      if (haveBarrel && haveEndcap) {
        //cluster spans barren and endcap
        float eBarrel=getSamplVar(barrelSample, samplingPattern, e_sampl);  //Check for errorcode? Should not happen...
        float eEndcap=getSamplVar(endcapSample, samplingPattern, e_sampl);
        float eSum=eBarrel+eEndcap;
        float phiBarrel=getSamplVar(barrelSample, samplingPattern, phi_sampl);
        float phiEndcap=getSamplVar(endcapSample, samplingPattern, phi_sampl);
        if (eSum != 0.0) {
          if ((eBarrel > 0 && eEndcap > 0) || (eBarrel < 0 && eEndcap < 0)) {
            float phiSum = eSum * phiBarrel + eEndcap * CaloPhiRange::diff(phiEndcap, phiBarrel);
            return CaloPhiRange::fix(phiSum / eSum);
          } else if (eBarrel > 0)
            return phiBarrel;
          else
            return phiEndcap;
        }
        // energy==0 case, should never happen
        return CaloPhiRange::fix(0.5 * (phiBarrel + phiEndcap));
      }
      if  (haveBarrel) {
        return getSamplVar(barrelSample, samplingPattern, phi_sampl);
      }
      if (haveEndcap) {
        return getSamplVar(endcapSample, samplingPattern, phi_sampl);
      }
 
      //Should never reach this point ...
      return defaultErrorValue;
    }
  }
}
 
#endif