Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

MTD simulation: update Track ID stored in MTD SimHits #41318

Merged
merged 12 commits into from
Apr 14, 2023
Merged
2 changes: 2 additions & 0 deletions Geometry/MTDNumberingBuilder/src/CmsMTDStringToEnum.cc
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,8 @@ const CmsMTDStringToEnum::Impl CmsMTDStringToEnum::m_impl;

CmsMTDStringToEnum::Impl::Impl() {
map_.emplace("FastTimerRegion", GeometricTimingDet::MTD);
map_.emplace("FastTimerRegionBTL", GeometricTimingDet::MTD);
map_.emplace("FastTimerRegionETL", GeometricTimingDet::MTD);
map_.emplace("BarrelTimingLayer", GeometricTimingDet::BTL);
map_.emplace("Layer1", GeometricTimingDet::BTLLayer);
map_.emplace("Layer1Timing", GeometricTimingDet::BTLLayer);
Expand Down
11 changes: 9 additions & 2 deletions Geometry/MTDSimData/data/v2/mtdProdCuts.xml
Original file line number Diff line number Diff line change
Expand Up @@ -3,10 +3,17 @@
<DDDefinition xmlns="http://www.cern.ch/cms/DDL" xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance" xsi:schemaLocation="http://www.cern.ch/cms/DDL ../../../DetectorDescription/Schema/DDLSchema.xsd">

<SpecParSection label="mtdProdCuts.xml" eval="true">
<SpecPar name="fasttime">
<SpecPar name="fasttime-btl">
<PartSelector path="//BarrelTimingLayer"/>
<Parameter name="CMSCutsRegion" value="FastTimerRegionBTL" eval="false"/>
<Parameter name="ProdCutsForElectrons" value="0.1*mm"/>
<Parameter name="ProdCutsForPositrons" value="0.1*mm"/>
<Parameter name="ProdCutsForGamma" value="0.1*mm"/>
<Parameter name="ProdCutsForProtons" value="0.1*mm"/>
</SpecPar>
<SpecPar name="fasttime-etl">
<PartSelector path="//EndcapTimingLayer"/>
<Parameter name="CMSCutsRegion" value="FastTimerRegion" eval="false"/>
<Parameter name="CMSCutsRegion" value="FastTimerRegionETL" eval="false"/>
<Parameter name="ProdCutsForElectrons" value="0.1*mm"/>
<Parameter name="ProdCutsForPositrons" value="0.1*mm"/>
<Parameter name="ProdCutsForGamma" value="0.1*mm"/>
Expand Down
11 changes: 9 additions & 2 deletions Geometry/MTDSimData/data/v4/mtdProdCuts.xml
Original file line number Diff line number Diff line change
Expand Up @@ -3,11 +3,18 @@
<DDDefinition>

<SpecParSection label="mtdProdCuts.xml" eval="true">
<SpecPar name="fasttime">
<SpecPar name="fasttime-btl">
<PartSelector path="//BarrelTimingLayer"/>
<Parameter name="CMSCutsRegion" value="FastTimerRegionBTL" eval="false"/>
<Parameter name="ProdCutsForElectrons" value="0.1*mm"/>
<Parameter name="ProdCutsForPositrons" value="0.1*mm"/>
<Parameter name="ProdCutsForGamma" value="0.1*mm"/>
<Parameter name="ProdCutsForProtons" value="0.1*mm"/>
</SpecPar>
<SpecPar name="fasttime-etl">
<PartSelector path="//EndcapTimingLayer"/>
<PartSelector path="//CALOECTSFront"/>
<Parameter name="CMSCutsRegion" value="FastTimerRegion" eval="false"/>
<Parameter name="CMSCutsRegion" value="FastTimerRegionETL" eval="false"/>
<Parameter name="ProdCutsForElectrons" value="0.1*mm"/>
<Parameter name="ProdCutsForPositrons" value="0.1*mm"/>
<Parameter name="ProdCutsForGamma" value="0.1*mm"/>
Expand Down
6 changes: 6 additions & 0 deletions SimG4CMS/Forward/interface/MtdSD.h
Original file line number Diff line number Diff line change
Expand Up @@ -24,7 +24,13 @@ class MtdSD : public TimingSD {

uint32_t setDetUnitId(const G4Step *) override;

protected:
int getTrackID(const G4Track *) override;

private:
double energyCut;
double energyHistoryCut;

void setNumberingScheme(MTDNumberingScheme *);
void getBaseNumber(const G4Step *);

Expand Down
3 changes: 3 additions & 0 deletions SimG4CMS/Forward/interface/TimingSD.h
Original file line number Diff line number Diff line change
Expand Up @@ -62,6 +62,9 @@ class TimingSD : public SensitiveTkDetector, public Observer<const BeginOfEvent*
const G4ThreeVector& getLocalEntryPoint() const { return hitPointLocal; };
const G4ThreeVector& getGlobalEntryPoint() const { return hitPoint; };

// general method to assign track ID to be stored in hits
virtual int getTrackID(const G4Track*);

private:
void getStepInfo(const G4Step*);
bool hitExists(const G4Step*);
Expand Down
39 changes: 38 additions & 1 deletion SimG4CMS/Forward/src/MtdSD.cc
Original file line number Diff line number Diff line change
@@ -1,3 +1,5 @@
//#define EDM_ML_DEBUG

#include "SimG4CMS/Forward/interface/MtdSD.h"

#include "FWCore/ParameterSet/interface/ParameterSet.h"
Expand All @@ -13,7 +15,6 @@

#include <iostream>

//#define EDM_ML_DEBUG
//-------------------------------------------------------------------
MtdSD::MtdSD(const std::string& name,
const SensitiveDetectorCatalog& clg,
Expand All @@ -40,6 +41,11 @@ MtdSD::MtdSD(const std::string& name,
if (scheme)
setNumberingScheme(scheme);

energyCut = m_p.getParameter<double>("EnergyThresholdForPersistencyInGeV") * CLHEP::GeV; //default must be 0.5
energyHistoryCut = m_p.getParameter<double>("EnergyThresholdForHistoryInGeV") * CLHEP::GeV; //default must be 0.05

setCuts(energyCut, energyHistoryCut);
Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

this statement would treat tracks in BTL and in ETL in the same way, which is probably something we do not want. I suggest we possibly discuss this for a further development.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@fabiocos , yes, in the next rounds I would prefer remove these lines from all SD classes - it is not their job. For now this PR is fine. @srimanob, would you agree?


double newTimeFactor = 1. / m_p.getParameter<double>("TimeSliceUnit");
edm::LogVerbatim("MtdSim") << "New time factor = " << newTimeFactor;
setTimeFactor(newTimeFactor);
Expand Down Expand Up @@ -90,3 +96,34 @@ void MtdSD::getBaseNumber(const G4Step* aStep) {
}
}
}

int MtdSD::getTrackID(const G4Track* aTrack) {
int theID = aTrack->GetTrackID();
TrackInformation* trkInfo = cmsTrackInformation(aTrack);
const G4String& rname = aTrack->GetVolume()->GetLogicalVolume()->GetRegion()->GetName();
if (trkInfo != nullptr) {
#ifdef EDM_ML_DEBUG
trkInfo->Print();
#endif
if (rname == "FastTimerRegionSensBTL") {
if (trkInfo->isBTLdaughter()) {
theID = trkInfo->idAtBTLentrance();
}
#ifdef EDM_ML_DEBUG
edm::LogVerbatim("MtdSim") << "BTL Track ID: " << trkInfo->idAtBTLentrance() << ":" << theID;
#endif
} else if (rname == "FastTimerRegionSensETL") {
theID = trkInfo->getIDonCaloSurface();
#ifdef EDM_ML_DEBUG
edm::LogVerbatim("MtdSim") << "ETL Track ID: " << trkInfo->getIDonCaloSurface() << ":" << theID;
#endif
} else {
throw cms::Exception("MtdSDError") << "MtdSD called in incorrect region " << rname;
}
} else {
#ifdef EDM_ML_DEBUG
edm::LogWarning("MtdSim") << "MtdSD: Problem with primaryID **** set by force to TkID **** " << theID;
#endif
}
return theID;
}
19 changes: 15 additions & 4 deletions SimG4CMS/Forward/src/TimingSD.cc
Original file line number Diff line number Diff line change
@@ -1,3 +1,5 @@
//#define EDM_ML_DEBUG

///////////////////////////////////////////////////////////////////////////////
// File: TimingSD.cc
// Date: 02.2006
Expand Down Expand Up @@ -27,8 +29,6 @@
#include <vector>
#include <iostream>

//#define EDM_ML_DEBUG

static const float invgev = 1.0 / CLHEP::GeV;
static const double invns = 1.0 / CLHEP::nanosecond;
static const double invdeg = 1.0 / CLHEP::deg;
Expand Down Expand Up @@ -161,6 +161,12 @@ void TimingSD::getStepInfo(const G4Step* aStep) {
}
info->putInHistory();
}
#ifdef EDM_ML_DEBUG
if (info != nullptr) {
LogDebug("TimingSim") << "TrackInformation for ID = " << theTrack->GetTrackID();
info->Print();
}
#endif
}

edeposit *= invgev;
Expand All @@ -176,7 +182,7 @@ void TimingSD::getStepInfo(const G4Step* aStep) {
tSliceID = (int)tSlice;

unitID = setDetUnitId(aStep);
primaryID = theTrack->GetTrackID();
primaryID = getTrackID(theTrack);
}

bool TimingSD::hitExists(const G4Step* aStep) {
Expand Down Expand Up @@ -255,7 +261,7 @@ void TimingSD::createNewHit(const G4Step* aStep) {
edm::LogVerbatim("TimingSim") << "TimingSD CreateNewHit for " << GetName() << " PV " << currentPV->GetName()
<< " PVid = " << currentPV->GetCopyNo() << " Unit " << unitID << "\n primary "
<< primaryID << " Tof(ns)= " << tof << " time slice " << tSliceID
<< " E(GeV)= " << incidentEnergy << " trackID " << theTrack->GetTrackID() << " "
<< " E(MeV)= " << incidentEnergy << " trackID " << theTrack->GetTrackID() << " "
<< theTrack->GetDefinition()->GetParticleName() << " parentID "
<< theTrack->GetParentID();

Expand Down Expand Up @@ -365,3 +371,8 @@ void TimingSD::update(const BeginOfEvent* i) {
}

void TimingSD::clearHits() { slave->Initialize(); }

int TimingSD::getTrackID(const G4Track* aTrack) {
LogDebug("TimingSim") << "primary ID: " << aTrack->GetTrackID();
return aTrack->GetTrackID();
}
1 change: 1 addition & 0 deletions SimG4Core/Application/interface/SteppingAction.h
Original file line number Diff line number Diff line change
Expand Up @@ -48,6 +48,7 @@ class SteppingAction : public G4UserSteppingAction {

const G4VPhysicalVolume* tracker{nullptr};
const G4VPhysicalVolume* calo{nullptr};
const G4VPhysicalVolume* btl{nullptr};
const CMSSteppingVerbose* steppingVerbose;
double theCriticalEnergyForVacuum;
double theCriticalDensity;
Expand Down
2 changes: 2 additions & 0 deletions SimG4Core/Application/python/g4SimHits_cfi.py
Original file line number Diff line number Diff line change
Expand Up @@ -518,6 +518,8 @@
IgnoreTrackID = cms.bool(False),
EminHit = cms.double(0.0),
CheckID = cms.untracked.bool(True),
EnergyThresholdForPersistencyInGeV = cms.double(1e9), # keep temporarily old behaviour
EnergyThresholdForHistoryInGeV = cms.double(1e9) # keep temporarily old behaviour)
),
HGCSD = cms.PSet(
Verbosity = cms.untracked.int32(0),
Expand Down
3 changes: 2 additions & 1 deletion SimG4Core/Application/src/StackingAction.cc
Original file line number Diff line number Diff line change
Expand Up @@ -418,7 +418,8 @@ void StackingAction::initPointer() {
if (savePDandCinTracker &&
(rname == "BeamPipe" || rname == "BeamPipeVacuum" || rname == "TrackerPixelSensRegion" ||
rname == "TrackerPixelDeadRegion" || rname == "TrackerDeadRegion" || rname == "TrackerSensRegion" ||
rname == "FastTimerRegion" || rname == "FastTimerRegionSensBTL" || rname == "FastTimerRegionSensETL")) {
rname == "FastTimerRegionBTL" || rname == "FastTimerRegionETL" || rname == "FastTimerRegionSensBTL" ||
rname == "FastTimerRegionSensETL")) {
trackerRegions.push_back(reg);
}
if (savePDandCinCalo && (rname == "HcalRegion" || rname == "EcalRegion" || rname == "PreshowerSensRegion" ||
Expand Down
41 changes: 37 additions & 4 deletions SimG4Core/Application/src/SteppingAction.cc
Original file line number Diff line number Diff line change
@@ -1,3 +1,5 @@
//#define EDM_ML_DEBUG

#include "SimG4Core/Application/interface/SteppingAction.h"

#include "SimG4Core/Notification/interface/TrackInformation.h"
Expand Down Expand Up @@ -158,10 +160,39 @@ void SteppingAction::UserSteppingAction(const G4Step* aStep) {
}
}
}
// check transition tracker/calo
// check transition tracker/btl and tracker/calo
bool isKilled = false;
if (sAlive == tstat || sVeryForward == tstat) {
if (preStep->GetPhysicalVolume() == tracker && postStep->GetPhysicalVolume() == calo) {
// store TrackInformation about transition from one envelope to another
if (preStep->GetPhysicalVolume() == tracker && postStep->GetPhysicalVolume() == btl) {
// store transition tracker -> BTL only for tracks entering BTL for the first time
TrackInformation* trkinfo = static_cast<TrackInformation*>(theTrack->GetUserInformation());
if (!trkinfo->isFromTtoBTL() && !trkinfo->isFromBTLtoT()) {
trkinfo->setFromTtoBTL();
trkinfo->setIdAtBTLentrance(theTrack->GetTrackID());
#ifdef DebugLog
LogDebug("SimG4CoreApplication") << "Setting flag for Tracker -> BTL " << trkinfo->isFromTtoBTL()
<< " IdAtBTLentrance = " << trkinfo->idAtBTLentrance();
#endif
} else {
trkinfo->setBTLlooper();
trkinfo->setIdAtBTLentrance(theTrack->GetTrackID());
#ifdef DebugLog
LogDebug("SimG4CoreApplication") << "Setting flag for BTL looper " << trkinfo->isBTLlooper();
trkinfo->Print();
#endif
}
} else if (preStep->GetPhysicalVolume() == btl && postStep->GetPhysicalVolume() == tracker) {
// store transition BTL -> tracker
TrackInformation* trkinfo = static_cast<TrackInformation*>(theTrack->GetUserInformation());
if (!trkinfo->isFromBTLtoT()) {
trkinfo->setFromBTLtoT();
#ifdef DebugLog
LogDebug("SimG4CoreApplication") << "Setting flag for BTL -> Tracker " << trkinfo->isFromBTLtoT();
#endif
}
} else if (preStep->GetPhysicalVolume() == tracker && postStep->GetPhysicalVolume() == calo) {
// store transition tracker -> calo
TrackInformation* trkinfo = static_cast<TrackInformation*>(theTrack->GetUserInformation());
if (!trkinfo->crossedBoundary()) {
trkinfo->setCrossedBoundary(theTrack);
Expand Down Expand Up @@ -204,12 +235,14 @@ bool SteppingAction::initPointer() {
tracker = pvcite;
} else if (pvname == "CALO" || pvname == "caloBase:CALO_1") {
calo = pvcite;
} else if (pvname == "BarrelTimingLayer" || pvname == "btl:BarrelTimingLayer_1") {
btl = pvcite;
}
if (tracker && calo)
if (tracker && calo && btl)
break;
}
edm::LogVerbatim("SimG4CoreApplication")
<< "SteppingAction: pointer for Tracker " << tracker << " and for Calo " << calo;
<< "SteppingAction: pointer for Tracker " << tracker << " and for Calo " << calo << " and for BTL " << btl;

const G4LogicalVolumeStore* lvs = G4LogicalVolumeStore::GetInstance();
if (numberEkins > 0) {
Expand Down
4 changes: 2 additions & 2 deletions SimG4Core/Application/src/TrackingAction.cc
Original file line number Diff line number Diff line change
@@ -1,3 +1,5 @@
//#define EDM_ML_DEBUG

#include "SimG4Core/Application/interface/TrackingAction.h"

#include "SimG4Core/Notification/interface/NewTrackAction.h"
Expand All @@ -15,8 +17,6 @@
#include "G4TrackingManager.hh"
#include "G4SystemOfUnits.hh"

//#define EDM_ML_DEBUG

TrackingAction::TrackingAction(SimTrackManager* stm, CMSSteppingVerbose* sv, const edm::ParameterSet& p)
: trackManager_(stm),
steppingVerbose_(sv),
Expand Down
2 changes: 2 additions & 0 deletions SimG4Core/Notification/interface/NewTrackAction.h
Original file line number Diff line number Diff line change
Expand Up @@ -21,6 +21,8 @@ class NewTrackAction {
private:
void addUserInfoToPrimary(G4Track* aTrack) const;
void addUserInfoToSecondary(G4Track* aTrack, const TrackInformation& motherInfo, int) const;

bool isInBTL(const G4Track*) const;
};

#endif
21 changes: 20 additions & 1 deletion SimG4Core/Notification/interface/TrackInformation.h
Original file line number Diff line number Diff line change
Expand Up @@ -80,6 +80,20 @@ class TrackInformation : public G4VUserTrackInformation {
}
int getCastorHitPID() const { return castorHitPID_; }

// methods for MTD info management
//
inline void setFromTtoBTL() { mtdStatus_ |= 1 << 0; } // 1st bit
inline bool isFromTtoBTL() const { return (mtdStatus_ >> 0) & 1; }
inline void setFromBTLtoT() { mtdStatus_ |= 1 << 1; } // 2nd bit
inline bool isFromBTLtoT() const { return (mtdStatus_ >> 1) & 1; }
inline void setBTLdaughter() { mtdStatus_ |= 1 << 2; } // 3rd bit
inline bool isBTLdaughter() const { return (mtdStatus_ >> 2) & 1; }
inline void setBTLlooper() { mtdStatus_ |= 1 << 3; } // 4th bit
inline bool isBTLlooper() const { return (mtdStatus_ >> 3) & 1; }

int idAtBTLentrance() const { return idAtBTLentrance_; }
void setIdAtBTLentrance(int id) { idAtBTLentrance_ = id; }

void Print() const override;

private:
Expand All @@ -105,6 +119,9 @@ class TrackInformation : public G4VUserTrackInformation {
bool hasCastorHit_;
int castorHitPID_;

uint8_t mtdStatus_;
int idAtBTLentrance_;

// Restrict construction to friends
TrackInformation()
: G4VUserTrackInformation(),
Expand All @@ -126,7 +143,9 @@ class TrackInformation : public G4VUserTrackInformation {
genParticleP_(0),
caloSurfaceParticleP_(0),
hasCastorHit_(false),
castorHitPID_(0) {}
castorHitPID_(0),
mtdStatus_(0),
idAtBTLentrance_(0) {}
friend class NewTrackAction;
};

Expand Down
24 changes: 24 additions & 0 deletions SimG4Core/Notification/src/NewTrackAction.cc
Original file line number Diff line number Diff line change
@@ -1,3 +1,4 @@
//#define EDM_ML_DEBUG

#include "SimG4Core/Notification/interface/NewTrackAction.h"
#include "SimG4Core/Notification/interface/TrackInformation.h"
Expand Down Expand Up @@ -51,18 +52,41 @@ void NewTrackAction::addUserInfoToSecondary(G4Track *aTrack, const TrackInformat
motherInfo.getIDLastVolume(),
aTrack->GetDefinition()->GetPDGEncoding(),
aTrack->GetMomentum().mag());
LogDebug("SimG4CoreApplication") << "NewTrackAction: Id on calo surface " << trkInfo->getIDonCaloSurface()
<< " to be stored " << trkInfo->storeTrack();
} else {
// transfer calo ID from mother (to be checked in TrackingAction)
trkInfo->setIDonCaloSurface(motherInfo.getIDonCaloSurface(),
motherInfo.getIDCaloVolume(),
motherInfo.getIDLastVolume(),
motherInfo.caloSurfaceParticlePID(),
motherInfo.caloSurfaceParticleP());
LogDebug("SimG4CoreApplication") << "NewTrackAction: Id on calo surface " << trkInfo->getIDonCaloSurface();
}

if (motherInfo.hasCastorHit()) {
trkInfo->setCastorHitPID(motherInfo.getCastorHitPID());
}

// manage ID of tracks in BTL to map them to SimTracks to be stored
if (isInBTL(aTrack)) {
if ((motherInfo.storeTrack() && motherInfo.isFromTtoBTL()) || motherInfo.isBTLdaughter()) {
trkInfo->setBTLdaughter();
trkInfo->setIdAtBTLentrance(motherInfo.idAtBTLentrance());
LogDebug("SimG4CoreApplication") << "NewTrackAction: secondary in BTL " << trkInfo->isBTLdaughter()
<< " from mother ID " << trkInfo->idAtBTLentrance();
}
}

aTrack->SetUserInformation(trkInfo);
}

bool NewTrackAction::isInBTL(const G4Track *aTrack) const {
bool out = false;
G4String tName(aTrack->GetVolume()->GetLogicalVolume()->GetRegion()->GetName());
if (tName == "FastTimerRegionBTL" || tName == "FastTimerRegionSensBTL") {
out = true;
}

return out;
}
Loading