Geant4 Cross Reference |
1 // 2 // ******************************************************************** 3 // * License and Disclaimer * 4 // * * 5 // * The Geant4 software is copyright of the Copyright Holders of * 6 // * the Geant4 Collaboration. It is provided under the terms and * 7 // * conditions of the Geant4 Software License, included in the file * 8 // * LICENSE and available at http://cern.ch/geant4/license . These * 9 // * include a list of copyright holders. * 10 // * * 11 // * Neither the authors of this software system, nor their employing * 12 // * institutes,nor the agencies providing financial support for this * 13 // * work make any representation or warranty, express or implied, * 14 // * regarding this software system or assume any liability for its * 15 // * use. Please see the license in the file LICENSE and URL above * 16 // * for the full disclaimer and the limitation of liability. * 17 // * * 18 // * This code implementation is the result of the scientific and * 19 // * technical work of the GEANT4 collaboration. * 20 // * By using, copying, modifying or distributing the software (or * 21 // * any work based on the software) you agree to acknowledge its * 22 // * use in resulting scientific publications, and indicate your * 23 // * acceptance of all terms of the Geant4 Software license. * 24 // ******************************************************************** 25 // 26 #include "G4LEPTSIonisationModel.hh" 27 #include "CLHEP/Units/PhysicalConstants.h" 28 29 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 30 G4LEPTSIonisationModel::G4LEPTSIonisationModel(const G4String& modelName) 31 : G4VLEPTSModel( modelName ) 32 { 33 SetDeexcitationFlag(true); 34 fParticleChangeForGamma = nullptr; 35 theXSType = XSIonisation; 36 37 } // constructor 38 39 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 40 G4LEPTSIonisationModel::~G4LEPTSIonisationModel() 41 = default; 42 43 44 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 45 void G4LEPTSIonisationModel::Initialise(const G4ParticleDefinition* aParticle, 46 const G4DataVector&) 47 { 48 Init(); 49 BuildPhysicsTable( *aParticle ); 50 fParticleChangeForGamma = GetParticleChangeForGamma(); 51 52 } 53 54 55 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 56 G4double G4LEPTSIonisationModel::CrossSectionPerVolume(const G4Material* mate, 57 const G4ParticleDefinition* aParticle, 58 G4double kineticEnergy, 59 G4double, 60 G4double) 61 { 62 return 1./GetMeanFreePath( mate, aParticle, kineticEnergy ); 63 64 } 65 66 67 void G4LEPTSIonisationModel::SampleSecondaries(std::vector<G4DynamicParticle*>* fvect, 68 const G4MaterialCutsCouple* mateCuts, 69 const G4DynamicParticle* aDynamicParticle, 70 G4double, 71 G4double) 72 { 73 G4double P0KinEn = aDynamicParticle->GetKineticEnergy(); 74 75 G4double Edep=0; 76 G4double Energylost=0; 77 G4ThreeVector P0Dir = aDynamicParticle->GetMomentumDirection(); 78 79 const G4Material* aMaterial = mateCuts->GetMaterial(); 80 if(P0KinEn < theIonisPot[aMaterial]) { 81 theIonisPot[aMaterial] = P0KinEn; 82 } 83 Energylost = SampleEnergyLoss(aMaterial, theIonisPot[aMaterial], P0KinEn); 84 G4ThreeVector P1Dir = SampleNewDirection(aMaterial, P0Dir, P0KinEn/CLHEP::eV, Energylost/CLHEP::eV); 85 G4double P1KinEn = std::max(0., P0KinEn - Energylost); 86 fParticleChangeForGamma->ProposeMomentumDirection( P1Dir); 87 fParticleChangeForGamma->SetProposedKineticEnergy( P1KinEn); 88 #ifdef DEBUG_LEPTS 89 G4cout << " G4LEPTSIonisationModel::SampleSecondaries( SetProposedKineticEnergy " << P1KinEn << " " << P0KinEn << " - " << Energylost << G4endl; 90 #endif 91 92 G4double P2KinEn; 93 94 if( Energylost < theIonisPotInt[aMaterial]) { // External Ionisation 95 //- SetModelName("Ionisation"); 96 Edep = theIonisPot[aMaterial]; 97 P2KinEn = std::max(0.001*CLHEP::eV, (Energylost - theIonisPot[aMaterial]) ); 98 } 99 else { // Auger 100 //- SetModelName("IonisAuger"); 101 Edep = 35*CLHEP::eV; 102 P2KinEn = std::max(0.0, (Energylost - theIonisPotInt[aMaterial]) ); 103 G4double P3KinEn = std::max(0.0, theIonisPotInt[aMaterial] - Edep); 104 105 G4ThreeVector P3Dir; 106 P3Dir.setX( G4UniformRand() ); 107 P3Dir.setY( G4UniformRand() ); 108 P3Dir.setZ( G4UniformRand() ); 109 P3Dir /= P3Dir.mag(); 110 111 auto e3 = new G4DynamicParticle(G4Electron::Electron(), P3Dir, P3KinEn); 112 fvect->push_back(e3); 113 } 114 115 fParticleChangeForGamma->ProposeLocalEnergyDeposit(Edep); 116 117 if( P2KinEn > theLowestEnergyLimit) { 118 G4double cp0 = std::sqrt(P0KinEn*(P0KinEn + 2.*CLHEP::electron_mass_c2)); 119 G4double cp1 = std::sqrt(P1KinEn*(P1KinEn + 2.*CLHEP::electron_mass_c2)); 120 G4ThreeVector P2Momentum = cp0*P0Dir -cp1*P1Dir; 121 G4ThreeVector P2Dir = P2Momentum / P2Momentum.mag(); 122 P2Dir.rotateUz(P0Dir); 123 auto e2 = new G4DynamicParticle(G4Electron::Electron(), P2Dir, P2KinEn); 124 fvect->push_back(e2); 125 } 126 127 } 128