Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-10-25 10:04:15

0001 #include "SimG4CMS/HGCalTestBeam/interface/AHCalSD.h"
0002 #include "SimG4CMS/HGCalTestBeam/interface/AHCalDetId.h"
0003 #include "SimG4Core/Notification/interface/TrackInformation.h"
0004 #include "Geometry/HGCalCommonData/interface/AHCalParameters.h"
0005 
0006 #include "DataFormats/DetId/interface/DetId.h"
0007 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
0008 
0009 #include "G4LogicalVolume.hh"
0010 #include "G4LogicalVolumeStore.hh"
0011 #include "G4ParticleTable.hh"
0012 #include "G4Track.hh"
0013 #include "G4VProcess.hh"
0014 
0015 #include "G4PhysicalConstants.hh"
0016 #include "G4SystemOfUnits.hh"
0017 
0018 #include <iomanip>
0019 #include <map>
0020 
0021 //#define EDM_ML_DEBUG
0022 
0023 AHCalSD::AHCalSD(const std::string& name,
0024                  const SensitiveDetectorCatalog& clg,
0025                  edm::ParameterSet const& p,
0026                  const SimTrackManager* manager)
0027     : CaloSD(name,
0028              clg,
0029              p,
0030              manager,
0031              (float)(p.getParameter<edm::ParameterSet>("AHCalSD").getParameter<double>("TimeSliceUnit")),
0032              p.getParameter<edm::ParameterSet>("AHCalSD").getParameter<bool>("IgnoreTrackID")) {
0033   edm::ParameterSet m_HC = p.getParameter<edm::ParameterSet>("AHCalSD");
0034   useBirk = m_HC.getParameter<bool>("UseBirkLaw");
0035   birk1 = m_HC.getParameter<double>("BirkC1") * (CLHEP::g / (CLHEP::MeV * CLHEP::cm2));
0036   birk2 = m_HC.getParameter<double>("BirkC2");
0037   birk3 = m_HC.getParameter<double>("BirkC3");
0038   eminHit = m_HC.getParameter<double>("EminHit") * CLHEP::MeV;
0039 
0040   edm::LogVerbatim("HcalSim") << "AHCalSD::  Use of Birks law is set to      " << useBirk
0041                               << "  with three constants kB = " << birk1 << ", C1 = " << birk2 << ", C2 = " << birk3
0042                               << "\nAHCalSD:: Threshold for storing"
0043                               << " hits: " << eminHit;
0044 }
0045 
0046 double AHCalSD::getEnergyDeposit(const G4Step* aStep) {
0047   double destep = aStep->GetTotalEnergyDeposit();
0048   double wt2 = aStep->GetTrack()->GetWeight();
0049   double weight = (wt2 > 0.0) ? wt2 : 1.0;
0050 #ifdef EDM_ML_DEBUG
0051   double weight0 = weight;
0052 #endif
0053   if (useBirk)
0054     weight *= getAttenuation(aStep, birk1, birk2, birk3);
0055 #ifdef EDM_ML_DEBUG
0056   edm::LogVerbatim("HcalSim") << "AHCalSD: weight " << weight0 << " " << weight;
0057 #endif
0058   double edep = weight * destep;
0059   return edep;
0060 }
0061 
0062 uint32_t AHCalSD::setDetUnitId(const G4Step* aStep) {
0063   const G4StepPoint* preStepPoint = aStep->GetPreStepPoint();
0064   const G4VTouchable* touch = preStepPoint->GetTouchable();
0065 
0066   int depth = (touch->GetReplicaNumber(1));
0067   int incol = ((touch->GetReplicaNumber(0)) % AHCalParameters::kColumn_);
0068   int inrow = ((touch->GetReplicaNumber(0)) / AHCalParameters::kColumn_) % AHCalParameters::kRow_;
0069   int jncol = ((touch->GetReplicaNumber(0)) / AHCalParameters::kRowColumn_) % AHCalParameters::kSign_;
0070   int jnrow = ((touch->GetReplicaNumber(0)) / AHCalParameters::kSignRowColumn_) % AHCalParameters::kSign_;
0071   int col = (jncol == 0) ? incol : -incol;
0072   int row = (jnrow == 0) ? inrow : -inrow;
0073   uint32_t index = AHCalDetId(row, col, depth).rawId();
0074 #ifdef EDM_ML_DEBUG
0075   edm::LogVerbatim("HcalSim") << "AHCalSD: det = " << HcalOther << " depth = " << depth << " row = " << row
0076                               << " column = " << col << " packed index = 0x" << std::hex << index << std::dec;
0077   bool flag = unpackIndex(index, row, col, depth);
0078   edm::LogVerbatim("HcalSim") << "Results from unpacker for 0x" << std::hex << index << std::dec << " Flag " << flag
0079                               << " Row " << row << " Col " << col << " Depth " << depth;
0080 #endif
0081   return index;
0082 }
0083 
0084 bool AHCalSD::unpackIndex(const uint32_t& idx, int& row, int& col, int& depth) {
0085   DetId gen(idx);
0086   HcalSubdetector subdet = (HcalSubdetector(gen.subdetId()));
0087   bool rcode = (gen.det() == DetId::Hcal && subdet != HcalOther);
0088   row = col = depth = 0;
0089   if (rcode) {
0090     row = AHCalDetId(idx).irow();
0091     col = AHCalDetId(idx).icol();
0092     depth = AHCalDetId(idx).depth();
0093   }
0094 #ifdef EDM_ML_DEBUG
0095   edm::LogVerbatim("HcalSim") << "AHCalSD: packed index = 0x" << std::hex << idx << std::dec << " Row " << row
0096                               << " Column " << col << " Depth " << depth << " OK " << rcode;
0097 #endif
0098   return rcode;
0099 }
0100 
0101 bool AHCalSD::filterHit(CaloG4Hit* aHit, double time) {
0102   return ((time <= tmaxHit) && (aHit->getEnergyDeposit() > eminHit));
0103 }