Geant4 Cross Reference |
1 // 1 // 2 // ******************************************* 2 // ******************************************************************** 3 // * License and Disclaimer 3 // * License and Disclaimer * 4 // * 4 // * * 5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of * 6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and * 7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file * 8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These * 9 // * include a list of copyright holders. 9 // * include a list of copyright holders. * 10 // * 10 // * * 11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing * 12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this * 13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, * 14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its * 15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above * 16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. * 17 // * 17 // * * 18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and * 19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. * 20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or * 21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its * 22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your * 23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. * 24 // ******************************************* 24 // ******************************************************************** 25 // 25 // 26 // 26 // 27 /// \file GB05BOptnSplitAndKillByCrossSection. 27 /// \file GB05BOptnSplitAndKillByCrossSection.cc 28 /// \brief Implementation of the GB05BOptnSpli 28 /// \brief Implementation of the GB05BOptnSplitAndKillByCrossSection class 29 29 30 #include "GB05BOptnSplitAndKillByCrossSection. << 31 << 32 #include "Randomize.hh" 30 #include "Randomize.hh" >> 31 #include "GB05BOptnSplitAndKillByCrossSection.hh" 33 32 34 //....oooOO0OOooo........oooOO0OOooo........oo 33 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 35 34 36 GB05BOptnSplitAndKillByCrossSection::GB05BOptn 35 GB05BOptnSplitAndKillByCrossSection::GB05BOptnSplitAndKillByCrossSection(G4String name) 37 : G4VBiasingOperation(name), fParticleChange << 36 : G4VBiasingOperation(name), >> 37 fParticleChange(), >> 38 fInteractionLength(-1.0) 38 {} 39 {} 39 40 40 //....oooOO0OOooo........oooOO0OOooo........oo 41 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 41 42 42 GB05BOptnSplitAndKillByCrossSection::~GB05BOpt << 43 GB05BOptnSplitAndKillByCrossSection::~GB05BOptnSplitAndKillByCrossSection() >> 44 {} 43 45 44 //....oooOO0OOooo........oooOO0OOooo........oo 46 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 45 47 46 G4double GB05BOptnSplitAndKillByCrossSection:: << 48 G4double GB05BOptnSplitAndKillByCrossSection:: 47 << 49 DistanceToApplyOperation( const G4Track*, >> 50 G4double, >> 51 G4ForceCondition* condition) 48 { 52 { 49 *condition = NotForced; 53 *condition = NotForced; 50 54 51 // -- Sample the exponential law using the t 55 // -- Sample the exponential law using the total interaction length of processes 52 // -- to couterbalance for: 56 // -- to couterbalance for: 53 G4double proposedStepLength = -std::log(G4Un << 57 G4double proposedStepLength = -std::log( G4UniformRand() ) * fInteractionLength; 54 << 58 55 return proposedStepLength; 59 return proposedStepLength; 56 } 60 } 57 61 58 //....oooOO0OOooo........oooOO0OOooo........oo 62 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 59 63 60 G4VParticleChange* << 64 G4VParticleChange* 61 GB05BOptnSplitAndKillByCrossSection::GenerateB << 65 GB05BOptnSplitAndKillByCrossSection::GenerateBiasingFinalState( const G4Track* track, >> 66 const G4Step* ) 62 { 67 { >> 68 63 // -- This method is called if we have limit 69 // -- This method is called if we have limited the step. 64 // -- We hence make the splitting or killing 70 // -- We hence make the splitting or killing. 65 71 66 // Get track weight: 72 // Get track weight: 67 G4double initialWeight = track->GetWeight(); 73 G4double initialWeight = track->GetWeight(); 68 << 74 69 // The "particle change" is the object to be 75 // The "particle change" is the object to be used to communicate to 70 // the tracking the update of the primary st 76 // the tracking the update of the primary state and/or creation 71 // secondary tracks. 77 // secondary tracks. 72 fParticleChange.Initialize(*track); 78 fParticleChange.Initialize(*track); 73 79 74 // -- Splitting and killing factors. 80 // -- Splitting and killing factors. 75 // -- They are taken the same, but the killi 81 // -- They are taken the same, but the killing factor can be make bigger. 76 G4double splittingFactor = 2.0; << 82 G4double splittingFactor = 2.0; 77 G4double killingFactor = 2.0; << 83 G4double killingFactor = 2.0; 78 84 79 if (track->GetMomentumDirection().z() > 0) { << 80 // -- We split if the track is moving forw << 81 85 82 // Define the tracks weight: << 86 if ( track->GetMomentumDirection().z() > 0 ) 83 G4double weightOfTrack = initialWeight / s << 87 { 84 << 88 // -- We split if the track is moving forward: 85 // Ask currect track weight to be changed << 89 86 fParticleChange.ProposeParentWeight(weight << 90 // Define the tracks weight: 87 // Now make clones of this track (this is << 91 G4double weightOfTrack = initialWeight/splittingFactor; 88 // we will then have the primary and clone << 92 89 // splitting by a factor 2: << 93 // Ask currect track weight to be changed to new value: 90 G4Track* clone = new G4Track(*track); << 94 fParticleChange.ProposeParentWeight( weightOfTrack ); 91 clone->SetWeight(weightOfTrack); << 95 // Now make clones of this track (this is the actual splitting): 92 fParticleChange.AddSecondary(clone); << 96 // we will then have the primary and clone of it, hence the 93 // -- Below's call added for safety & illu << 97 // splitting by a factor 2: 94 // -- modify the clone (ie : daughter) wei << 98 G4Track* clone = new G4Track( *track ); 95 // -- primary. Here call not mandatory and << 99 clone->SetWeight( weightOfTrack ); 96 fParticleChange.SetSecondaryWeightByProces << 100 fParticleChange.AddSecondary( clone ); 97 } << 101 // -- Below's call added for safety & illustration : inform particle change to not 98 else { << 102 // -- modify the clone (ie : daughter) weight to male it that of the 99 // -- We apply Russian roulette if the tra << 103 // -- primary. Here call not mandatory and both tracks have same weights. 100 << 104 fParticleChange.SetSecondaryWeightByProcess(true); 101 // Shoot a random number (in ]0,1[ segment << 102 G4double random = G4UniformRand(); << 103 G4double killingProbability = 1.0 - 1.0 / << 104 if (random < killingProbability) { << 105 // We ask for the the track to be killed << 106 fParticleChange.ProposeTrackStatus(fStop << 107 } 105 } 108 else { << 106 else 109 // In this case, the track survives. We << 107 { 110 // to conserve weight among killed and s << 108 // -- We apply Russian roulette if the track is moving backward: 111 fParticleChange.ProposeParentWeight(init << 109 >> 110 // Shoot a random number (in ]0,1[ segment): >> 111 G4double random = G4UniformRand(); >> 112 G4double killingProbability = 1.0 - 1.0/killingFactor; >> 113 if ( random < killingProbability ) >> 114 { >> 115 // We ask for the the track to be killed: >> 116 fParticleChange.ProposeTrackStatus(fStopAndKill); >> 117 } >> 118 else >> 119 { >> 120 // In this case, the track survives. We change its weight >> 121 // to conserve weight among killed and survival tracks: >> 122 fParticleChange.ProposeParentWeight( initialWeight*killingFactor ); >> 123 } 112 } 124 } 113 } << 125 114 << 115 return &fParticleChange; 126 return &fParticleChange; >> 127 116 } 128 } 117 129 118 //....oooOO0OOooo........oooOO0OOooo........oo 130 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 119 131