Geant4  10.02.p03
G4ParameterisationPolyconeZ Class Reference

#include <G4ParameterisationPolycone.hh>

Inheritance diagram for G4ParameterisationPolyconeZ:
Collaboration diagram for G4ParameterisationPolyconeZ:

Public Member Functions

 G4ParameterisationPolyconeZ (EAxis axis, G4int nCopies, G4double offset, G4double step, G4VSolid *motherSolid, DivisionType divType)
 
 ~G4ParameterisationPolyconeZ ()
 
void CheckParametersValidity ()
 
G4double GetMaxParameter () const
 
void ComputeTransformation (const G4int copyNo, G4VPhysicalVolume *physVol) const
 
void ComputeDimensions (G4Polycone &pcone, const G4int copyNo, const G4VPhysicalVolume *physVol) const
 
- Public Member Functions inherited from G4VParameterisationPolycone
 G4VParameterisationPolycone (EAxis axis, G4int nCopies, G4double offset, G4double step, G4VSolid *msolid, DivisionType divType)
 
virtual ~G4VParameterisationPolycone ()
 
- Public Member Functions inherited from G4VDivisionParameterisation
 G4VDivisionParameterisation (EAxis axis, G4int nDiv, G4double width, G4double offset, DivisionType divType, G4VSolid *motherSolid=0)
 
virtual ~G4VDivisionParameterisation ()
 
virtual G4VSolidComputeSolid (const G4int, G4VPhysicalVolume *)
 
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 ()
 
virtual ~G4VPVParameterisation ()
 
virtual G4MaterialComputeMaterial (const G4int repNo, G4VPhysicalVolume *currentVol, const G4VTouchable *parentTouch=0)
 
virtual G4bool IsNested () const
 
virtual G4VVolumeMaterialScannerGetMaterialScanner ()
 

Private Member Functions

G4double GetR (G4double z, G4double z1, G4double r1, G4double z2, G4double r2) const
 
G4double GetRmin (G4double z, G4int nsegment) const
 
G4double GetRmax (G4double z, G4int nsegment) const
 
void ComputeDimensions (G4Trd &, const G4int, const G4VPhysicalVolume *) const
 
void ComputeDimensions (G4Trap &, const G4int, const G4VPhysicalVolume *) const
 
void ComputeDimensions (G4Box &, const G4int, const G4VPhysicalVolume *) const
 
void ComputeDimensions (G4Orb &, const G4int, const G4VPhysicalVolume *) const
 
void ComputeDimensions (G4Sphere &, const G4int, const G4VPhysicalVolume *) const
 
void ComputeDimensions (G4Ellipsoid &, const G4int, const G4VPhysicalVolume *) const
 
void ComputeDimensions (G4Torus &, const G4int, const G4VPhysicalVolume *) const
 
void ComputeDimensions (G4Para &, const G4int, const G4VPhysicalVolume *) const
 
void ComputeDimensions (G4Hype &, const G4int, const G4VPhysicalVolume *) const
 
void ComputeDimensions (G4Tubs &, const G4int, const G4VPhysicalVolume *) const
 
void ComputeDimensions (G4Cons &, const G4int, const G4VPhysicalVolume *) const
 
void ComputeDimensions (G4Polyhedra &, const G4int, const G4VPhysicalVolume *) const
 

Private Attributes

G4int fNSegment
 
G4PolyconeHistoricalfOrigParamMother
 

Additional Inherited Members

- Protected Member Functions inherited from G4VDivisionParameterisation
void ChangeRotMatrix (G4VPhysicalVolume *physVol, G4double rotZ=0.) const
 
G4int CalculateNDiv (G4double motherDim, G4double width, G4double offset) const
 
G4double CalculateWidth (G4double motherDim, G4int nDiv, G4double offset) const
 
void CheckOffset (G4double maxPar)
 
void CheckNDivAndWidth (G4double maxPar)
 
G4double OffsetZ () const
 
- Protected Attributes inherited from G4VDivisionParameterisation
G4String ftype
 
EAxis faxis
 
G4int fnDiv
 
G4double fwidth
 
G4double foffset
 
DivisionType fDivisionType
 
G4VSolidfmotherSolid
 
G4bool fReflectedSolid
 
G4bool fDeleteSolid
 
G4int theVoluFirstCopyNo
 
G4double kCarTolerance
 
G4double fhgap
 
- Static Protected Attributes inherited from G4VDivisionParameterisation
static G4ThreadLocal G4RotationMatrixfRot = 0
 
static const G4int verbose = 5
 

Detailed Description

Definition at line 180 of file G4ParameterisationPolycone.hh.

Constructor & Destructor Documentation

◆ G4ParameterisationPolyconeZ()

G4ParameterisationPolyconeZ::G4ParameterisationPolyconeZ ( EAxis  axis,
G4int  nCopies,
G4double  offset,
G4double  step,
G4VSolid motherSolid,
DivisionType  divType 
)

Definition at line 352 of file G4ParameterisationPolycone.cc.

355  : G4VParameterisationPolycone( axis, nDiv, width, offset, msolid, divType ),
356  fNSegment(0),
357  fOrigParamMother(((G4Polycone*)fmotherSolid)->GetOriginalParameters())
358 {
359 
361  SetType( "DivisionPolyconeZ" );
362 
363  if( divType == DivWIDTH )
364  {
365  fnDiv =
367  - fOrigParamMother->Z_values[0] , width, offset );
368  }
369  else if( divType == DivNDIV )
370  {
371  fwidth =
373  - fOrigParamMother->Z_values[0] , nDiv, offset );
374  }
375 
376 #ifdef G4DIVDEBUG
377  if( verbose >= 1 )
378  {
379  G4cout << " G4ParameterisationPolyconeZ - # divisions " << fnDiv << " = "
380  << nDiv << G4endl
381  << " Offset " << foffset << " = " << offset << G4endl
382  << " Width " << fwidth << " = " << width << G4endl;
383  }
384 #endif
385 }
G4int CalculateNDiv(G4double motherDim, G4double width, G4double offset) const
void SetType(const G4String &type)
#define width
G4GLOB_DLL std::ostream G4cout
G4VParameterisationPolycone(EAxis axis, G4int nCopies, G4double offset, G4double step, G4VSolid *msolid, DivisionType divType)
#define G4endl
Definition: G4ios.hh:61
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ~G4ParameterisationPolyconeZ()

G4ParameterisationPolyconeZ::~G4ParameterisationPolyconeZ ( )

Definition at line 388 of file G4ParameterisationPolycone.cc.

389 {
390 }

Member Function Documentation

◆ CheckParametersValidity()

void G4ParameterisationPolyconeZ::CheckParametersValidity ( )
virtual

Reimplemented from G4VDivisionParameterisation.

Definition at line 437 of file G4ParameterisationPolycone.cc.

438 {
440 
441  // Division will be following the mother polycone segments
442  if( fDivisionType == DivNDIV ) {
443  if( fnDiv > fOrigParamMother->Num_z_planes-1 ) {
444  std::ostringstream error;
445  error << "Configuration not supported." << G4endl
446  << "Division along Z will be done by splitting in the defined"
447  << G4endl
448  << "Z planes, i.e, the number of division would be: "
450  << ", instead of: " << fnDiv << " !";
451  G4Exception("G4ParameterisationPolyconeZ::CheckParametersValidity()",
452  "GeomDiv0001", FatalException, error);
453  }
454  }
455 
456  // Division will be done within one polycone segment
457  // with applying given width and offset
459  // Check if divided region does not span over more
460  // than one z segment
461 
462  G4int isegstart = -1; // number of the segment containing start position
463  G4int isegend = -1; // number of the segment containing end position
464 
465  if ( ! fReflectedSolid ) {
466  // The start/end position of the divided region
467  G4double zstart
469  G4double zend
471 
472  G4int counter = 0;
473  while ( isegend < 0 && counter < fOrigParamMother->Num_z_planes-1 ) {
474  // first segment
475  if ( zstart >= fOrigParamMother->Z_values[counter] &&
476  zstart < fOrigParamMother->Z_values[counter+1] ) {
477  isegstart = counter;
478  }
479  // last segment
480  if ( zend > fOrigParamMother->Z_values[counter] &&
481  zend <= fOrigParamMother->Z_values[counter+1] ) {
482  isegend = counter;
483  }
484  ++counter;
485  } // Loop checking, 06.08.2015, G.Cosmo
486  }
487  else {
488  // The start/end position of the divided region
489  G4double zstart
491  G4double zend
493 
494  G4int counter = 0;
495  while ( isegend < 0 && counter < fOrigParamMother->Num_z_planes-1 ) {
496  // first segment
497  if ( zstart <= fOrigParamMother->Z_values[counter] &&
498  zstart > fOrigParamMother->Z_values[counter+1] ) {
499  isegstart = counter;
500  }
501  // last segment
502  if ( zend < fOrigParamMother->Z_values[counter] &&
503  zend >= fOrigParamMother->Z_values[counter+1] ) {
504  isegend = counter;
505  }
506  ++counter;
507  } // Loop checking, 06.08.2015, G.Cosmo
508  }
509 
510 
511  if ( isegstart != isegend ) {
512  std::ostringstream message;
513  message << "Condiguration not supported." << G4endl
514  << "Division with user defined width." << G4endl
515  << "Solid " << fmotherSolid->GetName() << G4endl
516  << "Divided region is not between two z planes.";
517  G4Exception("G4ParameterisationPolyconeZ::CheckParametersValidity()",
518  "GeomDiv0001", FatalException, message);
519  }
520 
521  fNSegment = isegstart;
522  }
523 }
int G4int
Definition: G4Types.hh:78
G4String GetName() const
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
#define G4endl
Definition: G4ios.hh:61
static PROLOG_HANDLER error
Definition: xmlrole.cc:112
double G4double
Definition: G4Types.hh:76
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ComputeDimensions() [1/13]

void G4ParameterisationPolyconeZ::ComputeDimensions ( G4Polycone pcone,
const G4int  copyNo,
const G4VPhysicalVolume physVol 
) const
virtual

Reimplemented from G4VPVParameterisation.

Definition at line 577 of file G4ParameterisationPolycone.cc.

579 {
580 
581  // Define division solid
582  G4PolyconeHistorical origparam;
583  G4int nz = 2;
584  origparam.Num_z_planes = nz;
587 
588  // Define division solid z sections
589  origparam.Z_values = new G4double[nz];
590  origparam.Rmin = new G4double[nz];
591  origparam.Rmax = new G4double[nz];
592 
593  if ( fDivisionType == DivNDIV ) {
594  // The position of the centre of copyNo-th mother polycone segment
595  G4double posi = (fOrigParamMother->Z_values[copyNo]
596  + fOrigParamMother->Z_values[copyNo+1])/2;
597 
598  origparam.Z_values[0] = fOrigParamMother->Z_values[copyNo] - posi;
599  origparam.Z_values[1] = fOrigParamMother->Z_values[copyNo+1] - posi;
600  origparam.Rmin[0] = fOrigParamMother->Rmin[copyNo];
601  origparam.Rmin[1] = fOrigParamMother->Rmin[copyNo+1];
602  origparam.Rmax[0] = fOrigParamMother->Rmax[copyNo];
603  origparam.Rmax[1] = fOrigParamMother->Rmax[copyNo+1];
604  }
605 
607  if ( ! fReflectedSolid ) {
608  origparam.Z_values[0] = - fwidth/2.;
609  origparam.Z_values[1] = fwidth/2.;
610 
611  // The position of the centre of copyNo-th division
612  G4double posi
613  = fOrigParamMother->Z_values[0] + foffset + (2*copyNo + 1) * fwidth/2.;
614 
615  // The first and last z sides z values
616  G4double zstart = posi - fwidth/2.;
617  G4double zend = posi + fwidth/2.;
618  origparam.Rmin[0] = GetRmin(zstart, fNSegment);
619  origparam.Rmax[0] = GetRmax(zstart, fNSegment);
620  origparam.Rmin[1] = GetRmin(zend, fNSegment);
621  origparam.Rmax[1] = GetRmax(zend, fNSegment);
622  }
623  else {
624  origparam.Z_values[0] = fwidth/2.;
625  origparam.Z_values[1] = - fwidth/2.;
626 
627  // The position of the centre of copyNo-th division
628  G4double posi
629  = fOrigParamMother->Z_values[0] - ( foffset + (2*copyNo + 1) * fwidth/2.);
630 
631  // The first and last z sides z values
632  G4double zstart = posi + fwidth/2.;
633  G4double zend = posi - fwidth/2.;
634  origparam.Rmin[0] = GetRmin(zstart, fNSegment);
635  origparam.Rmax[0] = GetRmax(zstart, fNSegment);
636  origparam.Rmin[1] = GetRmin(zend, fNSegment);
637  origparam.Rmax[1] = GetRmax(zend, fNSegment);
638  }
639 
640  // It can happen due to rounding errors
641  if ( origparam.Rmin[0] < 0.0 ) origparam.Rmin[0] = 0.0;
642  if ( origparam.Rmin[nz-1] < 0.0 ) origparam.Rmin[1] = 0.0;
643  }
644 
645  pcone.SetOriginalParameters(&origparam); // copy values & transfer pointers
646  pcone.Reset(); // reset to new solid parameters
647 
648 #ifdef G4DIVDEBUG
649  if( verbose >= 2 )
650  {
651  G4cout << "G4ParameterisationPolyconeZ::ComputeDimensions()" << G4endl
652  << "-- Parametrised pcone copy-number: " << copyNo << G4endl;
653  pcone.DumpInfo();
654  }
655 #endif
656 }
G4double GetRmax(G4double z, G4int nsegment) const
int G4int
Definition: G4Types.hh:78
void DumpInfo() const
G4double GetRmin(G4double z, G4int nsegment) const
G4GLOB_DLL std::ostream G4cout
void SetOriginalParameters(G4PolyconeHistorical *pars)
G4bool Reset()
Definition: G4Polycone.cc:442
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ComputeDimensions() [2/13]

void G4ParameterisationPolyconeZ::ComputeDimensions ( G4Trd ,
const G4int  ,
const G4VPhysicalVolume  
) const
inlineprivatevirtual

Reimplemented from G4VPVParameterisation.

Definition at line 208 of file G4ParameterisationPolycone.hh.

209  {}

◆ ComputeDimensions() [3/13]

void G4ParameterisationPolyconeZ::ComputeDimensions ( G4Trap ,
const G4int  ,
const G4VPhysicalVolume  
) const
inlineprivatevirtual

Reimplemented from G4VPVParameterisation.

Definition at line 210 of file G4ParameterisationPolycone.hh.

211  {}

◆ ComputeDimensions() [4/13]

void G4ParameterisationPolyconeZ::ComputeDimensions ( G4Box ,
const G4int  ,
const G4VPhysicalVolume  
) const
inlineprivatevirtual

Reimplemented from G4VPVParameterisation.

Definition at line 212 of file G4ParameterisationPolycone.hh.

213  {}

◆ ComputeDimensions() [5/13]

void G4ParameterisationPolyconeZ::ComputeDimensions ( G4Orb ,
const G4int  ,
const G4VPhysicalVolume  
) const
inlineprivatevirtual

Reimplemented from G4VPVParameterisation.

Definition at line 214 of file G4ParameterisationPolycone.hh.

215  {}

◆ ComputeDimensions() [6/13]

void G4ParameterisationPolyconeZ::ComputeDimensions ( G4Sphere ,
const G4int  ,
const G4VPhysicalVolume  
) const
inlineprivatevirtual

Reimplemented from G4VPVParameterisation.

Definition at line 216 of file G4ParameterisationPolycone.hh.

217  {}

◆ ComputeDimensions() [7/13]

void G4ParameterisationPolyconeZ::ComputeDimensions ( G4Ellipsoid ,
const G4int  ,
const G4VPhysicalVolume  
) const
inlineprivatevirtual

Reimplemented from G4VPVParameterisation.

Definition at line 218 of file G4ParameterisationPolycone.hh.

219  {}

◆ ComputeDimensions() [8/13]

void G4ParameterisationPolyconeZ::ComputeDimensions ( G4Torus ,
const G4int  ,
const G4VPhysicalVolume  
) const
inlineprivatevirtual

Reimplemented from G4VPVParameterisation.

Definition at line 220 of file G4ParameterisationPolycone.hh.

221  {}

◆ ComputeDimensions() [9/13]

void G4ParameterisationPolyconeZ::ComputeDimensions ( G4Para ,
const G4int  ,
const G4VPhysicalVolume  
) const
inlineprivatevirtual

Reimplemented from G4VPVParameterisation.

Definition at line 222 of file G4ParameterisationPolycone.hh.

223  {}

◆ ComputeDimensions() [10/13]

void G4ParameterisationPolyconeZ::ComputeDimensions ( G4Hype ,
const G4int  ,
const G4VPhysicalVolume  
) const
inlineprivatevirtual

Reimplemented from G4VPVParameterisation.

Definition at line 224 of file G4ParameterisationPolycone.hh.

225  {}

◆ ComputeDimensions() [11/13]

void G4ParameterisationPolyconeZ::ComputeDimensions ( G4Tubs ,
const G4int  ,
const G4VPhysicalVolume  
) const
inlineprivatevirtual

Reimplemented from G4VPVParameterisation.

Definition at line 226 of file G4ParameterisationPolycone.hh.

227  {}

◆ ComputeDimensions() [12/13]

void G4ParameterisationPolyconeZ::ComputeDimensions ( G4Cons ,
const G4int  ,
const G4VPhysicalVolume  
) const
inlineprivatevirtual

Reimplemented from G4VPVParameterisation.

Definition at line 228 of file G4ParameterisationPolycone.hh.

229  {}

◆ ComputeDimensions() [13/13]

void G4ParameterisationPolyconeZ::ComputeDimensions ( G4Polyhedra ,
const G4int  ,
const G4VPhysicalVolume  
) const
inlineprivatevirtual

Reimplemented from G4VPVParameterisation.

Definition at line 230 of file G4ParameterisationPolycone.hh.

231  {}

◆ ComputeTransformation()

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

Implements G4VDivisionParameterisation.

Definition at line 528 of file G4ParameterisationPolycone.cc.

529 {
530  if ( fDivisionType == DivNDIV ) {
531  // The position of the centre of copyNo-th mother polycone segment
532  G4double posi
533  = ( fOrigParamMother->Z_values[copyNo]
534  + fOrigParamMother->Z_values[copyNo+1])/2;
535  physVol->SetTranslation( G4ThreeVector(0, 0, posi) );
536  }
537 
539  // The position of the centre of copyNo-th division
541 
542  if ( ! fReflectedSolid )
543  posi += foffset + (2*copyNo + 1) * fwidth/2.;
544  else
545  posi -= foffset + (2*copyNo + 1) * fwidth/2.;
546 
547  physVol->SetTranslation( G4ThreeVector(0, 0, posi) );
548  }
549 
550  //----- calculate rotation matrix: unit
551 
552 #ifdef G4DIVDEBUG
553  if( verbose >= 2 )
554  {
555  G4cout << " G4ParameterisationPolyconeZ - position: " << posi << G4endl
556  << " copyNo: " << copyNo << " - foffset: " << foffset/deg
557  << " - fwidth: " << fwidth/deg << G4endl;
558  }
559 #endif
560 
561  ChangeRotMatrix( physVol );
562 
563 #ifdef G4DIVDEBUG
564  if( verbose >= 2 )
565  {
566  G4cout << std::setprecision(8) << " G4ParameterisationPolyconeZ "
567  << copyNo << G4endl
568  << " Position: " << origin << " - Width: " << fwidth
569  << " - Axis: " << faxis << G4endl;
570  }
571 #endif
572 }
CLHEP::Hep3Vector G4ThreeVector
void ChangeRotMatrix(G4VPhysicalVolume *physVol, G4double rotZ=0.) const
G4GLOB_DLL std::ostream G4cout
static const double deg
Definition: G4SIunits.hh:151
void SetTranslation(const G4ThreeVector &v)
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
Here is the call graph for this function:
Here is the caller graph for this function:

◆ GetMaxParameter()

G4double G4ParameterisationPolyconeZ::GetMaxParameter ( ) const
virtual

◆ GetR()

G4double G4ParameterisationPolyconeZ::GetR ( G4double  z,
G4double  z1,
G4double  r1,
G4double  z2,
G4double  r2 
) const
private

Definition at line 393 of file G4ParameterisationPolycone.cc.

396 {
397  // Linear parameterisation:
398  // r = az + b
399  // a = (r1 - r2)/(z1-z2)
400  // b = r1 - a*z1
401 
402  return (r1-r2)/(z1-z2)*z + ( r1 - (r1-r2)/(z1-z2)*z1 ) ;
403 }
Here is the caller graph for this function:

◆ GetRmax()

G4double G4ParameterisationPolyconeZ::GetRmax ( G4double  z,
G4int  nsegment 
) const
private

Definition at line 418 of file G4ParameterisationPolycone.cc.

419 {
420 // Get Rmax in the given z position for the given polycone segment
421 
422  return GetR(z,
423  fOrigParamMother->Z_values[nseg],
424  fOrigParamMother->Rmax[nseg],
425  fOrigParamMother->Z_values[nseg+1],
426  fOrigParamMother->Rmax[nseg+1]);
427 }
G4double GetR(G4double z, G4double z1, G4double r1, G4double z2, G4double r2) const
Here is the call graph for this function:
Here is the caller graph for this function:

◆ GetRmin()

G4double G4ParameterisationPolyconeZ::GetRmin ( G4double  z,
G4int  nsegment 
) const
private

Definition at line 406 of file G4ParameterisationPolycone.cc.

407 {
408 // Get Rmin in the given z position for the given polycone segment
409 
410  return GetR(z,
411  fOrigParamMother->Z_values[nseg],
412  fOrigParamMother->Rmin[nseg],
413  fOrigParamMother->Z_values[nseg+1],
414  fOrigParamMother->Rmin[nseg+1]);
415 }
G4double GetR(G4double z, G4double z1, G4double r1, G4double z2, G4double r2) const
Here is the call graph for this function:
Here is the caller graph for this function:

Member Data Documentation

◆ fNSegment

G4int G4ParameterisationPolyconeZ::fNSegment
private

Definition at line 234 of file G4ParameterisationPolycone.hh.

◆ fOrigParamMother

G4PolyconeHistorical* G4ParameterisationPolyconeZ::fOrigParamMother
private

Definition at line 235 of file G4ParameterisationPolycone.hh.


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