| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer << 3 // * DISCLAIMER *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th << 5 // * The following disclaimer summarizes all the specific disclaimers *
6 // * the Geant4 Collaboration. It is provided << 6 // * of contributors to this software. The specific disclaimers,which *
7 // * conditions of the Geant4 Software License << 7 // * govern, are listed with their locations in: *
8 // * LICENSE and available at http://cern.ch/ << 8 // * http://cern.ch/geant4/license *
9 // * include a list of copyright holders. <<
10 // * 9 // * *
11 // * Neither the authors of this software syst 10 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 11 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 12 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 13 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file << 14 // * use. *
16 // * for the full disclaimer and the limitatio <<
17 // * 15 // * *
18 // * This code implementation is the result << 16 // * This code implementation is the intellectual property of the *
19 // * technical work of the GEANT4 collaboratio << 17 // * GEANT4 collaboration. *
20 // * By using, copying, modifying or distri << 18 // * By copying, distributing or modifying the Program (or any work *
21 // * any work based on the software) you ag << 19 // * based on the Program) you indicate your acceptance of this *
22 // * use in resulting scientific publicati << 20 // * statement, and all its terms. *
23 // * acceptance of all terms of the Geant4 Sof <<
24 // ******************************************* 21 // ********************************************************************
25 // 22 //
26 // << 23 // $Id: HistoManager.cc,v 1.10 2005/06/01 10:48:40 maire Exp $
>> 24 // GEANT4 tag $Name: geant4-07-01 $
>> 25 //
27 //....oooOO0OOooo........oooOO0OOooo........oo 26 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
28 //....oooOO0OOooo........oooOO0OOooo........oo 27 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
29 28
30 #include "HistoManager.hh" 29 #include "HistoManager.hh"
31 << 30 #include "HistoMessenger.hh"
32 #include "DetectorConstruction.hh" <<
33 <<
34 #include "G4UnitsTable.hh" 31 #include "G4UnitsTable.hh"
35 32
>> 33 #ifdef G4ANALYSIS_USE
>> 34 #include "AIDA/AIDA.h"
>> 35 #endif
>> 36
36 //....oooOO0OOooo........oooOO0OOooo........oo 37 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
37 38
38 HistoManager::HistoManager() 39 HistoManager::HistoManager()
>> 40 :af(0),tree(0),factoryOn(false)
39 { 41 {
40 Book(); << 42 #ifdef G4ANALYSIS_USE
>> 43 // Creating the analysis factory
>> 44 af = AIDA_createAnalysisFactory();
>> 45 if(!af) {
>> 46 G4cout << "TestEm1::HistoManager::HistoManager :"
>> 47 << " problem creating the AIDA analysis factory."
>> 48 << G4endl;
>> 49 }
>> 50 #endif
>> 51
>> 52 fileName[0] = "testem3";
>> 53 fileType = "hbook";
>> 54 fileOption = "--noErrors uncompress";
>> 55 // histograms
>> 56 for (G4int k=0; k<MaxHisto; k++) {
>> 57 histo[k] = 0;
>> 58 exist[k] = false;
>> 59 Unit[k] = 1.0;
>> 60 Width[k] = 1.0;
>> 61 }
>> 62 histoMessenger = new HistoMessenger(this);
>> 63 }
>> 64
>> 65 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
>> 66
>> 67 HistoManager::~HistoManager()
>> 68 {
>> 69 delete histoMessenger;
>> 70
>> 71 #ifdef G4ANALYSIS_USE
>> 72 delete af;
>> 73 #endif
41 } 74 }
42 75
43 //....oooOO0OOooo........oooOO0OOooo........oo 76 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
44 77
45 void HistoManager::Book() << 78 void HistoManager::book()
>> 79 {
>> 80 #ifdef G4ANALYSIS_USE
>> 81 if(!af) return;
>> 82
>> 83 // Creating a tree mapped to an hbook file.
>> 84 fileName[1] = fileName[0] + "." + fileType;
>> 85 G4bool readOnly = false;
>> 86 G4bool createNew = true;
>> 87 AIDA::ITreeFactory* tf = af->createTreeFactory();
>> 88 tree = tf->create(fileName[1], fileType, readOnly, createNew, fileOption);
>> 89 delete tf;
>> 90 if(!tree) {
>> 91 G4cout << "TestEm1::HistoManager::book :"
>> 92 << " problem creating the AIDA tree with "
>> 93 << " storeName = " << fileName[1]
>> 94 << " storeType = " << fileType
>> 95 << " readOnly = " << readOnly
>> 96 << " createNew = " << createNew
>> 97 << " options = " << fileOption
>> 98 << G4endl;
>> 99 return;
>> 100 }
>> 101
>> 102 // Creating a histogram factory, whose histograms will be handled by the tree
>> 103 AIDA::IHistogramFactory* hf = af->createHistogramFactory(*tree);
>> 104
>> 105 // create selected histograms
>> 106 for (G4int k=1; k<MaxHisto; k++) {
>> 107 if (exist[k]) {
>> 108 histo[k] = hf->createHistogram1D( Label[k], Title[k],
>> 109 Nbins[k], Vmin[k], Vmax[k]);
>> 110 factoryOn = true;
>> 111 }
>> 112 }
>> 113 delete hf;
>> 114 if (factoryOn)
>> 115 G4cout << "\n----> Histogram Tree is opened in " << fileName[1] << G4endl;
>> 116 #endif
>> 117 }
>> 118
>> 119 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
>> 120
>> 121 void HistoManager::save()
>> 122 {
>> 123 #ifdef G4ANALYSIS_USE
>> 124 if (factoryOn) {
>> 125 tree->commit(); // Writing the histograms to the file
>> 126 tree->close(); // and closing the tree (and the file)
>> 127 G4cout << "\n----> Histogram Tree is saved in " << fileName[1] << G4endl;
>> 128
>> 129 delete tree;
>> 130 tree = 0;
>> 131 factoryOn = false;
>> 132 }
>> 133 #endif
>> 134 }
>> 135
>> 136 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
>> 137
>> 138 void HistoManager::FillHisto(G4int ih, G4double xbin, G4double weight)
46 { 139 {
47 // Create or get analysis manager << 140 if (ih >= MaxHisto) {
48 // The choice of analysis technology is done << 141 G4cout << "---> warning from HistoManager::FillHisto() : histo " << ih
49 // in HistoManager.hh << 142 << G4endl;
50 G4AnalysisManager* analysisManager = G4Analy << 143 return;
51 analysisManager->SetDefaultFileType("root"); << 144 }
52 analysisManager->SetFileName(fFileName); << 145 #ifdef G4ANALYSIS_USE
53 analysisManager->SetVerboseLevel(1); << 146 if(exist[ih]) histo[ih]->fill(xbin/Unit[ih], weight);
54 analysisManager->SetActivation(true); // en << 147 #endif
55 << 148 }
56 // Define histograms start values << 149
57 << 150 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
58 const G4String id[] = {"0", "1", "2", "3" << 151
59 "9", "10", "11", "12 << 152 void HistoManager::SetHisto(G4int ih,
60 "18", "19", "20", "21 << 153 G4int nbins, G4double valmin, G4double valmax, const G4String& unit)
>> 154 {
>> 155 if (ih < 1 || ih >= MaxHisto) {
>> 156 G4cout << "---> warning from HistoManager::SetHisto() : histo " << ih
>> 157 << "does not exist" << G4endl;
>> 158 return;
>> 159 }
>> 160
>> 161 // histo 1 : energy deposit in absorber 1
>> 162 // histo 2 : energy deposit in absorber 2
>> 163 // ...etc...........
>> 164 // MaxAbsor = 10 (-1)
>> 165 //
>> 166 // histo 11 : longitudinal profile of energy deposit in absorber 1 (MeV)
>> 167 // histo 12 : longitudinal profile of energy deposit in absorber 2 (MeV)
>> 168 // ...etc...........
>> 169 //
>> 170 // histo 21 : forward energy flow (MeV)
>> 171 // histo 22 : backward energy flow (MeV)
>> 172
>> 173 const G4String id[] = { "0", "1", "2", "3", "4", "5", "6", "7", "8", "9",
>> 174 "10","11","12","13","14","15","16","17","18","19",
>> 175 "20","21","22"};
>> 176
61 G4String title; 177 G4String title;
>> 178 G4double vmin = valmin, vmax = valmax;
>> 179
>> 180 if (ih < MaxAbsor) {
>> 181 title = "Edep in absorber " + id[ih] + " (" + unit + ")";
>> 182 Unit[ih] = G4UnitDefinition::GetValueOf(unit);
>> 183 vmin = valmin/Unit[ih]; vmax = valmax/Unit[ih];
>> 184 } else if (ih > MaxAbsor && ih < 2*MaxAbsor) {
>> 185 title = "longit. profile of Edep (MeV/event) in absorber "
>> 186 + id[ih-MaxAbsor];
>> 187 } else if (ih == 2*MaxAbsor+1) {
>> 188 title = "Forward energy flow (MeV/event)";
>> 189 } else if (ih == 2*MaxAbsor+2) {
>> 190 title = "Backward energy flow (MeV/event)";
>> 191 } else return;
>> 192
>> 193 exist[ih] = true;
>> 194 Label[ih] = id[ih];
>> 195 Title[ih] = title;
>> 196 Nbins[ih] = nbins;
>> 197 Vmin[ih] = vmin;
>> 198 Vmax[ih] = vmax;
>> 199 Width[ih] = (valmax-valmin)/nbins;
>> 200
>> 201 G4cout << "----> SetHisto " << ih << ": " << title << "; "
>> 202 << nbins << " bins from "
>> 203 << vmin << " " << unit << " to " << vmax << " " << unit << G4endl;
>> 204
>> 205 }
>> 206
>> 207 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
>> 208
>> 209 void HistoManager::Normalize(G4int ih, G4double fac)
>> 210 {
>> 211 if (ih >= MaxHisto) {
>> 212 G4cout << "---> warning from HistoManager::Normalize() : histo " << ih
>> 213 << G4endl;
>> 214 return;
>> 215 }
>> 216 #ifdef G4ANALYSIS_USE
>> 217 if(exist[ih]) histo[ih]->scale(fac);
>> 218 #endif
>> 219 }
>> 220
>> 221 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
62 222
63 // Default values (to be reset via /analysis << 223 void HistoManager::RemoveHisto(G4int ih)
64 G4int nbins = 100; << 224 {
65 G4double vmin = 0.; << 225 if (ih >= MaxHisto) {
66 G4double vmax = 100.; << 226 G4cout << "---> warning from HistoManager::RemoveHisto() : histo " << ih
67 << 227 << "does not exist" << G4endl;
68 // Create all histograms as inactivated << 228 return;
69 // as we have not yet set nbins, vmin, vmax <<
70 for (G4int k = 0; k < kMaxHisto; k++) { <<
71 if (k < kMaxAbsor) title = "Edep in absorb <<
72 if (k == kMaxAbsor) title = "total energy <<
73 if (k > kMaxAbsor) title = "Edep longit. p <<
74 if (k == 2 * kMaxAbsor + 1) title = "energ <<
75 if (k == 2 * kMaxAbsor + 2) title = "later <<
76 if (k == 2 * kMaxAbsor + 3) title = "total <<
77 if (k == 2 * kMaxAbsor + 4) title = "total <<
78 G4int ih = analysisManager->CreateH1(id[k] <<
79 analysisManager->SetH1Activation(ih, false <<
80 } 229 }
>> 230
>> 231 histo[ih] = 0; exist[ih] = false;
81 } 232 }
>> 233
>> 234 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
>> 235
>> 236
82 237