Application/src/Phase2TrackingAction.cc

0001 #include "SimG4Core/Application/interface/Phase2TrackingAction.h"
0002 #include "SimG4Core/Physics/interface/CMSG4TrackInterface.h"
0003
0004 #include "SimG4Core/Notification/interface/BeginOfTrack.h"
0005 #include "SimG4Core/Notification/interface/EndOfTrack.h"
0006 #include "SimG4Core/Notification/interface/TrackInformation.h"
0007 #include "SimG4Core/Notification/interface/TrackWithHistory.h"
0008 #include "SimG4Core/Notification/interface/SimTrackManager.h"
0009 #include "SimG4Core/Notification/interface/CMSSteppingVerbose.h"
0010
0011 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0012
0013 #include "G4UImanager.hh"
0014 #include "G4TrackingManager.hh"
0015 #include <CLHEP/Units/SystemOfUnits.h>
0016
0017 //#define EDM_ML_DEBUG
0018
0019 Phase2TrackingAction::Phase2TrackingAction(SimTrackManager* stm, CMSSteppingVerbose* sv, const edm::ParameterSet& p)
0020     : trackManager_(stm),
0021       steppingVerbose_(sv),
0022       endPrintTrackID_(p.getParameter<int>("EndPrintTrackID")),
0023       checkTrack_(p.getUntrackedParameter<bool>("CheckTrack", false)),
0024       doFineCalo_(p.getParameter<bool>("DoFineCalo")),
0025       saveCaloBoundaryInformation_(p.getParameter<bool>("SaveCaloBoundaryInformation")),
0026       ekinMin_(p.getParameter<double>("PersistencyEmin") * CLHEP::GeV),
0027       ekinMinRegion_(p.getParameter<std::vector<double>>("RegionEmin")) {
0028   trackInterface_ = CMSG4TrackInterface::instance();
0029   double eth = p.getParameter<double>("EminFineTrack") * CLHEP::MeV;
0030   if (doFineCalo_ && eth < ekinMin_) {
0031     ekinMin_ = eth;
0032   }
0033   edm::LogVerbatim("SimG4CoreApplication") << "Phase2TrackingAction: boundary: " << saveCaloBoundaryInformation_
0034                                            << "; DoFineCalo: " << doFineCalo_ << "; ekinMin(MeV)=" << ekinMin_;
0035   if (!ekinMinRegion_.empty()) {
0036     ptrRegion_.resize(ekinMinRegion_.size(), nullptr);
0037   }
0038 }
0039
0040 void Phase2TrackingAction::PreUserTrackingAction(const G4Track* aTrack) {
0041   g4Track_ = aTrack;
0042   currentTrack_ = new TrackWithHistory(aTrack, aTrack->GetParentID());
0043   trackInterface_->setCurrentTrack(aTrack);
0044
0045   BeginOfTrack bt(aTrack);
0046   m_beginOfTrackSignal(&bt);
0047
0048   trkInfo_ = static_cast<TrackInformation*>(aTrack->GetUserInformation());
0049
0050   // Always save primaries
0051   if (trkInfo_->isPrimary()) {
0052     trackManager_->cleanTracksWithHistory();
0053     currentTrack_->setToBeSaved();
0054   }
0055
0056   if (nullptr != steppingVerbose_) {
0057     steppingVerbose_->trackStarted(aTrack, false);
0058     if (aTrack->GetTrackID() == endPrintTrackID_) {
0059       steppingVerbose_->stopEventPrint();
0060     }
0061   }
0062   double ekin = aTrack->GetKineticEnergy();
0063
0064 #ifdef EDM_ML_DEBUG
0065   edm::LogVerbatim("DoFineCalo") << "PreUserPhase2TrackingAction: Start processing track " << aTrack->GetTrackID()
0066                                  << " pdgid=" << aTrack->GetDefinition()->GetPDGEncoding()
0067                                  << " ekin[GeV]=" << ekin / CLHEP::GeV << " vertex[cm]=("
0068                                  << aTrack->GetVertexPosition().x() / CLHEP::cm << ","
0069                                  << aTrack->GetVertexPosition().y() / CLHEP::cm << ","
0070                                  << aTrack->GetVertexPosition().z() / CLHEP::cm << ")"
0071                                  << " parentid=" << aTrack->GetParentID();
0072 #endif
0073   if (ekin > ekinMin_) {
0074     // Each track with energy above the threshold should be saved
0075     trkInfo_->putInHistory();
0076   }
0077 }
0078
0079 void Phase2TrackingAction::PostUserTrackingAction(const G4Track* aTrack) {
0080   // Tracks in history may be upgraded to stored secondary tracks,
0081   // which cross the boundary between Tracker and Calo
0082   int id = aTrack->GetTrackID();
0083   bool ok = (trkInfo_->storeTrack() || currentTrack_->saved());
0084   if (trkInfo_->crossedBoundary()) {
0085     currentTrack_->setCrossedBoundaryPosMom(id, trkInfo_->getPositionAtBoundary(), trkInfo_->getMomentumAtBoundary());
0086     ok = (ok || saveCaloBoundaryInformation_ || doFineCalo_);
0087   }
0088   if (ok) {
0089     currentTrack_->setToBeSaved();
0090   }
0091
0092   bool withAncestor = (trkInfo_->getIDonCaloSurface() == id || trkInfo_->isAncestor());
0093   bool isInHistory = trkInfo_->isInHistory();
0094
0095   trackManager_->addTrack(currentTrack_, aTrack, isInHistory, withAncestor);
0096
0097 #ifdef EDM_ML_DEBUG
0098   edm::LogVerbatim("Phase2TrackingAction")
0099       << "Phase2TrackingAction end track=" << id << "  " << aTrack->GetDefinition()->GetParticleName()
0100       << " proposed to be saved= " << ok << " end point " << aTrack->GetPosition();
0101 #endif
0102
0103   if (!isInHistory) {
0104     delete currentTrack_;
0105   }
0106
0107   EndOfTrack et(aTrack);
0108   m_endOfTrackSignal(&et);
0109 }