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

#include <G4tgbPlaceParamCircle.hh>

+ Inheritance diagram for G4tgbPlaceParamCircle:

Public Member Functions

 G4tgbPlaceParamCircle (G4tgrPlaceParameterisation *)
 
 ~G4tgbPlaceParamCircle ()
 
void ComputeTransformation (const G4int copyNo, G4VPhysicalVolume *physVol) const
 
- Public Member Functions inherited from G4tgbPlaceParameterisation
 G4tgbPlaceParameterisation (G4tgrPlaceParameterisation *tgrParam)
 
virtual ~G4tgbPlaceParameterisation ()
 
void CheckNExtraData (G4tgrPlaceParameterisation *tgrParam, G4int nWcheck, WLSIZEtype st, const G4String &methodName)
 
G4int GetNCopies () const
 
EAxis GetAxis () const
 
- Public Member Functions inherited from G4VPVParameterisation
 G4VPVParameterisation ()=default
 
virtual ~G4VPVParameterisation ()=default
 
virtual G4VSolidComputeSolid (const G4int, G4VPhysicalVolume *)
 
virtual G4MaterialComputeMaterial (const G4int repNo, G4VPhysicalVolume *currentVol, const G4VTouchable *parentTouch=nullptr)
 
virtual G4bool IsNested () const
 
virtual G4VVolumeMaterialScannerGetMaterialScanner ()
 
virtual void ComputeDimensions (G4Box &, const G4int, const G4VPhysicalVolume *) const
 
virtual void ComputeDimensions (G4Tubs &, const G4int, const G4VPhysicalVolume *) const
 
virtual void ComputeDimensions (G4Trd &, const G4int, const G4VPhysicalVolume *) const
 
virtual void ComputeDimensions (G4Trap &, const G4int, const G4VPhysicalVolume *) const
 
virtual void ComputeDimensions (G4Cons &, const G4int, const G4VPhysicalVolume *) const
 
virtual void ComputeDimensions (G4Sphere &, const G4int, const G4VPhysicalVolume *) const
 
virtual void ComputeDimensions (G4Orb &, const G4int, const G4VPhysicalVolume *) const
 
virtual void ComputeDimensions (G4Ellipsoid &, const G4int, const G4VPhysicalVolume *) const
 
virtual void ComputeDimensions (G4Torus &, const G4int, const G4VPhysicalVolume *) const
 
virtual void ComputeDimensions (G4Para &, const G4int, const G4VPhysicalVolume *) const
 
virtual void ComputeDimensions (G4Polycone &, const G4int, const G4VPhysicalVolume *) const
 
virtual void ComputeDimensions (G4Polyhedra &, const G4int, const G4VPhysicalVolume *) const
 
virtual void ComputeDimensions (G4Hype &, const G4int, const G4VPhysicalVolume *) const
 

Additional Inherited Members

- Protected Attributes inherited from G4tgbPlaceParameterisation
G4int theNCopies = 0
 
EAxis theAxis = kUndefined
 
G4ThreeVector theTranslation
 
G4RotationMatrixtheRotationMatrix = nullptr
 

Detailed Description

Definition at line 46 of file G4tgbPlaceParamCircle.hh.

Constructor & Destructor Documentation

◆ G4tgbPlaceParamCircle()

G4tgbPlaceParamCircle::G4tgbPlaceParamCircle ( G4tgrPlaceParameterisation * tgrParam)

Definition at line 45 of file G4tgbPlaceParamCircle.cc.

48{
49 //---- Get translation and rotation
50 if(tgrParam->GetParamType() == "CIRCLE")
51 {
52 CheckNExtraData(tgrParam, 7, WLSIZE_EQ, "G4tgbPlaceParamCircle:");
53 theCircleAxis =
54 G4ThreeVector(tgrParam->GetExtraData()[4], tgrParam->GetExtraData()[5],
55 tgrParam->GetExtraData()[6]);
56
57 G4ThreeVector zaxis(0., 0., -1.);
58 if(zaxis.cross(theCircleAxis).mag() > 1.E-6)
59 {
60 theDirInPlane = zaxis.cross(theCircleAxis);
61 }
62 else
63 {
64 theDirInPlane = theCircleAxis.cross(G4ThreeVector(0., -1., 0.));
65 }
67 }
68 else
69 {
70 CheckNExtraData(tgrParam, 4, WLSIZE_EQ, "G4tgbPlaceParamCircle:");
71 if(tgrParam->GetParamType() == "CIRCLE_XY")
72 {
73 theCircleAxis = G4ThreeVector(0., 0., 1.);
74 theDirInPlane = G4ThreeVector(1., 0., 0.);
76 }
77 else if(tgrParam->GetParamType() == "CIRCLE_XZ")
78 {
79 theCircleAxis = G4ThreeVector(0., 1., 0.);
80 theDirInPlane = G4ThreeVector(1., 0., 0.);
82 }
83 else if(tgrParam->GetParamType() == "CIRCLE_YZ")
84 {
85 theCircleAxis = G4ThreeVector(1., 0., 0.);
86 theDirInPlane = G4ThreeVector(0., 1., 0.);
88 }
89 }
90
91 if(theCircleAxis.mag() == 0.)
92 {
93 G4Exception("G4tgbPlaceParamCircle::G4tgbPlaceParamCircle()",
94 "InvalidSetup", FatalException, "Circle axis is zero !");
95 }
96 theCircleAxis /= theCircleAxis.mag();
97
99
100 theNCopies = G4int(tgrParam->GetExtraData()[0]);
101 theStep = tgrParam->GetExtraData()[1];
102 theOffset = tgrParam->GetExtraData()[2];
103 theRadius = tgrParam->GetExtraData()[3];
104
105#ifdef G4VERBOSE
107 {
108 G4cout << " G4tgbPlaceParamCircle::G4tgbPlaceParamCircle():" << G4endl
109 << " param type " << tgrParam->GetParamType() << G4endl
110 << " no copies - " << theNCopies << G4endl << " step - "
111 << theStep << G4endl << " offset - " << theOffset << G4endl
112 << " radius - " << theRadius << G4endl << " circle axis - "
113 << theCircleAxis << G4endl << " dir in plane - " << theDirInPlane
114 << G4endl;
115 }
116#endif
117}
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
CLHEP::Hep3Vector G4ThreeVector
int G4int
Definition G4Types.hh:85
#define G4endl
Definition G4ios.hh:67
G4GLOB_DLL std::ostream G4cout
@ WLSIZE_EQ
Definition G4tgrUtils.hh:47
Hep3Vector cross(const Hep3Vector &) const
double mag() const
G4tgbPlaceParameterisation(G4tgrPlaceParameterisation *tgrParam)
void CheckNExtraData(G4tgrPlaceParameterisation *tgrParam, G4int nWcheck, WLSIZEtype st, const G4String &methodName)
static G4int GetVerboseLevel()
std::vector< G4double > GetExtraData() const
const G4String & GetParamType() const
@ kYAxis
Definition geomdefs.hh:56
@ kXAxis
Definition geomdefs.hh:55
@ kZAxis
Definition geomdefs.hh:57

◆ ~G4tgbPlaceParamCircle()

G4tgbPlaceParamCircle::~G4tgbPlaceParamCircle ( )

Definition at line 40 of file G4tgbPlaceParamCircle.cc.

41{
42}

Member Function Documentation

◆ ComputeTransformation()

void G4tgbPlaceParamCircle::ComputeTransformation ( const G4int copyNo,
G4VPhysicalVolume * physVol ) const
virtual

Reimplemented from G4tgbPlaceParameterisation.

Definition at line 120 of file G4tgbPlaceParamCircle.cc.

122{
123 G4double posi = theOffset + copyNo * theStep;
124 G4ThreeVector origin = theDirInPlane * theRadius;
125 origin.rotate(posi, theCircleAxis);
126
127 //----- Calculate rotation matrix (so that all volumes point to the centre)
129 rm.rotate(-posi, theCircleAxis);
130
131 //----- Set translation and rotation
132 physVol->SetTranslation(origin);
133 G4RotationMatrix* pvRm = physVol->GetRotation();
134 if(pvRm == nullptr)
135 {
136 pvRm = new G4RotationMatrix;
137 }
138 *pvRm = *theRotationMatrix * rm;
139 physVol->SetRotation(pvRm);
140 physVol->SetCopyNo(copyNo);
141
142#ifdef G4VERBOSE
144 {
145 G4cout << " G4tgbPlaceParamCircle::ComputeTransformation():"
146 << physVol->GetName() << G4endl << " no copies - " << theNCopies
147 << G4endl << " centre - " << origin << G4endl
148 << " rotation-matrix - " << *pvRm << G4endl;
149 }
150#endif
151}
CLHEP::HepRotation G4RotationMatrix
double G4double
Definition G4Types.hh:83
Hep3Vector & rotate(double, const Hep3Vector &)
HepRotation & rotate(double delta, const Hep3Vector &axis)
Definition Rotation.cc:42
const G4RotationMatrix * GetRotation() const
virtual void SetCopyNo(G4int CopyNo)=0
const G4String & GetName() const
void SetTranslation(const G4ThreeVector &v)
void SetRotation(G4RotationMatrix *)

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