source: trunk/source/processes/hadronic/models/binary_cascade/src/G4GeneratorPrecompoundInterface.cc@ 1014

Last change on this file since 1014 was 819, checked in by garnier, 17 years ago

import all except CVS
File size: 7.7 KB
Line 
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 "G4GeneratorPrecompoundInterface.hh"
27#include "G4DynamicParticleVector.hh"
28#include "G4IonTable.hh"
29
30//
31// HPW, 10DEC 98, the decay part originally written by Gunter Folger in his FTF-test-program.
32//
33
34 G4HadFinalState* G4GeneratorPrecompoundInterface::
35 ApplyYourself(const G4HadProjectile &, G4Nucleus & )
36 {
37 std::cout << "G4GeneratorPrecompoundInterface: ApplyYourself interface called stand-allone."<< G4endl;
38 std::cout << "This class is only a mediator between generator and precompound"<<G4endl;
39 std::cout << "Please remove from your physics list."<<G4endl;
40 throw G4HadronicException(__FILE__, __LINE__, "SEVERE: G4GeneratorPrecompoundInterface model interface called stand-allone.");
41 return new G4HadFinalState;
42 }
43
44 G4ReactionProductVector* G4GeneratorPrecompoundInterface::
45 Propagate(G4KineticTrackVector* theSecondaries, G4V3DNucleus* theNucleus)
46 {
47 G4ReactionProductVector * theTotalResult = new G4ReactionProductVector;
48
49 // decay the strong resonances
50 G4KineticTrackVector *result1, *secondaries, *result;
51 result1=theSecondaries;
52 result=new G4KineticTrackVector();
53
54 for (unsigned int aResult=0; aResult < result1->size(); aResult++)
55 {
56 G4ParticleDefinition * pdef;
57 pdef=result1->operator[](aResult)->GetDefinition();
58 secondaries=NULL;
59 if ( pdef->IsShortLived() )
60 {
61 secondaries = result1->operator[](aResult)->Decay();
62 }
63 if ( secondaries == NULL )
64 {
65 result->push_back(result1->operator[](aResult));
66 result1->operator[](aResult)=NULL; //protect for clearAndDestroy
67 }
68 else
69 {
70 for (unsigned int aSecondary=0; aSecondary<secondaries->size(); aSecondary++)
71 {
72 result1->push_back(secondaries->operator[](aSecondary));
73 }
74 delete secondaries;
75 }
76 }
77 std::for_each(result1->begin(), result1->end(), DeleteKineticTrack());
78 delete result1;
79
80
81
82 // prepare the fragment
83 G4Fragment anInitialState;
84 G4int anA=theNucleus->GetMassNumber();
85 G4int aZ=theNucleus->GetCharge();
86 G4int numberOfEx = 0;
87 G4int numberOfCh = 0;
88 G4int numberOfHoles = 0;
89 G4double exEnergy = 0;
90 G4ThreeVector exciton3Momentum(0,0,0);
91 // loop over secondaries
92 for(unsigned int list=0; list < result->size(); list++)
93 {
94 G4KineticTrack *aTrack = result->operator[](list);
95 if(aTrack->GetDefinition() != G4Proton::Proton() &&
96 aTrack->GetDefinition() != G4Neutron::Neutron())
97 {
98 G4ReactionProduct * theNew = new G4ReactionProduct(aTrack->GetDefinition());
99 theNew->SetMomentum(aTrack->Get4Momentum().vect());
100 theNew->SetTotalEnergy(aTrack->Get4Momentum().e());
101 theTotalResult->push_back(theNew);
102 }
103 else if(aTrack->Get4Momentum().t() - aTrack->Get4Momentum().mag()>80*MeV)
104 {
105 G4ReactionProduct * theNew = new G4ReactionProduct(aTrack->GetDefinition());
106 theNew->SetMomentum(aTrack->Get4Momentum().vect());
107 theNew->SetTotalEnergy(aTrack->Get4Momentum().e());
108 theTotalResult->push_back(theNew);
109 }
110 else if(aTrack->GetPosition().mag() > theNucleus->GetNuclearRadius())
111 {
112 G4ReactionProduct * theNew = new G4ReactionProduct(aTrack->GetDefinition());
113 theNew->SetMomentum(aTrack->Get4Momentum().vect());
114 theNew->SetTotalEnergy(aTrack->Get4Momentum().e());
115 theTotalResult->push_back(theNew);
116 }
117 else
118 {
119 // within the nucleus, neutron or proton
120 // now calculate A, Z of the fragment, momentum, number of exciton states
121 anA++;;
122 numberOfEx++;
123 aZ += G4int(aTrack->GetDefinition()->GetPDGCharge());
124 numberOfCh += G4int(aTrack->GetDefinition()->GetPDGCharge());
125 exciton3Momentum += aTrack->Get4Momentum().vect();
126 exEnergy += (aTrack->Get4Momentum().t()-aTrack->Get4Momentum().m());
127 }
128 }
129
130 // loop over wounded nucleus
131 G4Nucleon * theCurrentNucleon = theNucleus->StartLoop() ? theNucleus->GetNextNucleon() : NULL;
132 while(theCurrentNucleon != NULL)
133 {
134 if(theCurrentNucleon->AreYouHit())
135 {
136 numberOfHoles++;
137 numberOfEx++;
138 anA--;
139 aZ -= G4int(theCurrentNucleon->GetDefinition()->GetPDGCharge());
140 exciton3Momentum -= theCurrentNucleon->Get4Momentum().vect();
141 exEnergy+=theCurrentNucleon->GetBindingEnergy();
142 }
143 theCurrentNucleon = theNucleus->GetNextNucleon();
144 }
145
146 if(!theDeExcitation)
147 {
148 // throw G4HadronicException(__FILE__, __LINE__, "Please register an evaporation phase with G4GeneratorPrecompoundInterface.");
149 }
150 else if(0!=anA && 0!=aZ)
151 {
152 G4double residualMass =
153 G4ParticleTable::GetParticleTable()->GetIonTable()->GetIonMass(aZ ,anA);
154 residualMass += exEnergy;
155 G4LorentzVector exciton4Momentum(exciton3Momentum,
156 std::sqrt(exciton3Momentum.mag2()+residualMass*residualMass));
157
158 anInitialState.SetA(anA);
159 anInitialState.SetZ(aZ);
160 anInitialState.SetNumberOfParticles(numberOfEx-numberOfHoles);
161 anInitialState.SetNumberOfCharged(numberOfCh);
162 anInitialState.SetNumberOfHoles(numberOfHoles);
163 anInitialState.SetMomentum(exciton4Momentum);
164// anInitialState.SetExcitationEnergy(exEnergy); // now a redundant call.
165
166 // call pre-compound
167 const G4Fragment aFragment(anInitialState);
168 G4ReactionProductVector * aPreResult = theDeExcitation->DeExcite(aFragment);
169// G4ReactionProductVector * aPreResult = new G4ReactionProductVector;
170 // fill pre-compound part into the result, and return
171 for(unsigned int ll=0; ll<aPreResult->size(); ll++)
172 {
173 theTotalResult->push_back(aPreResult->operator[](ll));
174 }
175 delete aPreResult;
176 }
177 else
178 {
179 // throw G4HadronicException(__FILE__, __LINE__, "Please register an evaporation phase with G4GeneratorPrecompoundInterface.");
180 }
181 // now return
182
183 std::for_each(result->begin(), result->end(), DeleteKineticTrack());
184 delete result;
185 return theTotalResult;
186 }
187
Note: See TracBrowser for help on using the repository browser.