Geant4  10.02.p01
G4P1ToolsManager.cc
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 // $Id$
27 
28 // Author: Ivana Hrivnacova, 24/07/2014 (ivana@ipno.in2p3.fr)
29 
30 #include "G4P1ToolsManager.hh"
32 #include "G4BaseHistoUtilities.hh"
33 #include "G4AnalysisUtilities.hh"
34 
35 #include "tools/histo/p1d"
36 
37 #include <fstream>
38 
39 using namespace G4Analysis;
40 
41 // static definitions
43 
44 //
45 // Constructors, destructor
46 //
47 
48 //_____________________________________________________________________________
50  : G4VP1Manager(),
51  G4THnManager<tools::histo::p1d>(state, "P1")
52 {}
53 
54 //_____________________________________________________________________________
56 {}
57 
58 //
59 // Utility functions
60 //
61 
62 namespace {
63 
64 //_____________________________________________________________________________
65 void UpdateP1Information(G4HnInformation* hnInformation,
66  const G4String& xunitName,
67  const G4String& yunitName,
68  const G4String& xfcnName,
69  const G4String& yfcnName,
70  G4BinScheme xbinScheme)
71 {
72  hnInformation->SetDimension(kX, xunitName, xfcnName, xbinScheme);
73  hnInformation->SetDimension(kY, yunitName, yfcnName, G4BinScheme::kLinear);
74 }
75 
76 //_____________________________________________________________________________
77 void AddP1Annotation(tools::histo::p1d* p1d,
78  const G4String& xunitName,
79  const G4String& yunitName,
80  const G4String& xfcnName,
81  const G4String& yfcnName)
82 {
83  G4String xaxisTitle;
84  G4String yaxisTitle;
85  UpdateTitle(xaxisTitle, xunitName, xfcnName);
86  UpdateTitle(yaxisTitle, yunitName, yfcnName);
87  p1d->add_annotation(tools::histo::key_axis_x_title(), xaxisTitle);
88  p1d->add_annotation(tools::histo::key_axis_y_title(), yaxisTitle);
89 }
90 
91 //_____________________________________________________________________________
92 tools::histo::p1d* CreateToolsP1(const G4String& title,
93  G4int nbins, G4double xmin, G4double xmax,
94  G4double ymin, G4double ymax,
95  const G4String& xunitName,
96  const G4String& yunitName,
97  const G4String& xfcnName,
98  const G4String& yfcnName,
99  const G4String& xbinSchemeName)
100 {
101  auto xunit = GetUnitValue(xunitName);
102  auto yunit = GetUnitValue(yunitName);
103  auto xfcn = GetFunction(xfcnName);
104  auto yfcn = GetFunction(yfcnName);
105  auto xbinScheme = GetBinScheme(xbinSchemeName);
106 
107  if ( xbinScheme != G4BinScheme::kLog ) {
108  if ( xbinScheme == G4BinScheme::kUser ) {
109  // This should never happen, but let's make sure about it
110  // by issuing a warning
111  G4ExceptionDescription description;
112  description
113  << " User binning scheme setting was ignored." << G4endl
114  << " Linear binning will be applied with given (nbins, xmin, xmax) values";
115  G4Exception("G4P1ToolsManager::CreateP1",
116  "Analysis_W013", JustWarning, description);
117  }
118  return new tools::histo::p1d(title,
119  nbins, xfcn(xmin/xunit), xfcn(xmax/xunit),
120  yfcn(ymin/yunit), yfcn(ymax/yunit));
121  }
122  else {
123  // Compute edges
124  std::vector<G4double> edges;
125  ComputeEdges(nbins, xmin, xmax, xunit, xfcn, xbinScheme, edges);
126  return new tools::histo::p1d(title, edges, yfcn(ymin/yunit), yfcn(ymax/yunit));
127  }
128 }
129 
130 //_____________________________________________________________________________
131 tools::histo::p1d* CreateToolsP1(const G4String& title,
132  const std::vector<G4double>& edges,
133  G4double ymin, G4double ymax,
134  const G4String& xunitName,
135  const G4String& yunitName,
136  const G4String& xfcnName,
137  const G4String& yfcnName)
138 {
139  auto xunit = GetUnitValue(xunitName);
140  auto yunit = GetUnitValue(yunitName);
141  auto xfcn = GetFunction(xfcnName);
142  auto yfcn = GetFunction(yfcnName);
143 
144  // Apply function
145  std::vector<G4double> newEdges;
146  ComputeEdges(edges, xunit, xfcn, newEdges);
147 
148  return new tools::histo::p1d(title, newEdges, yfcn(ymin/yunit), yfcn(ymax/yunit));
149 }
150 
151 //_____________________________________________________________________________
152 void ConfigureToolsP1(tools::histo::p1d* p1d,
153  G4int nbins, G4double xmin, G4double xmax,
154  G4double ymin, G4double ymax,
155  const G4String& xunitName,
156  const G4String& yunitName,
157  const G4String& xfcnName,
158  const G4String& yfcnName,
159  const G4String& xbinSchemeName)
160 {
161  auto xunit = GetUnitValue(xunitName);
162  auto yunit = GetUnitValue(yunitName);
163  auto xfcn = GetFunction(xfcnName);
164  auto yfcn = GetFunction(yfcnName);
165  auto xbinScheme = GetBinScheme(xbinSchemeName);
166 
167  if ( xbinScheme != G4BinScheme::kLog ) {
168  if ( xbinScheme == G4BinScheme::kUser ) {
169  // This should never happen, but let's make sure about it
170  // by issuing a warning
171  G4ExceptionDescription description;
172  description
173  << " User binning scheme setting was ignored." << G4endl
174  << " Linear binning will be applied with given (nbins, xmin, xmax) values";
175  G4Exception("G4P1ToolsManager::SetP1",
176  "Analysis_W013", JustWarning, description);
177  }
178  p1d->configure(nbins, xfcn(xmin/xunit), xfcn(xmax/xunit),
179  yfcn(ymin/yunit), yfcn(ymax/yunit));
180  }
181  else {
182  // Compute bins
183  std::vector<G4double> edges;
184  ComputeEdges(nbins, xmin, xmax, xunit, xfcn, xbinScheme, edges);
185  p1d->configure(edges, yfcn(ymin/yunit), yfcn(ymax/yunit));
186  }
187 }
188 
189 //_____________________________________________________________________________
190 void ConfigureToolsP1(tools::histo::p1d* p1d,
191  const std::vector<G4double>& edges,
192  G4double ymin, G4double ymax,
193  const G4String& xunitName,
194  const G4String& yunitName,
195  const G4String& xfcnName,
196  const G4String& yfcnName)
197 {
198  // Apply function to edges
199  auto xunit = GetUnitValue(xunitName);
200  auto yunit = GetUnitValue(yunitName);
201  auto xfcn = GetFunction(xfcnName);
202  auto yfcn = GetFunction(yfcnName);
203  std::vector<G4double> newEdges;
204  ComputeEdges(edges, xunit, xfcn, newEdges);
205 
206  p1d->configure(newEdges, yfcn(ymin/yunit), yfcn(ymax/yunit));
207 }
208 }
209 
210 //
211 // private methods
212 //
213 
214 //_____________________________________________________________________________
216  const G4String& xunitName,
217  const G4String& yunitName,
218  const G4String& xfcnName,
219  const G4String& yfcnName,
220  G4BinScheme xbinScheme) const
221 {
222  auto hnInformation = fHnManager->AddHnInformation(name, 2);
223  hnInformation->AddDimension(xunitName, xfcnName, xbinScheme);
224  hnInformation->AddDimension(yunitName, yfcnName, G4BinScheme::kLinear);
225 }
226 
227 //
228 // protected methods
229 //
230 
231 //_____________________________________________________________________________
233  G4int nbins, G4double xmin, G4double xmax,
234  G4double ymin, G4double ymax,
235  const G4String& xunitName, const G4String& yunitName,
236  const G4String& xfcnName, const G4String& yfcnName,
237  const G4String& xbinSchemeName)
238 {
239 #ifdef G4VERBOSE
240  if ( fState.GetVerboseL4() )
241  fState.GetVerboseL4()->Message("create", "P1", name);
242 #endif
243  tools::histo::p1d* p1d
244  = CreateToolsP1(title, nbins, xmin, xmax, ymin, ymax,
245  xunitName, yunitName, xfcnName, yfcnName,
246  xbinSchemeName);
247 
248  // Add annotation
249  AddP1Annotation(p1d, xunitName, yunitName, xfcnName, yfcnName);
250 
251  // Save P1 information
252  auto xbinScheme = GetBinScheme(xbinSchemeName);
254  name, xunitName, yunitName, xfcnName, yfcnName, xbinScheme);
255 
256  // Register profile
257  G4int id = RegisterT(p1d, name);
258 
259 #ifdef G4VERBOSE
260  if ( fState.GetVerboseL2() )
261  fState.GetVerboseL2()->Message("create", "P1", name);
262 #endif
263  return id;
264 }
265 
266 //_____________________________________________________________________________
268  const std::vector<G4double>& edges,
269  G4double ymin, G4double ymax,
270  const G4String& xunitName, const G4String& yunitName,
271  const G4String& xfcnName, const G4String& yfcnName)
272 {
273 #ifdef G4VERBOSE
274  if ( fState.GetVerboseL4() )
275  fState.GetVerboseL4()->Message("create", "P1", name);
276 #endif
277  tools::histo::p1d* p1d
278  = CreateToolsP1(title, edges, ymin, ymax,
279  xunitName, yunitName, xfcnName, yfcnName);
280 
281  // Add annotation
282  AddP1Annotation(p1d, xunitName, yunitName, xfcnName, yfcnName);
283 
284  // Save P1 information
286  name, xunitName, yunitName, xfcnName, yfcnName, G4BinScheme::kUser);
287 
288  // Register profile
289  G4int id = RegisterT(p1d, name);
290 
291 #ifdef G4VERBOSE
292  if ( fState.GetVerboseL2() )
293  fState.GetVerboseL2()->Message("create", "P1", name);
294 #endif
295  return id;
296 }
297 
298 //_____________________________________________________________________________
300  G4int nbins, G4double xmin, G4double xmax,
301  G4double ymin, G4double ymax,
302  const G4String& xunitName, const G4String& yunitName,
303  const G4String& xfcnName, const G4String& yfcnName,
304  const G4String& xbinSchemeName)
305 {
306  auto p1d = GetTInFunction(id, "SetP1", false, false);
307  if ( ! p1d ) return false;
308 
309  auto info = fHnManager->GetHnInformation(id,"SetP1");
310 #ifdef G4VERBOSE
311  if ( fState.GetVerboseL4() )
312  fState.GetVerboseL4()->Message("configure", "P1", info->GetName());
313 #endif
314 
315  // Configure tools p1
316  ConfigureToolsP1(
317  p1d, nbins, xmin, xmax, ymin, ymax,
318  xunitName, yunitName, xfcnName, yfcnName, xbinSchemeName);
319 
320  // Add annotation
321  AddP1Annotation(p1d, xunitName, yunitName, xfcnName, yfcnName);
322 
323  // Update information
324  auto xbinScheme = GetBinScheme(xbinSchemeName);
325  UpdateP1Information(
326  info, xunitName, yunitName, xfcnName, yfcnName, xbinScheme);
327 
328  // Set activation
329  fHnManager->SetActivation(id, true);
330 
331  return true;
332 }
333 
334 //_____________________________________________________________________________
336  const std::vector<G4double>& edges,
337  G4double ymin, G4double ymax,
338  const G4String& xunitName, const G4String& yunitName,
339  const G4String& xfcnName, const G4String& yfcnName)
340 {
341  auto p1d = GetTInFunction(id, "SetP1", false, false);
342  if ( ! p1d ) return false;
343 
344  auto info = fHnManager->GetHnInformation(id,"SetP1");
345 #ifdef G4VERBOSE
346  if ( fState.GetVerboseL4() )
347  fState.GetVerboseL4()->Message("configure", "P1", info->GetName());
348 #endif
349 
350  // Configure tools p1
351  ConfigureToolsP1(p1d, edges, ymin, ymax,
352  xunitName, yunitName, xfcnName, yfcnName);
353 
354  // Add annotation
355  AddP1Annotation(p1d, xunitName, yunitName, xfcnName, yfcnName);
356 
357  // Update information
358  UpdateP1Information(
359  info, xunitName, yunitName, xfcnName, yfcnName, G4BinScheme::kUser);
360 
361  // Set activation
362  fHnManager->SetActivation(id, true);
363 
364  return true;
365 }
366 
367 
368 //_____________________________________________________________________________
370 {
371  auto p1d = GetTInFunction(id, "ScaleP1", false, false);
372  if ( ! p1d ) return false;
373 
374  return p1d->scale(factor);
375 }
376 
377 //_____________________________________________________________________________
379  G4double weight)
380 {
381  auto p1d = GetTInFunction(id, "FillP1", true, false);
382  if ( ! p1d ) return false;
383 
384  if ( fState.GetIsActivation() && ( ! fHnManager->GetActivation(id) ) ) {
385  //G4cout << "Skipping FillP1 for " << id << G4endl;
386  return false;
387  }
388 
389  auto xInfo
390  = fHnManager->GetHnDimensionInformation(id, kX, "FillP1");
391  auto yInfo
392  = fHnManager->GetHnDimensionInformation(id, kY, "FillP1");
393 
394  p1d->fill(xInfo->fFcn(xvalue/xInfo->fUnit),
395  yInfo->fFcn(yvalue/yInfo->fUnit), weight);
396 
397 #ifdef G4VERBOSE
398  if ( fState.GetVerboseL4() ) {
399  G4ExceptionDescription description;
400  description << " id " << id
401  << " xvalue " << xvalue
402  << " xfcn(xvalue/xunit) " << xInfo->fFcn(xvalue/xInfo->fUnit)
403  << " yvalue " << yvalue
404  << " yfcn(yvalue/yunit) " << yInfo->fFcn(yvalue/yInfo->fUnit)
405  << " weight " << weight;
406  fState.GetVerboseL4()->Message("fill", "P1", description);
407  }
408 #endif
409  return true;
410 }
411 
412 //_____________________________________________________________________________
414 {
415  return GetTId(name, warn);
416 }
417 
418 //_____________________________________________________________________________
420 {
421  auto p1d = GetTInFunction(id, "GetP1Nbins");
422  if ( ! p1d ) return 0;
423 
424  return GetNbins(*p1d, kX);
425 }
426 
427 //_____________________________________________________________________________
429 {
430 // Returns xmin value with applied unit and profile function
431 
432  auto p1d = GetTInFunction(id, "GetP1Xmin");
433  if ( ! p1d ) return 0;
434 
435  return GetMin(*p1d, kX);
436 }
437 
438 //_____________________________________________________________________________
440 {
441  auto p1d = GetTInFunction(id, "GetP1Xmax");
442  if ( ! p1d ) return 0;
443 
444  return GetMax(*p1d, kX);
445 }
446 
447 //_____________________________________________________________________________
449 {
450  auto p1d = GetTInFunction(id, "GetP1XWidth", true, false);
451  if ( ! p1d ) return 0;
452 
453  return GetWidth(*p1d, kX, fHnManager->GetHnType());
454 }
455 
456 //_____________________________________________________________________________
458 {
459 // Returns xmin value with applied unit and profile function
460 
461  auto p1d = GetTInFunction(id, "GetP1Ymin");
462  if ( ! p1d ) return 0;
463 
464  return p1d->min_v();
465 }
466 
467 //_____________________________________________________________________________
469 {
470  auto p1d = GetTInFunction(id, "GetP1Ymax");
471  if ( ! p1d ) return 0;
472 
473  return p1d->max_v();
474 }
475 
476 //_____________________________________________________________________________
478 {
479  auto p1d = GetTInFunction(id, "SetP1Title");
480  if ( ! p1d ) return false;
481 
482  return SetTitle(*p1d, title);
483 }
484 
485 //_____________________________________________________________________________
487 {
488  auto p1d = GetTInFunction(id, "SetP1XAxisTitle");
489  if ( ! p1d ) return false;
490 
491  return SetAxisTitle(*p1d, kX, title);
492 }
493 
494 //_____________________________________________________________________________
496 {
497  auto p1d = GetTInFunction(id, "SetP1YAxisTitle");
498  if ( ! p1d ) return false;
499 
500  return SetAxisTitle(*p1d, kY, title);
501 }
502 
503 //_____________________________________________________________________________
505 {
506  auto p1d = GetTInFunction(id, "GetP1Title");
507  if ( ! p1d ) return "";
508 
509  return GetTitle(*p1d);
510 }
511 
512 
513 //_____________________________________________________________________________
515 {
516  auto p1d = GetTInFunction(id, "GetP1XAxisTitle");
517  if ( ! p1d ) return "";
518 
519  return GetAxisTitle(*p1d, kX, fHnManager->GetHnType());
520 }
521 
522 //_____________________________________________________________________________
524 {
525  auto p1d = GetTInFunction(id, "GetP1YAxisTitle");
526  if ( ! p1d ) return "";
527 
528  return GetAxisTitle(*p1d, kY, fHnManager->GetHnType());
529 }
530 
531 //
532 // public methods
533 //
534 
535 //_____________________________________________________________________________
536 G4int G4P1ToolsManager::AddP1(const G4String& name, tools::histo::p1d* p1d)
537 {
538 #ifdef G4VERBOSE
539  if ( fState.GetVerboseL4() )
540  fState.GetVerboseL4()->Message("add", "P1", name);
541 #endif
542 
543  // Add annotation
544  AddP1Annotation(p1d, "none", "none", "none", "none");
545  // Add information
546  AddP1Information(name, "none", "none", "none", "none", G4BinScheme::kLinear);
547 
548  // Register profile
549  auto id = RegisterT(p1d, name);
550 
551 #ifdef G4VERBOSE
552  if ( fState.GetVerboseL2() )
553  fState.GetVerboseL2()->Message("add", "P1", name);
554 #endif
555  return id;
556 }
557 
558 //_____________________________________________________________________________
560  const std::vector<tools::histo::p1d*>& p1Vector)
561 {
562  AddTVector(p1Vector);
563 }
564 
565 //_____________________________________________________________________________
566 tools::histo::p1d* G4P1ToolsManager::GetP1(G4int id, G4bool warn,
567  G4bool onlyIfActive) const
568 {
569  return GetTInFunction(id, "GetP1", warn, onlyIfActive);
570 }
571 
void AddP1Information(const G4String &name, const G4String &xunitName, const G4String &yunitName, const G4String &xfcnName, const G4String &yfcnName, G4BinScheme xbinScheme) const
G4double GetWidth(const G4ToolsBaseHisto &baseHisto, G4int dimension, const G4String &hnType)
virtual ~G4P1ToolsManager()
void Message(const G4String &action, const G4String &object, const G4String &objectName, G4bool success=true) const
G4String GetTitle(const G4ToolsBaseHisto &baseHisto)
G4String GetAxisTitle(const G4ToolsBaseHisto &baseHisto, G4int dimension, const G4String &hnType)
virtual G4double GetP1Xmin(G4int id) const final
G4double GetMin(const G4ToolsBaseHisto &baseHisto, G4int dimension)
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
void ComputeEdges(G4int nbins, G4double xmin, G4double xmax, G4double unit, G4Fcn fcn, G4BinScheme, std::vector< G4double > &edges)
Definition: G4BinScheme.cc:56
virtual G4int GetP1Id(const G4String &name, G4bool warn=true) const final
virtual G4String GetP1Title(G4int id) const final
G4String name
Definition: TRTMaterials.hh:40
virtual G4bool ScaleP1(G4int id, G4double factor) final
void AddP1Vector(const std::vector< tools::histo::p1d * > &p1Vector)
virtual G4bool FillP1(G4int id, G4double xvalue, G4double yvalue, G4double weight=1.0) final
virtual G4bool SetP1XAxisTitle(G4int id, const G4String &title) final
G4bool SetAxisTitle(G4ToolsBaseHisto &baseHisto, G4int dimension, const G4String &title)
virtual G4double GetP1Xmax(G4int id) const final
int G4int
Definition: G4Types.hh:78
std::shared_ptr< G4HnManager > fHnManager
Definition: G4THnManager.hh:83
virtual G4double GetP1Ymax(G4int id) const final
G4int RegisterT(tools::histo::p1d *t, const G4String &name)
const G4AnalysisVerbose * GetVerboseL2() const
void UpdateTitle(G4String &title, const G4String &unitName, const G4String &fcnName)
virtual G4bool SetP1YAxisTitle(G4int id, const G4String &title) final
const G4AnalysisManagerState & fState
Definition: G4THnManager.hh:80
G4int GetTId(const G4String &name, G4bool warn=true) const
const G4AnalysisVerbose * GetVerboseL4() const
bool G4bool
Definition: G4Types.hh:79
const G4int kX
tools::histo::p1d * GetTInFunction(G4int id, G4String functionName, G4bool warn=true, G4bool onlyIfActive=true) const
virtual G4bool SetP1Title(G4int id, const G4String &title) final
G4double GetUnitValue(const G4String &unit)
virtual G4double GetP1XWidth(G4int id) const final
G4bool SetTitle(G4ToolsBaseHisto &baseHisto, const G4String &title)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4int AddP1(const G4String &name, tools::histo::p1d *p1d)
static const G4double factor
virtual G4bool SetP1(G4int id, G4int nbins, G4double xmin, G4double xmax, G4double ymin=0, G4double ymax=0, const G4String &xunitName="none", const G4String &yunitName="none", const G4String &xfcnName="none", const G4String &yfcnName="none", const G4String &xbinScheme="linear") final
G4Fcn GetFunction(const G4String &fcnName)
Definition: G4Fcn.cc:36
virtual G4int GetP1Nbins(G4int id) const final
void AddTVector(const std::vector< tools::histo::p1d * > &tVector)
G4BinScheme GetBinScheme(const G4String &binSchemeName)
Definition: G4BinScheme.cc:36
virtual G4String GetP1XAxisTitle(G4int id) const final
#define G4endl
Definition: G4ios.hh:61
G4double GetMax(const G4ToolsBaseHisto &baseHisto, G4int dimension)
G4BinScheme
Definition: G4BinScheme.hh:40
static const G4int kDimension
G4int GetNbins(const G4ToolsBaseHisto &baseHisto, G4int dimension)
double G4double
Definition: G4Types.hh:76
void AddDimension(const G4String &unitName, const G4String &fcnName, G4BinScheme binScheme)
const G4int kY
virtual G4String GetP1YAxisTitle(G4int id) const final
G4P1ToolsManager(const G4AnalysisManagerState &state)
virtual G4double GetP1Ymin(G4int id) const final
virtual G4int CreateP1(const G4String &name, const G4String &title, G4int nbins, G4double xmin, G4double xmax, G4double ymin=0, G4double ymax=0, const G4String &xunitName="none", const G4String &yunitName="none", const G4String &xfcnName="none", const G4String &yfcnName="none", const G4String &xbinScheme="linear") final
void SetDimension(G4int dimension, const G4String &unitName, const G4String &fcnName, G4BinScheme binScheme)
tools::histo::p1d * GetP1(G4int id, G4bool warn=true, G4bool onlyIfActive=true) const