MuonPhaseII/MuonDigitization/sTgcDigitizationR4/sTgcDigitizationR4/sTgcDigitMaker.h

Go to the documentation of this file.

/*
  Copyright (C) 2002-2025 CERN for the benefit of the ATLAS collaboration
*/

 
#ifndef sTGCDigitizationR4_STGCDIGITMAKER_H
#define sTGCDigitizationR4_STGCDIGITMAKER_H
 
#include "AthenaBaseComps/AthMessaging.h"
#include "MuonCondData/DigitEffiData.h"
#include "MuonCondData/NswCalibDbThresholdData.h"
#include <MuonReadoutGeometryR4/sTgcReadoutElement.h>
#include "MuonReadoutGeometryR4/MuonDetectorManager.h"
#include "MuonDigitContainer/sTgcDigitContainer.h"
#include "CxxUtils/ArrayHelper.h"
#include "xAODMuonSimHit/MuonSimHitContainer.h"
#include "HitManagement/TimedHitPtr.h"
 
namespace CLHEP {
  class HepRandomEngine;
  class HepRandom;
}
 
namespace MuonR4{
class sTgcDigitMaker : public AthMessaging {
 public:
  enum class digitMode: std::uint8_t  {
    StripsOnly = 1,
    StripsAndPads = 2,
    AllChType = 3,
  };
  sTgcDigitMaker(const MuonGMR4::MuonDetectorManager* detMgr,
                 digitMode mode,
                 double meanGasGain, 
                 bool doPadChargeSharing);

  virtual ~sTgcDigitMaker();

  StatusCode initialize();

  struct DigiConditions {
     const Muon::DigitEffiData* efficiencies{nullptr};
     const NswCalibDbThresholdData* thresholdData{nullptr};
     CLHEP::HepRandomEngine* rndEngine{nullptr};
  };
  
using TimedHit = TimedHitPtr<xAOD::MuonSimHit>;
using ReadoutChannelType = sTgcIdHelper::sTgcChannelTypes;

  using sTgcDigitVec = std::vector<std::unique_ptr<sTgcDigit>>;
  sTgcDigitVec executeDigi(const DigiConditions& condContainers, const TimedHit& hit) const;
 
 private:
  struct GammaParameter {
    double lowEdge{0.};
    double kParameter{0.};
    double thetaParameter{0.};
  };

  struct Ionization {
    double distance{-9.99}; //smallest distance bet the wire and particle trajectory
    double time{0.};        // time of arrival
    Amg::Vector3D posOnSegment{Amg::Vector3D::Zero()}; // Point of closest approach
    Amg::Vector3D posOnWire{Amg::Vector3D::Zero()}; // Position on the wire
  };

  struct DigiInput {
      Identifier hitId{};
      Amg::Vector3D posOnSurf{Amg::Vector3D::Zero()};
      Amg::Vector3D hitDir{Amg::Vector3D::Zero()};
      double totalCharge{0.};
      double time{0.};
      const MuonGMR4::sTgcReadoutElement* reEle{nullptr};
  };

  bool getIonizationPoint(const TimedHit& hit, const DigiConditions& condContainers, Ionization& ionization) const;

  double calculateTotalCharge(double energyDeposit, CLHEP::HepRandomEngine* rndEngine) const;

  sTgcDigitVec processStripDigitization(const DigiConditions& condContainers,
                                        const DigiInput& digiInput) const;

  sTgcDigitVec processPadDigitization(const DigiInput& digiInput) const;

  sTgcDigitVec processWireDigitization(const DigiInput& digiInput) const;

  sTgcDigitVec processStripChargeSharing(const DigiInput& digiInput,
                                         const double peak_position,
                                         const int stripNumber) const;

  sTgcDigitVec processPadChargeSharing(const DigiInput& digiInput,
                                       const int padEta,
                                       const int padPhi) const;

  static void addDigit(sTgcDigitVec& digits, 
                       const Identifier& id, 
                       double digittime, 
                       double charge);

  StatusCode readFileOfTimeArrival();
  
  StatusCode readFileOfTimeOffsetStrip();
  
  Ionization pointClosestApproach(const MuonGMR4::StripLayer& stripLayer,
                                  int wireNumber, 
                                  const Amg::Vector3D& locHitPos,
                                  const Amg::Vector3D& locHitDir,
                                  const double stepLength) const;

  double getTimeOffsetStrip(size_t neighbor_index) const;
  
  static double getPadChargeFraction(double distance);

  GammaParameter getGammaParameter(double distance) const;
  
  double getMostProbableArrivalTime(double distance) const;
 
  // Parameters of the gamma probability distribution function
  std::vector<GammaParameter> m_gammaParameter;
  std::array<double, 5> m_mostProbableArrivalTime{make_array<double, 5>(0.)};
  std::array<double, 6> m_timeOffsetStrip{make_array<double, 6>(0.)};
 
  const MuonGMR4::MuonDetectorManager* m_detMgr{nullptr};
  const Muon::IMuonIdHelperSvc* m_idHelperSvc{m_detMgr->idHelperSvc()};
  const sTgcIdHelper& m_idHelper{m_idHelperSvc->stgcIdHelper()};
  // Computes the charge on a strip given the limits of integral using the error function
  // In R3, strip cluster charge profile was defined by a double gaussian for every digit
  // which is computationally inefficient, here we use result of double gaussian integral
  // to find the charge on each strip of the cluster, where M, N are upper and lower 
  // limits respectively 
   // Strip cluster charge profile: [0] = sigma of inner Gaussian, [1] = sigma of outer Gaussian
  static constexpr std::array<double, 2> m_clusterParams{0.573, 1.092};
  double chargeIntegral(double N, double M) const;

  // Digitization parameters
  digitMode m_digitMode{digitMode::AllChType};
  double m_theta{10}; // theta=10 value best matches the PDF
  double m_meanGasGain{5.e4};  // mean gain estimated from ATLAS note "ATL-MUON-PUB-2014-001"
  bool m_doPadSharing{false};
  
  // Flag to enable strip time offset
  bool m_doTimeOffsetStrip{false};
  // Angular strip resolution parameter
  double m_StripResolution{0.0949};
  double m_posResIncident{1.};
  double m_posResAngular{0.305/m_StripResolution};
 
  // Dependence of energy deposited on incident angle
  double m_chargeAngularFactor{4.0};
};
}
#endif