| 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 * model )
75 { 74 {
76 productionModel = model; 75 productionModel = model;
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