| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 // 25 //
26 // << 26 // $Id: GammaRayTelAnalysis.cc 82268 2014-06-13 13:47:30Z gcosmo $
27 // ------------------------------------------- 27 // ------------------------------------------------------------
28 // GEANT 4 class implementation file 28 // GEANT 4 class implementation file
29 // CERN Geneva Switzerland 29 // CERN Geneva Switzerland
30 // 30 //
31 // 31 //
32 // ------------ GammaRayAnalysisManager 32 // ------------ GammaRayAnalysisManager ------
33 // by R.Giannitrapani, F.Longo & G.S 33 // by R.Giannitrapani, F.Longo & G.Santin (03 dic 2000)
34 // 34 //
35 // 03.04.2013 F.Longo/L.Pandola 35 // 03.04.2013 F.Longo/L.Pandola
36 // - migrated to G4tools 36 // - migrated to G4tools
37 // 37 //
38 // 29.05.2003 F.Longo 38 // 29.05.2003 F.Longo
39 // - Anaphe 5.0.5 compliant << 39 // - anaphe 5.0.5 compliant
40 // 40 //
41 // 18.06.2002 R.Giannitrapani, F.Longo & G.San 41 // 18.06.2002 R.Giannitrapani, F.Longo & G.Santin
42 // - new release for Anaphe 4.0.3 42 // - new release for Anaphe 4.0.3
43 // 43 //
44 // 07.12.2001 A.Pfeiffer 44 // 07.12.2001 A.Pfeiffer
45 // - integrated Guy's addition of the ntuple 45 // - integrated Guy's addition of the ntuple
46 // 46 //
47 // 06.12.2001 A.Pfeiffer 47 // 06.12.2001 A.Pfeiffer
48 // - updating to new design (singleton) 48 // - updating to new design (singleton)
49 // 49 //
50 // 22.11.2001 G.Barrand 50 // 22.11.2001 G.Barrand
51 // - Adaptation to AIDA 51 // - Adaptation to AIDA
52 // 52 //
53 // ******************************************* 53 // ************************************************************
54 <<
55 #include <fstream> 54 #include <fstream>
56 #include <iomanip> 55 #include <iomanip>
57 56
>> 57 #include "G4RunManager.hh"
>> 58
58 #include "GammaRayTelAnalysis.hh" 59 #include "GammaRayTelAnalysis.hh"
59 #include "GammaRayTelAnalysisMessenger.hh" <<
60 #include "GammaRayTelDetectorConstruction.hh" 60 #include "GammaRayTelDetectorConstruction.hh"
>> 61 #include "GammaRayTelAnalysisMessenger.hh"
61 62
62 #include "G4RunManager.hh" << 63 GammaRayTelAnalysis* GammaRayTelAnalysis::instance = 0;
63 <<
64 //....oooOO0OOooo........oooOO0OOooo........oo <<
65 <<
66 GammaRayTelAnalysis *GammaRayTelAnalysis::inst <<
67 64
>> 65 //--------------------------------------------------------------------------------
68 //....oooOO0OOooo........oooOO0OOooo........oo 66 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
69 67
70 GammaRayTelAnalysis::GammaRayTelAnalysis() : d << 68 GammaRayTelAnalysis::GammaRayTelAnalysis()
71 detector = static_cast<const GammaRayTelDe << 69 :GammaRayTelDetector(0),histo2DMode("strip")
>> 70 {
>> 71 GammaRayTelDetector =
>> 72 static_cast<const GammaRayTelDetectorConstruction*>
>> 73 (G4RunManager::GetRunManager()->GetUserDetectorConstruction());
72 74
73 // Define the messenger and the analysis s << 75 // Define the messenger and the analysis system
74 analysisMessenger = new GammaRayTelAnalysi << 76 analysisMessenger = new GammaRayTelAnalysisMessenger(this);
75 histogramFileName = "gammaraytel"; << 77 histoFileName = "gammaraytel";
76 } 78 }
77 79
78 //....oooOO0OOooo........oooOO0OOooo........oo <<
79 80
80 GammaRayTelAnalysis::~GammaRayTelAnalysis() { 81 GammaRayTelAnalysis::~GammaRayTelAnalysis() {
81 Finish(); << 82 Finish();
82 } << 83 // Complete clean-up
83 << 84 delete G4AnalysisManager::Instance();
84 //....oooOO0OOooo........oooOO0OOooo........oo <<
85 <<
86 void GammaRayTelAnalysis::Init() { <<
87 } 85 }
88 86
89 //....oooOO0OOooo........oooOO0OOooo........oo <<
90 87
91 void GammaRayTelAnalysis::Finish() { << 88 void GammaRayTelAnalysis::Init()
92 delete analysisMessenger; << 89 {;}
93 analysisMessenger = nullptr; <<
94 } <<
95 90
96 //....oooOO0OOooo........oooOO0OOooo........oo << 91 void GammaRayTelAnalysis::Finish()
>> 92 {
>> 93 delete analysisMessenger;
>> 94 analysisMessenger = 0;
>> 95 }
97 96
98 auto GammaRayTelAnalysis::getInstance() -> Gam << 97 GammaRayTelAnalysis* GammaRayTelAnalysis::getInstance()
99 if (instance == nullptr) { << 98 {
100 instance = new GammaRayTelAnalysis(); << 99 if (instance == 0) instance = new GammaRayTelAnalysis();
101 } << 100 return instance;
102 return instance; <<
103 } 101 }
104 102
105 //....oooOO0OOooo........oooOO0OOooo........oo <<
106 <<
107 // This function fill the 2d histogram of the 103 // This function fill the 2d histogram of the XZ positions
108 void GammaRayTelAnalysis::InsertPositionXZ(G4d << 104 void GammaRayTelAnalysis::InsertPositionXZ(double x, double z)
109 auto *manager = G4AnalysisManager::Instance( << 105 {
110 manager->FillH2(1, x, z); << 106 G4AnalysisManager* man = G4AnalysisManager::Instance();
>> 107 man->FillH2(1,x,z);
111 } 108 }
112 109
113 //....oooOO0OOooo........oooOO0OOooo........oo <<
114 <<
115 // This function fill the 2d histogram of the 110 // This function fill the 2d histogram of the YZ positions
116 void GammaRayTelAnalysis::InsertPositionYZ(G4d << 111 void GammaRayTelAnalysis::InsertPositionYZ(double y, double z)
117 auto *manager = G4AnalysisManager::Instance( << 112 {
118 manager->FillH2(2, y, z); << 113 G4AnalysisManager* man = G4AnalysisManager::Instance();
>> 114 man->FillH2(2,y,z);
119 } 115 }
120 116
121 //....oooOO0OOooo........oooOO0OOooo........oo <<
122 <<
123 // This function fill the 1d histogram of the 117 // This function fill the 1d histogram of the energy released in the last Si plane
124 void GammaRayTelAnalysis::InsertEnergy(G4doubl << 118 void GammaRayTelAnalysis::InsertEnergy(double en)
125 auto *manager = G4AnalysisManager::Instance( << 119 {
126 manager->FillH1(1, energy); << 120 G4AnalysisManager* man = G4AnalysisManager::Instance();
>> 121 man->FillH1(1,en);
127 } 122 }
128 123
129 //....oooOO0OOooo........oooOO0OOooo........oo <<
130 <<
131 // This function fill the 1d histogram of the 124 // This function fill the 1d histogram of the hits distribution along the TKR planes
132 void GammaRayTelAnalysis::InsertHits(G4int pla << 125 void GammaRayTelAnalysis::InsertHits(int nplane)
133 auto *manager = G4AnalysisManager::Instance( << 126 {
134 manager->FillH1(2, planeNumber); << 127 G4AnalysisManager* man = G4AnalysisManager::Instance();
>> 128 man->FillH1(2,nplane);
>> 129 }
>> 130
>> 131 void GammaRayTelAnalysis::setNtuple(float E, float p, float x,
>> 132 float y, float z)
>> 133 {
>> 134 G4AnalysisManager* man = G4AnalysisManager::Instance();
>> 135 man->FillNtupleDColumn(0,E);
>> 136 man->FillNtupleDColumn(1,p);
>> 137 man->FillNtupleDColumn(2,x);
>> 138 man->FillNtupleDColumn(3,y);
>> 139 man->FillNtupleDColumn(4,z);
>> 140 man->AddNtupleRow();
135 } 141 }
136 142
137 //....oooOO0OOooo........oooOO0OOooo........oo 143 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
>> 144 /*
>> 145 This member reset the histograms and it is called at the begin
>> 146 of each run; here we put the inizialization so that the histograms have
>> 147 always the right dimensions depending from the detector geometry
>> 148 */
>> 149
>> 150 void GammaRayTelAnalysis::BeginOfRun()
>> 151 {
>> 152 G4AnalysisManager* man = G4AnalysisManager::Instance();
>> 153
>> 154 // Open an output file
>> 155
>> 156 G4cout << "Opening output file " << histoFileName << " ... ";
>> 157 man->OpenFile(histoFileName);
>> 158 man->SetFirstHistoId(1);
>> 159 G4cout << " done" << G4endl;
>> 160
>> 161
>> 162 int Nplane = GammaRayTelDetector->GetNbOfTKRLayers();
>> 163 int Nstrip = GammaRayTelDetector->GetNbOfTKRStrips();
>> 164 int Ntile = GammaRayTelDetector->GetNbOfTKRTiles();
>> 165 double sizexy = GammaRayTelDetector->GetTKRSizeXY();
>> 166 double sizez = GammaRayTelDetector->GetTKRSizeZ();
>> 167 int N = Nstrip*Ntile;
>> 168
>> 169 // Book1D histograms
>> 170 //------------------
>> 171
>> 172 // 1D histogram that store the energy deposition of the
>> 173 // particle in the last (number 0) TKR X-plane
>> 174 man->CreateH1("1","Edep in the last X plane (keV)", 100, 50, 200);
>> 175
>> 176 // 1D histogram that store the hits distribution along the TKR X-planes
>> 177 man->CreateH1("2","Hits dist in TKR X planes",Nplane, 0, Nplane-1);
>> 178
>> 179 // Book 2D histograms
>> 180 //-------------------
>> 181
>> 182 // 2D histogram that store the position (mm) of the hits (XZ projection)
>> 183
>> 184 if (histo2DMode == "strip")
>> 185 {
>> 186 man->CreateH2("1","Tracker Hits XZ (strip,plane)",
>> 187 N, 0, N-1,
>> 188 2*Nplane, 0, Nplane-1);
>> 189 }
>> 190 else
>> 191 {
>> 192 man->CreateH2("1","Tracker Hits XZ (x,z) in mm",
>> 193 int(sizexy/5), -sizexy/2, sizexy/2,
>> 194 int(sizez/5), -sizez/2, sizez/2);
>> 195 }
>> 196
>> 197 // 2D histogram that store the position (mm) of the hits (YZ projection)
>> 198 if (histo2DMode == "strip")
>> 199 {
>> 200 man->CreateH2("2","Tracker Hits YZ (strip,plane)",
>> 201 N, 0, N-1,
>> 202 2*Nplane, 0, Nplane-1);
>> 203 }
>> 204 else
>> 205 {
>> 206 man->CreateH2("2","Tracker Hits YZ (x,z) in mm",
>> 207 int(sizexy/5), -sizexy/2, sizexy/2,
>> 208 int(sizez/5), -sizez/2, sizez/2);
>> 209 }
>> 210
>> 211
>> 212 // Book Ntuples (energy / plane/ x / y / z)
>> 213 //------------------------------------------
>> 214 man->CreateNtuple("1","Track ntuple");
>> 215 man->CreateNtupleDColumn("energy");
>> 216 man->CreateNtupleDColumn("plane"); // can I use Int values?
>> 217 man->CreateNtupleDColumn("x"); // can I use Int values?
>> 218 man->CreateNtupleDColumn("y"); // can I use Int values?
>> 219 man->CreateNtupleDColumn("z"); // can I use Int values?
>> 220 man->FinishNtuple();
138 221
139 void GammaRayTelAnalysis::setNtuple(G4double e <<
140 auto *manager = G4AnalysisManager::Instance( <<
141 manager->FillNtupleDColumn(0, energy); <<
142 manager->FillNtupleDColumn(1, planeNumber); <<
143 manager->FillNtupleDColumn(2, x); <<
144 manager->FillNtupleDColumn(3, y); <<
145 manager->FillNtupleDColumn(4, z); <<
146 manager->AddNtupleRow(); <<
147 } 222 }
148 223
149 //....oooOO0OOooo........oooOO0OOooo........oo 224 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
150 225
151 /* 226 /*
152 This member reset the histograms and it is ca << 227 This member is called at the end of each run
153 here we put the inizialization so that the hi << 228 */
154 always the right dimensions depending from th << 229 void GammaRayTelAnalysis::EndOfRun()
155 */ << 230 {
156 void GammaRayTelAnalysis::BeginOfRun() { << 231 //Save histograms
157 auto *manager = G4AnalysisManager::Instance( << 232 G4AnalysisManager* man = G4AnalysisManager::Instance();
158 manager->SetDefaultFileType("root"); << 233 man->Write();
159 << 234 man->CloseFile();
160 // Open an output file << 235 }
161 << 236
162 G4cout << "Opening output file " << histogra << 237 /* This member is called at the end of every event */
163 manager->OpenFile(histogramFileName); << 238
164 manager->SetFirstHistoId(1); << 239 void GammaRayTelAnalysis::EndOfEvent(G4int /* flag */ )
165 G4cout << " done" << G4endl; << 240 {;}
166 << 241
167 auto Nplane = detector->GetNbOfTKRLayers(); << 242
168 auto numberOfStrips = detector->GetNbOfTKRSt <<
169 auto numberOfTiles = detector->GetNbOfTKRTil <<
170 auto sizeXY = detector->GetTKRSizeXY(); <<
171 auto sizeZ = detector->GetTKRSizeZ(); <<
172 auto N = numberOfStrips * numberOfTiles; <<
173 <<
174 // Book 1D histograms <<
175 //------------------- <<
176 <<
177 constexpr auto NUMBER_OF_BINS{100}; <<
178 constexpr auto LOWER_BOUND{50}; <<
179 constexpr auto UPPER_BOUND{200}; <<
180 <<
181 // 1D histogram that store the energy deposi <<
182 manager->CreateH1("1", "Deposited energy in <<
183 <<
184 // 1D histogram that store the hits distribu <<
185 manager->CreateH1("2", "Hits distribution in <<
186 <<
187 // Book 2D histograms <<
188 //------------------- <<
189 <<
190 // 2D histogram that store the position (mm) <<
191 <<
192 if (histo2DMode == "strip") { <<
193 manager->CreateH2("1", "Tracker hits X <<
194 } else { <<
195 manager->CreateH2("1", "Tracker hits X <<
196 } <<
197 243
198 // 2D histogram that store the position (m <<
199 if (histo2DMode == "strip") { <<
200 manager->CreateH2("2", "Tracker hits Y <<
201 } else { <<
202 manager->CreateH2("2", "Tracker hits Y <<
203 } <<
204 244
205 // Book n-tuple (energy, plane, x, y, z) <<
206 //------------------------------------------ <<
207 manager->CreateNtuple("1", "Track n-tuple"); <<
208 manager->CreateNtupleDColumn("energy"); <<
209 manager->CreateNtupleDColumn("plane"); <<
210 manager->CreateNtupleDColumn("x"); <<
211 manager->CreateNtupleDColumn("y"); <<
212 manager->CreateNtupleDColumn("z"); <<
213 manager->FinishNtuple(); <<
214 } <<
215 245
216 /* <<
217 This member is called at the end of each run <<
218 */ <<
219 void GammaRayTelAnalysis::EndOfRun() { <<
220 // Save histograms <<
221 auto *manager = G4AnalysisManager::Instance( <<
222 manager->Write(); <<
223 manager->CloseFile(); <<
224 } <<
225 246
226 //....oooOO0OOooo........oooOO0OOooo........oo <<
227 247
228 // This member is called at the end of every e <<
229 void GammaRayTelAnalysis::EndOfEvent(G4int /* <<
230 } <<
231 248