| Geant4 Cross Reference |
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 /// \file G4H1Wrapper.cc
27 /// \brief This allows statistics compatible with FLUKA.
28 /// It is fully relying on G4VAnalysisManager's G4H1.
29 //
30 // Author: G.Hugo, 08 December 2022
31 //
32 // ***************************************************************************
33 //
34 // G4H1Wrapper
35 //
36 /// This allows event-level statistics.
37 /// This is necessary IN THIS SPECIFIC CASE ONLY, to be able to compare errors on results
38 /// with e.g. FLUKA, MCNP, PENELOPE, SERPENT, which all use an event-based error approach.
39 ///
40 /// It is fully relying on G4VAnalysisManager's G4H1.
41 /// Collects info within an event, then flush it to G4VAnalysisManager's G4H1.
42 //
43 // ***************************************************************************
44
45 #include "G4H1Wrapper.hh"
46
47 #include "G4RootAnalysisManager.hh"
48
49 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
50
51 G4H1Wrapper::G4H1Wrapper(G4VAnalysisManager* const analysisManager, const G4int histoIndex)
52 : fAnalysisManager(analysisManager),
53 fHistoIndex(histoIndex),
54 fHisto(dynamic_cast<G4RootAnalysisManager*>(analysisManager)->GetH1(histoIndex)),
55 // fAxis(fHisto->axis()),
56 fNumBins(fHisto->axis().bins()),
57 fIsSet(false),
58 fUnderflow(0.),
59 fOverflow(0.),
60 fSumSquaredEventTotals(0.),
61 fSumSquaredEventInRangeTotals(0.)
62 {
63 fEventData.assign(fNumBins, 0.);
64 }
65
66 // ***************************************************************************
67 // Fill data WITHIN THE EVENT. Can be called an arbitrary number of times.
68 // ***************************************************************************
69 void G4H1Wrapper::Fill(const G4double abscissaValue, const G4double weight)
70 {
71 if (!fIsSet) {
72 fIsSet = true;
73 }
74
75 const G4int binIndex = fHisto->axis().coord_to_index(abscissaValue);
76
77 // underflow
78 // See external/g4tools/include/tools/histo/axis.
79 if (binIndex == tools::histo::axis_UNDERFLOW_BIN) {
80 fUnderflow += weight;
81 }
82 // overflow
83 // See external/g4tools/include/tools/histo/axis.
84 else if (binIndex == tools::histo::axis_OVERFLOW_BIN) {
85 fOverflow += weight;
86 }
87 // in range
88 else {
89 fEventData[binIndex] += weight;
90 }
91 }
92
93 // ***************************************************************************
94 // End of event: flush all data collected during the event
95 // to G4VAnalysisManager's G4H1, then reset event data.
96 // ***************************************************************************
97 void G4H1Wrapper::EndOfEvent()
98 {
99 // Only do the flush if there was at least one fill during the event.
100 if (fIsSet) {
101 G4double eventTotal = 0.;
102 // FLUSH IN RANGE DATA.
103 for (G4int binIndex = 0; binIndex < fNumBins; ++binIndex) {
104 if (fEventData[binIndex] > 0.) {
105 const G4double binKineticEnergy = fHisto->axis().bin_center(binIndex);
106 const G4double binEventWeight = fEventData[binIndex];
107 fAnalysisManager->FillH1(fHistoIndex, binKineticEnergy, binEventWeight);
108 eventTotal += binEventWeight;
109
110 fEventData[binIndex] = 0.;
111 }
112 }
113 // Also keep track of the event's total (INTEGRAL over abscissa range).
114 const G4double squaredEventTotal = eventTotal * eventTotal;
115 fSumSquaredEventInRangeTotals += squaredEventTotal;
116 fSumSquaredEventTotals += squaredEventTotal;
117
118 // FLUSH UNDERFLOW DATA.
119 if (fUnderflow > 0.) {
120 fAnalysisManager->FillH1(fHistoIndex, fHisto->axis().lower_edge() - 1., fUnderflow);
121 fSumSquaredEventTotals += fUnderflow;
122 fUnderflow = 0.;
123 }
124
125 // FLUSH OVERFLOW DATA.
126 if (fOverflow > 0.) {
127 fAnalysisManager->FillH1(fHistoIndex, fHisto->axis().upper_edge() + 1., fOverflow);
128 fSumSquaredEventTotals += fOverflow;
129 fOverflow = 0.;
130 }
131
132 fIsSet = false;
133 }
134 }
135
136 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
137