Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/hadronic/models/de_excitation/photon_evaporation/src/G4GammaTransition.cc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

  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 //
 27 // -------------------------------------------------------------------
 28 //      GEANT 4 class file 
 29 //
 30 //      CERN, Geneva, Switzerland
 31 //
 32 //      File name:     G4GammaTransition
 33 //
 34 //      Author V.Ivanchenko 27 February 2015
 35 //
 36 // -------------------------------------------------------------------
 37 
 38 #include "G4GammaTransition.hh"
 39 #include "G4AtomicShells.hh"
 40 #include "Randomize.hh"
 41 #include "G4RandomDirection.hh"
 42 #include "G4Gamma.hh"
 43 #include "G4Electron.hh"
 44 #include "G4LorentzVector.hh"
 45 #include "G4SystemOfUnits.hh"
 46 #include "G4PhysicalConstants.hh"
 47 
 48 G4GammaTransition::G4GammaTransition() 
 49   : polarFlag(false), fDirection(0.,0.,0.), fTwoJMAX(10), fVerbose(0)
 50 {}
 51 
 52 G4GammaTransition::~G4GammaTransition() 
 53 {}
 54   
 55 G4Fragment* 
 56 G4GammaTransition::SampleTransition(G4Fragment* nucleus,
 57             G4double newExcEnergy,
 58             G4double mpRatio,
 59             G4int  JP1,
 60             G4int  JP2,
 61             G4int  MP,
 62             G4int  shell,
 63             G4bool isDiscrete,
 64             G4bool isGamma)
 65 {
 66   G4Fragment* result = nullptr;
 67   G4double bond_energy = 0.0;
 68 
 69   if (!isGamma) { 
 70     if(0 <= shell) {
 71       G4int Z = nucleus->GetZ_asInt();
 72       if(Z <= 104) {
 73   G4int idx = std::min(shell, G4AtomicShells::GetNumberOfShells(Z)-1);
 74   bond_energy = G4AtomicShells::GetBindingEnergy(Z, idx);
 75       }
 76     }
 77   }
 78   G4double etrans = nucleus->GetExcitationEnergy() - newExcEnergy 
 79     - bond_energy;
 80   if(fVerbose > 2) {
 81     G4cout << "G4GammaTransition::GenerateGamma - Etrans(MeV)= " 
 82      << etrans << "  Eexnew= " << newExcEnergy 
 83      << " Ebond= " << bond_energy << G4endl;
 84   } 
 85   if(etrans <= 0.0) { 
 86     etrans += bond_energy;
 87     bond_energy = 0.0;
 88   }
 89   
 90   // Do complete Lorentz computation 
 91   G4LorentzVector lv = nucleus->GetMomentum();
 92   G4double mass = nucleus->GetGroundStateMass() + newExcEnergy;
 93 
 94   // select secondary
 95   G4ParticleDefinition* part;
 96 
 97   if(isGamma) { part = G4Gamma::Gamma(); }
 98   else {
 99     part = G4Electron::Electron();
100     G4int ne = std::max(nucleus->GetNumberOfElectrons() - 1, 0);
101     nucleus->SetNumberOfElectrons(ne);
102   }
103 
104   if(polarFlag && isDiscrete && JP1 <= fTwoJMAX) {
105     SampleDirection(nucleus, mpRatio, JP1, JP2, MP);
106   } else {
107     fDirection = G4RandomDirection();
108   }
109 
110   G4double emass = part->GetPDGMass();
111 
112   // 2-body decay in rest frame
113   G4double ecm       = lv.mag();
114   G4ThreeVector bst  = lv.boostVector();
115   if(!isGamma) { ecm += (CLHEP::electron_mass_c2 - bond_energy); }
116 
117   //G4cout << "Ecm= " << ecm << " mass= " << mass << " emass= " << emass << G4endl;
118 
119   ecm = std::max(ecm, mass + emass);
120   G4double energy = 0.5*((ecm - mass)*(ecm + mass) + emass*emass)/ecm;
121   G4double mom = (emass > 0.0) ? std::sqrt((energy - emass)*(energy + emass))
122     : energy;
123 
124   // emitted gamma or e-
125   G4LorentzVector res4mom(mom * fDirection.x(),
126         mom * fDirection.y(),
127         mom * fDirection.z(), energy);
128   // residual
129   energy = std::max(ecm - energy, mass);
130   lv.set(-mom*fDirection.x(), -mom*fDirection.y(), -mom*fDirection.z(), energy);
131 
132   // Lab system transform for short lived level
133   lv.boost(bst);
134 
135   // modified primary fragment 
136   nucleus->SetExcEnergyAndMomentum(newExcEnergy, lv);
137 
138   // gamma or e- are produced
139   res4mom.boost(bst); 
140   result = new G4Fragment(res4mom, part);
141 
142   //G4cout << " DeltaE= " << e0 - lv.e() - res4mom.e() + emass
143   //   << "   Emass= " << emass << G4endl;
144   if(fVerbose > 2) {
145     G4cout << "G4GammaTransition::SampleTransition : " << *result << G4endl;
146     G4cout << "       Left nucleus: " << *nucleus << G4endl;
147   }
148   return result;
149 }
150 
151 void G4GammaTransition::SampleDirection(G4Fragment* nuc, G4double ratio, 
152           G4int twoJ1, G4int twoJ2, G4int mp)
153 {
154   G4double cosTheta, phi;
155   G4NuclearPolarization* np = nuc->GetNuclearPolarization(); 
156   if(fVerbose > 2) {
157     G4cout << "G4GammaTransition::SampleDirection : 2J1= " << twoJ1 
158      << " 2J2= " << twoJ2 << " ratio= " << ratio 
159      << " mp= " << mp << G4endl;
160     G4cout << "  Nucleus: " << *nuc << G4endl;
161   }
162   if(nullptr == np) {
163     cosTheta = 2*G4UniformRand() - 1.0;
164     phi = CLHEP::twopi*G4UniformRand();
165 
166   } else {
167     // PhotonEvaporation dataset:
168     // The multipolarity number with 1,2,3,4,5,6,7 representing E0,E1,M1,E2,M2,E3,M3
169     // monopole transition and 100*Nx+Ny representing multipolarity transition with
170     // Ny and Ny taking the value 1,2,3,4,5,6,7 referring to E0,E1,M1,E2,M2,E3,M3,..
171     // For example a M1+E2 transition would be written 304.
172     // M1 is the primary transition (L) and E2 is the secondary (L')
173 
174     G4double mpRatio = ratio;
175 
176     G4int L0 = 0, Lp = 0;
177     if (mp > 99) {
178       L0 = mp/200;
179       Lp = (mp%100)/2;
180     } else {
181       L0 = mp/2;
182       Lp = 0;
183       mpRatio = 0.;
184     } 
185     fPolTrans.SampleGammaTransition(np, twoJ1, twoJ2, L0, Lp, mpRatio, cosTheta, phi);
186   }
187 
188   G4double sinTheta = std::sqrt((1.-cosTheta)*(1.+cosTheta));
189   fDirection.set(sinTheta*std::cos(phi),sinTheta*std::sin(phi),cosTheta);
190   if(fVerbose > 3) {
191     G4cout << "G4GammaTransition::SampleDirection done: " << fDirection << G4endl;
192     if(np) { G4cout << *np << G4endl; }
193   }
194 }
195