Ignore:
Timestamp:
Apr 6, 2009, 12:21:12 PM (15 years ago)
Author:
garnier
Message:

update processes

Location:
trunk/source/processes/electromagnetic/highenergy/src
Files:
12 edited

Legend:

Unmodified
Added
Removed

  • trunk/source/processes/electromagnetic/highenergy/src/G4AnnihiToMuPair.cc

    r819 r961  
    2525//
    2626//
    27 // $Id: G4AnnihiToMuPair.cc,v 1.3 2006/06/29 19:32:34 gunter Exp $
    28 // GEANT4 tag $Name: geant4-09-01-patch-02 $
     27// $Id: G4AnnihiToMuPair.cc,v 1.5 2008/10/16 14:29:48 vnivanch Exp $
     28// GEANT4 tag $Name: geant4-09-02-ref-02 $
    2929//
    3030//         ------------ G4AnnihiToMuPair physics process ------
    … …  
    6666 
    6767 CrossSecFactor = 1.;
     68 SetProcessSubType(6);
     69
    6870}
    6971
    … …  
    244246void G4AnnihiToMuPair::PrintInfoDefinition()
    245247{
    246   G4String comments ="e+e->mu+mu- annihilation, atomic e- at rest.\n";
    247   G4cout << G4endl << GetProcessName() << ":  " << comments
    248          << "        threshold at " << LowestEnergyLimit/GeV << " GeV"
     248  G4String comments ="e+e->mu+mu- annihilation, atomic e- at rest, SubType=.";
     249  G4cout << G4endl << GetProcessName() << ":  " << comments
     250         << GetProcessSubType() << G4endl;
     251  G4cout << "        threshold at " << LowestEnergyLimit/GeV << " GeV"
    249252         << " good description up to "
    250253         << HighestEnergyLimit/TeV << " TeV for all Z." << G4endl;

  • trunk/source/processes/electromagnetic/highenergy/src/G4BetheBlochNoDeltaModel.cc

    r819 r961  
    2424// ********************************************************************
    2525//
    26 // $Id: G4BetheBlochNoDeltaModel.cc,v 1.3 2006/06/29 19:32:36 gunter Exp $
    27 // GEANT4 tag $Name: $
     26// $Id: G4BetheBlochNoDeltaModel.cc,v 1.4 2009/02/20 16:38:33 vnivanch Exp $
     27// GEANT4 tag $Name: geant4-09-02-ref-02 $
    2828//
    2929// -------------------------------------------------------------------
    … …  
    5858{}
    5959
     60//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     61
    6062G4BetheBlochNoDeltaModel::~G4BetheBlochNoDeltaModel()
    6163{}
    … …  
    6365//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
    6466
     67G4double G4BetheBlochNoDeltaModel::ComputeDEDXPerVolume(
     68                            const G4Material* material,
     69                            const G4ParticleDefinition* pd,
     70                            G4double kinEnergy, G4double)
     71{
     72  return
     73    G4BetheBlochModel::ComputeDEDXPerVolume(material, pd, kinEnergy, DBL_MAX);
     74}
    6575
     76//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     77
     78G4double G4BetheBlochNoDeltaModel::CrossSectionPerVolume(
     79                    const G4Material*,const G4ParticleDefinition*,
     80                    G4double, G4double, G4double)
     81{
     82  return 0.0;
     83}
     84
     85//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     86
     87

  • trunk/source/processes/electromagnetic/highenergy/src/G4BraggNoDeltaModel.cc

    r819 r961  
    2424// ********************************************************************
    2525//
    26 // $Id: G4BraggNoDeltaModel.cc,v 1.3 2006/06/29 19:32:38 gunter Exp $
    27 // GEANT4 tag $Name: $
     26// $Id: G4BraggNoDeltaModel.cc,v 1.4 2009/02/20 16:38:33 vnivanch Exp $
     27// GEANT4 tag $Name: geant4-09-02-ref-02 $
    2828//
    2929// -------------------------------------------------------------------
    … …  
    5353//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
    5454
    55 
    5655G4BraggNoDeltaModel::G4BraggNoDeltaModel(const G4ParticleDefinition*p,
    5756                                         const G4String& nam) :
    5857  G4BraggIonModel(p, nam)
    5958{}
     59
     60//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
    6061
    6162G4BraggNoDeltaModel::~G4BraggNoDeltaModel()
    … …  
    6465//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
    6566
     67G4double G4BraggNoDeltaModel::ComputeDEDXPerVolume(
     68                            const G4Material* material,
     69                            const G4ParticleDefinition* pd,
     70                            G4double kinEnergy, G4double)
     71{
     72  return
     73    G4BraggIonModel::ComputeDEDXPerVolume(material, pd, kinEnergy, DBL_MAX);
     74}
    6675
     76//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     77
     78G4double G4BraggNoDeltaModel::CrossSectionPerVolume(
     79                            const G4Material*,
     80                            const G4ParticleDefinition*,
     81                            G4double, G4double, G4double)
     82{
     83  return 0.0;
     84}
     85
     86//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     87
     88

  • trunk/source/processes/electromagnetic/highenergy/src/G4GammaConversionToMuons.cc

    r819 r961  
    2525//
    2626//
    27 // $Id: G4GammaConversionToMuons.cc,v 1.4 2006/06/29 19:32:40 gunter Exp $
    28 // GEANT4 tag $Name: geant4-09-01-patch-02 $
     27// $Id: G4GammaConversionToMuons.cc,v 1.7 2008/10/16 14:29:48 vnivanch Exp $
     28// GEANT4 tag $Name: geant4-09-02-ref-02 $
    2929//
    3030//         ------------ G4GammaConversionToMuons physics process ------
    … …  
    5151    HighestEnergyLimit(1e21*eV), // ok to 1e21eV=1e12GeV, then LPM suppression
    5252    CrossSecFactor(1.)
    53 { }
     53{
     54  SetProcessSubType(15);
     55}
    5456
    5557//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
    … …  
    261263    G4double xxp=1.-4./3.*xPM; // the main xPlus dependence
    262264    result=xxp*log(W)*LogWmaxInv;
    263     if(result>1.)
    264     { G4cout << "error in dSigxPlusGen, result=" << result << " is >1" << '\n';
    265       exit(10);
     265    if(result>1.) {
     266      G4cout << "G4GammaConversionToMuons::PostStepDoIt WARNING:"
     267             << " in dSigxPlusGen, result=" << result << " > 1" << G4endl;
    266268    }
    267269  }
    … …  
    285287      f1=(1.-2.*xPM+4.*xPM*t*(1.-t)) / (1.+C1/(t*t));
    286288      if(f1<0 || f1> f1_max) // should never happend
    287       { G4cout << "outside allowed range f1=" << f1 << G4endl;
    288         exit(1);
    289       }
     289        {
     290          G4cout << "G4GammaConversionToMuons::PostStepDoIt WARNING:"
     291                 << "outside allowed range f1=" << f1 << " is set to zero"
     292                 << G4endl;
     293          f1 = 0.0;
     294        }
    290295    }
    291296    while ( G4UniformRand()*f1_max > f1);
    … …  
    299304      f2=1.-2.*xPM+4.*xPM*t*(1.-t)*(1.+cos(2.*psi));
    300305      if(f2<0 || f2> f2_max) // should never happend
    301       { G4cout << "outside allowed range f2=" << f2 << G4endl;
    302         exit(1);
    303       }
     306        {
     307          G4cout << "G4GammaConversionToMuons::PostStepDoIt WARNING:"
     308                 << "outside allowed range f2=" << f2 << " is set to zero"
     309                 << G4endl;
     310          f2 = 0.0;
     311        }
    304312    }
    305313    while ( G4UniformRand()*f2_max > f2);
    … …  
    387395void G4GammaConversionToMuons::PrintInfoDefinition()
    388396{
    389   G4String comments ="gamma->mu+mu- Bethe Heitler process.\n";
     397  G4String comments ="gamma->mu+mu- Bethe Heitler process, SubType= ";
    390398  G4cout << G4endl << GetProcessName() << ":  " << comments
    391          << "        good cross section parametrization from "
     399         << GetProcessSubType() << G4endl;
     400  G4cout << "        good cross section parametrization from "
    392401         << G4BestUnit(LowestEnergyLimit,"Energy")
    393402         << " to " << HighestEnergyLimit/GeV << " GeV for all Z." << G4endl;

  • trunk/source/processes/electromagnetic/highenergy/src/G4eeCrossSections.cc

    r819 r961  
    2424// ********************************************************************
    2525//
    26 // $Id: G4eeCrossSections.cc,v 1.6 2006/06/29 19:32:42 gunter Exp $
    27 // GEANT4 tag $Name: $
     26// $Id: G4eeCrossSections.cc,v 1.7 2008/07/10 18:06:39 vnivanch Exp $
     27// GEANT4 tag $Name: geant4-09-02-ref-02 $
    2828//
    2929// -------------------------------------------------------------------
    … …  
    3939//
    4040// Modifications:
    41 //
     41// 10.07.2008 Updated for PDG Jour. Physics, G33, 1 (2006)
    4242//
    4343// -------------------------------------------------------------------
    … …  
    5252#include "G4PionMinus.hh"
    5353#include "G4PionZero.hh"
     54#include "G4Eta.hh"
    5455#include "G4KaonPlus.hh"
    5556#include "G4KaonMinus.hh"
    … …  
    8081  MsPi = G4PionPlus::PionPlus()->GetPDGMass();
    8182  MsPi0= G4PionZero::PionZero()->GetPDGMass();
    82   MsEta= 547.30*MeV;
     83  MsEta= G4Eta::Eta()->GetPDGMass();
    8384  MsEtap=957.78*MeV;
    8485  MsKs = G4KaonZeroLong::KaonZeroLong()->GetPDGMass();
    85   MsKc=G4KaonPlus::KaonPlus()->GetPDGMass();
    86   MsRho= 770.0*MeV;
    87   MsOm = 781.94*MeV;
     86  MsKc = G4KaonPlus::KaonPlus()->GetPDGMass();
     87  MsRho= 775.5*MeV;
     88  MsOm = 782.62*MeV;
    8889  MsF0 = 980.0*MeV;
    89   MsA0 = 983.4*MeV;
    90   MsPhi= 1019.413*MeV;
     90  MsA0 = 984.7*MeV;
     91  MsPhi= 1019.46*MeV;
    9192  MsK892 = 891.66*MeV;
    92   MsK0892 = 896.10*MeV;
    93   GRho = 150.7*MeV;
    94   GOm = 8.41*MeV;
    95   GPhi = 4.43*MeV;
     93  MsK0892 = 896.0*MeV;
     94  GRho = 149.4*MeV;
     95  GOm = 8.49*MeV;
     96  GPhi = 4.26*MeV;
    9697  GK892 = 50.8*MeV;
    97   GK0892 = 50.5*MeV;
     98  GK0892 = 50.3*MeV;
    9899  PhRho = 0.0;
    99100  PhOm = 0.0;
    … …  
    103104  BrRhoPiG = 4.5e-4;
    104105  BrRhoPi0G= 6.8e-4;
    105   BrRhoEtaG= 2.4e-4;
    106   BrRhoEe = 4.49e-5;
    107   BrOm3Pi = 0.888;
    108   BrOmPi0G= 0.085;
    109   BrOmEtaG= 6.5e-4;
    110   BrOm2Pi = 0.0221;
     106  BrRhoEtaG= 2.95e-4;
     107  BrRhoEe = 4.7e-5;
     108  BrOm3Pi = 0.891;
     109  BrOmPi0G= 0.089;
     110  BrOmEtaG= 4.9e-4;
     111  BrOm2Pi = 0.017;
    111112  PhOm2Pi = 90.0;
    112   BrOmEe = 7.07e-5;
    113   BrPhi2Kc = 0.491;
    114   BrPhiKsKl= 0.341;
    115   BrPhi3Pi = 0.155;
    116   BrPhiPi0G= 1.31e-3;
    117   BrPhiEtaG= 1.26e-2;
    118   BrPhi2Pi = 8.e-5;
     113  BrOmEe = 7.18e-5;
     114  BrPhi2Kc = 0.492;
     115  BrPhiKsKl= 0.34;
     116  BrPhi3Pi = 0.153;
     117  BrPhiPi0G= 1.25e-3;
     118  BrPhiEtaG= 1.301e-2;
     119  BrPhi2Pi = 7.3e-5;
    119120  PhPhi2Pi = -20.0*degree;
    120   BrPhiEe = 2.99e-4;
     121  BrPhiEe = 2.97e-4;
    121122
    122123  MsRho3 = MsRho*MsRho*MsRho;
    … …  
    125126
    126127  MeVnb = 3.8938e+11*nanobarn;
    127   Alpha = 1.0/137.036;
     128  Alpha = fine_structure_const;
    128129
    129130  AOmRho = 3.0;
    … …  
    134135  brsigpipi = 1.;
    135136
    136   msrho1450 = 1465.*MeV;
    137   msrho1700 = 1700.*MeV;
    138   grho1450 = 310.*MeV;
    139   grho1700 = 240.*MeV;
     137  msrho1450 = 1459.*MeV;
     138  msrho1700 = 1688.8*MeV;
     139  grho1450 = 171.*MeV;
     140  grho1700 = 161.*MeV;
    140141  arhoompi0 = 1.;
    141142  arho1450ompi0 = 1.;
    … …  
    190191
    191192G4double G4eeCrossSections::CrossSection2pi(G4double e)
    192 {
    193  
     193{
    194194  complex<G4double> xr(cos(PhRho),sin(PhRho));
    195195  complex<G4double> xo(cos(PhOm2Pi),sin(PhOm2Pi));
    … …  
    205205    + sqrt(Width2p(s,MsOm,GOm,BrOm2Pi,MsPi)*MsOm3*BrOmEe*GOm)*xo/dom
    206206    + sqrt(Width2p(s,MsPhi,GPhi,BrPhi2Pi,MsPi)*MsPhi3*BrPhiEe*GPhi)*xf/dphi;
     207
     208  G4double cross = 12.0*pi*MeVnb*norm(amp)/(e*s);
     209
     210  return cross;
     211}
     212
     213//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     214
     215G4double G4eeCrossSections::CrossSection3pi(G4double e)
     216{
     217  complex<G4double> xf(cos(PhPhi2Pi),sin(PhPhi));
     218
     219  G4double s = e*e;
     220  complex<G4double> dom  = DpOm(e);
     221  complex<G4double> dphi = DpPhi(e);
     222
     223  complex<G4double> amp =
     224    sqrt(Width3p(s,MsOm,GOm,BrOm3Pi)*MsOm3*BrOmEe*GOm)/dom
     225    + sqrt(Width3p(s,MsPhi,GPhi,BrPhi3Pi)*MsPhi3*BrPhiEe*GPhi)*xf/dphi;
     226
     227  G4double cross = 12.0*pi*MeVnb*norm(amp)/(e*s);
     228
     229  return cross;
     230}
     231
     232//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     233
     234G4double G4eeCrossSections::CrossSectionPi0G(G4double e)
     235{
     236  complex<G4double> xf(cos(PhPhi),sin(PhPhi));
     237
     238  G4double s = e*e;
     239  complex<G4double> drho = DpRho(e);
     240  complex<G4double> dom  = DpOm(e);
     241  complex<G4double> dphi = DpPhi(e);
     242
     243  complex<G4double> amp =
     244      sqrt(WidthPg(s,MsRho,GRho,BrRhoPi0G,MsPi0)*MsRho3*BrRhoEe*GRho)/drho
     245    + sqrt(WidthPg(s,MsOm,GOm,BrOmPi0G,MsPi0)*MsOm3*BrOmEe*GOm)/dom
     246    + sqrt(WidthPg(s,MsPhi,GPhi,BrPhiPi0G,MsPi0)*MsPhi3*BrPhiEe*GPhi)*xf/dphi;
     247
     248  G4double cross = 12.0*pi*MeVnb*norm(amp)/(e*s);
     249
     250  return cross;
     251}
     252
     253//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     254
     255G4double G4eeCrossSections::CrossSectionEtaG(G4double e)
     256{
     257  complex<G4double> xf(cos(PhPhi),sin(PhPhi));
     258
     259  G4double s = e*e;
     260  complex<G4double> drho = DpRho(e);
     261  complex<G4double> dom  = DpOm(e);
     262  complex<G4double> dphi = DpPhi(e);
     263
     264  complex<G4double> amp =
     265      sqrt(WidthPg(s,MsRho,GRho,BrRhoEtaG,MsEta)*MsRho3*BrRhoEe*GRho)/drho
     266    + sqrt(WidthPg(s,MsOm,GOm,BrOmEtaG,MsEta)*MsOm3*BrOmEe*GOm)/dom
     267    + sqrt(WidthPg(s,MsPhi,GPhi,BrPhiEtaG,MsEta)*MsPhi3*BrPhiEe*GPhi)*xf/dphi;
     268
     269  G4double cross = 12.0*pi*MeVnb*norm(amp)/(e*s);
     270
     271  return cross;
     272}
     273
     274//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     275
     276G4double G4eeCrossSections::CrossSection2Kcharged(G4double e)
     277{
     278  G4double s = e*e;
     279  complex<G4double> dphi = DpPhi(e);
     280
     281  complex<G4double> amp =
     282    sqrt(Width2p(s,MsPhi,GPhi,BrPhi2Kc,MsKc)*MsPhi3*BrPhiEe*GPhi)/dphi;
     283
     284  G4double cross = 12.0*pi*MeVnb*norm(amp)/(e*s);
     285
     286  return cross;
     287}
     288
     289//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     290
     291G4double G4eeCrossSections::CrossSection2Kneutral(G4double e)
     292{
     293  G4double s = e*e;
     294  complex<G4double> dphi = DpPhi(e);
     295
     296  complex<G4double> amp =
     297    sqrt(Width2p(s,MsPhi,GPhi,BrPhiKsKl,MsKs)*MsPhi3*BrPhiEe*GPhi)/dphi;
    207298
    208299  G4double cross = 12.0*pi*MeVnb*norm(amp)/(e*s);
    … …  
    238329G4double G4eeCrossSections::PhaseSpace3p(G4double e)
    239330{
    240  
     331  //     E.A.Kuraev, Z.K.Silagadze.
     332  //  Once more about the omega->3 pi contact term.
     333  //  Yadernaya Phisica, 1995, V58, N9, p.1678-1694. 
     334
    241335  //  G4bool b;
    242336  //  G4double x = ph3p->GetValue(e, b);

  • trunk/source/processes/electromagnetic/highenergy/src/G4eeToHadrons.cc

    r819 r961  
    2424// ********************************************************************
    2525//
    26 // $Id: G4eeToHadrons.cc,v 1.7 2006/06/29 19:32:44 gunter Exp $
    27 // GEANT4 tag $Name: $
     26// $Id: G4eeToHadrons.cc,v 1.9 2009/02/20 16:38:33 vnivanch Exp $
     27// GEANT4 tag $Name: geant4-09-02-ref-02 $
    2828//
    2929// -------------------------------------------------------------------
    … …  
    6464    isInitialised(false)
    6565{
    66     SetVerboseLevel(1);
     66  SetVerboseLevel(1);
     67  SetProcessSubType(fAnnihilationToHadrons);
    6768}
    6869
    … …  
    7172G4eeToHadrons::~G4eeToHadrons()
    7273{}
     74
     75//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     76
     77G4bool G4eeToHadrons::IsApplicable(const G4ParticleDefinition& p)
     78{
     79  return (&p == G4Positron::Positron());
     80}
    7381
    7482//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....

  • trunk/source/processes/electromagnetic/highenergy/src/G4eeToHadronsModel.cc

    r819 r961  
    2424// ********************************************************************
    2525//
    26 // $Id: G4eeToHadronsModel.cc,v 1.8 2007/05/22 17:37:30 vnivanch Exp $
    27 // GEANT4 tag $Name: $
     26// $Id: G4eeToHadronsModel.cc,v 1.9 2008/07/10 18:06:39 vnivanch Exp $
     27// GEANT4 tag $Name: geant4-09-02-ref-02 $
    2828//
    2929// -------------------------------------------------------------------
    … …  
    6565using namespace std;
    6666
    67 G4eeToHadronsModel::G4eeToHadronsModel(const G4Vee2hadrons* m,
    68                                              G4int ver,
     67G4eeToHadronsModel::G4eeToHadronsModel(G4Vee2hadrons* m, G4int ver,
    6968                                       const G4String& nam)
    7069  : G4VEmModel(nam),
    71   model(m),
    72   crossPerElectron(0),
    73   crossBornPerElectron(0),
    74   isInitialised(false),
    75   nbins(100),
    76   verbose(ver)
    77 {
    78   theGamma      = G4Gamma::Gamma();
     70    model(m),
     71    crossPerElectron(0),
     72    crossBornPerElectron(0),
     73    isInitialised(false),
     74    nbins(100),
     75    verbose(ver)
     76{
     77  theGamma = G4Gamma::Gamma();
    7978}
    8079
    … …  
    9695  isInitialised  = true;
    9796
     97  // Lab system
    9898  highKinEnergy = HighEnergyLimit();
    9999  lowKinEnergy  = LowEnergyLimit();
    100100
    101   emin  = model->ThresholdEnergy();
    102   emax = 2.0*electron_mass_c2*sqrt(1.0 + 0.5*highKinEnergy/electron_mass_c2);
    103   if(emin > emax) emin = emax;
    104 
    105   lowKinEnergy  = 0.5*emin*emin/electron_mass_c2 - 2.0*electron_mass_c2;
    106 
    107   epeak = min(model->PeakEnergy(), emax);
     101  // CM system
     102  emin  = model->LowEnergy();
     103  emax  = model->HighEnergy();
     104
     105  G4double emin0 =
     106    2.0*electron_mass_c2*sqrt(1.0 + 0.5*lowKinEnergy/electron_mass_c2);
     107  G4double emax0 =
     108    2.0*electron_mass_c2*sqrt(1.0 + 0.5*highKinEnergy/electron_mass_c2);
     109
     110  // recompute low energy
     111  if(emin0 > emax) {
     112    emin0 = emax;
     113    model->SetLowEnergy(emin0);
     114  }
     115  if(emin > emin0) {
     116    emin0 = emin;
     117    lowKinEnergy  = 0.5*emin*emin/electron_mass_c2 - 2.0*electron_mass_c2;
     118    SetLowEnergyLimit(lowKinEnergy);
     119  }
     120
     121  // recompute high energy
     122  if(emax < emax0) {
     123    emax0 = emax;
     124    highKinEnergy = 0.5*emax*emax/electron_mass_c2 - 2.0*electron_mass_c2;
     125    SetHighEnergyLimit(highKinEnergy);
     126  }
     127
     128  // peak energy
     129  epeak = std::min(model->PeakEnergy(), emax);
    108130  peakKinEnergy  = 0.5*epeak*epeak/electron_mass_c2 - 2.0*electron_mass_c2;
    109131
    … …  
    148170    }
    149171  }
     172}
     173
     174//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     175
     176G4double G4eeToHadronsModel::CrossSectionPerVolume(
     177                                      const G4Material* mat,
     178                                      const G4ParticleDefinition* p,
     179                                      G4double kineticEnergy,
     180                                      G4double, G4double)
     181{
     182  return mat->GetElectronDensity()*
     183    ComputeCrossSectionPerElectron(p, kineticEnergy);
     184}
     185
     186//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     187
     188G4double G4eeToHadronsModel::ComputeCrossSectionPerAtom(
     189                                      const G4ParticleDefinition* p,
     190                                      G4double kineticEnergy,
     191                                      G4double Z, G4double,
     192                                      G4double, G4double)
     193{
     194  return Z*ComputeCrossSectionPerElectron(p, kineticEnergy);
    150195}
    151196

  • trunk/source/processes/electromagnetic/highenergy/src/G4eeToHadronsMultiModel.cc

    r819 r961  
    2424// ********************************************************************
    2525//
    26 // $Id: G4eeToHadronsMultiModel.cc,v 1.4 2007/05/23 08:50:41 vnivanch Exp $
    27 // GEANT4 tag $Name: $
     26// $Id: G4eeToHadronsMultiModel.cc,v 1.6 2008/07/11 17:49:11 vnivanch Exp $
     27// GEANT4 tag $Name: geant4-09-02-ref-02 $
    2828//
    2929// -------------------------------------------------------------------
    … …  
    3434// File name:     G4eeToHadronsMultiModel
    3535//
    36 // Author:        Vladimir Ivanchenko on base of Michel Maire code
     36// Author:        Vladimir Ivanchenko
    3737//
    3838// Creation date: 02.08.2004
    … …  
    5151#include "G4eeToHadronsMultiModel.hh"
    5252#include "G4eeToTwoPiModel.hh"
     53#include "G4eeTo3PiModel.hh"
     54#include "G4eeToPGammaModel.hh"
     55#include "G4ee2KNeutralModel.hh"
     56#include "G4ee2KChargedModel.hh"
    5357#include "G4eeCrossSections.hh"
     58#include "G4Vee2hadrons.hh"
    5459
    5560//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
    … …  
    8085//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
    8186
    82 void G4eeToHadronsMultiModel::Initialise(const G4ParticleDefinition* p, const G4DataVector& v)
     87void G4eeToHadronsMultiModel::Initialise(const G4ParticleDefinition*,
     88                                         const G4DataVector&)
    8389{
    8490  if(!isInitialised) {
    8591    isInitialised = true;
    8692
    87     thKineticEnergy = DBL_MAX;
    88     maxKineticEnergy = HighEnergyLimit();
     93    thKineticEnergy  = DBL_MAX;
     94    maxKineticEnergy = 1.2*GeV;
    8995
    9096    cross = new G4eeCrossSections();
    91     G4eeToHadronsModel* model =
    92       new G4eeToHadronsModel(new G4eeToTwoPiModel(cross), verbose);
    93     models.push_back(model);
    94     model->SetHighEnergyLimit(maxKineticEnergy);
    95     model->Initialise(p, v);
    96     G4double emin = model->LowEnergyLimit();
    97     if(emin < thKineticEnergy) thKineticEnergy = emin;
    98     ekinMin.push_back(emin);
    99     ekinMax.push_back(model->HighEnergyLimit());
    100     ekinPeak.push_back(model->PeakEnergy());
    101     cumSum.push_back(0.0);
    102     nModels = 1;
    103 
    104     if(pParticleChange)
     97
     98    G4eeToTwoPiModel* m2pi = new G4eeToTwoPiModel(cross);
     99    m2pi->SetHighEnergy(maxKineticEnergy);
     100    AddEEModel(m2pi);
     101
     102    G4eeTo3PiModel* m3pi1 = new G4eeTo3PiModel(cross);
     103    m3pi1->SetHighEnergy(0.95*GeV);
     104    AddEEModel(m3pi1);
     105
     106    G4eeTo3PiModel* m3pi2 = new G4eeTo3PiModel(cross);
     107    m3pi2->SetLowEnergy(0.95*GeV);
     108    m3pi2->SetHighEnergy(maxKineticEnergy);
     109    AddEEModel(m3pi2);
     110
     111    G4ee2KChargedModel* m2kc = new G4ee2KChargedModel(cross);
     112    m2kc->SetHighEnergy(maxKineticEnergy);
     113    AddEEModel(m2kc);
     114
     115    G4ee2KNeutralModel* m2kn = new G4ee2KNeutralModel(cross);
     116    m2kn->SetHighEnergy(maxKineticEnergy);
     117    AddEEModel(m2kn);
     118
     119    G4eeToPGammaModel* mpg1 = new G4eeToPGammaModel(cross,"pi0");
     120    mpg1->SetLowEnergy(0.7*GeV);
     121    mpg1->SetHighEnergy(maxKineticEnergy);
     122    AddEEModel(mpg1);
     123
     124    G4eeToPGammaModel* mpg2 = new G4eeToPGammaModel(cross,"eta");
     125    mpg2->SetLowEnergy(0.7*GeV);
     126    mpg2->SetHighEnergy(maxKineticEnergy);
     127    AddEEModel(mpg2);
     128
     129    nModels = models.size();
     130
     131    if(pParticleChange) {
    105132      fParticleChange =
    106133        reinterpret_cast<G4ParticleChangeForGamma*>(pParticleChange);
    107     else
     134    } else {
    108135      fParticleChange = new G4ParticleChangeForGamma();
     136    }
     137  }
     138}
     139
     140//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     141
     142void G4eeToHadronsMultiModel::AddEEModel(G4Vee2hadrons* mod)
     143{
     144  G4eeToHadronsModel* model = new G4eeToHadronsModel(mod, verbose);
     145  model->SetLowEnergyLimit(LowEnergyLimit());
     146  model->SetHighEnergyLimit(HighEnergyLimit());
     147  models.push_back(model);
     148  G4double elow = mod->ThresholdEnergy();
     149  ekinMin.push_back(elow);
     150  if(thKineticEnergy > elow) thKineticEnergy = elow;
     151  ekinMax.push_back(mod->HighEnergy());
     152  ekinPeak.push_back(mod->PeakEnergy());
     153  cumSum.push_back(0.0);
     154}
     155
     156//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     157
     158G4double G4eeToHadronsMultiModel::CrossSectionPerVolume(
     159                                      const G4Material* mat,
     160                                      const G4ParticleDefinition* p,
     161                                      G4double kineticEnergy,
     162                                      G4double, G4double)
     163{
     164  return mat->GetElectronDensity()*
     165    ComputeCrossSectionPerElectron(p, kineticEnergy);
     166}
     167
     168//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     169
     170G4double G4eeToHadronsMultiModel::ComputeCrossSectionPerAtom(
     171                                      const G4ParticleDefinition* p,
     172                                      G4double kineticEnergy,
     173                                      G4double Z, G4double,
     174                                      G4double, G4double)
     175{
     176  return Z*ComputeCrossSectionPerElectron(p, kineticEnergy);
     177}
     178
     179
     180//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     181
     182void G4eeToHadronsMultiModel::SampleSecondaries(std::vector<G4DynamicParticle*>* newp,
     183                                                const G4MaterialCutsCouple* couple,
     184                                                const G4DynamicParticle* dp,
     185                                                G4double, G4double)
     186{
     187  G4double kinEnergy = dp->GetKineticEnergy();
     188  if (kinEnergy > thKineticEnergy) {
     189    G4double q = cumSum[nModels-1]*G4UniformRand();
     190    for(G4int i=0; i<nModels; i++) {
     191      if(q <= cumSum[i]) {
     192        (models[i])->SampleSecondaries(newp, couple,dp);
     193        if(newp->size() > 0) fParticleChange->ProposeTrackStatus(fStopAndKill);
     194        break;
     195      }
     196    }
    109197  }
    110198}
    … …  
    115203{
    116204  if(verbose > 0) {
    117     G4cout << "      e+ annihilation into hadrons active above "
    118            << thKineticEnergy/GeV << " GeV"
     205    G4double e1 = 0.5*thKineticEnergy*thKineticEnergy/electron_mass_c2
     206      - 2.0*electron_mass_c2;
     207    G4double e2 = 0.5*maxKineticEnergy*maxKineticEnergy/electron_mass_c2
     208      - 2.0*electron_mass_c2;
     209    G4cout << "      e+ annihilation into hadrons active from "
     210           << e1/GeV << " GeV to " << e2/GeV << " GeV"
    119211           << G4endl;
    120212  }
    … …  
    128220    csFactor = fac;
    129221    if(verbose > 0)
    130       G4cout << "### G4eeToHadronsMultiModel: The cross section for G4eeToHadronsMultiModel is  "
    131              << "increased by the Factor= " << csFactor << G4endl;
    132   }
    133 }
    134 
    135 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     222      G4cout << "### G4eeToHadronsMultiModel: The cross section for G4eeToHadronsMultiModel "
     223             << " is increased by the Factor= " << csFactor << G4endl;
     224  }
     225}
     226
     227//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....

  • trunk/source/processes/electromagnetic/highenergy/src/G4eeToTwoPiModel.cc

    r819 r961  
    2424// ********************************************************************
    2525//
    26 // $Id: G4eeToTwoPiModel.cc,v 1.5 2007/05/22 17:37:30 vnivanch Exp $
    27 // GEANT4 tag $Name: $
     26// $Id: G4eeToTwoPiModel.cc,v 1.7 2009/02/20 16:38:33 vnivanch Exp $
     27// GEANT4 tag $Name: geant4-09-02-ref-02 $
    2828//
    2929// -------------------------------------------------------------------
    … …  
    6464  cross(cr)
    6565{
    66   Initialise();
     66  massPi = G4PionPlus::PionPlus()->GetPDGMass();
     67  massRho = 775.5*MeV;
    6768}
    6869
    … …  
    7475//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
    7576
    76 void G4eeToTwoPiModel::Initialise()
     77G4double G4eeToTwoPiModel::ThresholdEnergy() const
    7778{
    78   massPi = G4PionPlus::PionPlus()->GetPDGMass();
    79   massRho = 770.*MeV;
    80   highEnergy = 1.*GeV;
    81   cross = new G4eeCrossSections();
     79  return 2.0*massPi;
     80}
     81
     82//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     83
     84G4double G4eeToTwoPiModel::PeakEnergy() const
     85{
     86  return massRho;
     87}
     88
     89//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     90
     91G4double G4eeToTwoPiModel::ComputeCrossSection(G4double e) const
     92{
     93  G4double ee = std::min(HighEnergy(),e);
     94  return cross->CrossSection2pi(ee);
    8295}
    8396
    … …  
    87100                                                 G4double emax) const
    88101{
    89   G4double tmin = max(emin, 2.0*massPi);
    90   G4double tmax = max(tmin, emax);
     102  G4double tmin = std::max(emin, 2.0*massPi);
     103  G4double tmax = std::max(tmin, emax);
    91104  G4int nbins = (G4int)((tmax - tmin)/(5.*MeV));
    92105  G4PhysicsVector* v = new G4PhysicsLinearVector(emin,emax,nbins);
     106  v->SetSpline(true);
    93107  return v;
    94108}
    … …  
    97111
    98112void G4eeToTwoPiModel::SampleSecondaries(std::vector<G4DynamicParticle*>* newp,
    99             G4double e, const G4ThreeVector& direction) const
     113            G4double e, const G4ThreeVector& direction)
    100114{
    101115
    … …  
    113127  dir.rotateUz(direction);
    114128
    115   // create G4DynamicParticle object for delta ray
     129  // create G4DynamicParticle objects
    116130  G4DynamicParticle* pip =
    117131     new G4DynamicParticle(G4PionPlus::PionPlus(),dir,tkin);

  • trunk/source/processes/electromagnetic/highenergy/src/G4hhIonisation.cc

    r819 r961  
    2424// ********************************************************************
    2525//
    26 // $Id: G4hhIonisation.cc,v 1.6 2007/05/22 17:37:30 vnivanch Exp $
    27 // GEANT4 tag $Name: geant4-09-01-patch-02 $
     26// $Id: G4hhIonisation.cc,v 1.9 2009/02/20 16:38:33 vnivanch Exp $
     27// GEANT4 tag $Name: geant4-09-02-ref-02 $
    2828//
    2929// -------------------------------------------------------------------
    … …  
    6464    isInitialised(false)
    6565{
    66   minKinEnergy = 0.1*keV;
    67   maxKinEnergy = 100.*TeV;
    68   SetDEDXBinning(120);
    69   SetMinKinEnergy(minKinEnergy);
    70   SetMaxKinEnergy(maxKinEnergy);
    7166  SetStepFunction(0.1, 0.1*mm);
    7267  SetVerboseLevel(1);
     68  SetProcessSubType(fIonisation);
    7369  mass = 0.0;
    7470  ratio = 0.0;
    … …  
    7975G4hhIonisation::~G4hhIonisation()
    8076{}
     77
     78//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     79
     80G4bool G4hhIonisation::IsApplicable(const G4ParticleDefinition& p)
     81{
     82  return (p.GetPDGCharge() != 0.0 && p.GetPDGMass() > 100.0*MeV &&
     83         !p.IsShortLived());
     84}
     85
     86//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     87
     88G4double G4hhIonisation::MinPrimaryEnergy(const G4ParticleDefinition*,
     89                                          const G4Material*,
     90                                          G4double cut)
     91{
     92  G4double x = 0.5*cut/electron_mass_c2;
     93  G4double y = electron_mass_c2/mass;
     94  G4double g = x*y + std::sqrt((1. + x)*(1. + x*y*y));
     95  return mass*(g - 1.0);
     96}
    8197
    8298//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
    … …  
    101117  G4int nm = 1;
    102118
     119  minKinEnergy = MinKinEnergy();
     120
    103121  if(eth > minKinEnergy) {
    104122    G4VEmModel* em = new G4BraggNoDeltaModel();
    … …  
    109127  }
    110128
    111   if(eth < maxKinEnergy) {
     129  if(eth < MaxKinEnergy()) {
    112130    G4VEmModel* em1 = new G4BetheBlochNoDeltaModel();
    113131    em1->SetLowEnergyLimit(std::max(eth,minKinEnergy));
    114     em1->SetHighEnergyLimit(maxKinEnergy);
     132    em1->SetHighEnergyLimit(MaxKinEnergy());
    115133    AddEmModel(nm, em1, flucModel);
    116134  }
    117135
    118   if(verboseLevel>0)
     136  if(verboseLevel>1) {
    119137    G4cout << "G4hhIonisation is initialised: Nmodels= " << nm << G4endl;
    120 
     138  }
    121139  isInitialised = true;
    122140}
    … …  
    127145{
    128146  G4cout << "      Delta-ray will not be produced; "
    129          << "Bether-Bloch model for E > " << std::max(eth,minKinEnergy)
    130147         << G4endl;
    131   if(eth > minKinEnergy) G4cout
    132          << "      ICRU49 parametrisation scaled from protons below.";
    133   G4cout << G4endl;
    134148}
    135149

  • trunk/source/processes/electromagnetic/highenergy/src/G4mplIonisation.cc

    r819 r961  
    2424// ********************************************************************
    2525//
    26 // $Id: G4mplIonisation.cc,v 1.5 2007/05/31 11:13:31 vnivanch Exp $
    27 // GEANT4 tag $Name: geant4-09-01-patch-02 $
     26// $Id: G4mplIonisation.cc,v 1.8 2009/02/20 16:38:33 vnivanch Exp $
     27// GEANT4 tag $Name: geant4-09-02-ref-02 $
    2828//
    2929// -------------------------------------------------------------------
    … …  
    6464
    6565  SetVerboseLevel(0);
     66  SetProcessSubType(fIonisation);
     67  SetStepFunction(0.2, 1*mm);
    6668}
    6769
    … …  
    7072G4mplIonisation::~G4mplIonisation()
    7173{}
     74
     75//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     76
     77G4bool G4mplIonisation::IsApplicable(const G4ParticleDefinition& p)
     78{
     79  return (p.GetParticleName() == "monopole");
     80}
    7281
    7382//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
    … …  
    8291
    8392  G4mplIonisationModel* ion  = new G4mplIonisationModel(magneticCharge,"PAI");
    84   ion->SetLowEnergyLimit(0.1*keV);
    85   ion->SetHighEnergyLimit(100.*TeV);
     93  ion->SetLowEnergyLimit(MinKinEnergy());
     94  ion->SetHighEnergyLimit(MaxKinEnergy());
    8695  AddEmModel(0,ion,ion);
    87 
    88   SetStepFunction(0.2, 1*mm);
    8996
    9097  isInitialised = true;

  • trunk/source/processes/electromagnetic/highenergy/src/G4mplIonisationModel.cc

    r819 r961  
    2424// ********************************************************************
    2525//
    26 // $Id: G4mplIonisationModel.cc,v 1.5 2007/11/13 18:36:29 vnivanch Exp $
    27 // GEANT4 tag $Name: $
     26// $Id: G4mplIonisationModel.cc,v 1.6 2009/02/20 16:38:33 vnivanch Exp $
     27// GEANT4 tag $Name: geant4-09-02-ref-02 $
    2828//
    2929// -------------------------------------------------------------------
    … …  
    139139//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
    140140
    141 G4double G4mplIonisationModel::ComputeDEDXAhlen(const G4Material* material, G4double bg2)
     141G4double G4mplIonisationModel::ComputeDEDXAhlen(const G4Material* material,
     142                                                G4double bg2)
    142143{
    143144  G4double eDensity = material->GetElectronDensity();
    … …  
    176177  return dedx;
    177178}
     179
     180//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     181
     182void G4mplIonisationModel::SampleSecondaries(std::vector<G4DynamicParticle*>*,
     183                                             const G4MaterialCutsCouple*,
     184                                             const G4DynamicParticle*,
     185                                             G4double,
     186                                             G4double)
     187{}
    178188
    179189//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
    … …  
    204214  return loss;
    205215}
     216
     217//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     218
     219G4double G4mplIonisationModel::Dispersion(const G4Material* material,
     220                                          const G4DynamicParticle* dp,
     221                                          G4double& tmax,
     222                                          G4double& length)
     223{
     224  G4double siga = 0.0;
     225  G4double tau   = dp->GetKineticEnergy()/mass;
     226  if(tau > 0.0) {
     227    G4double electronDensity = material->GetElectronDensity();
     228    G4double gam   = tau + 1.0;
     229    G4double invbeta2 = (gam*gam)/(tau * (tau+2.0));
     230    siga  = (invbeta2 - 0.5) * twopi_mc2_rcl2 * tmax * length
     231      * electronDensity * chargeSquare;
     232  }
     233  return siga;
     234}
     235
     236//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
Note: See TracChangeset for help on using the changeset viewer.