#include <SteppingAction.h>

Inheritance diagram for SteppingAction:

Public Member Functions
	SteppingAction (EventAction *ea, const edm::ParameterSet &ps)

void	UserSteppingAction (const G4Step *aStep)

	~SteppingAction ()

Public Attributes
SimActivityRegistry::G4StepSignal	m_g4StepSignal

Private Member Functions
bool	catchLongLived (const G4Step *aStep)

bool	catchLowEnergyInVacuum (G4Track *theTrack, double theKenergy)

bool	initPointer ()

bool	isThisVolume (const G4VTouchable touch, G4VPhysicalVolume pv)

bool	killLowEnergy (const G4Step *aStep)

void	killTrack (const G4Step *aStep)

Private Attributes
G4VPhysicalVolume *	calo

std::vector< double >	ekinMins

std::vector< std::string >	ekinNames

std::vector< std::string >	ekinParticles

std::vector< int >	ekinPDG

std::vector< G4LogicalVolume * >	ekinVolumes

EventAction *	eventAction_

bool	initialized

bool	killBeamPipe

bool	killByEnergy

bool	killByTimeAtRegion

std::vector< std::string >	maxTimeNames

std::vector< G4Region * >	maxTimeRegions

double	maxTrackTime

std::vector< double >	maxTrackTimes

double	theCriticalDensity

double	theCriticalEnergyForVacuum

G4VPhysicalVolume *	tracker

int	verbose

Detailed Description

Definition at line 19 of file SteppingAction.h.

Constructor & Destructor Documentation

SteppingAction::SteppingAction	(	EventAction *	ea,
		const edm::ParameterSet &	ps
	)

Definition at line 14 of file SteppingAction.cc.

References ekinMins, ekinNames, ekinParticles, ekinVolumes, g, edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), i, killBeamPipe, killByEnergy, killByTimeAtRegion, maxTimeNames, maxTrackTime, maxTrackTimes, theCriticalDensity, and theCriticalEnergyForVacuum.

   : eventAction_(e), tracker(0), calo(0), initialized(false) {
 
   killBeamPipe = (p.getParameter<bool>("KillBeamPipe"));
   theCriticalEnergyForVacuum = (p.getParameter<double>("CriticalEnergyForVacuum")*MeV);
   theCriticalDensity = (p.getParameter<double>("CriticalDensity")*g/cm3);
   maxTrackTime  = (p.getParameter<double>("MaxTrackTime")*ns);
   maxTrackTimes = (p.getParameter<std::vector<double> >("MaxTrackTimes"));
   maxTimeNames  = (p.getParameter<std::vector<std::string> >("MaxTimeNames"));
   ekinMins      = (p.getParameter<std::vector<double> >("EkinThresholds"));
   ekinNames     = (p.getParameter<std::vector<std::string> >("EkinNames"));
   ekinParticles = (p.getParameter<std::vector<std::string> >("EkinParticles"));
   verbose       = (p.getUntrackedParameter<int>("Verbosity",0));
 
   edm::LogInfo("SimG4CoreApplication") << "SteppingAction:: KillBeamPipe = "
                        << killBeamPipe << " CriticalDensity = "
                        << theCriticalDensity << " g/cc;"
                        << " CriticalEnergyForVacuum = "
                        << theCriticalEnergyForVacuum << " Mev;"
                        << " MaxTrackTime = " << maxTrackTime 
                        << " ns";
 
   for (unsigned int i=0; i<maxTrackTimes.size(); i++) {
     maxTrackTimes[i] *= ns;
     edm::LogInfo("SimG4CoreApplication") << "SteppingAction::MaxTrackTime for "
                      << maxTimeNames[i] << " is " 
                      << maxTrackTimes[i];
   }
 
   killByTimeAtRegion = false;
   if(maxTrackTimes.size() > 0) { killByTimeAtRegion = true; }
 
   killByEnergy = false;
   if(ekinVolumes.size() > 0) {
 
     killByEnergy = true;
     edm::LogInfo("SimG4CoreApplication") << "SteppingAction::Kill following "
                      << ekinParticles.size() 
                      << " particles in " << ekinNames.size()
                      << " volumes";
     for (unsigned int i=0; i<ekinParticles.size(); i++) {
       ekinMins[i] *= GeV;
       edm::LogInfo("SimG4CoreApplication") << "SteppingAction::Particle " << i
                        << " " << ekinParticles[i]
                        << " (Threshold = " << ekinMins[i]
                        << " MeV)";
     }
     for (unsigned int i=0; i<ekinNames.size(); i++) 
       edm::LogInfo("SimG4CoreApplication") << "SteppingAction::Volume[" << i
                        << "] = " << ekinNames[i];
   }
 }

SteppingAction::~SteppingAction ( )

Definition at line 67 of file SteppingAction.cc.

67 {}

Member Function Documentation

bool SteppingAction::catchLongLived ( const G4Step * aStep )

private

Definition at line 153 of file SteppingAction.cc.

References i, killByTimeAtRegion, killTrack(), maxTimeRegions, maxTrackTime, maxTrackTimes, and cond::rpcobgas::time.

Referenced by UserSteppingAction().

 {
   bool flag   = true;
   double time = (aStep->GetPostStepPoint()->GetGlobalTime())/nanosecond;
   double tofM = maxTrackTime;
 
   if(killByTimeAtRegion) {
     G4Region* reg = 
       aStep->GetPreStepPoint()->GetPhysicalVolume()->GetLogicalVolume()->GetRegion();
     for (unsigned int i=0; i<maxTimeRegions.size(); i++) {
       if (reg == maxTimeRegions[i]) {
     tofM = maxTrackTimes[i];
     break;
       }
     }
   }
   if (time > tofM) {
     killTrack(aStep);
     flag = false;
   }
   return flag;
 }

bool SteppingAction::catchLowEnergyInVacuum	(	G4Track *	theTrack,
		double	theKenergy
	)

private

Definition at line 116 of file SteppingAction.cc.

References g, and theCriticalDensity.

Referenced by UserSteppingAction().

 {
   bool alive = true;
   if (theTrack->GetVolume()!=0) {
     double density = theTrack->GetVolume()->GetLogicalVolume()->GetMaterial()->GetDensity();
     if (density <= theCriticalDensity) {
       if (verbose>1) {
     edm::LogInfo("SimG4CoreApplication") 
       <<   " SteppingAction: LoopCatchSteppingAction:catchLowEnergyInVacuumHere: "
       << " Track from " << theTrack->GetDefinition()->GetParticleName()
       << " of kinetic energy " << theKenergy/MeV << " MeV "
       << " killed in " << theTrack->GetVolume()->GetLogicalVolume()->GetName()
       << " of density " << density/(g/cm3) << " g/cm3" ;
       }
       theTrack->SetTrackStatus(fStopAndKill);
       alive = false;
     }
     if(alive) {
       density = theTrack->GetNextVolume()->GetLogicalVolume()->GetMaterial()->GetDensity();
       if (density <= theCriticalDensity) {
     if (verbose>1) {
       edm::LogInfo("SimG4CoreApplication") 
         << " SteppingAction: LoopCatchSteppingAction::catchLowEnergyInVacuumNext: "
         << " Track from " << theTrack->GetDefinition()->GetParticleName()
         << " of kinetic energy " << theKenergy/MeV << " MeV "
         << " stopped in " << theTrack->GetVolume()->GetLogicalVolume()->GetName()
         << " before going into "<< theTrack->GetNextVolume()->GetLogicalVolume()->GetName()
         << " of density " << density/(g/cm3) << " g/cm3" ;
     }
     theTrack->SetTrackStatus(fStopButAlive);
     alive = false;
       }
     }
   }
   return alive;
 }

bool SteppingAction::initPointer ( )

private

Definition at line 206 of file SteppingAction.cc.

References calo, ekinMins, ekinNames, ekinParticles, ekinPDG, ekinVolumes, i, LogDebug, maxTimeNames, maxTimeRegions, maxTrackTimes, mergeVDriftHistosByStation::name, pileupDistInMC::num, and tracker.

Referenced by UserSteppingAction().

 {
   bool flag = true;
   const G4PhysicalVolumeStore * pvs = G4PhysicalVolumeStore::GetInstance();
   if (pvs) {
     std::vector<G4VPhysicalVolume*>::const_iterator pvcite;
     for (pvcite = pvs->begin(); pvcite != pvs->end(); pvcite++) {
       if ((*pvcite)->GetName() == "Tracker") tracker = (*pvcite);
       if ((*pvcite)->GetName() == "CALO")    calo    = (*pvcite);
       if (tracker && calo) break;
     }
     edm::LogInfo("SimG4CoreApplication") << "Pointer for Tracker " << tracker
                      << " and for Calo " << calo;
     if (tracker) LogDebug("SimG4CoreApplication") << "Tracker vol name "
                           << tracker->GetName();
     if (calo)    LogDebug("SimG4CoreApplication") << "Calorimeter vol name "
                           << calo->GetName();
   } else {
     flag = false;
   }
 
   const G4LogicalVolumeStore * lvs = G4LogicalVolumeStore::GetInstance();
   unsigned int num = ekinNames.size();
   if (num > 0) {
     if (lvs) {
       std::vector<G4LogicalVolume*>::const_iterator lvcite;
       for (lvcite = lvs->begin(); lvcite != lvs->end(); lvcite++) {
     for (unsigned int i=0; i<num; i++) {
       if ((*lvcite)->GetName() == (G4String)(ekinNames[i])) {
         ekinVolumes.push_back(*lvcite);
         break;
       }
     }
     if (ekinVolumes.size() == num) break;
       }
     }
     if (ekinVolumes.size() != num) flag = false;
     for (unsigned int i=0; i<ekinVolumes.size(); i++) {
       edm::LogInfo("SimG4CoreApplication") << ekinVolumes[i]->GetName()
                        <<" with pointer " <<ekinVolumes[i];
     }
   }
 
   G4ParticleTable * theParticleTable = G4ParticleTable::GetParticleTable();
   G4String particleName;
   for (unsigned int i=0; i<ekinParticles.size(); i++) {
     int pdg = theParticleTable->FindParticle(particleName=ekinParticles[i])->GetPDGEncoding();
     ekinPDG.push_back(pdg);
     edm::LogInfo("SimG4CoreApplication") << "Particle " << ekinParticles[i]
                      << " with code " << ekinPDG[i]
                      << " and KE cut off " << ekinMins[i];
   }
   if (!flag) edm::LogInfo("SimG4CoreApplication") << "SteppingAction fails to"
                           << " initialize some the "
                           << "LV pointers correctly";
 
   const G4RegionStore * rs = G4RegionStore::GetInstance();
   num = maxTimeNames.size();
   if (num > 0) {
     std::vector<double> tofs;
     if (rs) {
       std::vector<G4Region*>::const_iterator rcite;
       for (rcite = rs->begin(); rcite != rs->end(); rcite++) {
     for (unsigned int i=0; i<num; i++) {
       if ((*rcite)->GetName() == (G4String)(maxTimeNames[i])) {
         maxTimeRegions.push_back(*rcite);
         tofs.push_back(maxTrackTimes[i]);
         break;
       }
     }
     if (tofs.size() == num) break;
       }
     }
     for (unsigned int i=0; i<tofs.size(); i++) {
       maxTrackTimes[i] = tofs[i];
       G4String name = "Unknown";
       if (maxTimeRegions[i]) name = maxTimeRegions[i]->GetName();
       edm::LogInfo("SimG4CoreApplication") << name << " with pointer " 
                        << maxTimeRegions[i]<<" KE cut off "
                        << maxTrackTimes[i];
     }
     if (tofs.size() != num) 
       edm::LogInfo("SimG4CoreApplication") << "SteppingAction fails to "
                        << "initialize some the region "
                        << "pointers correctly";
   }
   return true;
 }

bool SteppingAction::isThisVolume	(	const G4VTouchable *	touch,
		G4VPhysicalVolume *	pv
	)

private

Definition at line 295 of file SteppingAction.cc.

References cuy::ii, and testEve_cfg::level.

Referenced by UserSteppingAction().

 {
 
   int level = ((touch->GetHistoryDepth())+1);
   if (level > 0 && level >= 3) {
     unsigned int ii = (unsigned int)(level - 3);
     return (touch->GetVolume(ii) == pv);
   }
   return false;
 }

bool SteppingAction::killLowEnergy ( const G4Step * aStep )

private

Definition at line 176 of file SteppingAction.cc.

References ekinMins, ekinPDG, ekinVolumes, i, killTrack(), and convertSQLiteXML::ok.

Referenced by UserSteppingAction().

 {
   bool ok = true;
   bool flag = false;
   G4LogicalVolume* lv = aStep->GetPreStepPoint()->GetPhysicalVolume()->GetLogicalVolume();
   for (unsigned int i=0; i<ekinVolumes.size(); i++) {
     if (lv == ekinVolumes[i]) {
       flag = true;
       break;
     }
   }
   if (flag) {
     G4Track * track = aStep->GetTrack();
     double    ekin  = track->GetKineticEnergy();
     double    ekinM = 0;
     int       pCode = track->GetDefinition()->GetPDGEncoding();
     for (unsigned int i=0; i<ekinPDG.size(); i++) {
       if (pCode == ekinPDG[i]) {
     ekinM = ekinMins[i];
     break;
       }
     }
     if (ekin < ekinM) {
       killTrack(aStep);
       ok = false;
     }
   }
   return ok;
 }

void SteppingAction::killTrack ( const G4Step * aStep )

private

Definition at line 307 of file SteppingAction.cc.

References GetVolume(), GetRecoTauVFromDQM_MC_cff::kk, and LogDebug.

Referenced by catchLongLived(), and killLowEnergy().

 {
   
   aStep->GetTrack()->SetTrackStatus(fStopAndKill);
   G4TrackVector tv = *(aStep->GetSecondary());
   for (unsigned int kk=0; kk<tv.size(); kk++) {
     if (tv[kk]->GetVolume() == aStep->GetPreStepPoint()->GetPhysicalVolume())
       tv[kk]->SetTrackStatus(fStopAndKill);
   }
   LogDebug("SimG4CoreApplication") 
     << "SteppingAction: Kills track " 
     << aStep->GetTrack()->GetTrackID() << " ("
     << aStep->GetTrack()->GetDefinition()->GetParticleName()
     << ") at " << aStep->GetPostStepPoint()->GetGlobalTime()/nanosecond << " ns in " 
     << aStep->GetPreStepPoint()->GetPhysicalVolume()->GetLogicalVolume()->GetName()
     << " from region " 
     << aStep->GetPreStepPoint()->GetPhysicalVolume()->GetLogicalVolume()->GetRegion()->GetName();
 }

void SteppingAction::UserSteppingAction ( const G4Step * aStep )

Definition at line 69 of file SteppingAction.cc.

References EventAction::addTkCaloStateInfo(), calo, catchLongLived(), catchLowEnergyInVacuum(), eventAction_, initialized, initPointer(), isThisVolume(), killBeamPipe, killByEnergy, killLowEnergy(), m_g4StepSignal, convertSQLiteXML::ok, AlCaHLTBitMon_ParallelJobs::p, theCriticalEnergyForVacuum, tracker, x, detailsBasic3DVector::y, and detailsBasic3DVector::z.

 {
   if (!initialized) { initialized = initPointer(); }
   m_g4StepSignal(aStep);
 
   G4Track * theTrack = aStep->GetTrack();
   bool ok = (theTrack->GetTrackStatus() == fAlive);
   G4double kinEnergy = theTrack->GetKineticEnergy();
 
   if (ok && killBeamPipe && kinEnergy < theCriticalEnergyForVacuum
       && theTrack->GetDefinition()->GetPDGCharge() != 0.0 && kinEnergy > 0.0) {
     ok = catchLowEnergyInVacuum(theTrack, kinEnergy);
 
     //ok = catchLowEnergyInVacuumHere(theTrack);
     //if(ok) { ok = catchLowEnergyInVacuumNext(aStep); }
   }
 
   if(ok && aStep->GetPostStepPoint()->GetPhysicalVolume() != 0) {
 
     ok = catchLongLived(aStep); 
 
     if(ok && killByEnergy) { ok = killLowEnergy(aStep); }
 
     if(ok) {
 
       G4StepPoint* preStep = aStep->GetPreStepPoint();
       if(isThisVolume(preStep->GetTouchable(),tracker) &&
      isThisVolume(aStep->GetPostStepPoint()->GetTouchable(),calo)) {
 
     math::XYZVectorD pos((preStep->GetPosition()).x(),
                  (preStep->GetPosition()).y(),
                  (preStep->GetPosition()).z());
       
     math::XYZTLorentzVectorD mom((preStep->GetMomentum()).x(),
                      (preStep->GetMomentum()).y(),
                      (preStep->GetMomentum()).z(),
                      preStep->GetTotalEnergy());
       
     uint32_t id = aStep->GetTrack()->GetTrackID();
       
     std::pair<math::XYZVectorD,math::XYZTLorentzVectorD> p(pos,mom);
     eventAction_->addTkCaloStateInfo(id,p);
       }
     }
   }
 }