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

#include <G4ParameterisationTrd.hh>

+ Inheritance diagram for G4ParameterisationTrdY:

Public Member Functions

 G4ParameterisationTrdY (EAxis axis, G4int nCopies, G4double width, G4double offset, G4VSolid *motherSolid, DivisionType divType)
 
 ~G4ParameterisationTrdY ()
 
G4double GetMaxParameter () const
 
void ComputeTransformation (const G4int copyNo, G4VPhysicalVolume *physVol) const
 
void ComputeDimensions (G4Trd &trd, const G4int copyNo, const G4VPhysicalVolume *pv) const
 
void ComputeDimensions (G4Trap &, const G4int, const G4VPhysicalVolume *) const
 
- Public Member Functions inherited from G4VParameterisationTrd
 G4VParameterisationTrd (EAxis axis, G4int nCopies, G4double offset, G4double step, G4VSolid *msolid, DivisionType divType)
 
virtual ~G4VParameterisationTrd ()
 
- Public Member Functions inherited from G4VDivisionParameterisation
 G4VDivisionParameterisation (EAxis axis, G4int nDiv, G4double width, G4double offset, DivisionType divType, G4VSolid *motherSolid=nullptr)
 
virtual ~G4VDivisionParameterisation ()
 
virtual G4VSolidComputeSolid (const G4int, G4VPhysicalVolume *)
 
virtual void ComputeTransformation (const G4int copyNo, G4VPhysicalVolume *physVol) const =0
 
const G4StringGetType () const
 
EAxis GetAxis () const
 
G4int GetNoDiv () const
 
G4double GetWidth () const
 
G4double GetOffset () const
 
G4VSolidGetMotherSolid () const
 
void SetType (const G4String &type)
 
G4int VolumeFirstCopyNo () const
 
void SetHalfGap (G4double hg)
 
G4double GetHalfGap () const
 
- Public Member Functions inherited from G4VPVParameterisation
 G4VPVParameterisation ()=default
 
virtual ~G4VPVParameterisation ()=default
 
virtual void ComputeTransformation (const G4int, G4VPhysicalVolume *) const =0
 
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 Member Functions inherited from G4VDivisionParameterisation
void ChangeRotMatrix (G4VPhysicalVolume *physVol, G4double rotZ=0.0) const
 
G4int CalculateNDiv (G4double motherDim, G4double width, G4double offset) const
 
G4double CalculateWidth (G4double motherDim, G4int nDiv, G4double offset) const
 
virtual void CheckParametersValidity ()
 
void CheckOffset (G4double maxPar)
 
void CheckNDivAndWidth (G4double maxPar)
 
virtual G4double GetMaxParameter () const =0
 
G4double OffsetZ () const
 
- Protected Attributes inherited from G4VParameterisationTrd
G4bool bDivInTrap = false
 
- Protected Attributes inherited from G4VDivisionParameterisation
G4String ftype
 
EAxis faxis
 
G4int fnDiv = 0
 
G4double fwidth = 0.0
 
G4double foffset = 0.0
 
DivisionType fDivisionType
 
G4VSolidfmotherSolid = nullptr
 
G4bool fReflectedSolid = false
 
G4bool fDeleteSolid = false
 
G4int theVoluFirstCopyNo = 1
 
G4double kCarTolerance
 
G4double fhgap = 0.0
 
- Static Protected Attributes inherited from G4VDivisionParameterisation
static G4ThreadLocal G4RotationMatrixfRot = nullptr
 
static const G4int verbose = 5
 

Detailed Description

Definition at line 122 of file G4ParameterisationTrd.hh.

Constructor & Destructor Documentation

◆ G4ParameterisationTrdY()

G4ParameterisationTrdY::G4ParameterisationTrdY ( EAxis  axis,
G4int  nCopies,
G4double  width,
G4double  offset,
G4VSolid motherSolid,
DivisionType  divType 
)

Definition at line 242 of file G4ParameterisationTrd.cc.

246 : G4VParameterisationTrd( axis, nDiv, width, offset, msolid, divType )
247{
249 SetType( "DivisionTrdY" );
250
251 G4Trd* msol = (G4Trd*)(fmotherSolid);
252 if( divType == DivWIDTH )
253 {
255 width, offset );
256 }
257 else if( divType == DivNDIV )
258 {
260 nDiv, offset );
261 }
262
263#ifdef G4DIVDEBUG
264 if( verbose >= 1 )
265 {
266 G4cout << " G4ParameterisationTrdY no divisions " << fnDiv
267 << " = " << nDiv << G4endl
268 << " Offset " << foffset << " = " << offset << G4endl
269 << " width " << fwidth << " = " << width << G4endl;
270 }
271#endif
272}
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
Definition: G4Trd.hh:63
G4double GetYHalfLength1() const
void SetType(const G4String &type)
G4double CalculateWidth(G4double motherDim, G4int nDiv, G4double offset) const
G4int CalculateNDiv(G4double motherDim, G4double width, G4double offset) const

◆ ~G4ParameterisationTrdY()

G4ParameterisationTrdY::~G4ParameterisationTrdY ( )

Definition at line 275 of file G4ParameterisationTrd.cc.

276{
277}

Member Function Documentation

◆ ComputeDimensions() [1/2]

void G4ParameterisationTrdY::ComputeDimensions ( G4Trap trap,
const G4int  copyNo,
const G4VPhysicalVolume  
) const
virtual

Reimplemented from G4VPVParameterisation.

Definition at line 350 of file G4ParameterisationTrd.cc.

352{
353 G4Trd* msol = (G4Trd*)(fmotherSolid);
354 G4double pDx1 = msol->GetXHalfLength1();
355 G4double pDx2 = msol->GetXHalfLength2();
356 G4double pDz = msol->GetZHalfLength();
357 //fwidth is at Z=0
358 G4double pDy1 = msol->GetYHalfLength1();
359 G4double pDy2 = msol->GetYHalfLength2();
360 // G4double pDxAVE = (pDy1+pDy2)/2.;
361 G4double yChangeRatio = (pDy2-pDy1) / (pDy2+pDy1);
362 G4double fWidChange = yChangeRatio*fwidth;
363 G4double fWid1 = fwidth - fWidChange;
364 G4double fWid2 = fwidth + fWidChange;
365 G4double fOffsetChange = yChangeRatio*foffset/2.;
366 G4double fOffset1 = foffset - fOffsetChange;
367 G4double fOffset2 = foffset + fOffsetChange;
368 G4double cyx1 = -pDy1+fOffset1 + (copyNo+0.5)*fWid1;// centre of the side at y=-pDy1
369 G4double cyx2 = -pDy2+fOffset2 + (copyNo+0.5)*fWid2;// centre of the side at y=+pDy1
370 G4double alp = std::atan( (cyx2-cyx1)/(pDz*2.) );
371
372 pDy1 = fwidth/2. - fWidChange/2.;
373 pDy2 = fwidth/2. + fWidChange/2.;
374
375
376 trap.SetAllParameters ( pDz,
377 alp,
378 90*CLHEP::deg,
379 pDy1,
380 pDx1,
381 pDx1,
382 0.,
383 pDy2,
384 pDx2 - fhgap,
385 pDx2 - fhgap * pDx2/pDx1,
386 0.);
387
388#ifdef G4DIVDEBUG
389 if( verbose >= 2 )
390 {
391 G4cout << " G4ParameterisationTrdY::ComputeDimensions():"
392 << copyNo << G4endl;
393 trap.DumpInfo();
394 }
395#endif
396}
double G4double
Definition: G4Types.hh:83
void SetAllParameters(G4double pDz, G4double pTheta, G4double pPhi, G4double pDy1, G4double pDx1, G4double pDx2, G4double pAlp1, G4double pDy2, G4double pDx3, G4double pDx4, G4double pAlp2)
Definition: G4Trap.cc:282
G4double GetXHalfLength2() const
G4double GetYHalfLength2() const
G4double GetXHalfLength1() const
G4double GetZHalfLength() const
void DumpInfo() const

◆ ComputeDimensions() [2/2]

void G4ParameterisationTrdY::ComputeDimensions ( G4Trd trd,
const G4int  copyNo,
const G4VPhysicalVolume pv 
) const
virtual

Reimplemented from G4VPVParameterisation.

Definition at line 325 of file G4ParameterisationTrd.cc.

327{
328 //---- The division along Y of a Trd will result a Trd, only
329 //--- if Y at -Z and +Z are equal, else use the G4Trap version
330 G4Trd* msol = (G4Trd*)(fmotherSolid);
331
332 G4double pDx1 = msol->GetXHalfLength1();
333 G4double pDx2 = msol->GetXHalfLength2();
334 G4double pDz = msol->GetZHalfLength();
335 G4double pDy = fwidth/2. - fhgap;
336
337 trd.SetAllParameters ( pDx1, pDx2, pDy, pDy, pDz );
338
339#ifdef G4DIVDEBUG
340 if( verbose >= 2 )
341 {
342 G4cout << " G4ParameterisationTrdY::ComputeDimensions():" << G4endl;
343 trd.DumpInfo();
344 }
345#endif
346}
void SetAllParameters(G4double pdx1, G4double pdx2, G4double pdy1, G4double pdy2, G4double pdz)
Definition: G4Trd.cc:128

◆ ComputeTransformation()

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

Implements G4VDivisionParameterisation.

Definition at line 288 of file G4ParameterisationTrd.cc.

290{
291 G4Trd* msol = (G4Trd*)(fmotherSolid );
292 G4double mdy = ( msol->GetYHalfLength1() + msol->GetYHalfLength2() ) / 2.;
293
294 //----- translation
295 G4ThreeVector origin(0.,0.,0.);
296 G4double posi = -mdy + foffset + (copyNo+0.5)*fwidth;
297
298 if( faxis == kYAxis )
299 {
300 origin.setY( posi );
301 }
302 else
303 {
304 std::ostringstream message;
305 message << "Only axes along Y are allowed ! Axis: " << faxis;
306 G4Exception("G4ParameterisationTrdY::ComputeTransformation()",
307 "GeomDiv0002", FatalException, message);
308 }
309
310#ifdef G4DIVDEBUG
311 if( verbose >= 2 )
312 {
313 G4cout << std::setprecision(8)
314 << " G4ParameterisationTrdY::ComputeTransformation " << copyNo
315 << " pos " << origin << " rot mat " << " axis " << faxis << G4endl;
316 }
317#endif
318
319 //----- set translation
320 physVol->SetTranslation( origin );
321}
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:59
void SetTranslation(const G4ThreeVector &v)
@ kYAxis
Definition: geomdefs.hh:56

◆ GetMaxParameter()

G4double G4ParameterisationTrdY::GetMaxParameter ( ) const
virtual

Implements G4VDivisionParameterisation.

Definition at line 280 of file G4ParameterisationTrd.cc.

281{
282 G4Trd* msol = (G4Trd*)(fmotherSolid);
283 return (msol->GetYHalfLength1()+msol->GetYHalfLength2());
284}

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