00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027 #include "G4LENDCapture.hh"
00028 #include "G4SystemOfUnits.hh"
00029 #include "G4Nucleus.hh"
00030 #include "G4ParticleTable.hh"
00031
00032 G4HadFinalState * G4LENDCapture::ApplyYourself(const G4HadProjectile& aTrack, G4Nucleus& aTarg )
00033 {
00034
00035 G4double temp = aTrack.GetMaterial()->GetTemperature();
00036
00037
00038
00039
00040 G4int iZ = aTarg.GetZ_asInt();
00041 G4int iA = aTarg.GetA_asInt();
00042
00043 G4double ke = aTrack.GetKineticEnergy();
00044
00045 G4HadFinalState* theResult = &theParticleChange;
00046 theResult->Clear();
00047
00048 G4GIDI_target* aTarget = usedTarget_map.find( lend_manager->GetNucleusEncoding( iZ , iA ) )->second->GetTarget();
00049 std::vector<G4GIDI_Product>* products = aTarget->getCaptureFinalState( ke*MeV, temp, NULL, NULL );
00050
00051
00052 if ( products != NULL )
00053 {
00054
00055 G4ThreeVector p(0,0,0);
00056 G4int totN = 0;
00057
00058 for ( G4int j = 0; j < int( products->size() ); j++ )
00059 {
00060 G4int jZ = (*products)[j].Z;
00061 G4int jA = (*products)[j].A;
00062
00063
00064
00065
00066
00067
00068
00069
00070 G4ThreeVector dp((*products)[j].px,(*products)[j].py,(*products)[j].pz);
00071 p += dp;
00072
00073 G4DynamicParticle* theSec = new G4DynamicParticle;
00074
00075 if ( jA == 1 && jZ == 1 )
00076 {
00077 theSec->SetDefinition( G4Proton::Proton() );
00078 totN += 1;
00079 }
00080 else if ( jA == 1 && jZ == 0 )
00081 {
00082 theSec->SetDefinition( G4Neutron::Neutron() );
00083 totN += 1;
00084 }
00085 else if ( jZ > 0 )
00086 {
00087 if ( jA != 0 )
00088 {
00089 theSec->SetDefinition( G4ParticleTable::GetParticleTable()->FindIon( jZ , jA , 0 , 0 ) );
00090 totN += jA;
00091 }
00092 else
00093 {
00094 theSec->SetDefinition( G4ParticleTable::GetParticleTable()->FindIon( jZ , iA+1-totN , 0 , 0 ) );
00095 }
00096 }
00097 else
00098 {
00099 theSec->SetDefinition( G4Gamma::Gamma() );
00100 }
00101
00102 theSec->SetMomentum( G4ThreeVector( (*products)[j].px*MeV , (*products)[j].py*MeV , (*products)[j].pz*MeV ) );
00103
00104 if ( dp.mag() == 0 )
00105 {
00106 theSec->SetMomentum( -p*MeV );
00107 }
00108
00109 theResult->AddSecondary( theSec );
00110 }
00111 }
00112 delete products;
00113
00114 theResult->SetStatusChange( stopAndKill );
00115
00116 return theResult;
00117
00118 }