File indexing completed on 2022-08-23 02:49:37
0001 #include "G4Version.hh"
0002 #if G4VERSION_NUMBER >= 1100
0003
0004 #include "SimG4Core/PhysicsLists/interface/CMSEmStandardPhysicsTrackingManager.h"
0005 #include "TrackingManagerHelper.h"
0006
0007 #include "G4CoulombScattering.hh"
0008 #include "G4UrbanMscModel.hh"
0009 #include "G4WentzelVIModel.hh"
0010 #include "G4eBremsstrahlung.hh"
0011 #include "G4eCoulombScatteringModel.hh"
0012 #include "G4eIonisation.hh"
0013 #include "G4eMultipleScattering.hh"
0014 #include "G4eplusAnnihilation.hh"
0015
0016 #include "G4ElectroVDNuclearModel.hh"
0017 #include "G4ElectronNuclearProcess.hh"
0018 #include "G4PositronNuclearProcess.hh"
0019
0020 #include "G4ComptonScattering.hh"
0021 #include "G4GammaConversion.hh"
0022 #include "G4PhotoElectricEffect.hh"
0023
0024 #include "G4CascadeInterface.hh"
0025 #include "G4CrossSectionDataSetRegistry.hh"
0026 #include "G4ExcitedStringDecay.hh"
0027 #include "G4GammaNuclearXS.hh"
0028 #include "G4GammaParticipants.hh"
0029 #include "G4GeneratorPrecompoundInterface.hh"
0030 #include "G4HadronInelasticProcess.hh"
0031 #include "G4HadronicParameters.hh"
0032 #include "G4LowEGammaNuclearModel.hh"
0033 #include "G4PhotoNuclearCrossSection.hh"
0034 #include "G4QGSMFragmentation.hh"
0035 #include "G4QGSModel.hh"
0036 #include "G4TheoFSGenerator.hh"
0037
0038 #include "G4GammaGeneralProcess.hh"
0039 #include "G4LossTableManager.hh"
0040
0041 #include "G4EmParameters.hh"
0042 #include "G4SystemOfUnits.hh"
0043
0044 #include "G4Electron.hh"
0045 #include "G4Gamma.hh"
0046 #include "G4Positron.hh"
0047
0048 #include "G4RegionStore.hh"
0049 #include "G4Region.hh"
0050 #include <string>
0051
0052 CMSEmStandardPhysicsTrackingManager *CMSEmStandardPhysicsTrackingManager::masterTrackingManager = nullptr;
0053
0054 CMSEmStandardPhysicsTrackingManager::CMSEmStandardPhysicsTrackingManager(const edm::ParameterSet &p) {
0055 fRangeFactor = p.getParameter<double>("G4MscRangeFactor");
0056 fGeomFactor = p.getParameter<double>("G4MscGeomFactor");
0057 fSafetyFactor = p.getParameter<double>("G4MscSafetyFactor");
0058 fLambdaLimit = p.getParameter<double>("G4MscLambdaLimit") * CLHEP::mm;
0059 std::string msc = p.getParameter<std::string>("G4MscStepLimit");
0060 fStepLimitType = fUseSafety;
0061 if (msc == "UseSafetyPlus") {
0062 fStepLimitType = fUseSafetyPlus;
0063 } else if (msc == "Minimal") {
0064 fStepLimitType = fMinimal;
0065 }
0066
0067 G4EmParameters *param = G4EmParameters::Instance();
0068 G4double highEnergyLimit = param->MscEnergyLimit();
0069
0070 const G4Region *aRegion = G4RegionStore::GetInstance()->GetRegion("HcalRegion", false);
0071 const G4Region *bRegion = G4RegionStore::GetInstance()->GetRegion("HGCalRegion", false);
0072
0073
0074 {
0075 G4eMultipleScattering *msc = new G4eMultipleScattering;
0076 G4UrbanMscModel *msc1 = new G4UrbanMscModel;
0077 G4WentzelVIModel *msc2 = new G4WentzelVIModel;
0078 msc1->SetHighEnergyLimit(highEnergyLimit);
0079 msc2->SetLowEnergyLimit(highEnergyLimit);
0080 msc->SetEmModel(msc1);
0081 msc->SetEmModel(msc2);
0082
0083
0084 if (nullptr != aRegion || nullptr != bRegion) {
0085 G4UrbanMscModel *msc3 = new G4UrbanMscModel();
0086 msc3->SetHighEnergyLimit(highEnergyLimit);
0087 msc3->SetRangeFactor(fRangeFactor);
0088 msc3->SetGeomFactor(fGeomFactor);
0089 msc3->SetSafetyFactor(fSafetyFactor);
0090 msc3->SetLambdaLimit(fLambdaLimit);
0091 msc3->SetStepLimitType(fStepLimitType);
0092 msc3->SetLocked(true);
0093 if (nullptr != aRegion) {
0094 msc->AddEmModel(-1, msc3, aRegion);
0095 }
0096 if (nullptr != bRegion) {
0097 msc->AddEmModel(-1, msc3, bRegion);
0098 }
0099 }
0100
0101 electron.msc = msc;
0102
0103 electron.ioni = new G4eIonisation;
0104 electron.brems = new G4eBremsstrahlung;
0105
0106 G4CoulombScattering *ss = new G4CoulombScattering;
0107 G4eCoulombScatteringModel *ssm = new G4eCoulombScatteringModel;
0108 ssm->SetLowEnergyLimit(highEnergyLimit);
0109 ssm->SetActivationLowEnergyLimit(highEnergyLimit);
0110 ss->SetEmModel(ssm);
0111 ss->SetMinKinEnergy(highEnergyLimit);
0112 electron.ss = ss;
0113 }
0114
0115
0116 {
0117 G4eMultipleScattering *msc = new G4eMultipleScattering;
0118 G4UrbanMscModel *msc1 = new G4UrbanMscModel;
0119 G4WentzelVIModel *msc2 = new G4WentzelVIModel;
0120 msc1->SetHighEnergyLimit(highEnergyLimit);
0121 msc2->SetLowEnergyLimit(highEnergyLimit);
0122 msc->SetEmModel(msc1);
0123 msc->SetEmModel(msc2);
0124
0125
0126 if (nullptr != aRegion || nullptr != bRegion) {
0127 G4UrbanMscModel *msc3 = new G4UrbanMscModel();
0128 msc3->SetHighEnergyLimit(highEnergyLimit);
0129 msc3->SetRangeFactor(fRangeFactor);
0130 msc3->SetGeomFactor(fGeomFactor);
0131 msc3->SetSafetyFactor(fSafetyFactor);
0132 msc3->SetLambdaLimit(fLambdaLimit);
0133 msc3->SetStepLimitType(fStepLimitType);
0134 msc3->SetLocked(true);
0135 if (nullptr != aRegion) {
0136 msc->AddEmModel(-1, msc3, aRegion);
0137 }
0138 if (nullptr != bRegion) {
0139 msc->AddEmModel(-1, msc3, bRegion);
0140 }
0141 }
0142
0143 positron.msc = msc;
0144
0145 positron.ioni = new G4eIonisation;
0146 positron.brems = new G4eBremsstrahlung;
0147 positron.annihilation = new G4eplusAnnihilation;
0148
0149 G4CoulombScattering *ss = new G4CoulombScattering;
0150 G4eCoulombScatteringModel *ssm = new G4eCoulombScatteringModel;
0151 ssm->SetLowEnergyLimit(highEnergyLimit);
0152 ssm->SetActivationLowEnergyLimit(highEnergyLimit);
0153 ss->SetEmModel(ssm);
0154 ss->SetMinKinEnergy(highEnergyLimit);
0155 positron.ss = ss;
0156 }
0157
0158
0159 {
0160 gammaProc = new G4GammaGeneralProcess;
0161
0162 gammaProc->AddEmProcess(new G4PhotoElectricEffect);
0163 gammaProc->AddEmProcess(new G4ComptonScattering);
0164 gammaProc->AddEmProcess(new G4GammaConversion);
0165
0166 G4HadronInelasticProcess *nuc = new G4HadronInelasticProcess("photonNuclear", G4Gamma::Definition());
0167 auto xsreg = G4CrossSectionDataSetRegistry::Instance();
0168 G4VCrossSectionDataSet *xs = nullptr;
0169 bool useGammaNuclearXS = true;
0170 if (useGammaNuclearXS) {
0171 xs = xsreg->GetCrossSectionDataSet("GammaNuclearXS");
0172 if (nullptr == xs)
0173 xs = new G4GammaNuclearXS;
0174 } else {
0175 xs = xsreg->GetCrossSectionDataSet("PhotoNuclearXS");
0176 if (nullptr == xs)
0177 xs = new G4PhotoNuclearCrossSection;
0178 }
0179 nuc->AddDataSet(xs);
0180
0181 G4QGSModel<G4GammaParticipants> *theStringModel = new G4QGSModel<G4GammaParticipants>;
0182 G4QGSMFragmentation *theFrag = new G4QGSMFragmentation;
0183 G4ExcitedStringDecay *theStringDecay = new G4ExcitedStringDecay(theFrag);
0184 theStringModel->SetFragmentationModel(theStringDecay);
0185
0186 G4GeneratorPrecompoundInterface *theCascade = new G4GeneratorPrecompoundInterface;
0187
0188 G4TheoFSGenerator *theModel = new G4TheoFSGenerator;
0189 theModel->SetTransport(theCascade);
0190 theModel->SetHighEnergyGenerator(theStringModel);
0191
0192 G4HadronicParameters *param = G4HadronicParameters::Instance();
0193
0194 G4CascadeInterface *cascade = new G4CascadeInterface;
0195
0196
0197 G4double gnLowEnergyLimit = 200 * CLHEP::MeV;
0198 if (gnLowEnergyLimit > 0.0) {
0199 G4LowEGammaNuclearModel *lemod = new G4LowEGammaNuclearModel;
0200 lemod->SetMaxEnergy(gnLowEnergyLimit);
0201 nuc->RegisterMe(lemod);
0202 cascade->SetMinEnergy(gnLowEnergyLimit - CLHEP::MeV);
0203 }
0204 cascade->SetMaxEnergy(param->GetMaxEnergyTransitionFTF_Cascade());
0205 nuc->RegisterMe(cascade);
0206 theModel->SetMinEnergy(param->GetMinEnergyTransitionFTF_Cascade());
0207 theModel->SetMaxEnergy(param->GetMaxEnergy());
0208 nuc->RegisterMe(theModel);
0209
0210 gammaProc->AddHadProcess(nuc);
0211 G4LossTableManager::Instance()->SetGammaGeneralProcess(gammaProc);
0212 }
0213
0214
0215
0216 G4ElectroVDNuclearModel *eModel = new G4ElectroVDNuclearModel;
0217
0218 {
0219 G4ElectronNuclearProcess *nuc = new G4ElectronNuclearProcess;
0220 nuc->RegisterMe(eModel);
0221 electron.nuc = nuc;
0222 }
0223 {
0224 G4PositronNuclearProcess *nuc = new G4PositronNuclearProcess;
0225 nuc->RegisterMe(eModel);
0226 positron.nuc = nuc;
0227 }
0228
0229 if (masterTrackingManager == nullptr) {
0230 masterTrackingManager = this;
0231 } else {
0232 electron.msc->SetMasterProcess(masterTrackingManager->electron.msc);
0233 electron.ss->SetMasterProcess(masterTrackingManager->electron.ss);
0234 electron.ioni->SetMasterProcess(masterTrackingManager->electron.ioni);
0235 electron.brems->SetMasterProcess(masterTrackingManager->electron.brems);
0236 electron.nuc->SetMasterProcess(masterTrackingManager->electron.nuc);
0237
0238 positron.msc->SetMasterProcess(masterTrackingManager->positron.msc);
0239 positron.ss->SetMasterProcess(masterTrackingManager->positron.ss);
0240 positron.ioni->SetMasterProcess(masterTrackingManager->positron.ioni);
0241 positron.brems->SetMasterProcess(masterTrackingManager->positron.brems);
0242 positron.annihilation->SetMasterProcess(masterTrackingManager->positron.annihilation);
0243 positron.nuc->SetMasterProcess(masterTrackingManager->positron.nuc);
0244
0245 gammaProc->SetMasterProcess(masterTrackingManager->gammaProc);
0246 }
0247 }
0248
0249 CMSEmStandardPhysicsTrackingManager::~CMSEmStandardPhysicsTrackingManager() {
0250 if (masterTrackingManager == this) {
0251 masterTrackingManager = nullptr;
0252 }
0253 }
0254
0255 void CMSEmStandardPhysicsTrackingManager::BuildPhysicsTable(const G4ParticleDefinition &part) {
0256 if (&part == G4Electron::Definition()) {
0257 electron.msc->BuildPhysicsTable(part);
0258 electron.ioni->BuildPhysicsTable(part);
0259 electron.brems->BuildPhysicsTable(part);
0260 electron.ss->BuildPhysicsTable(part);
0261 electron.nuc->BuildPhysicsTable(part);
0262 } else if (&part == G4Positron::Definition()) {
0263 positron.msc->BuildPhysicsTable(part);
0264 positron.ioni->BuildPhysicsTable(part);
0265 positron.brems->BuildPhysicsTable(part);
0266 positron.annihilation->BuildPhysicsTable(part);
0267 positron.ss->BuildPhysicsTable(part);
0268 positron.nuc->BuildPhysicsTable(part);
0269 } else if (&part == G4Gamma::Definition()) {
0270 gammaProc->BuildPhysicsTable(part);
0271 }
0272 }
0273
0274 void CMSEmStandardPhysicsTrackingManager::PreparePhysicsTable(const G4ParticleDefinition &part) {
0275 if (&part == G4Electron::Definition()) {
0276 electron.msc->PreparePhysicsTable(part);
0277 electron.ioni->PreparePhysicsTable(part);
0278 electron.brems->PreparePhysicsTable(part);
0279 electron.ss->PreparePhysicsTable(part);
0280 electron.nuc->PreparePhysicsTable(part);
0281 } else if (&part == G4Positron::Definition()) {
0282 positron.msc->PreparePhysicsTable(part);
0283 positron.ioni->PreparePhysicsTable(part);
0284 positron.brems->PreparePhysicsTable(part);
0285 positron.annihilation->PreparePhysicsTable(part);
0286 positron.ss->PreparePhysicsTable(part);
0287 positron.nuc->PreparePhysicsTable(part);
0288 } else if (&part == G4Gamma::Definition()) {
0289 gammaProc->PreparePhysicsTable(part);
0290 }
0291 }
0292
0293 void CMSEmStandardPhysicsTrackingManager::TrackElectron(G4Track *aTrack) {
0294 class ElectronPhysics final : public TrackingManagerHelper::Physics {
0295 public:
0296 ElectronPhysics(CMSEmStandardPhysicsTrackingManager &mgr) : fMgr(mgr) {}
0297
0298 void StartTracking(G4Track *aTrack) override {
0299 auto &electron = fMgr.electron;
0300
0301 electron.msc->StartTracking(aTrack);
0302 electron.ioni->StartTracking(aTrack);
0303 electron.brems->StartTracking(aTrack);
0304 electron.ss->StartTracking(aTrack);
0305 electron.nuc->StartTracking(aTrack);
0306
0307 fPreviousStepLength = 0;
0308 }
0309 void EndTracking() override {
0310 auto &electron = fMgr.electron;
0311
0312 electron.msc->EndTracking();
0313 electron.ioni->EndTracking();
0314 electron.brems->EndTracking();
0315 electron.ss->EndTracking();
0316 electron.nuc->EndTracking();
0317 }
0318
0319 G4double GetPhysicalInteractionLength(const G4Track &track) override {
0320 auto &electron = fMgr.electron;
0321 G4double physIntLength, proposedSafety = DBL_MAX;
0322 G4ForceCondition condition;
0323 G4GPILSelection selection;
0324
0325 fProposedStep = DBL_MAX;
0326 fSelected = -1;
0327
0328 physIntLength = electron.nuc->PostStepGPIL(track, fPreviousStepLength, &condition);
0329 if (physIntLength < fProposedStep) {
0330 fProposedStep = physIntLength;
0331 fSelected = 3;
0332 }
0333
0334 physIntLength = electron.ss->PostStepGPIL(track, fPreviousStepLength, &condition);
0335 if (physIntLength < fProposedStep) {
0336 fProposedStep = physIntLength;
0337 fSelected = 0;
0338 }
0339
0340 physIntLength = electron.brems->PostStepGPIL(track, fPreviousStepLength, &condition);
0341 if (physIntLength < fProposedStep) {
0342 fProposedStep = physIntLength;
0343 fSelected = 1;
0344 }
0345
0346 physIntLength = electron.ioni->PostStepGPIL(track, fPreviousStepLength, &condition);
0347 if (physIntLength < fProposedStep) {
0348 fProposedStep = physIntLength;
0349 fSelected = 2;
0350 }
0351
0352 physIntLength =
0353 electron.ioni->AlongStepGPIL(track, fPreviousStepLength, fProposedStep, proposedSafety, &selection);
0354 if (physIntLength < fProposedStep) {
0355 fProposedStep = physIntLength;
0356 fSelected = -1;
0357 }
0358
0359 physIntLength =
0360 electron.msc->AlongStepGPIL(track, fPreviousStepLength, fProposedStep, proposedSafety, &selection);
0361 if (physIntLength < fProposedStep) {
0362 fProposedStep = physIntLength;
0363
0364
0365 if (selection == CandidateForSelection) {
0366 fSelected = -1;
0367 }
0368 }
0369
0370 return fProposedStep;
0371 }
0372
0373 void AlongStepDoIt(G4Track &track, G4Step &step, G4TrackVector &) override {
0374 if (step.GetStepLength() == fProposedStep) {
0375 step.GetPostStepPoint()->SetStepStatus(fAlongStepDoItProc);
0376 } else {
0377
0378 fSelected = -1;
0379 }
0380 auto &electron = fMgr.electron;
0381 G4VParticleChange *particleChange;
0382
0383 particleChange = electron.msc->AlongStepDoIt(track, step);
0384 particleChange->UpdateStepForAlongStep(&step);
0385 track.SetTrackStatus(particleChange->GetTrackStatus());
0386 particleChange->Clear();
0387
0388 particleChange = electron.ioni->AlongStepDoIt(track, step);
0389 particleChange->UpdateStepForAlongStep(&step);
0390 track.SetTrackStatus(particleChange->GetTrackStatus());
0391 particleChange->Clear();
0392
0393 fPreviousStepLength = step.GetStepLength();
0394 }
0395
0396 void PostStepDoIt(G4Track &track, G4Step &step, G4TrackVector &secondaries) override {
0397 if (fSelected < 0) {
0398 return;
0399 }
0400 step.GetPostStepPoint()->SetStepStatus(fPostStepDoItProc);
0401
0402 auto &electron = fMgr.electron;
0403 G4VProcess *process = nullptr;
0404 G4VParticleChange *particleChange = nullptr;
0405
0406 switch (fSelected) {
0407 case 3:
0408 process = electron.nuc;
0409 particleChange = electron.nuc->PostStepDoIt(track, step);
0410 break;
0411 case 0:
0412 process = electron.ss;
0413 particleChange = electron.ss->PostStepDoIt(track, step);
0414 break;
0415 case 1:
0416 process = electron.brems;
0417 particleChange = electron.brems->PostStepDoIt(track, step);
0418 break;
0419 case 2:
0420 process = electron.ioni;
0421 particleChange = electron.ioni->PostStepDoIt(track, step);
0422 break;
0423 }
0424
0425 particleChange->UpdateStepForPostStep(&step);
0426 step.UpdateTrack();
0427
0428 int numSecondaries = particleChange->GetNumberOfSecondaries();
0429 for (int i = 0; i < numSecondaries; i++) {
0430 G4Track *secondary = particleChange->GetSecondary(i);
0431 secondary->SetParentID(track.GetTrackID());
0432 secondary->SetCreatorProcess(process);
0433 secondaries.push_back(secondary);
0434 }
0435
0436 track.SetTrackStatus(particleChange->GetTrackStatus());
0437 particleChange->Clear();
0438 }
0439
0440 private:
0441 CMSEmStandardPhysicsTrackingManager &fMgr;
0442 G4double fPreviousStepLength;
0443 G4double fProposedStep;
0444 G4int fSelected;
0445 };
0446
0447 ElectronPhysics physics(*this);
0448 TrackingManagerHelper::TrackChargedParticle(aTrack, physics);
0449 }
0450
0451 void CMSEmStandardPhysicsTrackingManager::TrackPositron(G4Track *aTrack) {
0452 class PositronPhysics final : public TrackingManagerHelper::Physics {
0453 public:
0454 PositronPhysics(CMSEmStandardPhysicsTrackingManager &mgr) : fMgr(mgr) {}
0455
0456 void StartTracking(G4Track *aTrack) override {
0457 auto &positron = fMgr.positron;
0458
0459 positron.msc->StartTracking(aTrack);
0460 positron.ioni->StartTracking(aTrack);
0461 positron.brems->StartTracking(aTrack);
0462 positron.annihilation->StartTracking(aTrack);
0463 positron.ss->StartTracking(aTrack);
0464 positron.nuc->StartTracking(aTrack);
0465
0466 fPreviousStepLength = 0;
0467 }
0468 void EndTracking() override {
0469 auto &positron = fMgr.positron;
0470
0471 positron.msc->EndTracking();
0472 positron.ioni->EndTracking();
0473 positron.brems->EndTracking();
0474 positron.annihilation->EndTracking();
0475 positron.ss->EndTracking();
0476 positron.nuc->EndTracking();
0477 }
0478
0479 G4double GetPhysicalInteractionLength(const G4Track &track) override {
0480 auto &positron = fMgr.positron;
0481 G4double physIntLength, proposedSafety = DBL_MAX;
0482 G4ForceCondition condition;
0483 G4GPILSelection selection;
0484
0485 fProposedStep = DBL_MAX;
0486 fSelected = -1;
0487
0488 physIntLength = positron.nuc->PostStepGPIL(track, fPreviousStepLength, &condition);
0489 if (physIntLength < fProposedStep) {
0490 fProposedStep = physIntLength;
0491 fSelected = 4;
0492 }
0493
0494 physIntLength = positron.ss->PostStepGPIL(track, fPreviousStepLength, &condition);
0495 if (physIntLength < fProposedStep) {
0496 fProposedStep = physIntLength;
0497 fSelected = 0;
0498 }
0499
0500 physIntLength = positron.annihilation->PostStepGPIL(track, fPreviousStepLength, &condition);
0501 if (physIntLength < fProposedStep) {
0502 fProposedStep = physIntLength;
0503 fSelected = 1;
0504 }
0505
0506 physIntLength = positron.brems->PostStepGPIL(track, fPreviousStepLength, &condition);
0507 if (physIntLength < fProposedStep) {
0508 fProposedStep = physIntLength;
0509 fSelected = 2;
0510 }
0511
0512 physIntLength = positron.ioni->PostStepGPIL(track, fPreviousStepLength, &condition);
0513 if (physIntLength < fProposedStep) {
0514 fProposedStep = physIntLength;
0515 fSelected = 3;
0516 }
0517
0518 physIntLength =
0519 positron.ioni->AlongStepGPIL(track, fPreviousStepLength, fProposedStep, proposedSafety, &selection);
0520 if (physIntLength < fProposedStep) {
0521 fProposedStep = physIntLength;
0522 fSelected = -1;
0523 }
0524
0525 physIntLength =
0526 positron.msc->AlongStepGPIL(track, fPreviousStepLength, fProposedStep, proposedSafety, &selection);
0527 if (physIntLength < fProposedStep) {
0528 fProposedStep = physIntLength;
0529
0530
0531 if (selection == CandidateForSelection) {
0532 fSelected = -1;
0533 }
0534 }
0535
0536 return fProposedStep;
0537 }
0538
0539 void AlongStepDoIt(G4Track &track, G4Step &step, G4TrackVector &) override {
0540 if (step.GetStepLength() == fProposedStep) {
0541 step.GetPostStepPoint()->SetStepStatus(fAlongStepDoItProc);
0542 } else {
0543
0544 fSelected = -1;
0545 }
0546 auto &positron = fMgr.positron;
0547 G4VParticleChange *particleChange;
0548
0549 particleChange = positron.msc->AlongStepDoIt(track, step);
0550 particleChange->UpdateStepForAlongStep(&step);
0551 track.SetTrackStatus(particleChange->GetTrackStatus());
0552 particleChange->Clear();
0553
0554 particleChange = positron.ioni->AlongStepDoIt(track, step);
0555 particleChange->UpdateStepForAlongStep(&step);
0556 track.SetTrackStatus(particleChange->GetTrackStatus());
0557 particleChange->Clear();
0558
0559 fPreviousStepLength = step.GetStepLength();
0560 }
0561
0562 void PostStepDoIt(G4Track &track, G4Step &step, G4TrackVector &secondaries) override {
0563 if (fSelected < 0) {
0564 return;
0565 }
0566 step.GetPostStepPoint()->SetStepStatus(fPostStepDoItProc);
0567
0568 auto &positron = fMgr.positron;
0569 G4VProcess *process;
0570 G4VParticleChange *particleChange = nullptr;
0571
0572 switch (fSelected) {
0573 case 4:
0574 process = positron.nuc;
0575 particleChange = positron.nuc->PostStepDoIt(track, step);
0576 break;
0577 case 0:
0578 process = positron.ss;
0579 particleChange = positron.ss->PostStepDoIt(track, step);
0580 break;
0581 case 1:
0582 process = positron.annihilation;
0583 particleChange = positron.annihilation->PostStepDoIt(track, step);
0584 break;
0585 case 2:
0586 process = positron.brems;
0587 particleChange = positron.brems->PostStepDoIt(track, step);
0588 break;
0589 case 3:
0590 process = positron.ioni;
0591 particleChange = positron.ioni->PostStepDoIt(track, step);
0592 break;
0593 }
0594
0595 particleChange->UpdateStepForPostStep(&step);
0596 step.UpdateTrack();
0597
0598 int numSecondaries = particleChange->GetNumberOfSecondaries();
0599 for (int i = 0; i < numSecondaries; i++) {
0600 G4Track *secondary = particleChange->GetSecondary(i);
0601 secondary->SetParentID(track.GetTrackID());
0602 secondary->SetCreatorProcess(process);
0603 secondaries.push_back(secondary);
0604 }
0605
0606 track.SetTrackStatus(particleChange->GetTrackStatus());
0607 particleChange->Clear();
0608 }
0609
0610 G4bool HasAtRestProcesses() override { return true; }
0611
0612 void AtRestDoIt(G4Track &track, G4Step &step, G4TrackVector &secondaries) override {
0613 auto &positron = fMgr.positron;
0614
0615 G4VParticleChange *particleChange = positron.annihilation->AtRestDoIt(track, step);
0616 particleChange->UpdateStepForAtRest(&step);
0617 step.UpdateTrack();
0618
0619 int numSecondaries = particleChange->GetNumberOfSecondaries();
0620 for (int i = 0; i < numSecondaries; i++) {
0621 G4Track *secondary = particleChange->GetSecondary(i);
0622 secondary->SetParentID(track.GetTrackID());
0623 secondary->SetCreatorProcess(positron.annihilation);
0624 secondaries.push_back(secondary);
0625 }
0626
0627 track.SetTrackStatus(particleChange->GetTrackStatus());
0628 particleChange->Clear();
0629 }
0630
0631 private:
0632 CMSEmStandardPhysicsTrackingManager &fMgr;
0633 G4double fPreviousStepLength;
0634 G4double fProposedStep;
0635 G4int fSelected;
0636 };
0637
0638 PositronPhysics physics(*this);
0639 TrackingManagerHelper::TrackChargedParticle(aTrack, physics);
0640 }
0641
0642 void CMSEmStandardPhysicsTrackingManager::TrackGamma(G4Track *aTrack) {
0643 class GammaPhysics final : public TrackingManagerHelper::Physics {
0644 public:
0645 GammaPhysics(CMSEmStandardPhysicsTrackingManager &mgr) : fMgr(mgr) {}
0646
0647 void StartTracking(G4Track *aTrack) override {
0648 fMgr.gammaProc->StartTracking(aTrack);
0649
0650 fPreviousStepLength = 0;
0651 }
0652 void EndTracking() override { fMgr.gammaProc->EndTracking(); }
0653
0654 G4double GetPhysicalInteractionLength(const G4Track &track) override {
0655 G4double physIntLength;
0656 G4ForceCondition condition;
0657
0658 fProposedStep = DBL_MAX;
0659 fSelected = -1;
0660
0661 physIntLength = fMgr.gammaProc->PostStepGPIL(track, fPreviousStepLength, &condition);
0662 if (physIntLength < fProposedStep) {
0663 fProposedStep = physIntLength;
0664 fSelected = 0;
0665 }
0666
0667 return fProposedStep;
0668 }
0669
0670 void AlongStepDoIt(G4Track &, G4Step &step, G4TrackVector &) override {
0671 if (step.GetStepLength() == fProposedStep) {
0672 step.GetPostStepPoint()->SetStepStatus(fAlongStepDoItProc);
0673 } else {
0674
0675 fSelected = -1;
0676 }
0677 fPreviousStepLength = step.GetStepLength();
0678 }
0679
0680 void PostStepDoIt(G4Track &track, G4Step &step, G4TrackVector &secondaries) override {
0681 if (fSelected < 0) {
0682 return;
0683 }
0684 step.GetPostStepPoint()->SetStepStatus(fPostStepDoItProc);
0685
0686 G4VProcess *process = fMgr.gammaProc;
0687 G4VParticleChange *particleChange = fMgr.gammaProc->PostStepDoIt(track, step);
0688
0689 particleChange->UpdateStepForPostStep(&step);
0690 step.UpdateTrack();
0691
0692 int numSecondaries = particleChange->GetNumberOfSecondaries();
0693 for (int i = 0; i < numSecondaries; i++) {
0694 G4Track *secondary = particleChange->GetSecondary(i);
0695 secondary->SetParentID(track.GetTrackID());
0696 secondary->SetCreatorProcess(process);
0697 secondaries.push_back(secondary);
0698 }
0699
0700 track.SetTrackStatus(particleChange->GetTrackStatus());
0701 particleChange->Clear();
0702 }
0703
0704 private:
0705 CMSEmStandardPhysicsTrackingManager &fMgr;
0706 G4double fPreviousStepLength;
0707 G4double fProposedStep;
0708 G4int fSelected;
0709 };
0710
0711 GammaPhysics physics(*this);
0712 TrackingManagerHelper::TrackNeutralParticle(aTrack, physics);
0713 }
0714
0715 void CMSEmStandardPhysicsTrackingManager::HandOverOneTrack(G4Track *aTrack) {
0716 const G4ParticleDefinition *part = aTrack->GetParticleDefinition();
0717
0718 if (part == G4Electron::Definition()) {
0719 TrackElectron(aTrack);
0720 } else if (part == G4Positron::Definition()) {
0721 TrackPositron(aTrack);
0722 } else if (part == G4Gamma::Definition()) {
0723 TrackGamma(aTrack);
0724 }
0725
0726 aTrack->SetTrackStatus(fStopAndKill);
0727 delete aTrack;
0728 }
0729
0730 #endif