Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
UltraFresnelLens Class Reference

#include <UltraFresnelLens.hh>

Public Member Functions

 UltraFresnelLens (G4double Diameter, G4int nGrooves, G4Material *Material, G4VPhysicalVolume *MotherPV, G4ThreeVector Pos)
 
 ~UltraFresnelLens ()
 
G4VPhysicalVolumeGetPhysicalVolume ()
 
G4MaterialGetMaterial ()
 
G4double GetDiameter ()
 
G4double GetThickness ()
 
G4double GetGrooveWidth ()
 
G4int GetNumberOfGrooves ()
 
G4double GetSagita (G4double)
 

Detailed Description

Definition at line 54 of file UltraFresnelLens.hh.

Constructor & Destructor Documentation

UltraFresnelLens::UltraFresnelLens ( G4double  Diameter,
G4int  nGrooves,
G4Material Material,
G4VPhysicalVolume MotherPV,
G4ThreeVector  Pos 
)

Definition at line 59 of file UltraFresnelLens.cc.

60 {
61 
62  LensMaterial = Material ;
63  LensDiameter = Diameter ;
64  NumberOfGrooves = nGrooves;
65  GrooveWidth = (Diameter/2.0)/nGrooves ;
66  LensPosition = Pos ;
67 
68  if( GrooveWidth <= 0 ){
69  G4Exception("UltraFresnelLens::UltraFresnelLens()","AirSh001",FatalException,
70  "UltraFresnelLens constructor: GrooveWidth<=0");
71  }
72 
73 
74 G4double Rmin1 = (NumberOfGrooves-1)*(GrooveWidth) ;
75 G4double Rmax1 = (NumberOfGrooves-0)*(GrooveWidth) ;
76 LensThickness = GetSagita(Rmax1)-GetSagita(Rmin1) ; // Height of the highest groove
77 
78 BuildLens(MotherPV) ;
79 
80 }
G4double GetSagita(G4double)
ush Pos
Definition: deflate.h:89
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

UltraFresnelLens::~UltraFresnelLens ( )

Definition at line 83 of file UltraFresnelLens.cc.

84 {;}

Member Function Documentation

G4double UltraFresnelLens::GetDiameter ( )
inline

Definition at line 64 of file UltraFresnelLens.hh.

64 {return LensDiameter ;}
G4double UltraFresnelLens::GetGrooveWidth ( )
inline

Definition at line 66 of file UltraFresnelLens.hh.

66 {return GrooveWidth ;}

Here is the caller graph for this function:

G4Material* UltraFresnelLens::GetMaterial ( )
inline

Definition at line 63 of file UltraFresnelLens.hh.

63 {return LensMaterial ;}
G4int UltraFresnelLens::GetNumberOfGrooves ( )
inline

Definition at line 67 of file UltraFresnelLens.hh.

67 {return NumberOfGrooves ;}

Here is the caller graph for this function:

G4VPhysicalVolume* UltraFresnelLens::GetPhysicalVolume ( )
inline

Definition at line 62 of file UltraFresnelLens.hh.

62 {return LensPhysicalVolume ;}
G4double UltraFresnelLens::GetSagita ( G4double  radius)

Definition at line 123 of file UltraFresnelLens.cc.

124 {
125 
126  G4double Conic = -1.0;
127  G4double Curvature = 0.00437636761488/mm ;
128  G4double Aspher[8] = { 4.206739256e-05/(mm),
129  9.6440152e-10/(mm3),
130  -1.4884317e-15/(mm2*mm3),
131  0.0/(mm*mm3*mm3),
132  0.0/(mm3*mm3*mm3),
133  0.0/(mm2*mm3*mm3*mm3),
134  0.0/(mm*mm3*mm3*mm3*mm3),
135  0.0/(mm*3*mm3*mm3*mm3*mm3)
136  };
137 
138  G4double TotAspher = 0.0*mm ;
139 
140  for(G4int k=1;k<9;k++){
141  TotAspher += Aspher[k-1]*std::pow(radius,2*k) ;
142  }
143 
144  G4double ArgSqrt = 1.0-(1.0+Conic)*std::pow(Curvature,2)*std::pow(radius,2) ;
145 
146  if (ArgSqrt < 0.0){
147  G4Exception("UltraFresnelLens::GetSagita()","AirSh002",
149  "UltraFresnelLensParameterisation::Sagita: Square Root of <0 !");
150  }
151  G4double Sagita_value = Curvature*std::pow(radius,2)/(1.0+std::sqrt(ArgSqrt)) + TotAspher;
152 
153  return Sagita_value ;
154 
155 
156 }
static constexpr double mm
Definition: G4SIunits.hh:115
int G4int
Definition: G4Types.hh:78
static constexpr double mm2
Definition: G4SIunits.hh:116
static constexpr double mm3
Definition: G4SIunits.hh:117
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4double UltraFresnelLens::GetThickness ( )
inline

Definition at line 65 of file UltraFresnelLens.hh.

65 {return LensThickness ;}

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