Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-02-23 09:22:21

0001 //
0002 // ********************************************************************
0003 // * License and Disclaimer                                           *
0004 // *                                                                  *
0005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
0006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
0007 // * conditions of the Geant4 Software License,  included in the file *
0008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
0009 // * include a list of copyright holders.                             *
0010 // *                                                                  *
0011 // * Neither the authors of this software system, nor their employing *
0012 // * institutes,nor the agencies providing financial support for this *
0013 // * work  make  any representation or  warranty, express or implied, *
0014 // * regarding  this  software system or assume any liability for its *
0015 // * use.  Please see the license in the file  LICENSE  and URL above *
0016 // * for the full disclaimer and the limitation of liability.         *
0017 // *                                                                  *
0018 // * This  code  implementation is the result of  the  scientific and *
0019 // * technical work of the GEANT4 collaboration.                      *
0020 // * By using,  copying,  modifying or  distributing the software (or *
0021 // * any work based  on the software)  you  agree  to acknowledge its *
0022 // * use  in  resulting  scientific  publications,  and indicate your *
0023 // * acceptance of all terms of the Geant4 Software license.          *
0024 // ********************************************************************
0025 //
0026 //
0027 /// \file optical/LXe/src/LXeDetectorConstruction.cc
0028 /// \brief Implementation of the LXeDetectorConstruction class
0029 //
0030 //
0031 #include "LXeDetectorConstruction.hh"
0032 
0033 #include "LXeDetectorMessenger.hh"
0034 #include "LXeMainVolume.hh"
0035 #include "LXePMTSD.hh"
0036 #include "LXeScintSD.hh"
0037 #include "LXeWLSSlab.hh"
0038 
0039 #include "G4Box.hh"
0040 #include "G4GeometryManager.hh"
0041 #include "G4LogicalBorderSurface.hh"
0042 #include "G4LogicalSkinSurface.hh"
0043 #include "G4LogicalVolume.hh"
0044 #include "G4LogicalVolumeStore.hh"
0045 #include "G4Material.hh"
0046 #include "G4MaterialTable.hh"
0047 #include "G4OpticalSurface.hh"
0048 #include "G4PVPlacement.hh"
0049 #include "G4PhysicalConstants.hh"
0050 #include "G4PhysicalVolumeStore.hh"
0051 #include "G4RunManager.hh"
0052 #include "G4SDManager.hh"
0053 #include "G4SolidStore.hh"
0054 #include "G4Sphere.hh"
0055 #include "G4SystemOfUnits.hh"
0056 #include "G4ThreeVector.hh"
0057 #include "G4Tubs.hh"
0058 #include "G4UImanager.hh"
0059 #include "G4VisAttributes.hh"
0060 #include "globals.hh"
0061 
0062 G4bool LXeDetectorConstruction::fSphereOn = true;
0063 
0064 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0065 
0066 LXeDetectorConstruction::LXeDetectorConstruction()
0067 {
0068   SetDefaults();
0069   DefineMaterials();
0070   fDetectorMessenger = new LXeDetectorMessenger(this);
0071 }
0072 
0073 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0074 
0075 LXeDetectorConstruction::~LXeDetectorConstruction()
0076 {
0077   delete fMainVolume;
0078   delete fLXe_mt;
0079   delete fDetectorMessenger;
0080   delete fMPTPStyrene;
0081 }
0082 
0083 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0084 
0085 void LXeDetectorConstruction::DefineMaterials()
0086 {
0087   G4double a;  // atomic mass
0088   G4double z;  // atomic number
0089   G4double density;
0090 
0091   G4int polyPMMA = 1;
0092   G4int nC_PMMA = 3 + 2 * polyPMMA;
0093   G4int nH_PMMA = 6 + 2 * polyPMMA;
0094 
0095   G4int polyeth = 1;
0096   G4int nC_eth = 2 * polyeth;
0097   G4int nH_eth = 4 * polyeth;
0098 
0099   //***Elements
0100   fH = new G4Element("H", "H", z = 1., a = 1.01 * g / mole);
0101   fC = new G4Element("C", "C", z = 6., a = 12.01 * g / mole);
0102   fN = new G4Element("N", "N", z = 7., a = 14.01 * g / mole);
0103   fO = new G4Element("O", "O", z = 8., a = 16.00 * g / mole);
0104 
0105   //***Materials
0106   // Liquid Xenon
0107   fLXe = new G4Material("LXe", z = 54., a = 131.29 * g / mole, density = 3.020 * g / cm3);
0108   // Aluminum
0109   fAl = new G4Material("Al", z = 13., a = 26.98 * g / mole, density = 2.7 * g / cm3);
0110   // Vacuum
0111   fVacuum = new G4Material("Vacuum", z = 1., a = 1.01 * g / mole, density = universe_mean_density,
0112                            kStateGas, 0.1 * kelvin, 1.e-19 * pascal);
0113   // Air
0114   fAir = new G4Material("Air", density = 1.29 * mg / cm3, 2);
0115   fAir->AddElement(fN, 70 * perCent);
0116   fAir->AddElement(fO, 30 * perCent);
0117   // Glass
0118   fGlass = new G4Material("Glass", density = 1.032 * g / cm3, 2);
0119   fGlass->AddElement(fC, 91.533 * perCent);
0120   fGlass->AddElement(fH, 8.467 * perCent);
0121   // Polystyrene
0122   fPstyrene = new G4Material("Polystyrene", density = 1.03 * g / cm3, 2);
0123   fPstyrene->AddElement(fC, 8);
0124   fPstyrene->AddElement(fH, 8);
0125   // Fiber(PMMA)
0126   fPMMA = new G4Material("PMMA", density = 1190. * kg / m3, 3);
0127   fPMMA->AddElement(fH, nH_PMMA);
0128   fPMMA->AddElement(fC, nC_PMMA);
0129   fPMMA->AddElement(fO, 2);
0130   // Cladding(polyethylene)
0131   fPethylene1 = new G4Material("Pethylene1", density = 1200. * kg / m3, 2);
0132   fPethylene1->AddElement(fH, nH_eth);
0133   fPethylene1->AddElement(fC, nC_eth);
0134   // Double cladding(flourinated polyethylene)
0135   fPethylene2 = new G4Material("Pethylene2", density = 1400. * kg / m3, 2);
0136   fPethylene2->AddElement(fH, nH_eth);
0137   fPethylene2->AddElement(fC, nC_eth);
0138 
0139   //***Material properties tables
0140 
0141   std::vector<G4double> lxe_Energy = {7.0 * eV, 7.07 * eV, 7.14 * eV};
0142 
0143   std::vector<G4double> lxe_SCINT = {0.1, 1.0, 0.1};
0144   std::vector<G4double> lxe_RIND = {1.59, 1.57, 1.54};
0145   std::vector<G4double> lxe_ABSL = {35. * cm, 35. * cm, 35. * cm};
0146   fLXe_mt = new G4MaterialPropertiesTable();
0147   fLXe_mt->AddProperty("SCINTILLATIONCOMPONENT1", lxe_Energy, lxe_SCINT);
0148   fLXe_mt->AddProperty("SCINTILLATIONCOMPONENT2", lxe_Energy, lxe_SCINT);
0149   fLXe_mt->AddProperty("RINDEX", lxe_Energy, lxe_RIND);
0150   fLXe_mt->AddProperty("ABSLENGTH", lxe_Energy, lxe_ABSL);
0151   fLXe_mt->AddConstProperty("SCINTILLATIONYIELD", 12000. / MeV);
0152   fLXe_mt->AddConstProperty("RESOLUTIONSCALE", 1.0);
0153   fLXe_mt->AddConstProperty("SCINTILLATIONTIMECONSTANT1", 20. * ns);
0154   fLXe_mt->AddConstProperty("SCINTILLATIONTIMECONSTANT2", 45. * ns);
0155   fLXe_mt->AddConstProperty("SCINTILLATIONYIELD1", 1.0);
0156   fLXe_mt->AddConstProperty("SCINTILLATIONYIELD2", 0.0);
0157   fLXe->SetMaterialPropertiesTable(fLXe_mt);
0158 
0159   // Set the Birks Constant for the LXe scintillator
0160   fLXe->GetIonisation()->SetBirksConstant(0.126 * mm / MeV);
0161 
0162   std::vector<G4double> glass_AbsLength = {420. * cm, 420. * cm, 420. * cm};
0163   auto glass_mt = new G4MaterialPropertiesTable();
0164   glass_mt->AddProperty("ABSLENGTH", lxe_Energy, glass_AbsLength);
0165   glass_mt->AddProperty("RINDEX", "Fused Silica");
0166   fGlass->SetMaterialPropertiesTable(glass_mt);
0167 
0168   auto vacuum_mt = new G4MaterialPropertiesTable();
0169   vacuum_mt->AddProperty("RINDEX", "Air");
0170   fVacuum->SetMaterialPropertiesTable(vacuum_mt);
0171   fAir->SetMaterialPropertiesTable(vacuum_mt);  // Give air the same rindex
0172 
0173   std::vector<G4double> wls_Energy = {2.00 * eV, 2.87 * eV, 2.90 * eV, 3.47 * eV};
0174 
0175   std::vector<G4double> rIndexPstyrene = {1.5, 1.5, 1.5, 1.5};
0176   std::vector<G4double> absorption1 = {2. * cm, 2. * cm, 2. * cm, 2. * cm};
0177   std::vector<G4double> scintilFast = {0.0, 0.0, 1.0, 1.0};
0178   fMPTPStyrene = new G4MaterialPropertiesTable();
0179   fMPTPStyrene->AddProperty("RINDEX", wls_Energy, rIndexPstyrene);
0180   fMPTPStyrene->AddProperty("ABSLENGTH", wls_Energy, absorption1);
0181   fMPTPStyrene->AddProperty("SCINTILLATIONCOMPONENT1", wls_Energy, scintilFast);
0182   fMPTPStyrene->AddConstProperty("SCINTILLATIONYIELD", 10. / keV);
0183   fMPTPStyrene->AddConstProperty("RESOLUTIONSCALE", 1.0);
0184   fMPTPStyrene->AddConstProperty("SCINTILLATIONTIMECONSTANT1", 10. * ns);
0185   fPstyrene->SetMaterialPropertiesTable(fMPTPStyrene);
0186 
0187   // Set the Birks Constant for the Polystyrene scintillator
0188   fPstyrene->GetIonisation()->SetBirksConstant(0.126 * mm / MeV);
0189 
0190   std::vector<G4double> AbsFiber = {9.0 * m, 9.0 * m, 0.1 * mm, 0.1 * mm};
0191   std::vector<G4double> EmissionFib = {1.0, 1.0, 0.0, 0.0};
0192   auto fiberProperty = new G4MaterialPropertiesTable();
0193   fiberProperty->AddProperty("RINDEX", "PMMA");
0194   fiberProperty->AddProperty("WLSABSLENGTH", wls_Energy, AbsFiber);
0195   fiberProperty->AddProperty("WLSCOMPONENT", wls_Energy, EmissionFib);
0196   fiberProperty->AddConstProperty("WLSTIMECONSTANT", 0.5 * ns);
0197   fPMMA->SetMaterialPropertiesTable(fiberProperty);
0198 
0199   std::vector<G4double> RefractiveIndexClad1 = {1.49, 1.49, 1.49, 1.49};
0200   auto clad1Property = new G4MaterialPropertiesTable();
0201   clad1Property->AddProperty("RINDEX", wls_Energy, RefractiveIndexClad1);
0202   clad1Property->AddProperty("ABSLENGTH", wls_Energy, AbsFiber);
0203   fPethylene1->SetMaterialPropertiesTable(clad1Property);
0204 
0205   std::vector<G4double> RefractiveIndexClad2 = {1.42, 1.42, 1.42, 1.42};
0206   auto clad2Property = new G4MaterialPropertiesTable();
0207   clad2Property->AddProperty("RINDEX", wls_Energy, RefractiveIndexClad2);
0208   clad2Property->AddProperty("ABSLENGTH", wls_Energy, AbsFiber);
0209   fPethylene2->SetMaterialPropertiesTable(clad2Property);
0210 }
0211 
0212 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0213 
0214 G4VPhysicalVolume* LXeDetectorConstruction::Construct()
0215 {
0216   // The experimental hall walls are all 1m away from housing walls
0217   G4double expHall_x = fScint_x + fD_mtl + 1. * m;
0218   G4double expHall_y = fScint_y + fD_mtl + 1. * m;
0219   G4double expHall_z = fScint_z + fD_mtl + 1. * m;
0220 
0221   // Create experimental hall
0222   fExperimentalHall_box = new G4Box("expHall_box", expHall_x, expHall_y, expHall_z);
0223   fExperimentalHall_log = new G4LogicalVolume(fExperimentalHall_box, fVacuum, "expHall_log");
0224   fExperimentalHall_phys = new G4PVPlacement(nullptr, G4ThreeVector(), fExperimentalHall_log,
0225                                              "expHall", nullptr, false, 0);
0226 
0227   fExperimentalHall_log->SetVisAttributes(G4VisAttributes::GetInvisible());
0228 
0229   // Place the main volume
0230   if (fMainVolumeOn) {
0231     fMainVolume =
0232       new LXeMainVolume(nullptr, G4ThreeVector(), fExperimentalHall_log, false, 0, this);
0233   }
0234 
0235   // Place the WLS slab
0236   if (fWLSslab) {
0237     G4VPhysicalVolume* slab =
0238       new LXeWLSSlab(nullptr, G4ThreeVector(0., 0., -fScint_z / 2. - fSlab_z - 1. * cm),
0239                      fExperimentalHall_log, false, 0, this);
0240 
0241     // Surface properties for the WLS slab
0242     auto scintWrap = new G4OpticalSurface("ScintWrap");
0243 
0244     new G4LogicalBorderSurface("ScintWrap", slab, fExperimentalHall_phys, scintWrap);
0245 
0246     scintWrap->SetType(dielectric_metal);
0247     scintWrap->SetFinish(polished);
0248     scintWrap->SetModel(glisur);
0249 
0250     std::vector<G4double> pp = {2.0 * eV, 3.5 * eV};
0251     std::vector<G4double> reflectivity = {1.0, 1.0};
0252     std::vector<G4double> efficiency = {0.0, 0.0};
0253 
0254     auto scintWrapProperty = new G4MaterialPropertiesTable();
0255 
0256     scintWrapProperty->AddProperty("REFLECTIVITY", pp, reflectivity);
0257     scintWrapProperty->AddProperty("EFFICIENCY", pp, efficiency);
0258     scintWrap->SetMaterialPropertiesTable(scintWrapProperty);
0259   }
0260 
0261   return fExperimentalHall_phys;
0262 }
0263 
0264 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0265 
0266 void LXeDetectorConstruction::ConstructSDandField()
0267 {
0268   if (!fMainVolume) return;
0269 
0270   // PMT SD
0271 
0272   LXePMTSD* pmt = fPmt_SD.Get();
0273   if (!pmt) {
0274     // Created here so it exists as pmts are being placed
0275     G4cout << "Construction /LXeDet/pmtSD" << G4endl;
0276     auto pmt_SD = new LXePMTSD("/LXeDet/pmtSD");
0277     fPmt_SD.Put(pmt_SD);
0278 
0279     pmt_SD->InitPMTs();
0280     pmt_SD->SetPmtPositions(fMainVolume->GetPmtPositions());
0281   }
0282   else {
0283     pmt->InitPMTs();
0284     pmt->SetPmtPositions(fMainVolume->GetPmtPositions());
0285   }
0286   G4SDManager::GetSDMpointer()->AddNewDetector(fPmt_SD.Get());
0287   // sensitive detector is not actually on the photocathode.
0288   // processHits gets done manually by the stepping action.
0289   // It is used to detect when photons hit and get absorbed & detected at the
0290   // boundary to the photocathode (which doesn't get done by attaching it to a
0291   // logical volume.
0292   // It does however need to be attached to something or else it doesn't get
0293   // reset at the begining of events
0294 
0295   SetSensitiveDetector(fMainVolume->GetLogPhotoCath(), fPmt_SD.Get());
0296 
0297   // Scint SD
0298 
0299   if (!fScint_SD.Get()) {
0300     G4cout << "Construction /LXeDet/scintSD" << G4endl;
0301     auto scint_SD = new LXeScintSD("/LXeDet/scintSD");
0302     fScint_SD.Put(scint_SD);
0303   }
0304   G4SDManager::GetSDMpointer()->AddNewDetector(fScint_SD.Get());
0305   SetSensitiveDetector(fMainVolume->GetLogScint(), fScint_SD.Get());
0306 }
0307 
0308 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0309 
0310 void LXeDetectorConstruction::SetDimensions(G4ThreeVector dims)
0311 {
0312   fScint_x = dims[0];
0313   fScint_y = dims[1];
0314   fScint_z = dims[2];
0315   G4RunManager::GetRunManager()->ReinitializeGeometry();
0316 }
0317 
0318 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0319 
0320 void LXeDetectorConstruction::SetHousingThickness(G4double d_mtl)
0321 {
0322   fD_mtl = d_mtl;
0323   G4RunManager::GetRunManager()->ReinitializeGeometry();
0324 }
0325 
0326 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0327 
0328 void LXeDetectorConstruction::SetNX(G4int nx)
0329 {
0330   fNx = nx;
0331   G4RunManager::GetRunManager()->ReinitializeGeometry();
0332 }
0333 
0334 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0335 
0336 void LXeDetectorConstruction::SetNY(G4int ny)
0337 {
0338   fNy = ny;
0339   G4RunManager::GetRunManager()->ReinitializeGeometry();
0340 }
0341 
0342 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0343 
0344 void LXeDetectorConstruction::SetNZ(G4int nz)
0345 {
0346   fNz = nz;
0347   G4RunManager::GetRunManager()->ReinitializeGeometry();
0348 }
0349 
0350 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0351 
0352 void LXeDetectorConstruction::SetPMTRadius(G4double outerRadius_pmt)
0353 {
0354   fOuterRadius_pmt = outerRadius_pmt;
0355   G4RunManager::GetRunManager()->ReinitializeGeometry();
0356 }
0357 
0358 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0359 
0360 void LXeDetectorConstruction::SetDefaults()
0361 {
0362   // Resets to default values
0363   fD_mtl = 0.0635 * cm;
0364 
0365   fScint_x = 17.8 * cm;
0366   fScint_y = 17.8 * cm;
0367   fScint_z = 22.6 * cm;
0368 
0369   fNx = 2;
0370   fNy = 2;
0371   fNz = 3;
0372 
0373   fOuterRadius_pmt = 2.3 * cm;
0374 
0375   fSphereOn = true;
0376   fRefl = 1.0;
0377 
0378   fNfibers = 15;
0379   fWLSslab = false;
0380   fMainVolumeOn = true;
0381   fMainVolume = nullptr;
0382   fSlab_z = 2.5 * mm;
0383 
0384   G4UImanager::GetUIpointer()->ApplyCommand("/LXe/detector/scintYieldFactor 1.");
0385 
0386   if (fLXe_mt) fLXe_mt->AddConstProperty("SCINTILLATIONYIELD", 12000. / MeV);
0387   if (fMPTPStyrene) fMPTPStyrene->AddConstProperty("SCINTILLATIONYIELD", 10. / keV);
0388 }
0389 
0390 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0391 
0392 void LXeDetectorConstruction::SetSphereOn(G4bool b)
0393 {
0394   fSphereOn = b;
0395   G4RunManager::GetRunManager()->ReinitializeGeometry();
0396 }
0397 
0398 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0399 
0400 void LXeDetectorConstruction::SetHousingReflectivity(G4double r)
0401 {
0402   fRefl = r;
0403   G4RunManager::GetRunManager()->ReinitializeGeometry();
0404 }
0405 
0406 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0407 
0408 void LXeDetectorConstruction::SetWLSSlabOn(G4bool b)
0409 {
0410   fWLSslab = b;
0411   G4RunManager::GetRunManager()->ReinitializeGeometry();
0412 }
0413 
0414 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0415 
0416 void LXeDetectorConstruction::SetMainVolumeOn(G4bool b)
0417 {
0418   fMainVolumeOn = b;
0419   G4RunManager::GetRunManager()->ReinitializeGeometry();
0420 }
0421 
0422 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0423 
0424 void LXeDetectorConstruction::SetNFibers(G4int n)
0425 {
0426   fNfibers = n;
0427   G4RunManager::GetRunManager()->ReinitializeGeometry();
0428 }
0429 
0430 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0431 
0432 void LXeDetectorConstruction::SetMainScintYield(G4double y)
0433 {
0434   fLXe_mt->AddConstProperty("SCINTILLATIONYIELD", y / MeV);
0435 }
0436 
0437 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0438 
0439 void LXeDetectorConstruction::SetWLSScintYield(G4double y)
0440 {
0441   fMPTPStyrene->AddConstProperty("SCINTILLATIONYIELD", y / MeV);
0442 }
0443 
0444 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0445 
0446 void LXeDetectorConstruction::SetSaveThreshold(G4int save)
0447 {
0448   // Sets the save threshold for the random number seed. If the number of
0449   // photons generated in an event is lower than this, then save the seed for
0450   // this event in a file called run###evt###.rndm
0451 
0452   fSaveThreshold = save;
0453   G4RunManager::GetRunManager()->SetRandomNumberStore(true);
0454 }
0455 
0456 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......