Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/particles/management/src/G4DalitzDecayChannel.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 ]

Diff markup

Differences between /particles/management/src/G4DalitzDecayChannel.cc (Version 11.3.0) and /particles/management/src/G4DalitzDecayChannel.cc (Version 10.6.p3)


  1 //                                                  1 //
  2 // *******************************************      2 // ********************************************************************
  3 // * License and Disclaimer                         3 // * License and Disclaimer                                           *
  4 // *                                                4 // *                                                                  *
  5 // * The  Geant4 software  is  copyright of th      5 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
  6 // * the Geant4 Collaboration.  It is provided      6 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
  7 // * conditions of the Geant4 Software License      7 // * conditions of the Geant4 Software License,  included in the file *
  8 // * LICENSE and available at  http://cern.ch/      8 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
  9 // * include a list of copyright holders.           9 // * include a list of copyright holders.                             *
 10 // *                                               10 // *                                                                  *
 11 // * Neither the authors of this software syst     11 // * Neither the authors of this software system, nor their employing *
 12 // * institutes,nor the agencies providing fin     12 // * institutes,nor the agencies providing financial support for this *
 13 // * work  make  any representation or  warran     13 // * work  make  any representation or  warranty, express or implied, *
 14 // * regarding  this  software system or assum     14 // * regarding  this  software system or assume any liability for its *
 15 // * use.  Please see the license in the file      15 // * use.  Please see the license in the file  LICENSE  and URL above *
 16 // * for the full disclaimer and the limitatio     16 // * for the full disclaimer and the limitation of liability.         *
 17 // *                                               17 // *                                                                  *
 18 // * This  code  implementation is the result      18 // * This  code  implementation is the result of  the  scientific and *
 19 // * technical work of the GEANT4 collaboratio     19 // * technical work of the GEANT4 collaboration.                      *
 20 // * By using,  copying,  modifying or  distri     20 // * By using,  copying,  modifying or  distributing the software (or *
 21 // * any work based  on the software)  you  ag     21 // * any work based  on the software)  you  agree  to acknowledge its *
 22 // * use  in  resulting  scientific  publicati     22 // * use  in  resulting  scientific  publications,  and indicate your *
 23 // * acceptance of all terms of the Geant4 Sof     23 // * acceptance of all terms of the Geant4 Software license.          *
 24 // *******************************************     24 // ********************************************************************
 25 //                                                 25 //
 26 // G4DalitzDecayChannel class implementation   << 
 27 //                                                 26 //
 28 // Author: H.Kurashige, 30 May 1997            <<  27 //
 29 // ------------------------------------------- <<  28 // 
 30                                                <<  29 // ------------------------------------------------------------
 31 #include "G4DalitzDecayChannel.hh"             <<  30 //      GEANT 4 class header file
                                                   >>  31 //
                                                   >>  32 //      History: first implementation, based on object model of
                                                   >>  33 //      30 May 1997 H.Kurashige
                                                   >>  34 // ------------------------------------------------------------
 32                                                    35 
 33 #include "G4DecayProducts.hh"                  << 
 34 #include "G4LorentzRotation.hh"                << 
 35 #include "G4LorentzVector.hh"                  << 
 36 #include "G4ParticleDefinition.hh"             << 
 37 #include "G4PhaseSpaceDecayChannel.hh"         << 
 38 #include "G4PhysicalConstants.hh"                  36 #include "G4PhysicalConstants.hh"
 39 #include "G4SystemOfUnits.hh"                      37 #include "G4SystemOfUnits.hh"
                                                   >>  38 #include "G4ParticleDefinition.hh"
                                                   >>  39 #include "G4DecayProducts.hh"
 40 #include "G4VDecayChannel.hh"                      40 #include "G4VDecayChannel.hh"
                                                   >>  41 #include "G4DalitzDecayChannel.hh"
                                                   >>  42 #include "G4PhaseSpaceDecayChannel.hh"
 41 #include "Randomize.hh"                            43 #include "Randomize.hh"
                                                   >>  44 #include "G4LorentzVector.hh"
                                                   >>  45 #include "G4LorentzRotation.hh"
                                                   >>  46 
                                                   >>  47 G4DalitzDecayChannel::G4DalitzDecayChannel()
                                                   >>  48   :G4VDecayChannel()
                                                   >>  49 {
                                                   >>  50 }
 42                                                    51 
 43 G4DalitzDecayChannel::G4DalitzDecayChannel(con <<  52 G4DalitzDecayChannel::G4DalitzDecayChannel(
 44                                            con <<  53          const G4String& theParentName,
 45                                            con <<  54          G4double        theBR,
 46   : G4VDecayChannel("Dalitz Decay", 1)         <<  55          const G4String& theLeptonName,
                                                   >>  56          const G4String& theAntiLeptonName)
                                                   >>  57                    :G4VDecayChannel("Dalitz Decay",1)
 47 {                                                  58 {
 48   // set names for daughter particles              59   // set names for daughter particles
 49   SetParent(theParentName);                        60   SetParent(theParentName);
 50   SetBR(theBR);                                    61   SetBR(theBR);
 51   SetNumberOfDaughters(3);                         62   SetNumberOfDaughters(3);
 52   G4String gammaName = "gamma";                    63   G4String gammaName = "gamma";
 53   SetDaughter(idGamma, gammaName);                 64   SetDaughter(idGamma, gammaName);
 54   SetDaughter(idLepton, theLeptonName);            65   SetDaughter(idLepton, theLeptonName);
 55   SetDaughter(idAntiLepton, theAntiLeptonName)     66   SetDaughter(idAntiLepton, theAntiLeptonName);
 56 }                                                  67 }
 57                                                    68 
 58 G4DalitzDecayChannel& G4DalitzDecayChannel::op <<  69 G4DalitzDecayChannel::~G4DalitzDecayChannel()
 59 {                                                  70 {
 60   if (this != &right) {                        <<  71 }
                                                   >>  72 
                                                   >>  73 G4DalitzDecayChannel::G4DalitzDecayChannel(const G4DalitzDecayChannel &right):
                                                   >>  74   G4VDecayChannel(right)
                                                   >>  75 {
                                                   >>  76 }
                                                   >>  77 
                                                   >>  78 G4DalitzDecayChannel & G4DalitzDecayChannel::operator=(const G4DalitzDecayChannel & right)
                                                   >>  79 {
                                                   >>  80   if (this != &right) { 
 61     kinematics_name = right.kinematics_name;       81     kinematics_name = right.kinematics_name;
 62     verboseLevel = right.verboseLevel;             82     verboseLevel = right.verboseLevel;
 63     rbranch = right.rbranch;                       83     rbranch = right.rbranch;
 64                                                    84 
 65     // copy parent name                            85     // copy parent name
 66     parent_name = new G4String(*right.parent_n     86     parent_name = new G4String(*right.parent_name);
 67                                                    87 
 68     // clear daughters_name array                  88     // clear daughters_name array
 69     ClearDaughtersName();                          89     ClearDaughtersName();
 70                                                    90 
 71     // recreate array                              91     // recreate array
 72     numberOfDaughters = right.numberOfDaughter     92     numberOfDaughters = right.numberOfDaughters;
 73     if (numberOfDaughters > 0) {               <<  93     if ( numberOfDaughters >0 ) {
 74       if (daughters_name != nullptr) ClearDaug <<  94       if (daughters_name !=0) ClearDaughtersName();
 75       daughters_name = new G4String*[numberOfD     95       daughters_name = new G4String*[numberOfDaughters];
 76       // copy daughters name                   <<  96       //copy daughters name
 77       for (G4int index = 0; index < numberOfDa <<  97       for (G4int index=0; index < numberOfDaughters; index++) {
 78         daughters_name[index] = new G4String(* <<  98           daughters_name[index] = new G4String(*right.daughters_name[index]);
 79       }                                            99       }
 80     }                                             100     }
 81   }                                               101   }
 82   return *this;                                   102   return *this;
 83 }                                                 103 }
 84                                                   104 
 85 G4DecayProducts* G4DalitzDecayChannel::DecayIt << 105 G4DecayProducts *G4DalitzDecayChannel::DecayIt(G4double) 
 86 {                                                 106 {
 87 #ifdef G4VERBOSE                                  107 #ifdef G4VERBOSE
 88   if (GetVerboseLevel() > 1) G4cout << "G4Dali << 108   if (GetVerboseLevel()>1) G4cout << "G4DalitzDecayChannel::DecayIt ";
 89 #endif                                         << 109 #endif 
 90   CheckAndFillParent();                           110   CheckAndFillParent();
 91   CheckAndFillDaughters();                        111   CheckAndFillDaughters();
 92                                                   112 
 93   // parent mass                                  113   // parent mass
 94   G4double parentmass = G4MT_parent->GetPDGMas    114   G4double parentmass = G4MT_parent->GetPDGMass();
 95                                                << 115  
 96   // create parent G4DynamicParticle at rest   << 116  //create parent G4DynamicParticle at rest
 97   G4ThreeVector dummy;                            117   G4ThreeVector dummy;
 98   auto parentparticle = new G4DynamicParticle( << 118   G4DynamicParticle * parentparticle = new G4DynamicParticle( G4MT_parent, dummy, 0.0);
 99                                                << 119  
100   // daughters' mass                           << 120   //daughters'mass
101   G4double leptonmass = G4MT_daughters[idLepto << 121   G4double leptonmass = G4MT_daughters[idLepton]->GetPDGMass(); 
102                                                << 122 
103   // Generate t ( = std::exp(x):mass Square of << 123  // Generate t ( = std::exp(x):mass Square of (l+ l-) system) 
104   G4double xmin = 2.0 * std::log(2.0 * leptonm << 124   G4double xmin  = 2.0*std::log(2.0*leptonmass);
105   G4double xmax = 2.0 * std::log(parentmass);  << 125   G4double xmax  = 2.0*std::log(parentmass);
106   G4double wmax = 1.5;                            126   G4double wmax = 1.5;
107   G4double x, w, ww, w1, w2, w3, t;               127   G4double x, w, ww, w1, w2, w3, t;
108   const std::size_t MAX_LOOP = 10000;          << 128   const size_t MAX_LOOP = 10000;
109   for (std::size_t loop_counter = 0; loop_coun << 129   for (size_t loop_counter=0; loop_counter <MAX_LOOP; ++loop_counter){
110     x = G4UniformRand() * (xmax - xmin) + xmin << 130     x = G4UniformRand()*(xmax-xmin) + xmin;
111     w = G4UniformRand() * wmax;                << 131     w = G4UniformRand()*wmax;
112     t = std::exp(x);                              132     t = std::exp(x);
113     w1 = (1.0 - 4.0 * leptonmass * leptonmass  << 133     w1 = (1.0-4.0*leptonmass*leptonmass/t);
114     if (w1 > 0.0) {                            << 134     if ( w1 > 0.0) {
115       w2 = (1.0 + 2.0 * leptonmass * leptonmas << 135       w2 = ( 1.0 + 2.0*leptonmass*leptonmass/t );
116       w3 = (1.0 - t / parentmass / parentmass) << 136       w3 = ( 1.0 - t/parentmass/parentmass );
117       w3 = w3 * w3 * w3;                          137       w3 = w3 * w3 * w3;
118       ww = w3 * w2 * std::sqrt(w1);               138       ww = w3 * w2 * std::sqrt(w1);
119     }                                          << 139     } else {
120     else {                                     << 
121       ww = 0.0;                                   140       ww = 0.0;
122     }                                             141     }
123     if (w <= ww) break;                        << 142     if (w <= ww) break;    
124   }                                               143   }
125                                                   144 
126   // calculate gamma momentum                     145   // calculate gamma momentum
127   G4double Pgamma = G4PhaseSpaceDecayChannel:: << 146   G4double Pgamma = 
128   G4double costheta = 2. * G4UniformRand() - 1 << 147       G4PhaseSpaceDecayChannel::Pmx(parentmass, 0.0, std::sqrt(t)); 
129   G4double sintheta = std::sqrt((1.0 - costhet << 148   G4double costheta = 2.*G4UniformRand()-1.0;
130   G4double phi = twopi * G4UniformRand() * rad << 149   G4double sintheta = std::sqrt((1.0 - costheta)*(1.0 + costheta));
131   G4ThreeVector gdirection(sintheta * std::cos << 150   G4double phi  = twopi*G4UniformRand()*rad;
                                                   >> 151   G4ThreeVector gdirection(sintheta*std::cos(phi),sintheta*std::sin(phi),costheta);
                                                   >> 152 
                                                   >> 153   //create G4DynamicParticle for gamma 
                                                   >> 154   G4DynamicParticle * gammaparticle
                                                   >> 155       = new G4DynamicParticle(G4MT_daughters[idGamma] , gdirection, Pgamma);
132                                                   156 
133   // create G4DynamicParticle for gamma        << 157   // calcurate beta of (l+ l-)system
134   auto gammaparticle = new G4DynamicParticle(G << 158   G4double beta = Pgamma/(parentmass-Pgamma);
135                                                << 
136   // calculate beta of (l+ l-)system           << 
137   G4double beta = Pgamma / (parentmass - Pgamm << 
138                                                   159 
139   // calculate lepton momentum in the rest fra    160   // calculate lepton momentum in the rest frame of (l+ l-)system
140   G4double Plepton = G4PhaseSpaceDecayChannel: << 161   G4double Plepton = 
141   G4double Elepton = std::sqrt(Plepton * Plept << 162       G4PhaseSpaceDecayChannel::Pmx(std::sqrt(t),leptonmass, leptonmass); 
142   costheta = 2. * G4UniformRand() - 1.0;       << 163   G4double Elepton = std::sqrt(Plepton*Plepton + leptonmass*leptonmass );
143   sintheta = std::sqrt((1.0 - costheta) * (1.0 << 164   costheta = 2.*G4UniformRand()-1.0;
144   phi = twopi * G4UniformRand() * rad;         << 165   sintheta = std::sqrt((1.0 - costheta)*(1.0 + costheta));
145   G4ThreeVector ldirection(sintheta * std::cos << 166   phi  = twopi*G4UniformRand()*rad;
146   // create G4DynamicParticle for leptons  in  << 167   G4ThreeVector ldirection(sintheta*std::cos(phi),sintheta*std::sin(phi),costheta);
147   auto leptonparticle =                        << 168   //create G4DynamicParticle for leptons  in the rest frame of (l+ l-)system
148     new G4DynamicParticle(G4MT_daughters[idLep << 169   G4DynamicParticle * leptonparticle 
149   auto antileptonparticle =                    << 170     = new G4DynamicParticle(G4MT_daughters[idLepton] , 
150     new G4DynamicParticle(G4MT_daughters[idAnt << 171           ldirection, Elepton-leptonmass );
151   // boost leptons in the rest frame of the pa << 172   G4DynamicParticle * antileptonparticle 
                                                   >> 173     = new G4DynamicParticle(G4MT_daughters[idAntiLepton] , 
                                                   >> 174           -1.0*ldirection, Elepton-leptonmass );
                                                   >> 175   //boost leptons in the rest frame of the parent 
152   G4LorentzVector p4 = leptonparticle->Get4Mom    176   G4LorentzVector p4 = leptonparticle->Get4Momentum();
153   p4.boost(-1.0 * gdirection.x() * beta, -1.0  << 177   p4.boost( -1.0*gdirection.x()*beta, -1.0*gdirection.y()*beta, -1.0*gdirection.z()*beta);
154            -1.0 * gdirection.z() * beta);      << 
155   leptonparticle->Set4Momentum(p4);               178   leptonparticle->Set4Momentum(p4);
156   p4 = antileptonparticle->Get4Momentum();        179   p4 = antileptonparticle->Get4Momentum();
157   p4.boost(-1.0 * gdirection.x() * beta, -1.0  << 180   p4.boost( -1.0*gdirection.x()*beta, -1.0*gdirection.y()*beta, -1.0*gdirection.z()*beta);
158            -1.0 * gdirection.z() * beta);      << 
159   antileptonparticle->Set4Momentum(p4);           181   antileptonparticle->Set4Momentum(p4);
160                                                   182 
161   // create G4Decayproducts                    << 183   //create G4Decayproducts
162   auto products = new G4DecayProducts(*parentp << 184   G4DecayProducts *products = new G4DecayProducts(*parentparticle);
163   delete parentparticle;                          185   delete parentparticle;
164   products->PushProducts(gammaparticle);          186   products->PushProducts(gammaparticle);
165   products->PushProducts(leptonparticle);         187   products->PushProducts(leptonparticle);
166   products->PushProducts(antileptonparticle);     188   products->PushProducts(antileptonparticle);
167                                                   189 
168 #ifdef G4VERBOSE                                  190 #ifdef G4VERBOSE
169   if (GetVerboseLevel() > 1) {                 << 191   if (GetVerboseLevel()>1) {
170     G4cout << "G4DalitzDecayChannel::DecayIt " << 192      G4cout << "G4DalitzDecayChannel::DecayIt ";
171     G4cout << "  create decay products in rest << 193      G4cout << "  create decay products in rest frame " <<G4endl;
172     products->DumpInfo();                      << 194      products->DumpInfo();
173   }                                               195   }
174 #endif                                            196 #endif
175   return products;                                197   return products;
176 }                                                 198 }
                                                   >> 199 
                                                   >> 200 
                                                   >> 201 
                                                   >> 202 
                                                   >> 203 
177                                                   204