Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4ScoringCylinder Class Reference

#include <G4ScoringCylinder.hh>

+ Inheritance diagram for G4ScoringCylinder:

Public Member Functions

 G4ScoringCylinder (G4String wName)
 
 ~G4ScoringCylinder ()
 
virtual void Construct (G4VPhysicalVolume *fWorldPhys)
 
virtual void List () const
 
virtual void Draw (std::map< G4int, G4double * > *map, G4VScoreColorMap *colorMap, G4int axflg=111)
 
virtual void DrawColumn (std::map< G4int, G4double * > *map, G4VScoreColorMap *colorMap, G4int idxProj, G4int idxColumn)
 
void SetRMax (G4double rMax)
 
void SetZSize (G4double zSize)
 
void RegisterPrimitives (std::vector< G4VPrimitiveScorer * > &vps)
 
void GetRZPhi (G4int index, G4int q[3]) const
 
- Public Member Functions inherited from G4VScoringMesh
 G4VScoringMesh (const G4String &wName)
 
virtual ~G4VScoringMesh ()
 
virtual void Construct (G4VPhysicalVolume *fWorldPhys)=0
 
virtual void List () const
 
const G4StringGetWorldName () const
 
G4bool IsActive () const
 
void Activate (G4bool vl=true)
 
MeshShape GetShape () const
 
void Accumulate (G4THitsMap< G4double > *map)
 
void Dump ()
 
void DrawMesh (const G4String &psName, G4VScoreColorMap *colorMap, G4int axflg=111)
 
void DrawMesh (const G4String &psName, G4int idxPlane, G4int iColumn, G4VScoreColorMap *colorMap)
 
virtual void Draw (std::map< G4int, G4double * > *map, G4VScoreColorMap *colorMap, G4int axflg=111)=0
 
virtual void DrawColumn (std::map< G4int, G4double * > *map, G4VScoreColorMap *colorMap, G4int idxProj, G4int idxColumn)=0
 
void ResetScore ()
 
void SetSize (G4double size[3])
 
G4ThreeVector GetSize () const
 
void SetCenterPosition (G4double centerPosition[3])
 
G4ThreeVector GetTranslation () const
 
void RotateX (G4double delta)
 
void RotateY (G4double delta)
 
void RotateZ (G4double delta)
 
G4RotationMatrix GetRotationMatrix () const
 
void SetNumberOfSegments (G4int nSegment[3])
 
void GetNumberOfSegments (G4int nSegment[3])
 
void SetPrimitiveScorer (G4VPrimitiveScorer *ps)
 
void SetFilter (G4VSDFilter *filter)
 
void SetCurrentPrimitiveScorer (const G4String &name)
 
G4bool FindPrimitiveScorer (const G4String &psname)
 
G4bool IsCurrentPrimitiveScorerNull ()
 
G4String GetPSUnit (const G4String &psname)
 
G4String GetCurrentPSUnit ()
 
void SetCurrentPSUnit (const G4String &unit)
 
G4double GetPSUnitValue (const G4String &psname)
 
void SetDrawPSName (const G4String &psname)
 
void GetDivisionAxisNames (G4String divisionAxisNames[3])
 
void SetNullToCurrentPrimitiveScorer ()
 
void SetVerboseLevel (G4int vl)
 
MeshScoreMap GetScoreMap ()
 
G4bool ReadyForQuantity () const
 

Additional Inherited Members

- Protected Member Functions inherited from G4VScoringMesh
G4VPrimitiveScorerGetPrimitiveScorer (const G4String &name)
 
- Protected Attributes inherited from G4VScoringMesh
G4String fWorldName
 
G4VPrimitiveScorerfCurrentPS
 
G4bool fConstructed
 
G4bool fActive
 
MeshShape fShape
 
G4double fSize [3]
 
G4ThreeVector fCenterPosition
 
G4RotationMatrixfRotationMatrix
 
G4int fNSegment [3]
 
std::map< G4String, G4THitsMap< G4double > * > fMap
 
G4MultiFunctionalDetectorfMFD
 
G4int verboseLevel
 
G4bool sizeIsSet
 
G4bool nMeshIsSet
 
G4String fDrawUnit
 
G4double fDrawUnitValue
 
G4String fDrawPSName
 
G4String fDivisionAxisNames [3]
 

Detailed Description

Definition at line 41 of file G4ScoringCylinder.hh.

Constructor & Destructor Documentation

◆ G4ScoringCylinder()

G4ScoringCylinder::G4ScoringCylinder ( G4String  wName)

Definition at line 55 of file G4ScoringCylinder.cc.

56 :G4VScoringMesh(wName), fMeshElementLogical(0)
57{
59
60 fDivisionAxisNames[0] = "Z";
61 fDivisionAxisNames[1] = "PHI";
62 fDivisionAxisNames[2] = "R";
63}
@ cylinderMesh
G4String fDivisionAxisNames[3]

◆ ~G4ScoringCylinder()

G4ScoringCylinder::~G4ScoringCylinder ( )

Definition at line 65 of file G4ScoringCylinder.cc.

66{;}

Member Function Documentation

◆ Construct()

void G4ScoringCylinder::Construct ( G4VPhysicalVolume fWorldPhys)
virtual

Implements G4VScoringMesh.

Definition at line 68 of file G4ScoringCylinder.cc.

69{
70 if(fConstructed) {
71
72 if(verboseLevel > 0)
73 G4cout << fWorldPhys->GetName() << " --- All quantities are reset." << G4endl;
74 ResetScore();
75
76 } else {
77 fConstructed = true;
78 SetupGeometry(fWorldPhys);
79 }
80}
#define G4endl
Definition: G4ios.hh:52
G4DLLIMPORT std::ostream G4cout
const G4String & GetName() const

◆ Draw()

void G4ScoringCylinder::Draw ( std::map< G4int, G4double * > *  map,
G4VScoreColorMap colorMap,
G4int  axflg = 111 
)
virtual

Implements G4VScoringMesh.

Definition at line 225 of file G4ScoringCylinder.cc.

225 {
226
228 if(pVisManager) {
229
230 // cell vectors
231 std::vector<double> ephi;
232 for(int phi = 0; phi < fNSegment[IPHI]; phi++) ephi.push_back(0.);
233 //-
234 std::vector<std::vector<double> > zphicell; // zphicell[Z][PHI]
235 for(int z = 0; z < fNSegment[IZ]; z++) zphicell.push_back(ephi);
236 //-
237 std::vector<std::vector<double> > rphicell; // rphicell[R][PHI]
238 for(int r = 0; r < fNSegment[IR]; r++) rphicell.push_back(ephi);
239
240 // projections
241 G4int q[3];
242 std::map<G4int, G4double*>::iterator itr = map->begin();
243 for(; itr != map->end(); itr++) {
244 if(itr->first < 0) {
245 G4cout << itr->first << G4endl;
246 continue;
247 }
248 GetRZPhi(itr->first, q);
249
250 zphicell[q[IZ]][q[IPHI]] += *(itr->second)/fDrawUnitValue;
251 rphicell[q[IR]][q[IPHI]] += *(itr->second)/fDrawUnitValue;
252 }
253
254 // search min./max. values
255 G4double zphimin = DBL_MAX, rphimin = DBL_MAX;
256 G4double zphimax = 0., rphimax = 0.;
257 for(int iphi = 0; iphi < fNSegment[IPHI]; iphi++) {
258 for(int iz = 0; iz < fNSegment[IZ]; iz++) {
259 if(zphimin > zphicell[iz][iphi]) zphimin = zphicell[iz][iphi];
260 if(zphimax < zphicell[iz][iphi]) zphimax = zphicell[iz][iphi];
261 }
262 for(int ir = 0; ir < fNSegment[IR]; ir++) {
263 if(rphimin > rphicell[ir][iphi]) rphimin = rphicell[ir][iphi];
264 if(rphimax < rphicell[ir][iphi]) rphimax = rphicell[ir][iphi];
265 }
266 }
267
268 G4VisAttributes att;
269 att.SetForceSolid(true);
270 att.SetForceAuxEdgeVisible(true);
271
272
273 G4Scale3D scale;
274 if(axflg/100==1) {
275 // rz plane
276 }
277 axflg = axflg%100;
278 if(axflg/10==1) {
279 // z-phi plane
280 if(colorMap->IfFloatMinMax()) { colorMap->SetMinMax(zphimin, zphimax); }
281
282 G4double zhalf = fSize[1]/fNSegment[IZ];
283 for(int phi = 0; phi < fNSegment[IPHI]; phi++) {
284 for(int z = 0; z < fNSegment[IZ]; z++) {
285 //-
286 G4double angle = twopi/fNSegment[IPHI]*phi;
287 G4double dphi = twopi/fNSegment[IPHI];
288 G4Tubs cylinder("z-phi", // name
289 fSize[0]*0.99, fSize[0], // inner radius, outer radius
290 zhalf, // half length in z
291 angle, dphi*0.99999); // starting phi angle, delta angle
292 //-
293 G4ThreeVector zpos(0., 0., -fSize[1] + fSize[1]/fNSegment[IZ]*(1 + 2.*z));
294 G4Transform3D trans;
295 if(fRotationMatrix) {
296 trans = G4Rotate3D(*fRotationMatrix).inverse()*G4Translate3D(zpos);
297 trans = G4Translate3D(fCenterPosition)*trans;
298 } else {
300 }
301 G4double c[4];
302 colorMap->GetMapColor(zphicell[z][phi], c);
303 att.SetColour(c[0], c[1], c[2]);//, c[3]);
304 //-
305 pVisManager->Draw(cylinder, att, trans);
306 }
307 }
308 }
309 axflg = axflg%10;
310 if(axflg==1) {
311 // r-phi plane
312 if(colorMap->IfFloatMinMax()) { colorMap->SetMinMax(rphimin, rphimax); }
313
314 G4double rsize = fSize[0]/fNSegment[IR];
315 for(int phi = 0; phi < fNSegment[IPHI]; phi++) {
316 for(int r = 0; r < fNSegment[IR]; r++) {
317
318 G4double rs[2] = {rsize*r, rsize*(r+1)};
319 G4double angle = twopi/fNSegment[IPHI]*phi;
320 G4double dphi = twopi/fNSegment[IPHI];
321 G4Tubs cylinder("z-phi", rs[0], rs[1], 0.001,
322 angle, dphi*0.99999);
323 /*
324 G4cout << ">>>> "
325 << rs[0] << " - " << rs[1] << " : "
326 << angle << " - " << angle + dphi
327 << G4endl;
328 */
329
330 G4ThreeVector zposn(0., 0., -fSize[1]);
331 G4ThreeVector zposp(0., 0., fSize[1]);
332 G4Transform3D transn, transp;
333 if(fRotationMatrix) {
334 transn = G4Rotate3D(*fRotationMatrix).inverse()*G4Translate3D(zposn);
335 transn = G4Translate3D(fCenterPosition)*transn;
336 transp = G4Rotate3D(*fRotationMatrix).inverse()*G4Translate3D(zposp);
337 transp = G4Translate3D(fCenterPosition)*transp;
338 } else {
341 }
342 G4double c[4];
343 colorMap->GetMapColor(rphicell[r][phi], c);
344 att.SetColour(c[0], c[1], c[2]);//, c[3]);
345 /*
346 G4cout << " " << c[0] << ", "
347 << c[1] << ", " << c[2] << G4endl;
348 */
349 pVisManager->Draw(cylinder, att, transn);
350 pVisManager->Draw(cylinder, att, transp);
351 }
352 }
353 }
354
355 colorMap->SetPSUnit(fDrawUnit);
356 colorMap->SetPSName(fDrawPSName);
357 colorMap->DrawColorChart();
358
359 }
360}
HepGeom::Translate3D G4Translate3D
HepGeom::Rotate3D G4Rotate3D
double G4double
Definition: G4Types.hh:64
int G4int
Definition: G4Types.hh:66
void GetRZPhi(G4int index, G4int q[3]) const
Definition: G4Tubs.hh:77
virtual void DrawColorChart(G4int nPoint=5)
G4bool IfFloatMinMax() const
void SetMinMax(G4double minVal, G4double maxVal)
void SetPSUnit(G4String &unit)
virtual void GetMapColor(G4double val, G4double color[4])=0
void SetPSName(G4String &psName)
G4RotationMatrix * fRotationMatrix
G4double fDrawUnitValue
G4String fDrawPSName
G4double fSize[3]
G4ThreeVector fCenterPosition
static G4VVisManager * GetConcreteInstance()
virtual void Draw(const G4Circle &, const G4Transform3D &objectTransformation=G4Transform3D())=0
void SetForceAuxEdgeVisible(G4bool)
void SetColour(const G4Colour &)
void SetForceSolid(G4bool)
Transform3D inverse() const
Definition: Transform3D.cc:142
#define DBL_MAX
Definition: templates.hh:83

◆ DrawColumn()

void G4ScoringCylinder::DrawColumn ( std::map< G4int, G4double * > *  map,
G4VScoreColorMap colorMap,
G4int  idxProj,
G4int  idxColumn 
)
virtual

Implements G4VScoringMesh.

Definition at line 362 of file G4ScoringCylinder.cc.

364{
365 G4int projAxis = 0;
366 switch(idxProj) {
367 case 0:
368 projAxis = IR;
369 break;
370 case 1:
371 projAxis = IZ;
372 break;
373 case 2:
374 projAxis = IPHI;
375 break;
376 }
377
378 if(idxColumn<0 || idxColumn>=fNSegment[projAxis])
379 {
380 G4cerr << "Warning : Column number " << idxColumn << " is out of scoring mesh [0," << fNSegment[projAxis]-1 <<
381 "]. Method ignored." << G4endl;
382 return;
383 }
385 if(pVisManager) {
386
387 // cell vectors
388 std::vector<std::vector<std::vector<double> > > cell; // cell[R][Z][PHI]
389 std::vector<double> ephi;
390 for(int phi = 0; phi < fNSegment[IPHI]; phi++) ephi.push_back(0.);
391 std::vector<std::vector<double> > ezphi;
392 for(int z = 0; z < fNSegment[IZ]; z++) ezphi.push_back(ephi);
393 for(int r = 0; r < fNSegment[IR]; r++) cell.push_back(ezphi);
394
395 std::vector<std::vector<double> > rzcell; // rzcell[R][Z]
396 std::vector<double> ez;
397 for(int z = 0; z < fNSegment[IZ]; z++) ez.push_back(0.);
398 for(int r = 0; r < fNSegment[IR]; r++) rzcell.push_back(ez);
399
400 std::vector<std::vector<double> > zphicell; // zphicell[Z][PHI]
401 for(int z = 0; z < fNSegment[IZ]; z++) zphicell.push_back(ephi);
402
403 std::vector<std::vector<double> > rphicell; // rphicell[R][PHI]
404 for(int r = 0; r < fNSegment[IR]; r++) rphicell.push_back(ephi);
405
406 // projections
407 G4int q[3];
408 std::map<G4int, G4double*>::iterator itr = map->begin();
409 for(; itr != map->end(); itr++) {
410 if(itr->first < 0) {
411 G4cout << itr->first << G4endl;
412 continue;
413 }
414 GetRZPhi(itr->first, q);
415
416 if(projAxis == IR && q[IR] == idxColumn) { // zphi plane
417 zphicell[q[IZ]][q[IPHI]] += *(itr->second)/fDrawUnitValue;
418 }
419 if(projAxis == IZ && q[IZ] == idxColumn) { // rphi plane
420 rphicell[q[IR]][q[IPHI]] += *(itr->second)/fDrawUnitValue;
421 }
422 if(projAxis == IPHI && q[IPHI] == idxColumn) { // rz plane
423 rzcell[q[IR]][q[IZ]] += *(itr->second)/fDrawUnitValue;
424 }
425 }
426
427 // search min./max. values
428 G4double rzmin = DBL_MAX, zphimin = DBL_MAX, rphimin = DBL_MAX;
429 G4double rzmax = 0., zphimax = 0., rphimax = 0.;
430 for(int r = 0; r < fNSegment[IR]; r++) {
431 for(int phi = 0; phi < fNSegment[IPHI]; phi++) {
432 if(rphimin > rphicell[r][phi]) rphimin = rphicell[r][phi];
433 if(rphimax < rphicell[r][phi]) rphimax = rphicell[r][phi];
434 }
435 for(int z = 0; z < fNSegment[IZ]; z++) {
436 if(rzmin > rzcell[r][z]) rzmin = rzcell[r][z];
437 if(rzmax < rzcell[r][z]) rzmax = rzcell[r][z];
438 }
439 }
440 for(int z = 0; z < fNSegment[IZ]; z++) {
441 for(int phi = 0; phi < fNSegment[IPHI]; phi++) {
442 if(zphimin > zphicell[z][phi]) zphimin = zphicell[z][phi];
443 if(zphimax < zphicell[z][phi]) zphimax = zphicell[z][phi];
444 }
445 }
446
447
448 G4VisAttributes att;
449 att.SetForceSolid(true);
450 att.SetForceAuxEdgeVisible(true);
451
452
453 G4Scale3D scale;
454 // z-phi plane
455 if(projAxis == IR) {
456 if(colorMap->IfFloatMinMax()) { colorMap->SetMinMax(zphimin,zphimax); }
457
458 G4double zhalf = fSize[1]/fNSegment[IZ];
459 G4double rsize[2] = {fSize[0]/fNSegment[IR]*idxColumn,
460 fSize[0]/fNSegment[IR]*(idxColumn+1)};
461 for(int phi = 0; phi < fNSegment[IPHI]; phi++) {
462 for(int z = 0; z < fNSegment[IZ]; z++) {
463
464 G4double angle = twopi/fNSegment[IPHI]*phi*radian;
465 G4double dphi = twopi/fNSegment[IPHI]*radian;
466 G4Tubs cylinder("z-phi", rsize[0], rsize[1], zhalf,
467 angle, dphi*0.99999);
468
469 G4ThreeVector zpos(0., 0., -fSize[1] + fSize[1]/fNSegment[IZ]*(1 + 2.*z));
470 G4Transform3D trans;
471 if(fRotationMatrix) {
472 trans = G4Rotate3D(*fRotationMatrix).inverse()*G4Translate3D(zpos);
473 trans = G4Translate3D(fCenterPosition)*trans;
474 } else {
476 }
477 G4double c[4];
478 colorMap->GetMapColor(zphicell[z][phi], c);
479 att.SetColour(c[0], c[1], c[2]);//, c[3]);
480 pVisManager->Draw(cylinder, att, trans);
481 }
482 }
483
484 // r-phi plane
485 } else if(projAxis == IZ) {
486 if(colorMap->IfFloatMinMax()) { colorMap->SetMinMax(rphimin,rphimax); }
487
488 G4double rsize = fSize[0]/fNSegment[IR];
489 for(int phi = 0; phi < fNSegment[IPHI]; phi++) {
490 for(int r = 0; r < fNSegment[IR]; r++) {
491
492 G4double rs[2] = {rsize*r, rsize*(r+1)};
493 G4double angle = twopi/fNSegment[IPHI]*phi*radian;
494 G4double dz = fSize[1]/fNSegment[IZ];
495 G4double dphi = twopi/fNSegment[IPHI]*radian;
496 G4Tubs cylinder("r-phi", rs[0], rs[1], dz,
497 angle, dphi*0.99999);
498 G4ThreeVector zpos(0., 0.,
499 -fSize[1]+fSize[1]/fNSegment[IZ]*(idxColumn*2+1));
500 G4Transform3D trans;
501 if(fRotationMatrix) {
502 trans = G4Rotate3D(*fRotationMatrix).inverse()*G4Translate3D(zpos);
503 trans = G4Translate3D(fCenterPosition)*trans;
504 } else {
506 }
507 G4double c[4];
508 colorMap->GetMapColor(rphicell[r][phi], c);
509 att.SetColour(c[0], c[1], c[2]);//, c[3]);
510 pVisManager->Draw(cylinder, att, trans);
511 }
512 }
513
514 // r-z plane
515 } else if(projAxis == IPHI) {
516 if(colorMap->IfFloatMinMax()) { colorMap->SetMinMax(rzmin,rzmax); }
517
518 G4double rsize = fSize[0]/fNSegment[IR];
519 G4double zhalf = fSize[1]/fNSegment[IZ];
520 G4double angle = twopi/fNSegment[IPHI]*idxColumn*radian;
521 G4double dphi = twopi/fNSegment[IPHI]*radian;
522 for(int z = 0; z < fNSegment[IZ]; z++) {
523 for(int r = 0; r < fNSegment[IR]; r++) {
524
525 G4double rs[2] = {rsize*r, rsize*(r+1)};
526 G4Tubs cylinder("z-phi", rs[0], rs[1], zhalf,
527 angle, dphi);
528
529 G4ThreeVector zpos(0., 0.,
530 -fSize[1]+fSize[1]/fNSegment[IZ]*(2.*z+1));
531 G4Transform3D trans;
532 if(fRotationMatrix) {
533 trans = G4Rotate3D(*fRotationMatrix).inverse()*G4Translate3D(zpos);
534 trans = G4Translate3D(fCenterPosition)*trans;
535 } else {
537 }
538 G4double c[4];
539 colorMap->GetMapColor(rzcell[r][z], c);
540 att.SetColour(c[0], c[1], c[2]);//, c[3]);
541 pVisManager->Draw(cylinder, att, trans);
542 }
543 }
544 }
545 }
546
547 colorMap->SetPSUnit(fDrawUnit);
548 colorMap->SetPSName(fDrawPSName);
549 colorMap->DrawColorChart();
550
551}
G4DLLIMPORT std::ostream G4cerr

◆ GetRZPhi()

void G4ScoringCylinder::GetRZPhi ( G4int  index,
G4int  q[3] 
) const

Definition at line 553 of file G4ScoringCylinder.cc.

553 {
554 // index = k + j * k-size + i * jk-plane-size
555
556 // nested : z -> phi -> r
557 G4int i = IZ;
558 G4int j = IPHI;
559 G4int k = IR;
560 G4int jk = fNSegment[j]*fNSegment[k];
561 q[i] = index/jk;
562 q[j] = (index - q[i]*jk)/fNSegment[k];
563 q[k] = index - q[j]*fNSegment[k] - q[i]*jk;
564}

Referenced by Draw(), and DrawColumn().

◆ List()

void G4ScoringCylinder::List ( ) const
virtual

Reimplemented from G4VScoringMesh.

Definition at line 213 of file G4ScoringCylinder.cc.

213 {
214 G4cout << "G4ScoringCylinder : " << fWorldName << " --- Shape: Cylindrical mesh" << G4endl;
215
216 G4cout << " Size (R, Dz): ("
217 << fSize[0]/cm << ", "
218 << fSize[1]/cm << ") [cm]"
219 << G4endl;
220
222}
virtual void List() const

◆ RegisterPrimitives()

void G4ScoringCylinder::RegisterPrimitives ( std::vector< G4VPrimitiveScorer * > &  vps)

◆ SetRMax()

void G4ScoringCylinder::SetRMax ( G4double  rMax)
inline

Definition at line 54 of file G4ScoringCylinder.hh.

54{fSize[0] = rMax;}

◆ SetZSize()

void G4ScoringCylinder::SetZSize ( G4double  zSize)
inline

Definition at line 55 of file G4ScoringCylinder.hh.

55{fSize[1] = zSize;} // half height

The documentation for this class was generated from the following files: