Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4VScoringMesh.cc
Go to the documentation of this file.
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// $Id$
28//
29// ---------------------------------------------------------------------
30// Modifications
31// 17-Apr-2012 T.Aso SetSize() and SetNumberOfSegments() is not allowed
32// to call twice in same geometrical mesh. Add warning
33// message to notify.
34//
35// ---------------------------------------------------------------------
36
37#include "G4VScoringMesh.hh"
38
39#include "G4SystemOfUnits.hh"
40#include "G4VPhysicalVolume.hh"
42#include "G4VPrimitiveScorer.hh"
43#include "G4VSDFilter.hh"
44#include "G4SDManager.hh"
45
47 : fWorldName(wName),fCurrentPS(0),fConstructed(false),fActive(true),
48 fRotationMatrix(0), fMFD(new G4MultiFunctionalDetector(wName)),
49 verboseLevel(0),sizeIsSet(false),nMeshIsSet(false),
50 fDrawUnit(""), fDrawUnitValue(1.)
51{
53
54 fSize[0] = fSize[1] = fSize[2] = 0.;
55 fNSegment[0] = fNSegment[1] = fNSegment[2] = 1;
57}
58
60{
61 ;
62}
63
65 if(verboseLevel > 9) G4cout << "G4VScoringMesh::ResetScore() is called." << G4endl;
66 std::map<G4String, G4THitsMap<G4double>* >::iterator itr = fMap.begin();
67 for(; itr != fMap.end(); itr++) {
68 if(verboseLevel > 9) G4cout << "G4VScoringMesh::ResetScore()" << itr->first << G4endl;
69 itr->second->clear();
70 }
71}
73 if ( !sizeIsSet ){
74 for(int i = 0; i < 3; i++) fSize[i] = size[i];
75 sizeIsSet = true;
76 }else{
77 G4String message = " The size of scoring mesh can not be changed.";
78 G4Exception("G4VScoringMesh::SetSize()",
79 "DigiHitsUtilsScoreVScoringMesh000", JustWarning,
80 message);
81 }
82}
84 if(sizeIsSet)
85 return G4ThreeVector(fSize[0], fSize[1], fSize[2]);
86 else
87 return G4ThreeVector(0., 0., 0.);
88}
90 fCenterPosition = G4ThreeVector(centerPosition[0], centerPosition[1], centerPosition[2]);
91}
93 if ( !nMeshIsSet ){
94 for(int i = 0; i < 3; i++) fNSegment[i] = nSegment[i];
95 nMeshIsSet = true;
96 } else {
97 G4String message = " The size of scoring segments can not be changed.";
98 G4Exception("G4VScoringMesh::SetNumberOfSegments()",
99 "DigiHitsUtilsScoreVScoringMesh000", JustWarning,
100 message);
101 }
102}
104 for(int i = 0; i < 3; i++) nSegment[i] = fNSegment[i];
105}
108 fRotationMatrix->rotateX(delta);
109}
110
113 fRotationMatrix->rotateY(delta);
114}
115
118 fRotationMatrix->rotateZ(delta);
119}
120
122
123 if(!ReadyForQuantity())
124 {
125 G4cerr << "ERROR : G4VScoringMesh::SetPrimitiveScorer() : " << ps->GetName()
126 << " does not yet have mesh size or number of bins. Set them first." << G4endl
127 << "This Method is ignored." << G4endl;
128 return;
129 }
130 if(verboseLevel > 0) G4cout << "G4VScoringMesh::SetPrimitiveScorer() : "
131 << ps->GetName() << " is registered."
132 << " 3D size: ("
133 << fNSegment[0] << ", "
134 << fNSegment[1] << ", "
135 << fNSegment[2] << ")" << G4endl;
136
137 ps->SetNijk(fNSegment[0], fNSegment[1], fNSegment[2]);
138 fCurrentPS = ps;
141 fMap[ps->GetName()] = map;
142}
143
145
146 if(fCurrentPS == 0) {
147 G4cerr << "ERROR : G4VScoringMesh::SetSDFilter() : a quantity must be defined first. This method is ignored." << G4endl;
148 return;
149 }
150 if(verboseLevel > 0) G4cout << "G4VScoringMesh::SetFilter() : "
151 << filter->GetName()
152 << " is set to "
153 << fCurrentPS->GetName() << G4endl;
154
155 G4VSDFilter* oldFilter = fCurrentPS->GetFilter();
156 if(oldFilter)
157 {
158 G4cout << "WARNING : G4VScoringMesh::SetFilter() : " << oldFilter->GetName()
159 << " is overwritten by " << filter->GetName() << G4endl;
160 }
161 fCurrentPS->SetFilter(filter);
162}
163
166 if(fCurrentPS == 0) {
167 G4cerr << "ERROR : G4VScoringMesh::SetCurrentPrimitiveScorer() : The primitive scorer <"
168 << name << "> does not found." << G4endl;
169 }
170}
171
173 std::map<G4String, G4THitsMap<G4double>* >::iterator itr = fMap.find(psname);;
174 if(itr == fMap.end()) return false;
175 return true;
176}
177
179 std::map<G4String, G4THitsMap<G4double>* >::iterator itr = fMap.find(psname);;
180 if(itr == fMap.end()) {
181 return G4String("");
182 } else {
183 return GetPrimitiveScorer(psname)->GetUnit();
184 }
185}
186
188 G4String unit = "";
189 if(fCurrentPS == 0) {
190 G4String msg = "ERROR : G4VScoringMesh::GetCurrentPSUnit() : ";
191 msg += " Current primitive scorer is null.";
192 G4cerr << msg << G4endl;
193 }else{
194 unit = fCurrentPS->GetUnit();
195 }
196 return unit;
197}
198
200 if(fCurrentPS == 0) {
201 G4String msg = "ERROR : G4VScoringMesh::GetCurrentPSUnit() : ";
202 msg += " Current primitive scorer is null.";
203 G4cerr << msg << G4endl;
204 }else{
205 fCurrentPS->SetUnit(unit);
206 }
207}
208
210 std::map<G4String, G4THitsMap<G4double>* >::iterator itr = fMap.find(psname);;
211 if(itr == fMap.end()) {
212 return 1.;
213 } else {
214 return GetPrimitiveScorer(psname)->GetUnitValue();
215 }
216}
217
219 for(int i = 0; i < 3; i++) divisionAxisNames[i] = fDivisionAxisNames[i];
220}
221
223 if(fMFD == 0) return 0;
224
226 for(G4int i = 0; i < nps; i++) {
228 if(name == ps->GetName()) return ps;
229 }
230
231 return 0;
232}
234
235 G4cout << " # of segments: ("
236 << fNSegment[0] << ", "
237 << fNSegment[1] << ", "
238 << fNSegment[2] << ")"
239 << G4endl;
240 G4cout << " displacement: ("
241 << fCenterPosition.x()/cm << ", "
242 << fCenterPosition.y()/cm << ", "
243 << fCenterPosition.z()/cm << ") [cm]"
244 << G4endl;
245 if(fRotationMatrix != 0) {
246 G4cout << " rotation matrix: "
247 << fRotationMatrix->xx() << " "
248 << fRotationMatrix->xy() << " "
249 << fRotationMatrix->xz() << G4endl
250 << " "
251 << fRotationMatrix->yx() << " "
252 << fRotationMatrix->yy() << " "
253 << fRotationMatrix->yz() << G4endl
254 << " "
255 << fRotationMatrix->zx() << " "
256 << fRotationMatrix->zy() << " "
257 << fRotationMatrix->zz() << G4endl;
258 }
259
260
261 G4cout << " registered primitve scorers : " << G4endl;
264 for(int i = 0; i < nps; i++) {
265 ps = fMFD->GetPrimitive(i);
266 G4cout << " " << i << " " << ps->GetName();
267 if(ps->GetFilter() != 0) G4cout << " with " << ps->GetFilter()->GetName();
268 G4cout << G4endl;
269 }
270
271
272}
273
275 G4cout << "scoring mesh name: " << fWorldName << G4endl;
276
277 G4cout << "# of G4THitsMap : " << fMap.size() << G4endl;
278 std::map<G4String, G4THitsMap<G4double>* >::iterator itr = fMap.begin();
279 for(; itr != fMap.end(); itr++) {
280 G4cout << "[" << itr->first << "]" << G4endl;
281 itr->second->PrintAllHits();
282 }
283 G4cout << G4endl;
284
285}
286
287
288void G4VScoringMesh::DrawMesh(const G4String& psName,G4VScoreColorMap* colorMap,G4int axflg)
289{
290 fDrawPSName = psName;
291 std::map<G4String, G4THitsMap<G4double>* >::const_iterator fMapItr = fMap.find(psName);
292 if(fMapItr!=fMap.end()) {
293 fDrawUnit = GetPSUnit(psName);
295 Draw(fMapItr->second->GetMap(), colorMap,axflg);
296 } else {
297 G4cerr << "Scorer <" << psName << "> is not defined. Method ignored." << G4endl;
298 }
299}
300
301void G4VScoringMesh::DrawMesh(const G4String& psName,G4int idxPlane,G4int iColumn,G4VScoreColorMap* colorMap)
302{
303 fDrawPSName = psName;
304 std::map<G4String, G4THitsMap<G4double>* >::const_iterator fMapItr = fMap.find(psName);
305 if(fMapItr!=fMap.end()) {
306 fDrawUnit = GetPSUnit(psName);
308 DrawColumn(fMapItr->second->GetMap(),colorMap,idxPlane,iColumn);
309 } else {
310 G4cerr << "Scorer <" << psName << "> is not defined. Method ignored." << G4endl;
311 }
312}
313
@ JustWarning
CLHEP::HepRotation G4RotationMatrix
CLHEP::Hep3Vector G4ThreeVector
double G4double
Definition: G4Types.hh:64
int G4int
Definition: G4Types.hh:66
bool G4bool
Definition: G4Types.hh:67
#define G4endl
Definition: G4ios.hh:52
G4DLLIMPORT std::ostream G4cerr
G4DLLIMPORT std::ostream G4cout
double z() const
double x() const
double y() const
double zz() const
double yz() const
double zx() const
double yx() const
double zy() const
double xx() const
HepRotation & rotateX(double delta)
Definition: Rotation.cc:66
HepRotation & rotateZ(double delta)
Definition: Rotation.cc:92
double yy() const
double xz() const
HepRotation & rotateY(double delta)
Definition: Rotation.cc:79
double xy() const
G4bool RegisterPrimitive(G4VPrimitiveScorer *)
G4VPrimitiveScorer * GetPrimitive(G4int id) const
static G4SDManager * GetSDMpointer()
Definition: G4SDManager.cc:40
void AddNewDetector(G4VSensitiveDetector *aSD)
Definition: G4SDManager.cc:67
void SetUnit(const G4String &unit)
G4VSDFilter * GetFilter() const
void SetNijk(G4int i, G4int j, G4int k)
void SetFilter(G4VSDFilter *f)
G4String GetName() const
const G4String & GetUnit() const
G4double GetUnitValue() const
G4String GetName() const
Definition: G4VSDFilter.hh:57
void SetFilter(G4VSDFilter *filter)
std::map< G4String, G4THitsMap< G4double > * > fMap
G4bool ReadyForQuantity() const
G4RotationMatrix * fRotationMatrix
G4ThreeVector GetSize() const
virtual void List() const
G4double fDrawUnitValue
G4String fDrawPSName
void RotateY(G4double delta)
void GetNumberOfSegments(G4int nSegment[3])
G4double GetPSUnitValue(const G4String &psname)
G4String GetPSUnit(const G4String &psname)
virtual void Draw(std::map< G4int, G4double * > *map, G4VScoreColorMap *colorMap, G4int axflg=111)=0
void SetCurrentPSUnit(const G4String &unit)
G4MultiFunctionalDetector * fMFD
void SetCurrentPrimitiveScorer(const G4String &name)
G4String fDivisionAxisNames[3]
void SetPrimitiveScorer(G4VPrimitiveScorer *ps)
G4VPrimitiveScorer * fCurrentPS
G4String GetCurrentPSUnit()
void GetDivisionAxisNames(G4String divisionAxisNames[3])
void DrawMesh(const G4String &psName, G4VScoreColorMap *colorMap, G4int axflg=111)
void SetNumberOfSegments(G4int nSegment[3])
G4double fSize[3]
G4VScoringMesh(const G4String &wName)
virtual ~G4VScoringMesh()
virtual void DrawColumn(std::map< G4int, G4double * > *map, G4VScoreColorMap *colorMap, G4int idxProj, G4int idxColumn)=0
void SetCenterPosition(G4double centerPosition[3])
void RotateX(G4double delta)
void SetSize(G4double size[3])
G4bool FindPrimitiveScorer(const G4String &psname)
void RotateZ(G4double delta)
G4VPrimitiveScorer * GetPrimitiveScorer(const G4String &name)
G4ThreeVector fCenterPosition
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41