| 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 // G4VScoringMesh
27 // ---------------------------------------------------------------------
28
29 #include "G4VScoringMesh.hh"
30 #include "G4THitsMap.hh"
31 #include "G4SystemOfUnits.hh"
32 #include "G4VPhysicalVolume.hh"
33 #include "G4MultiFunctionalDetector.hh"
34 #include "G4VPrimitiveScorer.hh"
35 #include "G4VSDFilter.hh"
36 #include "G4SDManager.hh"
37
38 G4VScoringMesh::G4VScoringMesh(const G4String& wName)
39 : fWorldName(wName)
40 , fCurrentPS(nullptr)
41 , fConstructed(false)
42 , fActive(true)
43 , fShape(MeshShape::undefined)
44 , fRotationMatrix(nullptr)
45 , fMFD(new G4MultiFunctionalDetector(wName))
46 , verboseLevel(0)
47 , sizeIsSet(false)
48 , nMeshIsSet(false)
49 , fDrawUnit("")
50 , fDrawUnitValue(1.)
51 , fMeshElementLogical(nullptr)
52 , fParallelWorldProcess(nullptr)
53 , fGeometryHasBeenDestroyed(false)
54 , copyNumberLevel(0)
55 , layeredMassFlg(false)
56 {
57 G4SDManager::GetSDMpointer()->AddNewDetector(fMFD);
58
59 fSize[0] = fSize[1] = fSize[2] = 0.;
60 fAngle[0] = 0.0;
61 fAngle[1] = CLHEP::twopi * rad;
62 fNSegment[0] = fNSegment[1] = fNSegment[2] = 1;
63 fDivisionAxisNames[0] = fDivisionAxisNames[1] = fDivisionAxisNames[2] = "";
64 }
65
66 void G4VScoringMesh::ResetScore()
67 {
68 if(verboseLevel > 9)
69 G4cout << "G4VScoringMesh::ResetScore() is called." << G4endl;
70 for(const auto& mp : fMap)
71 {
72 if(verboseLevel > 9)
73 G4cout << "G4VScoringMesh::ResetScore()" << mp.first << G4endl;
74 mp.second->clear();
75 }
76 }
77
78 void G4VScoringMesh::SetSize(G4double size[3])
79 {
80 if(!sizeIsSet)
81 {
82 sizeIsSet = true;
83 for(G4int i = 0; i < 3; ++i)
84 {
85 fSize[i] = size[i];
86 }
87 }
88 else
89 {
90 G4String message = " Mesh size has already been set and it cannot be changed.\n";
91 message += " This method is ignored.";
92 G4Exception("G4VScoringMesh::SetSize()",
93 "DigiHitsUtilsScoreVScoringMesh000", JustWarning, message);
94 }
95 }
96
97 G4ThreeVector G4VScoringMesh::GetSize() const
98 {
99 if(sizeIsSet)
100 return G4ThreeVector(fSize[0], fSize[1], fSize[2]);
101 return G4ThreeVector(0., 0., 0.);
102 }
103
104 void G4VScoringMesh::SetAngles(G4double startAngle, G4double spanAngle)
105 {
106 fAngle[0] = startAngle;
107 fAngle[1] = spanAngle;
108 }
109
110 void G4VScoringMesh::SetCenterPosition(G4double centerPosition[3])
111 {
112 fCenterPosition =
113 G4ThreeVector(centerPosition[0], centerPosition[1], centerPosition[2]);
114 }
115
116 void G4VScoringMesh::SetNumberOfSegments(G4int nSegment[3])
117 {
118 if(!nMeshIsSet || fShape == MeshShape::realWorldLogVol ||
119 fShape == MeshShape::probe)
120 {
121 for(G4int i = 0; i < 3; ++i)
122 fNSegment[i] = nSegment[i];
123 nMeshIsSet = true;
124 }
125 else
126 {
127 G4String message = " Number of bins has already been set and it cannot be changed.\n";
128 message += " This method is ignored.";
129 G4Exception("G4VScoringMesh::SetNumberOfSegments()",
130 "DigiHitsUtilsScoreVScoringMesh000", JustWarning, message);
131 }
132 }
133
134 void G4VScoringMesh::GetNumberOfSegments(G4int nSegment[3])
135 {
136 for(G4int i = 0; i < 3; ++i)
137 nSegment[i] = fNSegment[i];
138 }
139
140 void G4VScoringMesh::RotateX(G4double delta)
141 {
142 if(fRotationMatrix == nullptr)
143 fRotationMatrix = new G4RotationMatrix();
144 fRotationMatrix->rotateX(delta);
145 }
146
147 void G4VScoringMesh::RotateY(G4double delta)
148 {
149 if(fRotationMatrix == nullptr)
150 fRotationMatrix = new G4RotationMatrix();
151 fRotationMatrix->rotateY(delta);
152 }
153
154 void G4VScoringMesh::RotateZ(G4double delta)
155 {
156 if(fRotationMatrix == nullptr)
157 fRotationMatrix = new G4RotationMatrix();
158 fRotationMatrix->rotateZ(delta);
159 }
160
161 void G4VScoringMesh::SetPrimitiveScorer(G4VPrimitiveScorer* prs)
162 {
163 if(!ReadyForQuantity())
164 {
165 G4cerr << "ERROR : G4VScoringMesh::SetPrimitiveScorer() : "
166 << prs->GetName()
167 << " does not yet have mesh size or number of bins. Set them first."
168 << G4endl << "This Method is ignored." << G4endl;
169 return;
170 }
171 if(verboseLevel > 0)
172 G4cout << "G4VScoringMesh::SetPrimitiveScorer() : " << prs->GetName()
173 << " is registered."
174 << " 3D size: (" << fNSegment[0] << ", " << fNSegment[1] << ", "
175 << fNSegment[2] << ")" << G4endl;
176
177 prs->SetNijk(fNSegment[0], fNSegment[1], fNSegment[2]);
178 fCurrentPS = prs;
179 fMFD->RegisterPrimitive(prs);
180 auto map =
181 new G4THitsMap<G4StatDouble>(fWorldName, prs->GetName());
182 fMap[prs->GetName()] = map;
183 }
184
185 void G4VScoringMesh::SetFilter(G4VSDFilter* filter)
186 {
187 if(fCurrentPS == nullptr)
188 {
189 G4cerr << "ERROR : G4VScoringMesh::SetSDFilter() : a quantity must be "
190 "defined first. This method is ignored."
191 << G4endl;
192 return;
193 }
194 if(verboseLevel > 0)
195 G4cout << "G4VScoringMesh::SetFilter() : " << filter->GetName()
196 << " is set to " << fCurrentPS->GetName() << G4endl;
197
198 G4VSDFilter* oldFilter = fCurrentPS->GetFilter();
199 if(oldFilter != nullptr)
200 {
201 G4cout << "WARNING : G4VScoringMesh::SetFilter() : " << oldFilter->GetName()
202 << " is overwritten by " << filter->GetName() << G4endl;
203 }
204 fCurrentPS->SetFilter(filter);
205 }
206
207 void G4VScoringMesh::SetCurrentPrimitiveScorer(const G4String& name)
208 {
209 fCurrentPS = GetPrimitiveScorer(name);
210 if(fCurrentPS == nullptr)
211 {
212 G4cerr << "ERROR : G4VScoringMesh::SetCurrentPrimitiveScorer() : The "
213 "primitive scorer <"
214 << name << "> does not found." << G4endl;
215 }
216 }
217
218 G4bool G4VScoringMesh::FindPrimitiveScorer(const G4String& psname)
219 {
220 const auto itr = fMap.find(psname);
221 return itr != fMap.cend();
222 }
223
224 G4String G4VScoringMesh::GetPSUnit(const G4String& psname)
225 {
226 const auto itr = fMap.find(psname);
227 if(itr == fMap.cend())
228 {
229 return G4String("");
230 }
231
232 return GetPrimitiveScorer(psname)->GetUnit();
233 }
234
235 G4String G4VScoringMesh::GetCurrentPSUnit()
236 {
237 G4String unit = "";
238 if(fCurrentPS == nullptr)
239 {
240 G4String msg = "ERROR : G4VScoringMesh::GetCurrentPSUnit() : ";
241 msg += " Current primitive scorer is null.";
242 G4cerr << msg << G4endl;
243 }
244 else
245 {
246 unit = fCurrentPS->GetUnit();
247 }
248 return unit;
249 }
250
251 void G4VScoringMesh::SetCurrentPSUnit(const G4String& unit)
252 {
253 if(fCurrentPS == nullptr)
254 {
255 G4String msg = "ERROR : G4VScoringMesh::GetCurrentPSUnit() : ";
256 msg += " Current primitive scorer is null.";
257 G4cerr << msg << G4endl;
258 }
259 else
260 {
261 fCurrentPS->SetUnit(unit);
262 }
263 }
264
265 G4double G4VScoringMesh::GetPSUnitValue(const G4String& psname)
266 {
267 const auto itr = fMap.find(psname);
268 if(itr == fMap.cend())
269 {
270 return 1.;
271 }
272
273 return GetPrimitiveScorer(psname)->GetUnitValue();
274 }
275
276 void G4VScoringMesh::GetDivisionAxisNames(G4String divisionAxisNames[3])
277 {
278 for(G4int i = 0; i < 3; ++i)
279 divisionAxisNames[i] = fDivisionAxisNames[i];
280 }
281
282 G4VPrimitiveScorer* G4VScoringMesh::GetPrimitiveScorer(const G4String& name)
283 {
284 if(fMFD == nullptr)
285 return nullptr;
286
287 G4int nps = fMFD->GetNumberOfPrimitives();
288 for(G4int i = 0; i < nps; ++i)
289 {
290 G4VPrimitiveScorer* prs = fMFD->GetPrimitive(i);
291 if(name == prs->GetName())
292 return prs;
293 }
294
295 return nullptr;
296 }
297
298 void G4VScoringMesh::List() const
299 {
300 G4cout << " # of segments: (" << fNSegment[0] << ", " << fNSegment[1] << ", "
301 << fNSegment[2] << ")" << G4endl;
302 G4cout << " displacement: (" << fCenterPosition.x() / cm << ", "
303 << fCenterPosition.y() / cm << ", " << fCenterPosition.z() / cm
304 << ") [cm]" << G4endl;
305 if(fRotationMatrix != nullptr)
306 {
307 G4cout << " rotation matrix: " << fRotationMatrix->xx() << " "
308 << fRotationMatrix->xy() << " " << fRotationMatrix->xz() << G4endl
309 << " " << fRotationMatrix->yx() << " "
310 << fRotationMatrix->yy() << " " << fRotationMatrix->yz() << G4endl
311 << " " << fRotationMatrix->zx() << " "
312 << fRotationMatrix->zy() << " " << fRotationMatrix->zz() << G4endl;
313 }
314
315 G4cout << " registered primitve scorers : " << G4endl;
316 G4int nps = fMFD->GetNumberOfPrimitives();
317 G4VPrimitiveScorer* prs;
318 for(G4int i = 0; i < nps; ++i)
319 {
320 prs = fMFD->GetPrimitive(i);
321 G4cout << " " << i << " " << prs->GetName();
322 if(prs->GetFilter() != nullptr)
323 G4cout << " with " << prs->GetFilter()->GetName();
324 G4cout << G4endl;
325 }
326 }
327
328 void G4VScoringMesh::Dump()
329 {
330 G4cout << "scoring mesh name: " << fWorldName << G4endl;
331 G4cout << "# of G4THitsMap : " << fMap.size() << G4endl;
332 for(const auto& mp : fMap)
333 {
334 G4cout << "[" << mp.first << "]" << G4endl;
335 mp.second->PrintAllHits();
336 }
337 G4cout << G4endl;
338 }
339
340 void G4VScoringMesh::DrawMesh(const G4String& psName,
341 G4VScoreColorMap* colorMap, G4int axflg)
342 {
343 fDrawPSName = psName;
344 const auto fMapItr = fMap.find(psName);
345 if(fMapItr != fMap.cend())
346 {
347 fDrawUnit = GetPSUnit(psName);
348 fDrawUnitValue = GetPSUnitValue(psName);
349 Draw(fMapItr->second, colorMap, axflg);
350 }
351 else
352 {
353 G4cerr << "Scorer <" << psName << "> is not defined. Method ignored."
354 << G4endl;
355 }
356 }
357
358 void G4VScoringMesh::DrawMesh(const G4String& psName, G4int idxPlane,
359 G4int iColumn, G4VScoreColorMap* colorMap)
360 {
361 fDrawPSName = psName;
362 const auto fMapItr = fMap.find(psName);
363 if(fMapItr != fMap.cend())
364 {
365 fDrawUnit = GetPSUnit(psName);
366 fDrawUnitValue = GetPSUnitValue(psName);
367 DrawColumn(fMapItr->second, colorMap, idxPlane, iColumn);
368 }
369 else
370 {
371 G4cerr << "Scorer <" << psName << "> is not defined. Method ignored."
372 << G4endl;
373 }
374 }
375
376 void G4VScoringMesh::Accumulate(G4THitsMap<G4double>* map)
377 {
378 G4String psName = map->GetName();
379 const auto fMapItr = fMap.find(psName);
380 if (fMapItr != fMap.cend()) { *(fMapItr->second) += *map; }
381
382 if(verboseLevel > 9)
383 {
384 G4cout << G4endl;
385 G4cout << "G4VScoringMesh::Accumulate()" << G4endl;
386 G4cout << " PS name : " << psName << G4endl;
387 if(fMapItr == fMap.cend())
388 {
389 G4cout << " " << psName << " was not found." << G4endl;
390 }
391 else
392 {
393 G4cout << " map size : " << map->GetSize() << G4endl;
394 map->PrintAllHits();
395 }
396 G4cout << G4endl;
397 }
398 }
399
400 void G4VScoringMesh::Accumulate(G4THitsMap<G4StatDouble>* map)
401 {
402 G4String psName = map->GetName();
403 const auto fMapItr = fMap.find(psName);
404 if (fMapItr != fMap.cend()) { *(fMapItr->second) += *map; }
405
406 if(verboseLevel > 9)
407 {
408 G4cout << G4endl;
409 G4cout << "G4VScoringMesh::Accumulate()" << G4endl;
410 G4cout << " PS name : " << psName << G4endl;
411 if(fMapItr == fMap.cend())
412 {
413 G4cout << " " << psName << " was not found." << G4endl;
414 }
415 else
416 {
417 G4cout << " map size : " << map->GetSize() << G4endl;
418 map->PrintAllHits();
419 }
420 G4cout << G4endl;
421 }
422 }
423
424 void G4VScoringMesh::Construct(G4VPhysicalVolume* fWorldPhys)
425 {
426 if(fConstructed)
427 {
428 if(fGeometryHasBeenDestroyed)
429 {
430 SetupGeometry(fWorldPhys);
431 fGeometryHasBeenDestroyed = false;
432 }
433 if(verboseLevel > 0)
434 G4cout << fWorldName << " --- All quantities are reset." << G4endl;
435 ResetScore();
436 }
437 else
438 {
439 fConstructed = true;
440 SetupGeometry(fWorldPhys);
441 }
442 }
443
444 void G4VScoringMesh::WorkerConstruct(G4VPhysicalVolume* fWorldPhys)
445 {
446 if(fConstructed)
447 {
448 if(fGeometryHasBeenDestroyed)
449 {
450 fMeshElementLogical->SetSensitiveDetector(fMFD);
451 fGeometryHasBeenDestroyed = false;
452 }
453
454 if(verboseLevel > 0)
455 G4cout << fWorldPhys->GetName() << " --- All quantities are reset."
456 << G4endl;
457 ResetScore();
458 }
459 else
460 {
461 fConstructed = true;
462 fMeshElementLogical->SetSensitiveDetector(fMFD);
463 }
464 }
465
466 void G4VScoringMesh::Merge(const G4VScoringMesh* scMesh)
467 {
468 const MeshScoreMap scMap = scMesh->GetScoreMap();
469
470 auto fMapItr = fMap.cbegin();
471 auto mapItr = scMap.cbegin();
472 for(; fMapItr != fMap.cend(); ++fMapItr)
473 {
474 if(verboseLevel > 9)
475 G4cout << "G4VScoringMesh::Merge()" << fMapItr->first << G4endl;
476 *(fMapItr->second) += *(mapItr->second);
477 ++mapItr;
478 }
479 }
480