CustomPhysics/src/CMSSQLoopProcessDiscr.cc

0001
0002 #include "SimG4Core/CustomPhysics/interface/CMSSQLoopProcessDiscr.h"
0003 #include "G4SystemOfUnits.hh"
0004 #include "G4Step.hh"
0005 #include "G4ParticleDefinition.hh"
0006 #include "G4VParticleChange.hh"
0007
0008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0009
0010 CMSSQLoopProcessDiscr::CMSSQLoopProcessDiscr(double mass, const G4String& name, G4ProcessType type)
0011     : G4VDiscreteProcess(name, type) {
0012   fParticleChange = new G4ParticleChange();
0013   fParticleChange->ClearDebugFlag();
0014   GenMass = mass;
0015 }
0016
0017 CMSSQLoopProcessDiscr::~CMSSQLoopProcessDiscr() { delete fParticleChange; }
0018
0019 G4VParticleChange* CMSSQLoopProcessDiscr::PostStepDoIt(const G4Track& track, const G4Step& step) {
0020   G4Track* mytr = const_cast<G4Track*>(&track);
0021   mytr->SetPosition(posini);
0022   if (mytr->GetGlobalTime() / ns > 4990)
0023     edm::LogWarning("CMSSQLoopProcess::AlongStepDoIt")
0024         << "going to loose the particle because the GlobalTime is getting close to 5000" << std::endl;
0025
0026   fParticleChange->Clear();
0027   fParticleChange->Initialize(track);
0028
0029   //adding secondary antiS
0030   fParticleChange->SetNumberOfSecondaries(1);
0031   G4DynamicParticle* replacementParticle =
0032       new G4DynamicParticle(CMSAntiSQ::AntiSQ(GenMass), track.GetMomentumDirection(), track.GetKineticEnergy());
0033   fParticleChange->AddSecondary(replacementParticle, globaltimeini);
0034
0035   //killing original AntiS
0036   fParticleChange->ProposeTrackStatus(fStopAndKill);
0037
0038   // note SL: this way of working makes a very long history of the track,
0039   // which all get saved recursively in SimTracks. If the cross section
0040   // is too low such that 10's of thousands of iterations are needed, then
0041   // this becomes too heavy to swallow writing out this history.
0042   // So if we ever need very small cross sections, then we really need
0043   // to change this way of working such that we can throw away all original
0044   // tracks and only save the one that interacted.
0045
0046   return fParticleChange;
0047 }
0048
0049 G4double CMSSQLoopProcessDiscr::PostStepGetPhysicalInteractionLength(const G4Track& track,
0050                                                                      G4double previousStepSize,
0051                                                                      G4ForceCondition* condition) {
0052   *condition = NotForced;
0053   G4double intLength =
0054       DBL_MAX;  //by default the interaction length is super large, only when outside tracker make it 0 to be sure it will do the reset to the original position
0055   G4Track* mytr = const_cast<G4Track*>(&track);
0056   if (sqrt(pow(mytr->GetPosition().rho(), 2)) >
0057       2.45 *
0058           centimeter) {  //this is an important cut for the looping: if the radius of the particle is largher than 2.45cm its interaction length becomes 0 which means it will get killed
0059     // updated from 2.5 to 2.45 so that the Sbar does not start to hit the support of the new inner tracker which was added in 2018
0060     intLength = 0.0;  //0 interaction length means particle will directly interact.
0061   }
0062   return intLength;
0063 }
0064
0065 G4double CMSSQLoopProcessDiscr::GetMeanFreePath(const G4Track&, G4double, G4ForceCondition*) { return DBL_MAX; }
0066
0067 void CMSSQLoopProcessDiscr::StartTracking(G4Track* aTrack) {
0068   posini = aTrack->GetPosition();
0069   globaltimeini = aTrack->GetGlobalTime();
0070 }