File indexing completed on 2024-05-10 02:21:21
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0017 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0018
0019 #include "SimG4Core/Notification/interface/TrackInformation.h"
0020 #include "SimG4Core/Notification/interface/G4TrackToParticleID.h"
0021 #include "SimG4Core/Physics/interface/G4ProcessTypeEnumerator.h"
0022
0023 #include "SimDataFormats/TrackingHit/interface/UpdatablePSimHit.h"
0024 #include "SimDataFormats/SimHitMaker/interface/TrackingSlaveSD.h"
0025
0026 #include "SimG4CMS/PPS/interface/PPSPixelSD.h"
0027 #include "SimG4CMS/PPS/interface/PPSPixelOrganization.h"
0028
0029 #include "G4SDManager.hh"
0030 #include "G4Step.hh"
0031 #include "G4Track.hh"
0032 #include "G4VProcess.hh"
0033
0034 #include "G4PhysicalConstants.hh"
0035 #include <CLHEP/Units/SystemOfUnits.h>
0036
0037 PPSPixelSD::PPSPixelSD(const std::string& pname,
0038 const SensitiveDetectorCatalog& clg,
0039 edm::ParameterSet const& p,
0040 SimTrackManager const* manager)
0041 : SensitiveTkDetector(pname, clg), theManager_(manager) {
0042
0043 collectionName.insert(pname);
0044
0045
0046 edm::ParameterSet m_p = p.getParameter<edm::ParameterSet>("PPSPixelSD");
0047 int verbn = m_p.getUntrackedParameter<int>("Verbosity");
0048 SetVerboseLevel(verbn);
0049 slave_ = std::make_unique<TrackingSlaveSD>(pname);
0050 if (pname == "CTPPSPixelHits") {
0051 numberingScheme_ = std::make_unique<PPSPixelOrganization>();
0052 } else {
0053 edm::LogWarning("PPSSim") << "PPSPixelSD: ReadoutName " << pname << " not supported";
0054 }
0055 edm::LogVerbatim("PPSSim") << "PPSPixelSD: Instantiation completed for " << pname;
0056 }
0057
0058 PPSPixelSD::~PPSPixelSD() {}
0059
0060 bool PPSPixelSD::ProcessHits(G4Step* aStep, G4TouchableHistory*) {
0061 eloss_ = aStep->GetTotalEnergyDeposit();
0062 if (eloss_ > 0.0) {
0063 eloss_ /= CLHEP::GeV;
0064 stepInfo(aStep);
0065 LogDebug("PPSSim") << "PPSPixelSD: ProcessHits: edep=" << eloss_ << " "
0066 << theTrack_->GetDefinition()->GetParticleName();
0067 if (!hitExists())
0068 createNewHit();
0069 }
0070 return true;
0071 }
0072
0073 uint32_t PPSPixelSD::setDetUnitId(const G4Step* aStep) {
0074 return (numberingScheme_ == nullptr ? 0 : numberingScheme_->unitID(aStep));
0075 }
0076
0077 void PPSPixelSD::Initialize(G4HCofThisEvent* HCE) {
0078 LogDebug("PPSSim") << "PPSPixelSD : Initialize called for " << GetName();
0079
0080 theHC_ = new PPSPixelG4HitCollection(GetName(), collectionName[0]);
0081 G4SDManager::GetSDMpointer()->AddNewCollection(GetName(), collectionName[0]);
0082
0083 if (hcID_ < 0)
0084 hcID_ = G4SDManager::GetSDMpointer()->GetCollectionID(collectionName[0]);
0085 HCE->AddHitsCollection(hcID_, theHC_);
0086
0087 tsID_ = -2;
0088 edm::LogVerbatim("PPSSim") << "PPSPixelSD: is initialized " << GetName();
0089 }
0090
0091 void PPSPixelSD::EndOfEvent(G4HCofThisEvent*) {
0092
0093 for (unsigned int j = 0; j < (unsigned int)theHC_->entries(); ++j) {
0094 PPSPixelG4Hit* aHit = (*theHC_)[j];
0095 #ifdef debug
0096 LogDebug("PPSSim") << "HIT NUMERO " << j << "unit ID = " << aHit->unitID() << "\n"
0097 << " "
0098 << "enrty z " << aHit->entry().z() << "\n"
0099 << " "
0100 << "theta " << aHit->thetaAtEntry() << "\n";
0101 #endif
0102
0103 Local3DPoint Enter(aHit->entryPoint().x(), aHit->entryPoint().y(), aHit->entryPoint().z());
0104 Local3DPoint Exit(aHit->exitPoint().x(), aHit->exitPoint().y(), aHit->exitPoint().z());
0105 slave_->processHits(PSimHit(Enter,
0106 Exit,
0107 aHit->p(),
0108 aHit->tof(),
0109 aHit->energyLoss(),
0110 aHit->particleType(),
0111 aHit->unitID(),
0112 aHit->trackID(),
0113 aHit->thetaAtEntry(),
0114 aHit->phiAtEntry()));
0115 }
0116 }
0117
0118 void PPSPixelSD::PrintAll() {
0119 LogDebug("PPSSim") << "PPSPixelSD: Collection " << theHC_->GetName();
0120 theHC_->PrintAllHits();
0121 }
0122
0123 void PPSPixelSD::fillHits(edm::PSimHitContainer& c, const std::string& n) {
0124 if (slave_->name() == n) {
0125 c = slave_->hits();
0126 }
0127 }
0128
0129 void PPSPixelSD::update(const BeginOfEvent* i) {
0130 LogDebug("PPSSim") << "PPSPixelSD: dispatched BeginOfEvent for " << GetName();
0131 clearHits();
0132 eventno_ = (*i)()->GetEventID();
0133 }
0134
0135 void PPSPixelSD::update(const ::EndOfEvent*) {}
0136
0137 void PPSPixelSD::clearHits() { slave_->Initialize(); }
0138
0139 G4ThreeVector PPSPixelSD::setToLocal(const G4ThreeVector& global) {
0140 G4ThreeVector localPoint;
0141 const G4VTouchable* touch = preStepPoint_->GetTouchable();
0142 localPoint = touch->GetHistory()->GetTopTransform().TransformPoint(global);
0143 return localPoint;
0144 }
0145
0146 void PPSPixelSD::stepInfo(const G4Step* aStep) {
0147 preStepPoint_ = aStep->GetPreStepPoint();
0148 postStepPoint_ = aStep->GetPostStepPoint();
0149 theTrack_ = aStep->GetTrack();
0150
0151 hitPoint_ = preStepPoint_->GetPosition();
0152 exitPoint_ = postStepPoint_->GetPosition();
0153 currentPV_ = preStepPoint_->GetPhysicalVolume();
0154 theLocalEntryPoint_ = setToLocal(hitPoint_);
0155 theLocalExitPoint_ = setToLocal(exitPoint_);
0156
0157 tSlice_ = postStepPoint_->GetGlobalTime() / CLHEP::nanosecond;
0158 tSliceID_ = (int)tSlice_;
0159 unitID_ = setDetUnitId(aStep);
0160 #ifdef debug
0161 LogDebug("PPSSim") << "UNIT ID " << unitID_;
0162 #endif
0163 primaryID_ = theTrack_->GetTrackID();
0164 parentID_ = theTrack_->GetParentID();
0165
0166 incidentEnergy_ = theTrack_->GetTotalEnergy() / CLHEP::GeV;
0167
0168 pabs_ = aStep->GetPreStepPoint()->GetMomentum().mag() / CLHEP::GeV;
0169 tof_ = aStep->GetPostStepPoint()->GetGlobalTime() / CLHEP::nanosecond;
0170
0171 particleType_ = theTrack_->GetDefinition()->GetPDGEncoding();
0172
0173 thetaAtEntry_ = aStep->GetPreStepPoint()->GetPosition().theta();
0174 phiAtEntry_ = aStep->GetPreStepPoint()->GetPosition().phi();
0175
0176 vx_ = theTrack_->GetVertexPosition().x();
0177 vy_ = theTrack_->GetVertexPosition().y();
0178 vz_ = theTrack_->GetVertexPosition().z();
0179 }
0180
0181 bool PPSPixelSD::hitExists() {
0182
0183 if (tSliceID_ == tsID_ && unitID_ == previousUnitID_) {
0184 updateHit();
0185 return true;
0186 }
0187
0188
0189
0190
0191 int nhits = theHC_->entries();
0192 for (int j = 0; j < nhits; ++j) {
0193 PPSPixelG4Hit* aPreviousHit = (*theHC_)[j];
0194 if (aPreviousHit->timeSliceID() == tSliceID_ && aPreviousHit->unitID() == unitID_) {
0195 currentHit_ = aPreviousHit;
0196 updateHit();
0197 return true;
0198 }
0199 }
0200 return false;
0201 }
0202
0203 void PPSPixelSD::createNewHit() {
0204 #ifdef debug
0205 LogDebug("PPSSim") << "PPSPixelSD::CreateNewHit for"
0206 << " PV " << currentPV_->GetName() << " PVid = " << currentPV_->GetCopyNo()
0207 << " MVid = " << currentPV_->GetMother()->GetCopyNo() << " Unit " << unitID_ << "\n"
0208 << " primary " << primaryID_ << " time slice " << tSliceID_ << " For Track "
0209 << " which is a " << theTrack_->GetDefinition()->GetParticleName();
0210
0211 if (parentID == 0) {
0212 LogDebug("PPSSim") << "PPSPixelSD: primary of energy " << incidentEnergy_;
0213 } else {
0214 LogDebug("PPSSim") << " daughter of track: " << parentID;
0215 }
0216
0217 if (theTrack_->GetCreatorProcess() != nullptr)
0218 LogDebug("PPSSim") << theTrack_->GetCreatorProcess()->GetProcessName();
0219 else
0220 LogDebug("PPSSim") << "NO process";
0221 #endif
0222
0223 currentHit_ = new PPSPixelG4Hit;
0224 currentHit_->setTrackID(primaryID_);
0225 currentHit_->setTimeSlice(tSlice_);
0226 currentHit_->setUnitID(unitID_);
0227 currentHit_->setIncidentEnergy(incidentEnergy_);
0228
0229 currentHit_->setP(pabs_);
0230 currentHit_->setTof(tof_);
0231 currentHit_->setParticleType(particleType_);
0232 currentHit_->setThetaAtEntry(thetaAtEntry_);
0233 currentHit_->setPhiAtEntry(phiAtEntry_);
0234
0235 G4ThreeVector pos = 0.5 * (hitPoint_ + exitPoint_);
0236 currentHit_->setMeanPosition(pos);
0237 currentHit_->setEntryPoint(theLocalEntryPoint_);
0238 currentHit_->setExitPoint(theLocalExitPoint_);
0239
0240 currentHit_->setParentId(parentID_);
0241 currentHit_->setVx(vx_);
0242 currentHit_->setVy(vy_);
0243 currentHit_->setVz(vz_);
0244
0245 updateHit();
0246 storeHit(currentHit_);
0247 }
0248
0249 void PPSPixelSD::updateHit() {
0250 #ifdef debug
0251 LogDebug("PPSSim") << "G4PPSPixelSD updateHit: add eloss=" << eloss_ << "\nCurrentHit=" << currentHit_
0252 << ", PostStepPoint=" << postStepPoint_->GetPosition();
0253 #endif
0254 currentHit_->setEnergyLoss(eloss_);
0255
0256 tsID_ = tSliceID_;
0257 previousUnitID_ = unitID_;
0258 }
0259
0260 void PPSPixelSD::storeHit(PPSPixelG4Hit* hit) {
0261 if (hit == nullptr) {
0262 edm::LogWarning("PPSSim") << "PPSPixelSD: hit to be stored is NULL !!";
0263 return;
0264 }
0265
0266 theHC_->insert(hit);
0267 }