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 * =========================================== 27 * ============================================================================ 28 * 28 * 29 * Filename: CexmcHadronicProcess.cc 29 * Filename: CexmcHadronicProcess.cc 30 * 30 * 31 * Description: hadronic process with prod 31 * Description: hadronic process with production model 32 * 32 * 33 * Version: 1.0 33 * Version: 1.0 34 * Created: 31.10.2009 23:54:38 34 * Created: 31.10.2009 23:54:38 35 * Revision: none 35 * Revision: none 36 * Compiler: gcc 36 * Compiler: gcc 37 * 37 * 38 * Author: Alexey Radkov (), 38 * Author: Alexey Radkov (), 39 * Company: PNPI 39 * Company: PNPI 40 * 40 * 41 * =========================================== 41 * ============================================================================ 42 */ 42 */ 43 43 44 #include <G4ParticleChange.hh> 44 #include <G4ParticleChange.hh> 45 #include <G4ParticleDefinition.hh> 45 #include <G4ParticleDefinition.hh> 46 #include <G4HadronicInteraction.hh> 46 #include <G4HadronicInteraction.hh> 47 #include <G4Track.hh> 47 #include <G4Track.hh> 48 #include <G4Step.hh> 48 #include <G4Step.hh> 49 #include <G4Element.hh> 49 #include <G4Element.hh> 50 #include <G4StableIsotopes.hh> 50 #include <G4StableIsotopes.hh> 51 #include <G4TrackStatus.hh> 51 #include <G4TrackStatus.hh> 52 #include "CexmcHadronicProcess.hh" 52 #include "CexmcHadronicProcess.hh" 53 #include "CexmcProductionModel.hh" 53 #include "CexmcProductionModel.hh" 54 #include "CexmcIncidentParticleTrackInfo.hh" 54 #include "CexmcIncidentParticleTrackInfo.hh" 55 #include "CexmcException.hh" 55 #include "CexmcException.hh" 56 56 57 57 58 CexmcHadronicProcess::CexmcHadronicProcess( co 58 CexmcHadronicProcess::CexmcHadronicProcess( const G4String & name ) : 59 G4HadronicProcess( name ), productionModel 59 G4HadronicProcess( name ), productionModel( NULL ), interaction( NULL ), 60 theTotalResult( NULL ), isInitialized( fal 60 theTotalResult( NULL ), isInitialized( false ) 61 { 61 { 62 theTotalResult = new G4ParticleChange(); 62 theTotalResult = new G4ParticleChange(); 63 SetIntegral(false); << 64 } 63 } 65 64 66 65 67 CexmcHadronicProcess::~CexmcHadronicProcess() 66 CexmcHadronicProcess::~CexmcHadronicProcess() 68 { 67 { 69 delete theTotalResult; 68 delete theTotalResult; 70 } 69 } 71 70 72 71 73 void CexmcHadronicProcess::RegisterProduction 72 void CexmcHadronicProcess::RegisterProductionModel( 74 << 73 CexmcProductionModel * productionModel_ ) 75 { 74 { 76 productionModel = model; << 75 productionModel = productionModel_; 77 76 78 interaction = dynamic_cast< G4HadronicInte 77 interaction = dynamic_cast< G4HadronicInteraction * >( productionModel ); 79 78 80 if ( ! interaction ) 79 if ( ! interaction ) 81 throw CexmcException( CexmcIncompatibl 80 throw CexmcException( CexmcIncompatibleProductionModel ); 82 81 83 G4HadronicProcess::RegisterMe( interaction 82 G4HadronicProcess::RegisterMe( interaction ); 84 } 83 } 85 84 86 85 87 void CexmcHadronicProcess::CalculateTargetNuc 86 void CexmcHadronicProcess::CalculateTargetNucleus( 88 87 const G4Material * material ) 89 { 88 { 90 G4int numberOfElements( material->GetNumb 89 G4int numberOfElements( material->GetNumberOfElements() ); 91 if ( numberOfElements > 1 ) 90 if ( numberOfElements > 1 ) 92 { 91 { 93 G4cout << CEXMC_LINE_START "WARNING: N 92 G4cout << CEXMC_LINE_START "WARNING: Number of elements in target " 94 "material is more than 1.\n 93 "material is more than 1.\n Only the first " 95 "element will be chosen for 94 "element will be chosen for target nucleus" << G4endl; 96 } 95 } 97 96 98 const G4Element * element( material->GetE 97 const G4Element * element( material->GetElement( 0 ) ); 99 G4double ZZ( element->GetZ() ); 98 G4double ZZ( element->GetZ() ); 100 G4int Z( G4int( ZZ + 0.5 ) ); 99 G4int Z( G4int( ZZ + 0.5 ) ); 101 100 102 G4StableIsotopes stableIsotopes; 101 G4StableIsotopes stableIsotopes; 103 G4int index( stableIsotopes.Ge 102 G4int index( stableIsotopes.GetFirstIsotope( Z ) ); 104 G4double AA( stableIsotopes.GetIs 103 G4double AA( stableIsotopes.GetIsotopeNucleonCount( index ) ); 105 104 106 targetNucleus.SetParameters( AA, ZZ ); 105 targetNucleus.SetParameters( AA, ZZ ); 107 } 106 } 108 107 109 108 110 void CexmcHadronicProcess::FillTotalResult( G 109 void CexmcHadronicProcess::FillTotalResult( G4HadFinalState * hadFinalState, 111 c 110 const G4Track & track ) 112 { 111 { 113 G4int numberOfSecondaries( hadFinalState- 112 G4int numberOfSecondaries( hadFinalState->GetNumberOfSecondaries() ); 114 113 115 theTotalResult->Clear(); 114 theTotalResult->Clear(); 116 theTotalResult->Initialize( track ); 115 theTotalResult->Initialize( track ); 117 theTotalResult->SetSecondaryWeightByProces 116 theTotalResult->SetSecondaryWeightByProcess( true ); 118 theTotalResult->ProposeLocalEnergyDeposit( 117 theTotalResult->ProposeLocalEnergyDeposit( 119 hadFinalSt 118 hadFinalState->GetLocalEnergyDeposit() ); 120 theTotalResult->SetNumberOfSecondaries( nu 119 theTotalResult->SetNumberOfSecondaries( numberOfSecondaries ); 121 theTotalResult->ProposeEnergy( hadFinalSta 120 theTotalResult->ProposeEnergy( hadFinalState->GetEnergyChange() ); 122 theTotalResult->ProposeTrackStatus( fAlive 121 theTotalResult->ProposeTrackStatus( fAlive ); 123 if ( hadFinalState->GetStatusChange() == s 122 if ( hadFinalState->GetStatusChange() == stopAndKill ) 124 theTotalResult->ProposeTrackStatus( fS 123 theTotalResult->ProposeTrackStatus( fStopAndKill ); 125 124 126 for ( G4int i( 0 ); i < numberOfSecondari << 125 for( G4int i( 0 ); i < numberOfSecondaries; ++i ) 127 { 126 { 128 G4double time( hadFinalState->GetSec 127 G4double time( hadFinalState->GetSecondary( i )->GetTime() ); 129 if ( time < 0 ) 128 if ( time < 0 ) 130 time = track.GetGlobalTime(); 129 time = track.GetGlobalTime(); 131 130 132 G4Track * newTrack( new G4Track( 131 G4Track * newTrack( new G4Track( 133 hadFinalState->Ge 132 hadFinalState->GetSecondary( i )->GetParticle(), 134 time, track.GetPo 133 time, track.GetPosition() ) ); 135 134 136 G4double newWeight( track.GetWeight( 135 G4double newWeight( track.GetWeight() * 137 hadFinalState->G 136 hadFinalState->GetSecondary( i )->GetWeight() ); 138 newTrack->SetWeight( newWeight ); 137 newTrack->SetWeight( newWeight ); 139 newTrack->SetTouchableHandle( track.Ge 138 newTrack->SetTouchableHandle( track.GetTouchableHandle() ); 140 theTotalResult->AddSecondary( newTrack 139 theTotalResult->AddSecondary( newTrack ); 141 } 140 } 142 141 143 hadFinalState->Clear(); 142 hadFinalState->Clear(); 144 } 143 } 145 144 146 145 147 G4VParticleChange * CexmcHadronicProcess::Pos 146 G4VParticleChange * CexmcHadronicProcess::PostStepDoIt( const G4Track & track, 148 147 const G4Step & ) 149 { 148 { 150 G4TrackStatus trackStatus( track.GetTrack 149 G4TrackStatus trackStatus( track.GetTrackStatus() ); 151 150 152 if ( trackStatus != fAlive && trackStatus 151 if ( trackStatus != fAlive && trackStatus != fSuspend ) 153 { 152 { 154 theTotalResult->Clear(); 153 theTotalResult->Clear(); 155 theTotalResult->Initialize( track ); 154 theTotalResult->Initialize( track ); 156 155 157 return theTotalResult; 156 return theTotalResult; 158 } 157 } 159 158 160 /* NB: the target nucleus is chosen only o 159 /* NB: the target nucleus is chosen only once, it means that it will always 161 * have same Z and A, practically the firs 160 * have same Z and A, practically the first stable isotope of the first 162 * element in elements vector will be chos 161 * element in elements vector will be chosen. This simplification prompts 163 * the user to choose simple single-elemen 162 * the user to choose simple single-element material for the target, for 164 * example liquid hydrogen. On the other h 163 * example liquid hydrogen. On the other hand target nucleus is supposedly 165 * only needed if user decides to turn Fer 164 * only needed if user decides to turn Fermi motion on, so this 166 * simplification should not be very harmf 165 * simplification should not be very harmful */ 167 if ( ! isInitialized ) 166 if ( ! isInitialized ) 168 { 167 { 169 CalculateTargetNucleus( track.GetMater 168 CalculateTargetNucleus( track.GetMaterial() ); 170 isInitialized = true; 169 isInitialized = true; 171 } 170 } 172 171 173 G4HadProjectile projectile( track ); 172 G4HadProjectile projectile( track ); 174 G4HadFinalState * result( interaction->Ap 173 G4HadFinalState * result( interaction->ApplyYourself( projectile, 175 174 targetNucleus ) ); 176 FillTotalResult( result, track ); 175 FillTotalResult( result, track ); 177 176 178 if ( theTotalResult->GetTrackStatus() != f 177 if ( theTotalResult->GetTrackStatus() != fStopAndKill ) 179 { 178 { 180 CexmcTrackInfo * trackInfo( static_ca 179 CexmcTrackInfo * trackInfo( static_cast< CexmcTrackInfo * >( 181 180 track.GetUserInformation() ) ); 182 181 183 if ( trackInfo && 182 if ( trackInfo && 184 trackInfo->GetTypeInfo() == Cexmc 183 trackInfo->GetTypeInfo() == CexmcIncidentParticleTrackType ) 185 { 184 { 186 CexmcIncidentParticleTrackInfo * 185 CexmcIncidentParticleTrackInfo * theTrackInfo( 187 static_cast< CexmcIncidentPart 186 static_cast< CexmcIncidentParticleTrackInfo * >( trackInfo ) ); 188 theTrackInfo->SetNeedsTrackLengthR 187 theTrackInfo->SetNeedsTrackLengthResampling(); 189 } 188 } 190 } 189 } 191 190 192 return theTotalResult; 191 return theTotalResult; 193 } 192 } 194 193 195 194 196 G4bool CexmcHadronicProcess::IsApplicable( 195 G4bool CexmcHadronicProcess::IsApplicable( 197 const 196 const G4ParticleDefinition & particle ) 198 { 197 { 199 if ( ! productionModel ) 198 if ( ! productionModel ) 200 return false; 199 return false; 201 200 202 G4ParticleDefinition * incidentParticle( 201 G4ParticleDefinition * incidentParticle( 203 production 202 productionModel->GetIncidentParticle() ); 204 203 205 if ( ! incidentParticle ) 204 if ( ! incidentParticle ) 206 return false; 205 return false; 207 206 208 return particle == *incidentParticle; 207 return particle == *incidentParticle; 209 } 208 } 210 209 211 210