Geant4  10.03
XVCrystalCharacteristic.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
28 //
29 
31 
34 
35  fMaximum = DBL_MAX;
36  fMinimum = DBL_MAX;
37 }
38 
39 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
40 
42 }
43 
44 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
45 
48  return fLatticeManager->GetXPhysicalLattice(vVolume);
49 }
50 
51 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
52 
55  return GetXPhysicalLattice(vVolume)->GetXUnitCell();
56 }
57 
58 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
59 
62  return GetXPhysicalLattice(vVolume)->GetLogicalLattice();
63 }
64 
65 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
66 
69  if(fPhysicalLattice != vLattice){
70  fPhysicalLattice = vLattice;
72  }
73 }
74 
75 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
76 
78 GetEC(G4ThreeVector vPosition,XPhysicalLattice* vLattice){
79  if(IsInitialized(vLattice)){
80  return ComputeECFromVector(vPosition);
81  }
82  else{
83  return ComputeEC(vPosition,vLattice);
84  }
85 }
86 
87 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
88 
91 
92  G4double vTFSR = CLHEP::Bohr_radius * 0.88534;
93 
94  double vZ = vLattice->GetXUnitCell()->GetBase(0)->GetElement()->GetZ();
95  vTFSR /= (std::pow(vZ,0.333333333));
96 
97  return vTFSR;
98 }
99 
100 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
101 
104  return G4ThreeVector(-1.,-1.,-1.);
105 }
106 
107 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
108 
110  if(fMaximum == DBL_MAX){
111  fMaximum = ComputeMaximum(vLattice);
112  }
113  return fMaximum;
114 }
115 
116 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
117 
119  if(fMinimum == DBL_MAX){
120  fMinimum = ComputeMinimum(vLattice);
121  }
122  return fMinimum;
123 }
124 
125 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
126 
128  unsigned int vPrecisionX = 1024;
129  unsigned int vPrecisionY = 1024;
130  unsigned int vPrecisionZ = 1024;
131 
132  G4VPhysicalVolume* vVolume = vLattice->GetVolume();
133  G4double vStepX = GetXUnitCell(vVolume)->GetSize().x() / vPrecisionX;
134  G4double vStepY = GetXUnitCell(vVolume)->GetSize().y() / vPrecisionY;
135  G4double vStepZ = GetXUnitCell(vVolume)->GetSize().z() / vPrecisionZ;
136 
137  G4double vMaximum = -DBL_MAX;
138  G4double vValue;
139 
140  for(unsigned int i=0;i<vPrecisionX;i++){
141  for(unsigned int j=0;j<vPrecisionY;j++){
142  for(unsigned int k=0;k<vPrecisionZ;k++){
143  if( (vValue = GetEC(G4ThreeVector(vStepX * i,
144  vStepY * i,
145  vStepZ * i),
146  vLattice).mag() ) > vMaximum) {
147  vMaximum = vValue;
148  }
149  }
150  }
151  }
152 
153  return vMaximum;
154 }
155 
156 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
157 
159  unsigned int vPrecisionX = 1024;
160  unsigned int vPrecisionY = 1024;
161  unsigned int vPrecisionZ = 1024;
162 
163  G4VPhysicalVolume* vVolume = vLattice->GetVolume();
164  G4double vStepX = GetXUnitCell(vVolume)->GetSize().x() / vPrecisionX;
165  G4double vStepY = GetXUnitCell(vVolume)->GetSize().y() / vPrecisionY;
166  G4double vStepZ = GetXUnitCell(vVolume)->GetSize().z() / vPrecisionZ;
167 
168  G4double vMinimum = DBL_MAX;
169  G4double vValue;
170 
171  for(unsigned int i=0;i<vPrecisionX;i++){
172  for(unsigned int j=0;j<vPrecisionY;j++){
173  for(unsigned int k=0;k<vPrecisionZ;k++){
174  if( (vValue = GetEC(G4ThreeVector(vStepX * i,
175  vStepY * i,
176  vStepZ * i),
177  vLattice).mag() ) < vMinimum){
178  vMinimum = vValue;}
179  }
180  }
181  }
182 
183  return vMinimum;
184 }
185 
186 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
187 
189  if(vLattice == fPhysicalLattice){
190  return true;
191  }
192  else{
193  return false;
194  }
195 }
196 
197 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
static XLatticeManager3 * GetXLatticeManager()
G4Element * GetElement()
Definition: XLogicalBase.cc:56
CLHEP::Hep3Vector G4ThreeVector
XLogicalLattice * GetLogicalLattice(G4VPhysicalVolume *)
G4double GetZ() const
Definition: G4Element.hh:131
XPhysicalLattice * GetXPhysicalLattice(G4VPhysicalVolume *)
Definition of the XVCrystalCharacteristic class.
virtual G4ThreeVector ComputeECFromVector(G4ThreeVector)=0
virtual G4double GetMinimum(XPhysicalLattice *)
virtual G4double ComputeTFScreeningRadius(XPhysicalLattice *)
bool G4bool
Definition: G4Types.hh:79
G4ThreeVector GetSize()
Definition: XUnitCell.cc:51
XPhysicalLattice * fPhysicalLattice
virtual G4ThreeVector ComputePositionInUnitCell(G4ThreeVector, XPhysicalLattice *)
G4bool IsInitialized(XPhysicalLattice *)
XUnitCell * GetXUnitCell()
G4ThreeVector GetEC(G4ThreeVector, XPhysicalLattice *)
void InitializePhysicalLattice(XPhysicalLattice *)
XPhysicalLattice * GetXPhysicalLattice(G4VPhysicalVolume *)
virtual G4double GetMaximum(XPhysicalLattice *)
XLatticeManager3 * fLatticeManager
XLogicalLattice * GetLogicalLattice()
virtual G4ThreeVector ComputeEC(G4ThreeVector, XPhysicalLattice *)=0
virtual void InitializeVector()=0
XLogicalBase * GetBase(G4int)
Definition: XUnitCell.cc:63
XUnitCell * GetXUnitCell(G4VPhysicalVolume *)
double G4double
Definition: G4Types.hh:76
virtual G4double ComputeMaximum(XPhysicalLattice *)
virtual G4double ComputeMinimum(XPhysicalLattice *)
#define DBL_MAX
Definition: templates.hh:83
G4VPhysicalVolume * GetVolume()