Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/hadronic/util/src/G4ReactionProduct.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 /processes/hadronic/util/src/G4ReactionProduct.cc (Version 11.3.0) and /processes/hadronic/util/src/G4ReactionProduct.cc (Version 11.0.p4)


  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 // J.L. Chuma, TRIUMF, 31-Oct-1996                 26 // J.L. Chuma, TRIUMF, 31-Oct-1996
 27 // last modified: 19-Dec-1996                      27 // last modified: 19-Dec-1996
 28 // Modified by J.L.Chuma, 05-May-97                28 // Modified by J.L.Chuma, 05-May-97
 29 // M. Kelsey 29-Aug-2011 -- Use G4Allocator fo     29 // M. Kelsey 29-Aug-2011 -- Use G4Allocator for better memory management
 30                                                    30 
 31 #include "G4ReactionProduct.hh"                    31 #include "G4ReactionProduct.hh"
 32                                                    32 
 33 G4Allocator<G4ReactionProduct>*& aRPAllocator(     33 G4Allocator<G4ReactionProduct>*& aRPAllocator()
 34 {                                                  34 {
 35     G4ThreadLocalStatic G4Allocator<G4Reaction     35     G4ThreadLocalStatic G4Allocator<G4ReactionProduct>* _instance = nullptr;
 36     return _instance;                              36     return _instance;
 37 }                                                  37 }
 38                                                    38 
 39  G4ReactionProduct::G4ReactionProduct() :          39  G4ReactionProduct::G4ReactionProduct() :
 40     theParticleDefinition(NULL),                   40     theParticleDefinition(NULL),
 41     formationTime(0.0),                            41     formationTime(0.0),
 42     hasInitialStateParton(false),                  42     hasInitialStateParton(false),
 43     mass(0.0),                                     43     mass(0.0),
 44     totalEnergy(0.0),                              44     totalEnergy(0.0),
 45     kineticEnergy(0.0),                            45     kineticEnergy(0.0),
 46     timeOfFlight(0.0),                             46     timeOfFlight(0.0),
 47     side(0),                                       47     side(0),
 48     theCreatorModel(-1),                           48     theCreatorModel(-1),
 49     theParentResonanceDef(nullptr),            << 
 50     theParentResonanceID(0),                   << 
 51     NewlyAdded(false),                             49     NewlyAdded(false),
 52     MayBeKilled(true)                              50     MayBeKilled(true)
 53   {                                                51   {
 54     SetMomentum( 0.0, 0.0, 0.0 );                  52     SetMomentum( 0.0, 0.0, 0.0 );
 55     SetPositionInNucleus( 0.0, 0.0, 0.0 );         53     SetPositionInNucleus( 0.0, 0.0, 0.0 );
 56   }                                                54   }
 57                                                    55  
 58  G4ReactionProduct::G4ReactionProduct(             56  G4ReactionProduct::G4ReactionProduct(
 59   const G4ParticleDefinition *aParticleDefinit     57   const G4ParticleDefinition *aParticleDefinition )
 60   {                                                58   {
 61     SetMomentum( 0.0, 0.0, 0.0 );                  59     SetMomentum( 0.0, 0.0, 0.0 );
 62     SetPositionInNucleus( 0.0, 0.0, 0.0 );         60     SetPositionInNucleus( 0.0, 0.0, 0.0 );
 63     formationTime = 0.0;                           61     formationTime = 0.0;
 64     hasInitialStateParton = false;                 62     hasInitialStateParton = false;
 65     theParticleDefinition = aParticleDefinitio     63     theParticleDefinition = aParticleDefinition;
 66     mass = aParticleDefinition->GetPDGMass();      64     mass = aParticleDefinition->GetPDGMass();
 67     totalEnergy = mass;                            65     totalEnergy = mass;
 68     kineticEnergy = 0.0;                           66     kineticEnergy = 0.0;
 69     (aParticleDefinition->GetPDGEncoding()<0)      67     (aParticleDefinition->GetPDGEncoding()<0) ? timeOfFlight=-1.0 : timeOfFlight=1.0;
 70     side = 0;                                      68     side = 0;
 71     theCreatorModel = -1;                          69     theCreatorModel = -1;
 72     theParentResonanceDef = nullptr;           << 
 73     theParentResonanceID = 0;                  << 
 74     NewlyAdded = false;                            70     NewlyAdded = false;
 75     MayBeKilled = true;                            71     MayBeKilled = true;
 76   }                                                72   }
 77                                                    73  
 78  G4ReactionProduct::G4ReactionProduct(             74  G4ReactionProduct::G4ReactionProduct(
 79   const G4ReactionProduct &right )                 75   const G4ReactionProduct &right )
 80   {                                                76   {
 81     theParticleDefinition = right.theParticleD     77     theParticleDefinition = right.theParticleDefinition;
 82     positionInNucleus = right.positionInNucleu     78     positionInNucleus = right.positionInNucleus;
 83     formationTime = right.formationTime;           79     formationTime = right.formationTime;
 84     hasInitialStateParton = right.hasInitialSt     80     hasInitialStateParton = right.hasInitialStateParton;
 85     momentum = right.momentum;                     81     momentum = right.momentum;
 86     mass = right.mass;                             82     mass = right.mass;
 87     totalEnergy = right.totalEnergy;               83     totalEnergy = right.totalEnergy;
 88     kineticEnergy = right.kineticEnergy;           84     kineticEnergy = right.kineticEnergy;
 89     timeOfFlight = right.timeOfFlight;             85     timeOfFlight = right.timeOfFlight;
 90     side = right.side;                             86     side = right.side;
 91     theCreatorModel = right.theCreatorModel;       87     theCreatorModel = right.theCreatorModel;
 92     theParentResonanceDef = right.theParentRes << 
 93     theParentResonanceID = right.theParentReso << 
 94     NewlyAdded = right.NewlyAdded;                 88     NewlyAdded = right.NewlyAdded;
 95     MayBeKilled = right.MayBeKilled;               89     MayBeKilled = right.MayBeKilled;
 96   }                                                90   }
 97                                                    91  
 98  G4ReactionProduct &G4ReactionProduct::operato     92  G4ReactionProduct &G4ReactionProduct::operator=(
 99   const G4ReactionProduct &right )                 93   const G4ReactionProduct &right )
100   {                                                94   {
101     if( this != &right ) {                         95     if( this != &right ) {
102       theParticleDefinition = right.theParticl     96       theParticleDefinition = right.theParticleDefinition;
103       positionInNucleus = right.positionInNucl     97       positionInNucleus = right.positionInNucleus;
104       formationTime = right.formationTime;         98       formationTime = right.formationTime;
105       hasInitialStateParton = right.hasInitial     99       hasInitialStateParton = right.hasInitialStateParton;
106       momentum = right.momentum;                  100       momentum = right.momentum;
107       mass = right.mass;                          101       mass = right.mass;
108       totalEnergy = right.totalEnergy;            102       totalEnergy = right.totalEnergy;
109       kineticEnergy = right.kineticEnergy;        103       kineticEnergy = right.kineticEnergy;
110       timeOfFlight = right.timeOfFlight;          104       timeOfFlight = right.timeOfFlight;
111       side = right.side;                          105       side = right.side;
112       theCreatorModel = right.theCreatorModel;    106       theCreatorModel = right.theCreatorModel;
113       theParentResonanceDef = right.theParentR << 
114       theParentResonanceID = right.theParentRe << 
115       NewlyAdded = right.NewlyAdded;              107       NewlyAdded = right.NewlyAdded;
116       MayBeKilled = right.MayBeKilled;            108       MayBeKilled = right.MayBeKilled;
117     }                                             109     }
118     return *this;                                 110     return *this;
119   }                                               111   }
120                                                   112     
121  G4ReactionProduct &G4ReactionProduct::operato    113  G4ReactionProduct &G4ReactionProduct::operator=(
122   const G4DynamicParticle &right )                114   const G4DynamicParticle &right )
123   {                                               115   {
124     theParticleDefinition = right.GetDefinitio    116     theParticleDefinition = right.GetDefinition();
125     SetPositionInNucleus( 0.0, 0.0, 0.0 );        117     SetPositionInNucleus( 0.0, 0.0, 0.0 );
126     formationTime = 0.0;                          118     formationTime = 0.0;
127     hasInitialStateParton = false;                119     hasInitialStateParton = false;
128     momentum = right.GetMomentum();               120     momentum = right.GetMomentum();
129     mass = right.GetDefinition()->GetPDGMass()    121     mass = right.GetDefinition()->GetPDGMass();
130     totalEnergy = right.GetTotalEnergy();         122     totalEnergy = right.GetTotalEnergy();
131     kineticEnergy = right.GetKineticEnergy();     123     kineticEnergy = right.GetKineticEnergy();
132     (right.GetDefinition()->GetPDGEncoding()<0    124     (right.GetDefinition()->GetPDGEncoding()<0) ? timeOfFlight=-1.0 : timeOfFlight=1.0;
133     side = 0;                                     125     side = 0;
134     theCreatorModel = -1;                         126     theCreatorModel = -1;
135     theParentResonanceDef = nullptr;           << 
136     theParentResonanceID = 0;                  << 
137     NewlyAdded = false;                           127     NewlyAdded = false;
138     MayBeKilled = true;                           128     MayBeKilled = true;
139     return *this;                                 129     return *this;
140   }                                               130   }
141                                                   131  
142  G4ReactionProduct &G4ReactionProduct::operato    132  G4ReactionProduct &G4ReactionProduct::operator=(
143   const G4HadProjectile &right )                  133   const G4HadProjectile &right )
144   {                                               134   {
145     theParticleDefinition = right.GetDefinitio    135     theParticleDefinition = right.GetDefinition();
146     SetPositionInNucleus( 0.0, 0.0, 0.0 );        136     SetPositionInNucleus( 0.0, 0.0, 0.0 );
147     formationTime = 0.0;                          137     formationTime = 0.0;
148     hasInitialStateParton = false;                138     hasInitialStateParton = false;
149     momentum = right.Get4Momentum().vect();       139     momentum = right.Get4Momentum().vect();
150     mass = right.GetDefinition()->GetPDGMass()    140     mass = right.GetDefinition()->GetPDGMass();
151     totalEnergy = right.Get4Momentum().e();       141     totalEnergy = right.Get4Momentum().e();
152     kineticEnergy = right.GetKineticEnergy();     142     kineticEnergy = right.GetKineticEnergy();
153     (right.GetDefinition()->GetPDGEncoding()<0    143     (right.GetDefinition()->GetPDGEncoding()<0) ? timeOfFlight=-1.0 : timeOfFlight=1.0;
154     side = 0;                                     144     side = 0;
155     theCreatorModel = -1;                         145     theCreatorModel = -1;
156     theParentResonanceDef = nullptr;           << 
157     theParentResonanceID = 0;                  << 
158     NewlyAdded = false;                           146     NewlyAdded = false;
159     MayBeKilled = true;                           147     MayBeKilled = true;
160     return *this;                                 148     return *this;
161   }                                               149   }
162                                                   150  
163  void G4ReactionProduct::SetDefinitionAndUpdat    151  void G4ReactionProduct::SetDefinitionAndUpdateE(
164   const G4ParticleDefinition *aParticleDefinit    152   const G4ParticleDefinition *aParticleDefinition )
165   {    G4double aKineticEnergy = GetKineticEne    153   {    G4double aKineticEnergy = GetKineticEnergy();
166     G4double pp = GetMomentum().mag();            154     G4double pp = GetMomentum().mag();
167     G4ThreeVector aMomentum = GetMomentum();      155     G4ThreeVector aMomentum = GetMomentum();
168     SetDefinition( aParticleDefinition );         156     SetDefinition( aParticleDefinition );
169     SetKineticEnergy( aKineticEnergy );           157     SetKineticEnergy( aKineticEnergy );
170     if( pp > DBL_MIN )                            158     if( pp > DBL_MIN )
171       SetMomentum( aMomentum * (std::sqrt(aKin    159       SetMomentum( aMomentum * (std::sqrt(aKineticEnergy*aKineticEnergy +
172                                     2*aKinetic    160                                     2*aKineticEnergy*GetMass())/pp) );
173   }                                               161   }
174                                                   162 
175  void G4ReactionProduct::SetDefinition(           163  void G4ReactionProduct::SetDefinition(
176   const G4ParticleDefinition *aParticleDefinit    164   const G4ParticleDefinition *aParticleDefinition )
177   {                                               165   {
178     theParticleDefinition = aParticleDefinitio    166     theParticleDefinition = aParticleDefinition;
179     mass = aParticleDefinition->GetPDGMass();     167     mass = aParticleDefinition->GetPDGMass();
180     totalEnergy = mass;                           168     totalEnergy = mass;
181     kineticEnergy = 0.0;                          169     kineticEnergy = 0.0;
182     (aParticleDefinition->GetPDGEncoding()<0)     170     (aParticleDefinition->GetPDGEncoding()<0) ?
183       timeOfFlight=-1.0 : timeOfFlight=1.0;       171       timeOfFlight=-1.0 : timeOfFlight=1.0;
184   }                                               172   }
185                                                   173  
186  void G4ReactionProduct::SetMomentum(             174  void G4ReactionProduct::SetMomentum(
187   const G4double x, const G4double y, const G4    175   const G4double x, const G4double y, const G4double z )
188   {                                               176   {
189     momentum.setX( x );                           177     momentum.setX( x );
190     momentum.setY( y );                           178     momentum.setY( y );
191     momentum.setZ( z );                           179     momentum.setZ( z );
192   }                                               180   }
193                                                   181  
194  void G4ReactionProduct::SetMomentum(             182  void G4ReactionProduct::SetMomentum(
195   const G4double x, const G4double y )            183   const G4double x, const G4double y )
196   {                                               184   {
197     momentum.setX( x );                           185     momentum.setX( x );
198     momentum.setY( y );                           186     momentum.setY( y );
199   }                                               187   }
200                                                   188  
201  void G4ReactionProduct::SetMomentum( const G4    189  void G4ReactionProduct::SetMomentum( const G4double z )
202   {                                               190   {
203     momentum.setZ( z );                           191     momentum.setZ( z );
204   }                                               192   }
205                                                   193  
206  void G4ReactionProduct::SetZero()                194  void G4ReactionProduct::SetZero()
207   {                                               195   {
208     SetMomentum( 0.0, 0.0, 0.0 );                 196     SetMomentum( 0.0, 0.0, 0.0 );
209     totalEnergy = 0.0;                            197     totalEnergy = 0.0;
210     kineticEnergy = 0.0;                          198     kineticEnergy = 0.0;
211     mass = 0.0;                                   199     mass = 0.0;
212     timeOfFlight = 0.0;                           200     timeOfFlight = 0.0;
213     side = 0;                                     201     side = 0;
214     theCreatorModel = -1;                         202     theCreatorModel = -1;
215     theParentResonanceDef = nullptr;           << 
216     theParentResonanceID = 0;                  << 
217     NewlyAdded = false;                           203     NewlyAdded = false;
218     SetPositionInNucleus( 0.0, 0.0, 0.0 );        204     SetPositionInNucleus( 0.0, 0.0, 0.0 );
219     formationTime = 0.0;                          205     formationTime = 0.0;
220     hasInitialStateParton = false;                206     hasInitialStateParton = false;
221   }                                               207   }
222                                                   208  
223  void G4ReactionProduct::Lorentz(                 209  void G4ReactionProduct::Lorentz(
224    const G4ReactionProduct &p1, const G4Reacti    210    const G4ReactionProduct &p1, const G4ReactionProduct &p2 )
225   {                                               211   {
226     G4ThreeVector p1M = p1.momentum;              212     G4ThreeVector p1M = p1.momentum;
227     G4ThreeVector p2M = p2.momentum;              213     G4ThreeVector p2M = p2.momentum;
228     G4double p1x = p1M.x(); G4double p1y = p1M    214     G4double p1x = p1M.x(); G4double p1y = p1M.y(); G4double p1z = p1M.z();
229     G4double p2x = p2M.x(); G4double p2y = p2M    215     G4double p2x = p2M.x(); G4double p2y = p2M.y(); G4double p2z = p2M.z();
230     G4double a = ( (p1x*p2x+p1y*p2y+p1z*p2z)/(    216     G4double a = ( (p1x*p2x+p1y*p2y+p1z*p2z)/(p2.totalEnergy+p2.mass) -
231                    p1.totalEnergy ) / p2.mass;    217                    p1.totalEnergy ) / p2.mass;
232     G4double x = p1x+a*p2x;                       218     G4double x = p1x+a*p2x;
233     G4double y = p1y+a*p2y;                       219     G4double y = p1y+a*p2y;
234     G4double z = p1z+a*p2z;                       220     G4double z = p1z+a*p2z;
235     G4double p = std::sqrt(x*x+y*y+z*z);          221     G4double p = std::sqrt(x*x+y*y+z*z);
236     SetMass( p1.mass );                           222     SetMass( p1.mass );
237     SetTotalEnergy( std::sqrt( (p1.mass+p)*(p1    223     SetTotalEnergy( std::sqrt( (p1.mass+p)*(p1.mass+p) - 2.*p1.mass*p ) );
238     //SetTotalEnergy( std::sqrt( p1.mass*p1.ma    224     //SetTotalEnergy( std::sqrt( p1.mass*p1.mass + x*x + y*y + z*z ) );
239     SetMomentum( x, y, z );                       225     SetMomentum( x, y, z );
240   }                                               226   }
241                                                   227  
242  G4double G4ReactionProduct::Angle(               228  G4double G4ReactionProduct::Angle(
243   const G4ReactionProduct& p ) const              229   const G4ReactionProduct& p ) const
244   {                                               230   {
245     G4ThreeVector tM = momentum;                  231     G4ThreeVector tM = momentum;
246     G4ThreeVector pM = p.momentum;                232     G4ThreeVector pM = p.momentum;
247     G4double tx = tM.x(); G4double ty = tM.y()    233     G4double tx = tM.x(); G4double ty = tM.y(); G4double tz = tM.z();
248     G4double px = pM.x(); G4double py = pM.y()    234     G4double px = pM.x(); G4double py = pM.y(); G4double pz = pM.z();
249     G4double a = std::sqrt( ( px*px + py*py +     235     G4double a = std::sqrt( ( px*px + py*py + pz*pz ) * ( tx*tx + ty*ty + tz*tz ) );
250     if( a == 0.0 ) {                              236     if( a == 0.0 ) {
251       return 0.0;                                 237       return 0.0;
252     } else {                                      238     } else {
253       a = ( tx*px + ty*py + tz*pz ) / a;          239       a = ( tx*px + ty*py + tz*pz ) / a;
254       if( std::abs(a) > 1.0 ) { a<0.0 ? a=-1.0    240       if( std::abs(a) > 1.0 ) { a<0.0 ? a=-1.0 : a=1.0; }
255       return std::acos( a );                      241       return std::acos( a );
256     }                                             242     }
257   }                                               243   }
258                                                   244  
259  G4ReactionProduct operator+(                     245  G4ReactionProduct operator+(
260   const G4ReactionProduct& p1, const G4Reactio    246   const G4ReactionProduct& p1, const G4ReactionProduct& p2 )
261   {                                               247   {
262     G4double totEnergy = p1.totalEnergy + p2.t    248     G4double totEnergy = p1.totalEnergy + p2.totalEnergy;
263     G4double x = p1.momentum.x() + p2.momentum    249     G4double x = p1.momentum.x() + p2.momentum.x();
264     G4double y = p1.momentum.y() + p2.momentum    250     G4double y = p1.momentum.y() + p2.momentum.y();
265     G4double z = p1.momentum.z() + p2.momentum    251     G4double z = p1.momentum.z() + p2.momentum.z();
266     G4double newMass = totEnergy*totEnergy - (    252     G4double newMass = totEnergy*totEnergy - ( x*x + y*y + z*z );
267     if( newMass < 0.0 )                           253     if( newMass < 0.0 )
268       newMass = -1. * std::sqrt( -newMass );      254       newMass = -1. * std::sqrt( -newMass );
269     else                                          255     else
270       newMass = std::sqrt( newMass );             256       newMass = std::sqrt( newMass );
271     G4ReactionProduct result;                     257     G4ReactionProduct result;
272     result.SetMass( newMass );                    258     result.SetMass( newMass );
273     result.SetMomentum( x, y, z );                259     result.SetMomentum( x, y, z );
274     result.SetTotalEnergy( totEnergy );           260     result.SetTotalEnergy( totEnergy );
275     result.SetPositionInNucleus( 0.0, 0.0, 0.0    261     result.SetPositionInNucleus( 0.0, 0.0, 0.0 );
276     result.SetFormationTime(0.0);                 262     result.SetFormationTime(0.0);
277     result.HasInitialStateParton(false);          263     result.HasInitialStateParton(false);
278     return result;                                264     return result;
279   }                                               265   }
280                                                   266  
281  G4ReactionProduct operator-(                     267  G4ReactionProduct operator-(
282   const G4ReactionProduct& p1, const G4Reactio    268   const G4ReactionProduct& p1, const G4ReactionProduct& p2 )
283   {                                               269   {
284     G4double totEnergy = p1.totalEnergy - p2.t    270     G4double totEnergy = p1.totalEnergy - p2.totalEnergy;
285     G4double x = p1.momentum.x() - p2.momentum    271     G4double x = p1.momentum.x() - p2.momentum.x();
286     G4double y = p1.momentum.y() - p2.momentum    272     G4double y = p1.momentum.y() - p2.momentum.y();
287     G4double z = p1.momentum.z() - p2.momentum    273     G4double z = p1.momentum.z() - p2.momentum.z();
288     G4double newMass = totEnergy*totEnergy - (    274     G4double newMass = totEnergy*totEnergy - ( x*x + y*y + z*z );
289     if( newMass < 0.0 )                           275     if( newMass < 0.0 )
290       newMass = -1. * std::sqrt( -newMass );      276       newMass = -1. * std::sqrt( -newMass );
291     else                                          277     else
292       newMass = std::sqrt( newMass );             278       newMass = std::sqrt( newMass );
293     G4ReactionProduct result;                     279     G4ReactionProduct result;
294     result.SetMass( newMass );                    280     result.SetMass( newMass );
295     result.SetMomentum( x, y, z );                281     result.SetMomentum( x, y, z );
296     result.SetTotalEnergy( totEnergy );           282     result.SetTotalEnergy( totEnergy );
297     result.SetPositionInNucleus( 0.0, 0.0, 0.0    283     result.SetPositionInNucleus( 0.0, 0.0, 0.0 );
298     result.SetFormationTime(0.0);                 284     result.SetFormationTime(0.0);
299     result.HasInitialStateParton(false);          285     result.HasInitialStateParton(false);
300     return result;                                286     return result;
301   }                                               287   }
302  /* end of code */                                288  /* end of code */
303                                                   289  
304                                                   290