PPS/src/TotemRPSD.cc

0001 // File: TotemRPSD.cc
0002 // Date: 18.10.2005
0003 // Description: Sensitive Detector class for TOTEM RP Detectors
0004 // Modifications:
0005 #include "SimG4CMS/PPS/interface/TotemRPSD.h"
0006 #include "SimG4CMS/PPS/interface/PPSStripOrganization.h"
0007
0008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0009 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0010
0011 #include "SimG4Core/Notification/interface/TrackInformation.h"
0012 #include "SimG4Core/Notification/interface/G4TrackToParticleID.h"
0013 #include "SimG4Core/Geometry/interface/SensitiveDetectorCatalog.h"
0014 #include "SimG4Core/Physics/interface/G4ProcessTypeEnumerator.h"
0015
0016 #include "SimG4Core/SensitiveDetector/interface/SensitiveTkDetector.h"
0017
0018 #include "SimDataFormats/TrackingHit/interface/UpdatablePSimHit.h"
0019 #include "SimDataFormats/SimHitMaker/interface/TrackingSlaveSD.h"
0020
0021 #include "G4SDManager.hh"
0022 #include "G4Step.hh"
0023 #include "G4Track.hh"
0024 #include "G4VProcess.hh"
0025
0026 #include "G4PhysicalConstants.hh"
0027 #include <CLHEP/Units/SystemOfUnits.h>
0028
0029 #include <iostream>
0030 #include <vector>
0031 #include <string>
0032
0033 static constexpr double rp_garage_position_ = 40.0;
0034
0035 TotemRPSD::TotemRPSD(const std::string& pname,
0036                      const SensitiveDetectorCatalog& clg,
0037                      edm::ParameterSet const& p,
0038                      const SimTrackManager* manager)
0039     : SensitiveTkDetector(pname, clg) {
0040   collectionName.insert(pname);
0041
0042   edm::ParameterSet m_Anal = p.getParameter<edm::ParameterSet>("TotemRPSD");
0043   verbosity_ = m_Anal.getParameter<int>("Verbosity");
0044
0045   slave_ = std::make_unique<TrackingSlaveSD>(pname);
0046
0047   if (pname == "TotemHitsRP") {
0048     numberingScheme_ = std::make_unique<PPSStripOrganization>();
0049   } else {
0050     edm::LogWarning("TotemRP") << "TotemRPSD: ReadoutName " << pname << " not supported";
0051   }
0052
0053   edm::LogVerbatim("TotemRP") << "TotemRPSD: Instantiation completed for " << pname;
0054 }
0055
0056 TotemRPSD::~TotemRPSD() {}
0057
0058 void TotemRPSD::Initialize(G4HCofThisEvent* HCE) {
0059   LogDebug("TotemRP") << "TotemRPSD : Initialize called for " << GetName();
0060
0061   theHC_ = new TotemRPG4HitCollection(GetName(), collectionName[0]);
0062   G4SDManager::GetSDMpointer()->AddNewCollection(GetName(), collectionName[0]);
0063
0064   if (hcID_ < 0)
0065     hcID_ = G4SDManager::GetSDMpointer()->GetCollectionID(collectionName[0]);
0066   HCE->AddHitsCollection(hcID_, theHC_);
0067   LogDebug("TotemRP") << "TotemRPSD: is initialized for " << GetName();
0068 }
0069
0070 void TotemRPSD::printHitInfo() {
0071   LogDebug("TotemRP") << theTrack_->GetDefinition()->GetParticleName() << " TotemRPSD CreateNewHit for"
0072                       << " PV " << currentPV_->GetName() << " PVid = " << currentPV_->GetCopyNo() << " Unit "
0073                       << unitID_;
0074   LogDebug("TotemRP") << " primary " << primaryID_ << " time slice " << tSliceID_ << " of energy "
0075                       << theTrack_->GetTotalEnergy() << " Eloss_ " << eloss_ << " position pre: " << hitPoint_
0076                       << " post: " << exitPoint_;
0077
0078   if (parentID_ == 0) {
0079     LogDebug("TotemRP") << " primary particle ";
0080   } else {
0081     LogDebug("TotemRP") << " daughter of part. " << parentID_;
0082   }
0083
0084   LogDebug("TotemRP") << " and created by ";
0085
0086   if (theTrack_->GetCreatorProcess() != nullptr)
0087     LogDebug("TotemRP") << theTrack_->GetCreatorProcess()->GetProcessName();
0088   else
0089     LogDebug("TotemRP") << "NO process";
0090 }
0091
0092 bool TotemRPSD::ProcessHits(G4Step* aStep, G4TouchableHistory*) {
0093   eloss_ = aStep->GetTotalEnergyDeposit();
0094   if (eloss_ > 0.0) {
0095     eloss_ /= CLHEP::GeV;
0096     stepInfo(aStep);
0097     edm::LogVerbatim("TotemRP") << "TotemRPSD: ProcessHits 1: Eloss=" << eloss_ << " "
0098                                 << theTrack_->GetDefinition()->GetParticleName();
0099     createNewHit();
0100   }
0101   return true;
0102 }
0103
0104 void TotemRPSD::stepInfo(const G4Step* aStep) {
0105   preStepPoint_ = aStep->GetPreStepPoint();
0106   postStepPoint_ = aStep->GetPostStepPoint();
0107   theTrack_ = aStep->GetTrack();
0108   hitPoint_ = preStepPoint_->GetPosition();
0109   exitPoint_ = postStepPoint_->GetPosition();
0110   currentPV_ = preStepPoint_->GetPhysicalVolume();
0111   theLocalEntryPoint_ = setToLocal(hitPoint_);
0112   theLocalExitPoint_ = setToLocal(exitPoint_);
0113
0114   tSlice_ = (postStepPoint_->GetGlobalTime()) / CLHEP::nanosecond;
0115   tSliceID_ = (int)tSlice_;
0116   unitID_ = setDetUnitId(aStep);
0117
0118   if (verbosity_)
0119     LogDebug("TotemRP") << "UNIT " << unitID_;
0120
0121   primaryID_ = theTrack_->GetTrackID();
0122   parentID_ = theTrack_->GetParentID();
0123
0124   incidentEnergy_ = theTrack_->GetTotalEnergy() / CLHEP::GeV;
0125
0126   pabs_ = preStepPoint_->GetMomentum().mag() / CLHEP::GeV;
0127   thePx_ = preStepPoint_->GetMomentum().x() / CLHEP::GeV;
0128   thePy_ = preStepPoint_->GetMomentum().y() / CLHEP::GeV;
0129   thePz_ = preStepPoint_->GetMomentum().z() / CLHEP::GeV;
0130
0131   tof_ = postStepPoint_->GetGlobalTime() / CLHEP::nanosecond;
0132   particleType_ = theTrack_->GetDefinition()->GetPDGEncoding();
0133
0134   //corrected phi and theta treatment
0135   G4ThreeVector gmd = preStepPoint_->GetMomentumDirection();
0136   // convert it to local frame
0137   G4ThreeVector lmd =
0138       ((G4TouchableHistory*)(preStepPoint_->GetTouchable()))->GetHistory()->GetTopTransform().TransformAxis(gmd);
0139   Local3DPoint lnmd = ConvertToLocal3DPoint(lmd);
0140   thetaAtEntry_ = lnmd.theta();
0141   phiAtEntry_ = lnmd.phi();
0142
0143   vx_ = theTrack_->GetVertexPosition().x() / CLHEP::mm;
0144   vy_ = theTrack_->GetVertexPosition().y() / CLHEP::mm;
0145   vz_ = theTrack_->GetVertexPosition().z() / CLHEP::mm;
0146 }
0147
0148 uint32_t TotemRPSD::setDetUnitId(const G4Step* aStep) {
0149   return (numberingScheme_ == nullptr ? 0 : numberingScheme_->unitID(aStep));
0150 }
0151
0152 void TotemRPSD::storeHit(TotemRPG4Hit* hit) {
0153   if (hit == nullptr) {
0154     if (verbosity_)
0155       LogDebug("TotemRP") << "TotemRPSD: hit to be stored is NULL !!" << std::endl;
0156   } else {
0157     theHC_->insert(hit);
0158   }
0159 }
0160
0161 void TotemRPSD::createNewHit() {
0162   // Protect against creating hits in detectors not inserted
0163   double outrangeX = hitPoint_.x();
0164   double outrangeY = hitPoint_.y();
0165   if (std::abs(outrangeX) > rp_garage_position_ || std::abs(outrangeY) > rp_garage_position_)
0166     return;
0167   // end protection
0168
0169   currentHit_ = new TotemRPG4Hit;
0170   currentHit_->setTrackID(primaryID_);
0171   currentHit_->setTimeSlice(tSlice_);
0172   currentHit_->setUnitID(unitID_);
0173   currentHit_->setIncidentEnergy(incidentEnergy_);
0174
0175   currentHit_->setP(pabs_);
0176   currentHit_->setTof(tof_);
0177   currentHit_->setEnergyLoss(eloss_);
0178   currentHit_->setParticleType(particleType_);
0179   currentHit_->setThetaAtEntry(thetaAtEntry_);
0180   currentHit_->setPhiAtEntry(phiAtEntry_);
0181
0182   currentHit_->setEntry(hitPoint_);
0183   currentHit_->setExit(exitPoint_);
0184   currentHit_->setLocalEntry(theLocalEntryPoint_);
0185   currentHit_->setLocalExit(theLocalExitPoint_);
0186
0187   currentHit_->setParentId(parentID_);
0188   currentHit_->setVx(vx_);
0189   currentHit_->setVy(vy_);
0190   currentHit_->setVz(vz_);
0191
0192   currentHit_->setPx(thePx_);
0193   currentHit_->setPy(thePy_);
0194   currentHit_->setPz(thePz_);
0195
0196   storeHit(currentHit_);
0197 }
0198
0199 G4ThreeVector TotemRPSD::setToLocal(const G4ThreeVector& global) {
0200   G4ThreeVector localPoint;
0201   const G4VTouchable* touch = preStepPoint_->GetTouchable();
0202   localPoint = touch->GetHistory()->GetTopTransform().TransformPoint(global);
0203   return localPoint;
0204 }
0205
0206 void TotemRPSD::EndOfEvent(G4HCofThisEvent*) {
0207   // here we loop over transient hits and make them persistent
0208   for (unsigned int j = 0; j < (unsigned int)theHC_->entries(); ++j) {
0209     TotemRPG4Hit* aHit = (*theHC_)[j];
0210
0211     Local3DPoint entry(aHit->localEntry().x(), aHit->localEntry().y(), aHit->localEntry().z());
0212     Local3DPoint exit(aHit->localExit().x(), aHit->localExit().y(), aHit->localExit().z());
0213     slave_->processHits(PSimHit(entry,
0214                                 exit,
0215                                 aHit->p(),
0216                                 aHit->tof(),
0217                                 aHit->energyLoss(),
0218                                 aHit->particleType(),
0219                                 aHit->unitID(),
0220                                 aHit->trackID(),
0221                                 aHit->thetaAtEntry(),
0222                                 aHit->phiAtEntry()));
0223   }
0224 }
0225
0226 void TotemRPSD::PrintAll() {
0227   LogDebug("TotemRP") << "TotemRPSD: Collection " << theHC_->GetName() << std::endl;
0228   theHC_->PrintAllHits();
0229 }
0230
0231 void TotemRPSD::fillHits(edm::PSimHitContainer& c, const std::string& n) {
0232   if (slave_->name() == n) {
0233     c = slave_->hits();
0234   }
0235 }
0236
0237 void TotemRPSD::update(const BeginOfEvent* ptr) {
0238   clearHits();
0239   eventno_ = (*ptr)()->GetEventID();
0240 }
0241
0242 void TotemRPSD::update(const ::EndOfEvent*) {}
0243
0244 void TotemRPSD::clearHits() { slave_->Initialize(); }