ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4INCLGlobals.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4INCLGlobals.cc
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // INCL++ intra-nuclear cascade model
27 // Alain Boudard, CEA-Saclay, France
28 // Joseph Cugnon, University of Liege, Belgium
29 // Jean-Christophe David, CEA-Saclay, France
30 // Pekka Kaitaniemi, CEA-Saclay, France, and Helsinki Institute of Physics, Finland
31 // Sylvie Leray, CEA-Saclay, France
32 // Davide Mancusi, CEA-Saclay, France
33 //
34 #define INCLXX_IN_GEANT4_MODE 1
35 
36 #include "globals.hh"
37 
38 #include "G4INCLGlobals.hh"
39 #include "G4INCLParticle.hh"
40 #ifdef HAVE_WORDEXP
41 #include <wordexp.h>
42 #endif
43 
44 namespace G4INCL {
45  namespace Math {
46 
47  namespace {
48 
49  // constants for the Gaussian CDF approximation
50  const G4double gcdfa1 = 0.254829592;
51  const G4double gcdfa2 = -0.284496736;
52  const G4double gcdfa3 = 1.421413741;
53  const G4double gcdfa4 = -1.453152027;
54  const G4double gcdfa5 = 1.061405429;
55  const G4double gcdfp = 0.3275911;
56 
57  // constants for the inverse Gaussian CDF approximation
58  const G4double igcdfc1 = 2.515517;
59  const G4double igcdfc2 = 0.802853;
60  const G4double igcdfc3 = 0.010328;
61  const G4double igcdfd1 = 1.432788;
62  const G4double igcdfd2 = 0.189269;
63  const G4double igcdfd3 = 0.001308;
64 
65  G4double inverseGaussianCDFRational(const G4double t) {
66  // Abramowitz and Stegun formula 26.2.23.
67  // The absolute value of the error should be less than 4.5 e-4.
68  return t - ((igcdfc3*t + igcdfc2)*t + igcdfc1) /
69  (((igcdfd3*t + igcdfd2)*t + igcdfd1)*t + 1.0);
70  }
71 
72  }
73 
75  {
76  // Save the sign of x
77  const G4double sgn = sign(x);
78  const G4double z = std::fabs(x) * oneOverSqrtTwo;
79 
80  // A&S formula 7.1.26
81  G4double t = 1.0/(1.0 + gcdfp*z);
82  G4double y = 1.0 - (((((gcdfa5*t + gcdfa4)*t) + gcdfa3)*t + gcdfa2)*t + gcdfa1)*t*std::exp(-z*z);
83 
84  return 0.5*(1.0 + sgn*y);
85  }
86 
87  G4double gaussianCDF(const G4double x, const G4double x0, const G4double sigma) {
88  return gaussianCDF((x-x0)/sigma);
89  }
90 
92  if (x < 0.5)
93  return -inverseGaussianCDFRational( std::sqrt(-2.0*std::log(x)) );
94  else
95  return inverseGaussianCDFRational( std::sqrt(-2.0*std::log(1.-x)) );
96  }
97 
99 // assert(x>-1.000001 && x<1.000001);
100  return ((x > 1.) ? 0. : ((x<-1.) ? pi : std::asin(x)));
101  }
102 
104 // assert(x>-1.000001 && x<1.000001);
105  return ((x > 1.) ? 0. : ((x<-1.) ? pi : std::acos(x)));
106  }
107  }
108 
109  namespace ParticleConfig {
110  G4bool isPair(Particle const * const p1, Particle const * const p2, ParticleType t1, ParticleType t2) {
111  return ((p1->getType() == t1 && p2->getType() == t2) ||
112  (p1->getType() == t2 && p2->getType() == t1));
113  }
114  }
115 
116 #ifndef INCLXX_IN_GEANT4_MODE
117  namespace String {
118  void wrap(std::string &str, const size_t lineLength, const std::string &separators) {
119  const size_t len = str.size();
120  size_t startPos = 0;
121  while(len-startPos > lineLength) { /* Loop checking, 10.07.2015, D.Mancusi */
122  const size_t nextNewline = str.find('\n', startPos);
123  if(nextNewline!=std::string::npos && nextNewline-startPos<=lineLength)
124  startPos = nextNewline+1;
125  else {
126  size_t lastSeparator = str.find_last_of(separators, startPos+lineLength);
127  if(lastSeparator!=std::string::npos)
128  str[lastSeparator] = '\n';
129  startPos = lastSeparator+1;
130  }
131  }
132  }
133 
134  void replaceAll(std::string &str, const std::string &from, const std::string &to, const size_t maxPosition) {
135  if(from.empty())
136  return;
137  size_t start_pos = 0;
138  size_t cur_max_pos = maxPosition;
139  const size_t from_len = from.length();
140  const size_t to_len = to.length();
141  while((start_pos = str.find(from, start_pos)) != std::string::npos /* Loop checking, 10.07.2015, D.Mancusi */
142  && (cur_max_pos==std::string::npos || start_pos<cur_max_pos)) {
143  str.replace(start_pos, from_len, to);
144  start_pos += to_len; // In case 'to' contains 'from', like replacing 'x' with 'yx'
145  if(cur_max_pos!=std::string::npos)
146  cur_max_pos += to_len - from_len;
147  }
148  }
149 
150  std::vector<std::string> tokenize(std::string const &str, const std::string &delimiters) {
151  size_t startPos = 0, endPos;
152  std::vector<std::string> tokens;
153  do {
154  endPos = str.find_first_of(delimiters, startPos);
155  std::string token = str.substr(startPos, endPos-startPos);
156  tokens.push_back(token);
157  startPos = str.find_first_not_of(delimiters, endPos);
158  } while(endPos!=std::string::npos); /* Loop checking, 10.07.2015, D.Mancusi */
159 
160  return tokens;
161  }
162 
163  G4bool isInteger(std::string const &str) {
164  const size_t pos = str.find_first_not_of("0123456789");
165  return (pos==std::string::npos);
166  }
167 
168  std::string expandPath(std::string const &path) {
169 #ifdef HAVE_WORDEXP
170  wordexp_t expansionResult;
171  std::string result;
172  G4int err = wordexp(path.c_str(), &expansionResult, WRDE_NOCMD);
173  if(err)
174  result = path;
175  else
176  result = expansionResult.we_wordv[0];
177  wordfree(&expansionResult);
178  return result;
179 #else
180  // no-op if wordexp.h is not found
181  return path;
182 #endif
183  }
184  }
185 
186 #endif // INCLXX_IN_GEANT4_MODE
187 
188 }