When calculating a*a-b*b the form (a-b)*(a+b) is usually more numerically stable.

[u/mrichter/AliRoot.git] / EVGEN / AliGenGeVSim.cxx

diff --git a/EVGEN/AliGenGeVSim.cxx b/EVGEN/AliGenGeVSim.cxx
index ff744da6bb88166057aeb122099838149a90d6c7..3fdab521e37c4bd2a0b5b6717e64bf1b210101ce 100644 (file)
--- a/EVGEN/AliGenGeVSim.cxx
+++ b/EVGEN/AliGenGeVSim.cxx
@@ -265,11 +265,10 @@ static Double_t aPtYFormula2(Double_t *x, Double_t * par) {
   // mass -> [0] , temperature -> [1] , expansion velocity -> [2]
 
   Double_t aFormE = TMath::Sqrt(par[0]*par[0] + x[0]*x[0]) * TMath::CosH(x[1]);
-  Double_t aFormG = 1 / TMath::Sqrt( 1 - par[2]*par[2] );
+  Double_t aFormG = 1 / TMath::Sqrt((1.-par[2])*(1.+par[2]));
   Double_t aFormYp = par[2]*TMath::Sqrt( (par[0]*par[0] + x[0]*x[0]) 
-                                        * TMath::CosH(x[1])*TMath::CosH(x[1])
-                                        - par[0]*par[0] )
-    /( par[1]*TMath::Sqrt(1-par[2]*par[2]));
+                                        * (TMath::CosH(x[1])-par[0])*(TMath::CosH(x[1])+par[0]))
+    /( par[1]*TMath::Sqrt((1.-par[2])*(1.+par[2])));
 
   return x[0] * aFormE * TMath::Exp( - aFormG * aFormE / par[1])
     *( TMath::SinH(aFormYp)/aFormYp