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