G4VGammaDeexcitation.cc

Go to the documentation of this file.
00001 //
00002 // ********************************************************************
00003 // * License and Disclaimer                                           *
00004 // *                                                                  *
00005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
00006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
00007 // * conditions of the Geant4 Software License,  included in the file *
00008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
00009 // * include a list of copyright holders.                             *
00010 // *                                                                  *
00011 // * Neither the authors of this software system, nor their employing *
00012 // * institutes,nor the agencies providing financial support for this *
00013 // * work  make  any representation or  warranty, express or implied, *
00014 // * regarding  this  software system or assume any liability for its *
00015 // * use.  Please see the license in the file  LICENSE  and URL above *
00016 // * for the full disclaimer and the limitation of liability.         *
00017 // *                                                                  *
00018 // * This  code  implementation is the result of  the  scientific and *
00019 // * technical work of the GEANT4 collaboration.                      *
00020 // * By using,  copying,  modifying or  distributing the software (or *
00021 // * any work based  on the software)  you  agree  to acknowledge its *
00022 // * use  in  resulting  scientific  publications,  and indicate your *
00023 // * acceptance of all terms of the Geant4 Software license.          *
00024 // ********************************************************************
00025 //
00026 // $Id$
00027 //
00028 // -------------------------------------------------------------------
00029 //      GEANT 4 class file 
00030 //
00031 //      CERN, Geneva, Switzerland
00032 //
00033 //      File name:     G4VGammaDeexcitation
00034 //
00035 //      Author:        Maria Grazia Pia (pia@genova.infn.it)
00036 // 
00037 //      Creation date: 23 October 1998
00038 //
00039 //      Modifications: 
00040 //
00041 //        21 Nov 2001, Fan Lei (flei@space.qinetiq.com)
00042 //           Modified GenerateGamma() and UpdateUncleus() for implementation
00043 //           of Internal Conversion processs 
00044 //      
00045 //        15 April 1999, Alessandro Brunengo (Alessandro.Brunengo@ge.infn.it)
00046 //              Added creation time evaluation for products of evaporation
00047 //
00048 //        19 April 2010, J. M. Quesada calculations in CM system
00049 //              pending final boost to lab system  (not critical)
00050 //
00051 //        23 April 2010, V.Ivanchenko rewite kinematic part using PDG formula
00052 //                                    for 2-body decay
00053 //
00054 //        07 May   2011, V.Ivanchenko implement check ICM flag - produce or not e-
00055 //
00056 // -------------------------------------------------------------------
00057 
00058 #include "G4VGammaDeexcitation.hh"
00059 
00060 #include "globals.hh"
00061 #include "G4PhysicalConstants.hh"
00062 #include "Randomize.hh"
00063 #include "G4Gamma.hh"
00064 #include "G4Electron.hh"
00065 #include "G4LorentzVector.hh"
00066 #include "G4VGammaTransition.hh"
00067 #include "G4Fragment.hh"
00068 #include "G4FragmentVector.hh"
00069 
00070 #include "G4ParticleTable.hh"
00071 #include "G4IonTable.hh"
00072 
00073 #include "G4DiscreteGammaTransition.hh"
00074 
00075 
00076 G4VGammaDeexcitation::G4VGammaDeexcitation(): _transition(0), _verbose(0),
00077                                               _electronO (0), _vSN(-1)
00078 { 
00079   _nucleus = 0;
00080   fTimeLimit = DBL_MAX;
00081 }
00082 
00083 G4VGammaDeexcitation::~G4VGammaDeexcitation()
00084 { 
00085   if (_transition != 0) { delete _transition; }
00086 }
00087 
00088 G4FragmentVector* G4VGammaDeexcitation::DoTransition()
00089 {
00090   Initialize();
00091   G4FragmentVector* products = new G4FragmentVector();
00092  
00093   if (CanDoTransition())
00094     {
00095       G4Fragment* gamma = GenerateGamma();
00096       if (gamma != 0) { products->push_back(gamma); }
00097     }
00098  
00099   if (_verbose > 1) {
00100     G4cout << "G4VGammaDeexcitation::DoTransition - Transition deleted " << G4endl;
00101   }
00102  
00103   return products;
00104 }
00105 
00106 G4FragmentVector* G4VGammaDeexcitation::DoChain()
00107 {
00108   if (_verbose > 1) { G4cout << "G4VGammaDeexcitation::DoChain" << G4endl; }
00109   const G4double tolerance = CLHEP::keV;
00110 
00111   Initialize();
00112   G4FragmentVector* products = new G4FragmentVector();
00113   
00114   while (CanDoTransition())
00115     {      
00116       _transition->SetEnergyFrom(_nucleus->GetExcitationEnergy());
00117       G4Fragment* gamma = GenerateGamma();
00118       if (gamma != 0) 
00119         {
00120           products->push_back(gamma);
00121           //G4cout << "Eex(keV)= " << _nucleus->GetExcitationEnergy()/keV << G4endl;
00122           if(_nucleus->GetExcitationEnergy() <= tolerance) { break; }
00123           Update();
00124         }
00125     } 
00126   
00127   if (_verbose > 1) {
00128     G4cout << "G4VGammaDeexcitation::DoChain - Transition deleted, end of chain " << G4endl;
00129   }
00130   
00131   return products;
00132 }
00133 
00134 G4Fragment* G4VGammaDeexcitation::GenerateGamma()
00135 {
00136   // 23/04/10 V.Ivanchenko rewrite complitely
00137   G4double eGamma = 0.;
00138   
00139   if (_transition) {
00140     _transition->SelectGamma();  // it can be conversion electron too
00141     eGamma = _transition->GetGammaEnergy(); 
00142     //G4cout << "G4VGammaDeexcitation::GenerateGamma - Egam(MeV)= " 
00143     //     << eGamma << G4endl; 
00144     if(eGamma <= 0.0) { return 0; }
00145   } else { return 0; }
00146 
00147   G4double excitation = _nucleus->GetExcitationEnergy() - eGamma;
00148   if(excitation < 0.0) { excitation = 0.0; } 
00149   if (_verbose > 1) 
00150     {
00151       G4cout << "G4VGammaDeexcitation::GenerateGamma - Edeexc(MeV)= " << eGamma 
00152              << " ** left Eexc(MeV)= " << excitation
00153              << G4endl;
00154     }
00155 
00156   G4double gammaTime = _transition->GetGammaCreationTime();
00157   
00158   // Do complete Lorentz computation 
00159   G4LorentzVector lv = _nucleus->GetMomentum();
00160   G4double Mass = _nucleus->GetGroundStateMass() + excitation;
00161 
00162   // select secondary
00163   G4ParticleDefinition* gamma = G4Gamma::Gamma();
00164 
00165   G4DiscreteGammaTransition* dtransition = 
00166     dynamic_cast <G4DiscreteGammaTransition*> (_transition);
00167 
00168   G4bool eTransition = false;
00169   if (dtransition && !( dtransition->IsAGamma()) ) {
00170     eTransition = true; 
00171     gamma = G4Electron::Electron(); 
00172     _vSN = dtransition->GetOrbitNumber();   
00173     _electronO.RemoveElectron(_vSN);
00174     lv += G4LorentzVector(0.0,0.0,0.0,CLHEP::electron_mass_c2 - dtransition->GetBondEnergy());
00175   }
00176 
00177   G4double cosTheta = 1. - 2. * G4UniformRand(); 
00178   G4double sinTheta = std::sqrt(1. - cosTheta * cosTheta);
00179   G4double phi = twopi * G4UniformRand();
00180 
00181   G4double eMass = gamma->GetPDGMass();
00182   G4LorentzVector Gamma4P;
00183 
00184   //G4cout << "Egamma= " << eGamma << " Mass= " << eMass << " t= " << gammaTime
00185   //     << " tlim= " << fTimeLimit << G4endl;
00186 
00187   if(gammaTime > fTimeLimit) {
00188     // shortcut for long lived levels
00189     // not correct position of stopping ion gamma emission
00190     // 4-momentum balance is breaked
00191     G4double e = eGamma + eMass;
00192     G4double mom = std::sqrt(eGamma*(eGamma + 2*eMass));
00193     Gamma4P.set(mom * sinTheta * std::cos(phi),
00194                 mom * sinTheta * std::sin(phi),
00195                 mom * cosTheta, e); 
00196     lv -= Gamma4P;
00197     e = lv.e();
00198     if(e < Mass) { e = Mass; }
00199     mom = std::sqrt((e - Mass)*(e + Mass));
00200     G4ThreeVector v = lv.vect().unit();
00201     lv.set(mom*v.x(), mom*v.y(), mom*v.z(), e); 
00202 
00203   } else {
00204     // 2-body decay in rest frame
00205     G4double Ecm       = lv.mag();
00206     G4ThreeVector bst  = lv.boostVector();
00207 
00208     G4double GammaEnergy = 0.5*((Ecm - Mass)*(Ecm + Mass) + eMass*eMass)/Ecm;
00209     if(GammaEnergy <= eMass) { return 0; }
00210 
00211     G4double mom = std::sqrt((GammaEnergy - eMass)*(GammaEnergy + eMass));
00212     Gamma4P.set(mom * sinTheta * std::cos(phi),
00213                 mom * sinTheta * std::sin(phi),
00214                 mom * cosTheta,
00215                 GammaEnergy);
00216 
00217     Gamma4P.boost(bst);  
00218     lv -= Gamma4P;
00219   }
00220 
00221   // modified primary fragment 
00222   gammaTime += _nucleus->GetCreationTime();
00223 
00224   _nucleus->SetMomentum(lv);
00225   _nucleus->SetCreationTime(gammaTime);
00226 
00227   // e- is not produced
00228   if(eTransition && !dtransition->GetICM()) { return 0; }
00229 
00230   // gamma or e- are produced
00231   G4Fragment * thePhoton = new G4Fragment(Gamma4P,gamma);
00232   thePhoton->SetCreationTime(gammaTime);
00233 
00234   //G4cout << "G4VGammaDeexcitation::GenerateGamma : " << thePhoton << G4endl;
00235   //G4cout << "       Left nucleus: " << _nucleus << G4endl;
00236   return thePhoton;
00237 }
00238 
00239 void G4VGammaDeexcitation::Update()
00240 {
00241   if (_transition !=  0) 
00242     { 
00243       delete _transition;
00244       _transition = 0;
00245       if (_verbose > 1) {
00246         G4cout << "G4VGammaDeexcitation::Update - Transition deleted " << G4endl;
00247       }
00248     }
00249   
00250   _transition = CreateTransition();
00251   if (_transition != 0) 
00252     {
00253       _transition->SetEnergyFrom(_nucleus->GetExcitationEnergy());
00254       // if ( _vSN != -1) (dynamic_cast <G4DiscreteGammaTransition*> (_transition))->SetICM(false);
00255       // the above line is commented out for bug fix #952. It was intruduced for reason that
00256       // the k-shell electron is most likely one to be kicked out and there is no time for 
00257       // the atom to deexcite before the next IC. But this limitation is causing other problems as 
00258       // reported in #952
00259     }
00260   
00261   return;
00262 }

Generated on Mon May 27 17:50:13 2013 for Geant4 by  doxygen 1.4.7