Geant4  10.02.p03
UltraFresnelLens Class Reference

#include <UltraFresnelLens.hh>

Collaboration diagram for UltraFresnelLens:

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)
 

Private Member Functions

void BuildLens (G4VPhysicalVolume *)
 

Private Attributes

G4double LensDiameter
 
G4int NumberOfGrooves
 
G4MaterialLensMaterial
 
G4ThreeVector LensPosition
 
G4double GrooveWidth
 
G4double LensThickness
 
G4VPhysicalVolumeLensPhysicalVolume
 

Detailed Description

Definition at line 54 of file UltraFresnelLens.hh.

Constructor & Destructor Documentation

◆ UltraFresnelLens()

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 }
G4ThreeVector LensPosition
G4double GetSagita(G4double)
ush Pos
Definition: deflate.h:89
G4Material * LensMaterial
void BuildLens(G4VPhysicalVolume *)
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::~UltraFresnelLens ( )

Definition at line 83 of file UltraFresnelLens.cc.

84 {;}

Member Function Documentation

◆ BuildLens()

void UltraFresnelLens::BuildLens ( G4VPhysicalVolume MotherPV)
private

Definition at line 88 of file UltraFresnelLens.cc.

88  {
89 
90 G4double StartPhi = 0.0 ;
91 G4double DeltaPhi = twopi ;
92 
93 G4double LensMotherDz = LensThickness ;
94 
95 
96 G4Tubs *LensMotherCylinder
97  = new G4Tubs("LensMotherCylinder",0.0*mm,LensDiameter/2.0,LensMotherDz/2.0,StartPhi,DeltaPhi);
98 
99 G4Material *LensMotherMaterial = MotherPV->GetLogicalVolume()->GetMaterial() ;
100 G4LogicalVolume *LensMotherLV
101  = new G4LogicalVolume(LensMotherCylinder,LensMotherMaterial,"LensMotherLV",0,0,0);
102 G4VPhysicalVolume *LensMotherPV
103  = new G4PVPlacement(0,LensPosition,"LensMotherPV",LensMotherLV,MotherPV,false,0);
104 
106 
107 
108 G4Cons *solidGroove
109  = new G4Cons("Groove",40.0*mm,50.0*mm,40.0*mm,40.001*mm,1.0*mm,StartPhi,DeltaPhi);
110 
111 G4LogicalVolume *logicalGroove
112  = new G4LogicalVolume(solidGroove,LensMaterial,"Groove_log",0,0,0);
113 
114 G4VPVParameterisation *FresnelLensParam = new UltraFresnelLensParameterisation(this);
115 
117  new G4PVParameterised("LensPV",logicalGroove,LensMotherPV,kZAxis,NumberOfGrooves,FresnelLensParam) ;
118 
119 }
G4ThreeVector LensPosition
G4Material * GetMaterial() const
Definition: G4Tubs.hh:85
G4VPhysicalVolume * LensPhysicalVolume
G4Material * LensMaterial
Definition: G4Cons.hh:83
static const double twopi
Definition: G4SIunits.hh:75
static const G4VisAttributes Invisible
double G4double
Definition: G4Types.hh:76
G4LogicalVolume * GetLogicalVolume() const
static const double mm
Definition: G4SIunits.hh:114
void SetVisAttributes(const G4VisAttributes *pVA)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ GetDiameter()

G4double UltraFresnelLens::GetDiameter ( )
inline

Definition at line 64 of file UltraFresnelLens.hh.

64 {return LensDiameter ;}

◆ GetGrooveWidth()

G4double UltraFresnelLens::GetGrooveWidth ( )
inline

Definition at line 66 of file UltraFresnelLens.hh.

66 {return GrooveWidth ;}
Here is the caller graph for this function:

◆ GetMaterial()

G4Material* UltraFresnelLens::GetMaterial ( )
inline

Definition at line 63 of file UltraFresnelLens.hh.

63 {return LensMaterial ;}
G4Material * LensMaterial

◆ GetNumberOfGrooves()

G4int UltraFresnelLens::GetNumberOfGrooves ( )
inline

Definition at line 67 of file UltraFresnelLens.hh.

67 {return NumberOfGrooves ;}
Here is the call graph for this function:
Here is the caller graph for this function:

◆ GetPhysicalVolume()

G4VPhysicalVolume* UltraFresnelLens::GetPhysicalVolume ( )
inline

Definition at line 62 of file UltraFresnelLens.hh.

62 {return LensPhysicalVolume ;}
G4VPhysicalVolume * LensPhysicalVolume
Here is the caller graph for this function:

◆ GetSagita()

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 }
int G4int
Definition: G4Types.hh:78
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
static const double mm2
Definition: G4SIunits.hh:115
double G4double
Definition: G4Types.hh:76
static const double mm
Definition: G4SIunits.hh:114
static const double mm3
Definition: G4SIunits.hh:116
Here is the call graph for this function:
Here is the caller graph for this function:

◆ GetThickness()

G4double UltraFresnelLens::GetThickness ( )
inline

Definition at line 65 of file UltraFresnelLens.hh.

65 {return LensThickness ;}
Here is the caller graph for this function:

Member Data Documentation

◆ GrooveWidth

G4double UltraFresnelLens::GrooveWidth
private

Definition at line 81 of file UltraFresnelLens.hh.

◆ LensDiameter

G4double UltraFresnelLens::LensDiameter
private

Definition at line 76 of file UltraFresnelLens.hh.

◆ LensMaterial

G4Material* UltraFresnelLens::LensMaterial
private

Definition at line 78 of file UltraFresnelLens.hh.

◆ LensPhysicalVolume

G4VPhysicalVolume* UltraFresnelLens::LensPhysicalVolume
private

Definition at line 85 of file UltraFresnelLens.hh.

◆ LensPosition

G4ThreeVector UltraFresnelLens::LensPosition
private

Definition at line 79 of file UltraFresnelLens.hh.

◆ LensThickness

G4double UltraFresnelLens::LensThickness
private

Definition at line 82 of file UltraFresnelLens.hh.

◆ NumberOfGrooves

G4int UltraFresnelLens::NumberOfGrooves
private

Definition at line 77 of file UltraFresnelLens.hh.


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