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