FDsimG4Mirror.cc
Go to the documentation of this file.
1 #include "FDsimG4Mirror.hh"
2 #include "FDsimG4MirrorSegment.hh"
3 #include "FDsimG4MirrorSD.hh"
4 #include "FDsimG4DetectorConstruction.hh"
5 
6 #include "G4RunManager.hh"
7 #include "G4UnitsTable.hh"
8 #include "G4IntersectionSolid.hh"
9 #include "G4Box.hh"
10 
11 #include "G4SDManager.hh"
12 #include "G4Material.hh"
13 #include "G4MaterialTable.hh"
14 #include "G4Polyhedra.hh"
15 #include "G4BREPSolidPolyhedra.hh"
16 #include "G4LogicalVolume.hh"
17 #include "G4PVPlacement.hh"
18 #include "G4PVParameterised.hh"
19 #include "G4OpticalSurface.hh"
20 #include "G4LogicalBorderSurface.hh"
21 #include "G4ThreeVector.hh"
22 #include "G4VisAttributes.hh"
23 
24 #include "Randomize.hh" // Random Generators
25 
26 using namespace TelescopeSimulatorLX;
27 using namespace CLHEP;
28 
29 
30 FDsimG4Mirror::FDsimG4Mirror(const G4String& telescopeName, const G4String& mirrorDataFile, G4VPhysicalVolume* const mother)
31  fTelescopeName(telescopeName)
32  fMirrorDataFile(mirrorDataFile)
33 {
34  const G4String eyeName = fTelescopeName.substr(0, 2);
35 
36  if (eyeName == "LL")
37  fNumOfMirrors = 36;
38  else if (eyeName == "LM")
39  fNumOfMirrors = 36;
40  else if (eyeName == "LA")
41  fNumOfMirrors = 60;
42  else if (eyeName == "CO")
43  fNumOfMirrors = 60;
44  else
45  G4Exception("Error in FDsimG4Mirror !!! Unknown telescope.");
46 
47  fMotherPhysicalVolume = mother;
48 
49  m_MirrorParams = (PMIRRORPARAMS)new MIRRORPARAMS[fNumOfMirrors];
50 
51  const eTelescope eTelescopeName = GeteTelescope(fTelescopeName);
52  ReadMirrorSegmentPositions(eTelescopeName, fMirrorDataFile.c_str());
53 
54  DumpInfo();
55  CheckParameters();
56  BuildMirror();
57 }
58 
59 
60 FDsimG4Mirror::~FDsimG4Mirror()
61 {
62  delete[] m_MirrorParams;
63 }
64 
65 
66 void
67 FDsimG4Mirror::DumpInfo()
68 {
69  G4cerr << " ================================================================\n"
70  " = Mirror parameters\n"
71  " =\n"
72  " = Mirror parameters for telescope " << fTelescopeName << " read from:\n"
73  " = " << fMirrorDataFile << "\n"
74  " =\n"
75  " ================================================================" << G4endl;
76 }
77 
78 
79 void
80 FDsimG4Mirror::BuildMirror()
81 {
82  const G4RunManager* const runManager = G4RunManager::GetRunManager();
83  const FDsimG4DetectorConstruction* const detectorConstruction =
84  dynamic_cast<const FDsimG4DetectorConstruction*>(runManager->GetUserDetectorConstruction());
85 
86  const G4bool useSD = detectorConstruction->fUseSensitiveDetectors;
87  const G4int verbosityLevel = detectorConstruction->fVerbosityLevel;
88 
89  FDsimG4MirrorSD* mirrorSD = nullptr;
90  // Set sensitive detector
91  if (useSD) {
92  // Sensitive Detector Manager
93  G4SDManager* const sdManager = G4SDManager::GetSDMpointer();
94  mirrorSD = new FDsimG4MirrorSD("MirrorSD");
95  sdManager->AddNewDetector(mirrorSD);
96  }
97 
98  const G4int NUM = 2;
99  G4double fReflectivity[NUM] = { 0 };
100  G4double fPhotonEnergy[NUM] = { 0, 10*eV };
101 
102  const G4String& mirrorReflectionType = detectorConstruction->fMirrorReflectionType;
103  const SurfaceRoughness& mirrorRoughness = detectorConstruction->fMirrorRoughness;
104 
105  const G4double& sigma_alpha = mirrorRoughness.Sigma_alpha;
106  const vector<G4double>& wavelength = mirrorRoughness.Wavelength;
107  const vector<G4double>& specularSpike = mirrorRoughness.SpecularSpike;
108  const vector<G4double>& specularLobe = mirrorRoughness.SpecularLobe;
109  const vector<G4double>& backscatter = mirrorRoughness.Backscatter;
110 
111  const G4double& mirrorSegmentTiltSigma = detectorConstruction->fMirrorSegmentTiltSigma;
112 
113  // Initialize the "private" random number generator to be used to set the random orientation of the mirror segments.
114  // The seed is set according to the Eye and Telescope so that for each mirror the segments are always tilted in the same way.
115 
116  const eTelescope eTelescopeName = GeteTelescope(fTelescopeName);
117 
118  const G4int seed = eTelescopeName;
119  const G4int luxury = 3;
120  RanluxEngine myRanluxEngine(seed, luxury);
121  RandGauss GaussDist(myRanluxEngine);
122 
123  for (G4int i = 0; i < fNumOfMirrors; ++i) {
124  char szaux[100] = "";
125  sprintf(szaux, "MirrorOpticalSurface_%02i", i);
126  G4OpticalSurface* const opticalAirMirror = new G4OpticalSurface(szaux);
127 
128  opticalAirMirror->SetModel(unified);
129  opticalAirMirror->SetType(dielectric_metal);
130  opticalAirMirror->SetFinish(polished);
131 
132  for (G4int j = 0; j < NUM; ++j)
133  fReflectivity[j] = m_MirrorParams[i].fReflec375 / 100;
134 
135  G4MaterialPropertiesTable* const mirrorMPT = new G4MaterialPropertiesTable();
136  mirrorMPT->AddProperty("REFLECTIVITY", fPhotonEnergy, fReflectivity, NUM);
137 
138  if (MirrorReflectionType == "Diffusive") {
139  opticalAirMirror->SetFinish(ground);
140  opticalAirMirror->SetSigmaAlpha(sigma_alpha);
141 
142  mirrorMPT->AddProperty("SPECULARSPIKECONSTANT", new G4MaterialPropertyVector());
143  mirrorMPT->AddProperty("SPECULARLOBECONSTANT", new G4MaterialPropertyVector());
144  mirrorMPT->AddProperty("BACKSCATTERCONSTANT", new G4MaterialPropertyVector());
145 
146  for (unsigned iwl = 0, n = wavelength.size(); iwl < n; ++iwl) {
147  const G4double energy = h_Planck * c_light / wavelength[iwl];
148  mirrorMPT->AddEntry("SPECULARSPIKECONSTANT", energy, specularSpike[iwl]);
149  mirrorMPT->AddEntry("SPECULARLOBECONSTANT", energy, specularLobe[iwl]);
150  mirrorMPT->AddEntry("BACKSCATTERCONSTANT", energy, backscatter[iwl]);
151  }
152 
153  }
154 
155  opticalAirMirror->SetMaterialPropertiesTable(MirrorMPT);
156 
157  FDsimG4MirrorSegment* const mirrorSegment =
158  new FDsimG4MirrorSegment(m_MirrorParams[i].nType, m_MirrorParams[i].fRadius);
159  G4LogicalVolume* const mirrorSegment_log = mirrorSegment->GetLogicalVolume();
160  mirrorSegment_log->SetOptimisation(false);
161 
162  // Set the orientation of each mirror segment
163 
164  const G4double elevation = m_MirrorParams[i].fVertAnglePos*deg;
165  const G4double azimuth = m_MirrorParams[i].fHorzAnglePos*deg;
166 
167  const G4double xRotAngle = (elevation < 0) ? -90.0*deg : 90.0*deg;
168 
169  const G4double deltaElevation = GaussDist.fire(0, mirrorSegmentTiltSigma);
170  const G4double deltaAzimuth = GaussDist.fire(0, mirrorSegmentTiltSigma);
171 
172  G4RotationMatrix mirrorRot;
173  mirrorRot.rotateX(xRotAngle);
174  mirrorRot.rotateZ(-elevation + deltaElevation);
175  mirrorRot.rotateY(azimuth + deltaAzimuth);
176 
177  if (VerbosityLevel > 0 && mirrorSegmentTiltSigma > 0) {
178  G4cerr << "Mirror " << i << " tilted by (elevation, azimuth) "
179  << deltaElevation/deg << " (deg) "
180  << deltaAzimuth/deg << " (deg) "
181  << G4endl;
182  }
183 
184  // Set the position of each mirror segment
185 
186  G4RotationMatrix vectorRot;
187  vectorRot.rotateZ(-elevation);
188  vectorRot.rotateY(azimuth);
189  const G4double thickness = mirrorSegment->GetThickness();
190  G4ThreeVector mirrorPos(m_MirrorParams[i].fRadius + thickness/2, 0, 0);
191  mirrorPos = vectorRot * mirrorPos;
192 
193  sprintf(szaux, "MirrorSegment_%02i", i+1);
194  const G4Transform3D transform(mirrorRot, mirrorPos);
195 
196  sprintf(szaux,"MirrorLogicalBorderSurface_%02i", i);
197 
198 #if (G4VERSION_NUMBER >= 900)
199  G4int numOfCopy = 0;
200  G4PVPlacement* const mirrorSegment_phys =
201  new G4PVPlacement(transform, szaux, MirrorSegment_log, fMotherPhysicalVolume, false, ++numOfCopy);
202  if (VerbosityLevel > 2)
203  mirrorSegment_phys->CheckOverlaps(10000);
204 #endif
205 
206  if (useSD && mirrorSD)
207  mirrorSegment_log->SetSensitiveDetector(mirrorSD);
208 
209  }
210 
211  if (verbosityLevel > 0) {
212  G4cerr << "==== Mirror Surface properties ====\n"
213  "==== Mirror Reflection type : " << MirrorReflectionType << "====" << G4endl;
214  if (mirrorReflectionType == "Diffusive") {
215  G4cerr << "==== sigma_alpha = " << sigma_alpha/deg << " deg. ====\n"
216  "==== Parameters vs wavelength (Wl (nm) \t SpecularSpike \t SpecularLobe \t Backscatter) ====" << G4endl;
217  for (unsigned int i = 0, n = wavelength.size(); i < n; ++i) {
218  G4cerr << wavelength[i]/nm << '\t'
219  << specularSpike[i] << '\t'
220  << specularLobe[i] << '\t'
221  << backscatter[i]
222  << G4endl;
223  }
224  }
225  }
226 }
227 
228 
229 void
230 FDsimG4Mirror::ReadMirrorSegmentPositions(const eTelescope eTelescopeName, const char* const szFilename)
231 {
232  const int nMirrorCnt = fNumOfMirrors;
233  CMirrorParamsReaderFromFile reader(szFilename);
234 
235  stringstream errorstr;
236  if ((nMirrorCnt = reader.ReadMirrorParams(eTelescopeName)) != fNumOfMirrors) {
237  errorstr << " Reading of mirror params from file " << fMirrorDataFile << " failed !!!\n"
238  << "(info read for " << nMirrorCnt << " mirror segments, " << fNumOfMirrors << " required for this telescope.)";
239  G4Exception(errorstr.str());
240  } else
241  G4cerr << "read mirror params success - " << nMirrorCnt << " items readed" << G4endl;
242 
243  if (!reader.GetMirrorParams(m_MirrorParams, nMirrorCnt) || nMirrorCnt != fNumOfMirrors) {
244  G4cerr << "Get mirror params failed" << G4endl;
245  } else
246  G4cerr << "Get mirror params success, " << nMirrorCnt << " items added" << G4endl;
247 }
248 
249 
250 eTelescope
251 FDsimG4Mirror::GeteTelescope(const G4String& telescopeName)
252 {
253  eTelescope eTelescopeName = ET_None;
254 
255  if (telescopeName == "LL1")
256  eTelescopeName = ET_LL1;
257  else if (telescopeName == "LL2")
258  eTelescopeName = ET_LL2;
259  else if (telescopeName == "LL3")
260  eTelescopeName = ET_LL3;
261  else if (telescopeName == "LL4")
262  eTelescopeName = ET_LL4;
263  else if (telescopeName == "LL5")
264  eTelescopeName = ET_LL5;
265  else if (telescopeName == "LL6")
266  eTelescopeName = ET_LL6;
267  else if (telescopeName == "LM1")
268  eTelescopeName = ET_LM1;
269  else if (telescopeName == "LM2")
270  eTelescopeName = ET_LM2;
271  else if (telescopeName == "LM3")
272  eTelescopeName = ET_LM3;
273  else if (telescopeName == "LM4")
274  eTelescopeName = ET_LM4;
275  else if (telescopeName == "LM5")
276  eTelescopeName = ET_LM5;
277  else if (telescopeName == "LM6")
278  eTelescopeName = ET_LM6;
279  else if (telescopeName == "LA1")
280  eTelescopeName = ET_LA1;
281  else if (telescopeName == "LA2")
282  eTelescopeName = ET_LA2;
283  else if (telescopeName == "LA3")
284  eTelescopeName = ET_LA3;
285  else if (telescopeName == "LA4")
286  eTelescopeName = ET_LA4;
287  else if (telescopeName == "LA5")
288  eTelescopeName = ET_LA5;
289  else if (telescopeName == "LA6")
290  eTelescopeName = ET_LA6;
291  else if (telescopeName == "CO1")
292  eTelescopeName = ET_CO1;
293  else if (telescopeName == "CO2")
294  eTelescopeName = ET_CO2;
295  else if (telescopeName == "CO3")
296  eTelescopeName = ET_CO3;
297  else if (telescopeName == "CO4")
298  eTelescopeName = ET_CO4;
299  else if (telescopeName == "CO5")
300  eTelescopeName = ET_CO5;
301  else if (telescopeName == "CO6")
302  eTelescopeName = ET_CO6;
303 
304  return eTelescopeName;
305 }
galactic::eV
const double eV
Definition: GalacticUnits.h:35
TelescopeSimulatorLX::PMIRRORPARAMS
MIRRORPARAMS * PMIRRORPARAMS
Definition: MirrorParamsReader.hh:31
utl::deg
constexpr double deg
Definition: AugerUnits.h:140
TelescopeSimulatorLX::CMirrorParamsReaderFromFile::ReadMirrorParams
virtual G4int ReadMirrorParams(eTelescope eTelescopeName)
Definition: MirrorParamsReader.cc:56
TelescopeSimulatorLX::FDsimG4Mirror::ReadMirrorSegmentPositions
void ReadMirrorSegmentPositions(const eTelescope eTelescopeName, const char *const szFilename)
Definition: FDsimG4Mirror.cc:230
TelescopeSimulatorLX::CMirrorParamsReaderBase::GetMirrorParams
G4bool GetMirrorParams(PMIRRORPARAMS pparams, G4int &nParamsCnt)
Definition: MirrorParamsReader.cc:17
TelescopeSimulatorLX::FDsimG4Mirror::GeteTelescope
eTelescope GeteTelescope(const G4String &telescopeName)
Definition: FDsimG4Mirror.cc:251
, generated on Tue Sep 26 2023.