Geant4-11
G4StatAnalysis.icc
Go to the documentation of this file.
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// G4StatAnalysis inline methods implementation
27//
28// Author: J.Madsen, 25.10.2018
29// --------------------------------------------------------------------
30
31#include <cmath>
32#include <fstream>
33#include <iomanip>
34#include <iostream>
35#include <limits>
36
37#include "G4Timer.hh"
38#include "G4ios.hh"
39#include "globals.hh"
40#include "tls.hh"
41
42#include "G4Allocator.hh"
43#include "G4Types.hh"
44
45//----------------------------------------------------------------------------//
46
47G4StatAnalysis::G4StatAnalysis() {}
48
49//----------------------------------------------------------------------------//
50
51G4double G4StatAnalysis::GetMean() const
52{
53 return (fHits > 0) ? fSum1 / ((G4double) fHits) : 0.;
54}
55
56//----------------------------------------------------------------------------//
57
58const G4double& G4StatAnalysis::GetSum() const { return fSum1; }
59
60//----------------------------------------------------------------------------//
61
62const G4double& G4StatAnalysis::GetSumSquared() const { return fSum2; }
63
64//----------------------------------------------------------------------------//
65
66const G4double& G4StatAnalysis::GetSum1() const { return fSum1; }
67
68//----------------------------------------------------------------------------//
69
70const G4double& G4StatAnalysis::GetSum2() const { return fSum2; }
71
72//----------------------------------------------------------------------------//
73
74const G4int& G4StatAnalysis::GetHits() const { return fHits; }
75
76//----------------------------------------------------------------------------//
77
78G4int G4StatAnalysis::GetNumNonZero() const { return fHits - fZero; }
79
80//----------------------------------------------------------------------------//
81
82G4int G4StatAnalysis::GetNumZero() const { return fZero; }
83
84//----------------------------------------------------------------------------//
85
86void G4StatAnalysis::SetSum(const G4double& val) { fSum1 = val; }
87
88//----------------------------------------------------------------------------//
89
90void G4StatAnalysis::SetSumSquared(const G4double& val) { fSum2 = val; }
91
92//----------------------------------------------------------------------------//
93
94void G4StatAnalysis::SetSum1(const G4double& val) { fSum1 = val; }
95
96//----------------------------------------------------------------------------//
97
98void G4StatAnalysis::SetSum2(const G4double& val) { fSum2 = val; }
99
100//----------------------------------------------------------------------------//
101
102void G4StatAnalysis::SetHits(const G4int& val) { fHits = val; }
103
104//----------------------------------------------------------------------------//
105
106void G4StatAnalysis::SetZero(const G4int& val) { fZero = val; }
107
108//----------------------------------------------------------------------------//
109
110G4double G4StatAnalysis::GetFOM() const
111{
112 G4double elapsed_time = this->GetCpuTime();
113 G4double relative_err = this->GetRelativeError();
114 // lambda for equation clarity (will be inlined)
115 auto compute_figure_of_merit = [&]() {
116 return (1.0 / (relative_err * relative_err) / elapsed_time);
117 };
118 return (std::fabs(relative_err) > 0.0 && elapsed_time > 0.0)
119 ? compute_figure_of_merit()
120 : ((fHits > 0) ? 1.0 : 0.0);
121}
122
123//----------------------------------------------------------------------------//
124
125G4double G4StatAnalysis::GetRelativeError() const
126{
127 // lambda for equation clarity (will be inlined)
128 auto compute_relative_error = [&]() {
129 return (GetStdDev() / GetMean() / std::sqrt((G4double) fHits));
130 };
131 return (std::fabs(GetMean()) > 0 && fHits > 0) ? compute_relative_error()
132 : ((fHits > 0) ? 1.0 : 0.0);
133}
134
135//----------------------------------------------------------------------------//
136
137G4double G4StatAnalysis::GetStdDev() const
138{
139 return ::sqrt(std::fabs(GetVariance()));
140}
141
142//----------------------------------------------------------------------------//
143
144G4double G4StatAnalysis::GetVariance() const
145{
146 // lambda for equation clarity (will be inlined)
147 auto compute_variance = [&]() {
148 return ((fSum2 - (std::pow(fSum1, 2.0) / fHits)) /
149 (((G4double) fHits) - 1.0));
150 };
151 return (fHits > 1) ? compute_variance() : 0.0;
152}
153
154//----------------------------------------------------------------------------//
155
156G4double G4StatAnalysis::GetCoeffVariation() const
157{
158 // lambda for equation clarity (will be inlined)
159 auto coefficient_of_variation = [&]() {
160 G4double hits = fHits;
161 return ::sqrt((hits / (hits - 1.0)) *
162 ((fSum2 / (fSum1 * fSum1)) - (1.0 / hits)));
163 };
164 return (fHits > 1) ? coefficient_of_variation() : 0.0;
165 // return (fHits > 0 && fabs(fSum1) > 0.0)
166 // ? (100.0*GetStdDev()/GetMean()) : 0.0;
167}
168
169//----------------------------------------------------------------------------//
170
171G4double G4StatAnalysis::GetEfficiency() const
172{
173 G4double hits = fHits;
174 G4double nzero = fHits - fZero;
175 return (fHits > 0) ? (nzero / hits) : 0.0;
176}
177
178//----------------------------------------------------------------------------//
179
180G4double G4StatAnalysis::GetR2Int() const
181{
182 G4double hits = fHits;
183 return (fHits > 0)
184 ? (fSum2 / (fSum1 * fSum1)) - 1.0 / (GetEfficiency() * hits)
185 : 0.0;
186}
187
188//----------------------------------------------------------------------------//
189
190G4double G4StatAnalysis::GetR2Eff() const
191{
192 G4double hits = fHits;
193 return (fHits > 0) ? (1.0 - GetEfficiency()) / (GetEfficiency() * hits) : 0.0;
194}
195
196//----------------------------------------------------------------------------//
197
198G4StatAnalysis::operator G4double() const { return this->GetSum(); }
199
200//----------------------------------------------------------------------------//
201
202void G4StatAnalysis::Reset()
203{
204 fHits = 0;
205 fZero = 0;
206 fSum1 = 0.0;
207 fSum2 = 0.0;
208}
209
210//----------------------------------------------------------------------------//
211
212void G4StatAnalysis::Add(const G4double& val, const G4double& weight)
213{
214 fHits += 1;
215 fSum1 += val * weight;
216 fSum2 += val * val * weight;
217 if(std::fabs(val * weight) <
218 std::fabs(GetMean() * std::numeric_limits<double>::epsilon()))
219 fZero += 1;
220}
221
222//----------------------------------------------------------------------------//
223
224void G4StatAnalysis::Rescale(const G4double& factor)
225{
226 fSum1 *= factor;
227 fSum2 *= factor * factor;
228}
229
230//----------------------------------------------------------------------------//
231
232void G4StatAnalysis::PrintInfo(std::ostream& os, const std::string& tab) const
233{
234 G4int _hits = this->GetHits();
235 G4double _sum = this->GetSum();
236 G4double _sigma = this->GetStdDev();
237 G4double _coeff = this->GetCoeffVariation();
238 G4double _error = this->GetRelativeError();
239 G4double _eff = this->GetEfficiency();
240 G4double _fom = this->GetFOM();
241 G4double _r2int = this->GetR2Int();
242 G4double _r2eff = this->GetR2Eff();
243
244 using std::fixed;
245 using std::ios;
246 using std::left;
247 using std::right;
248 using std::scientific;
249 using std::setprecision;
250 using std::setw;
251
252 std::stringstream ss;
253 ss << tab //<< scientific
254 << setprecision(os.precision()) << right << _sum << left
255 << " [sigma: " << right << _sigma << left << " | error: " << right
256 << _error << left << " | coeff: " << right << _coeff << left
257 << " | eff: " << right << _eff << left << " | fom: " << right << _fom
258 << left << " | r2int: " << right << _r2int << left << " | r2eff: " << right
259 << _r2eff << left << " | hits: " << right << _hits << left << " ]";
260
261 os << ss.str();
262}
263
264//----------------------------------------------------------------------------//
265
266G4double G4StatAnalysis::GetCpuTime() const
267{
268 tms* startTime = GetCpuClock();
269 tms endTime;
270 times(&endTime);
271 return ((endTime.tms_stime - startTime->tms_stime) +
272 (endTime.tms_utime - startTime->tms_utime)) /
273 sysconf(_SC_CLK_TCK);
274}
275
276//----------------------------------------------------------------------------//
277
278G4StatAnalysis& G4StatAnalysis::operator+=(const G4double& _val)
279{
280 this->Add(_val);
281 return *this;
282}
283
284//----------------------------------------------------------------------------//
285
286G4StatAnalysis& G4StatAnalysis::operator/=(const G4double& _val)
287{
288 fSum1 /= _val;
289 fSum2 /= (_val * _val);
290 return *this;
291}
292
293//----------------------------------------------------------------------------//
294
295G4StatAnalysis& G4StatAnalysis::operator+=(const G4StatAnalysis& rhs)
296{
297 fHits += rhs.fHits;
298 fSum1 += rhs.fSum1;
299 fSum2 += rhs.fSum2;
300 fZero += rhs.fZero;
301 return *this;
302}
303
304//----------------------------------------------------------------------------//
305
306G4StatAnalysis& G4StatAnalysis::operator-=(const G4StatAnalysis& rhs)
307{
308 fHits -= rhs.fHits;
309 fSum1 -= rhs.fSum1;
310 fSum2 -= rhs.fSum2;
311 fZero -= rhs.fZero;
312 return *this;
313}
314
315//----------------------------------------------------------------------------//
316
317inline G4Allocator<G4StatAnalysis>*& _aStatAnalysisAllocator_G4MT_TLS_()
318{
319 G4ThreadLocalStatic G4Allocator<G4StatAnalysis>* _instance =
320 new G4Allocator<G4StatAnalysis>();
321 return _instance;
322}
323
324//----------------------------------------------------------------------------//
325
326void* G4StatAnalysis::operator new(std::size_t)
327{
328 G4Allocator<G4StatAnalysis>& _allocator =
329 *_aStatAnalysisAllocator_G4MT_TLS_();
330 return (void*) _allocator.MallocSingle();
331}
332
333//----------------------------------------------------------------------------//
334
335void G4StatAnalysis::operator delete(void* _ptr)
336{
337 G4Allocator<G4StatAnalysis>& _allocator =
338 *_aStatAnalysisAllocator_G4MT_TLS_();
339 _allocator.FreeSingle((G4StatAnalysis*) _ptr);
340}
341
342//----------------------------------------------------------------------------//
Generated by doxygen 1.9.3