Par04/src/Par04DetectorConstruction.cc

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 /// \file Par04DetectorConstruction.cc
0027 /// \brief Implementation of the Par04DetectorConstruction class
0028
0029 #include "Par04DetectorConstruction.hh"
0030
0031 #include "Par04DefineMeshModel.hh"  // for Par04DefineMeshModel
0032 #include "Par04DetectorMessenger.hh"  // for Par04DetectorMessenger
0033 #include "Par04SensitiveDetector.hh"  // for Par04SensitiveDetector
0034
0035 #include "G4Box.hh"  // for G4Box
0036 #include "G4LogicalVolume.hh"  // for G4LogicalVolume
0037 #include "G4Material.hh"  // for G4Material
0038 #include "G4NistManager.hh"  // for G4NistManager
0039 #include "G4PVPlacement.hh"  // for G4PVPlacement
0040 #include "G4Region.hh"  // for G4Region
0041 #include "G4RegionStore.hh"  // for G4RegionStore
0042 #include "G4RunManager.hh"  // for G4RunManager
0043 #include "G4SDManager.hh"  // for G4SDManager
0044 #include "G4Tubs.hh"  // for G4Tubs
0045 #include "G4UnitsTable.hh"  // for operator<<, G4BestUnit
0046 #include "G4VisAttributes.hh"  // for G4VisAttributes
0047
0048 #include <G4Colour.hh>  // for G4Colour
0049 #include <G4Exception.hh>  // for G4Exception
0050 #include <G4ExceptionSeverity.hh>  // for FatalException
0051 #include <G4SystemOfUnits.hh>  // for rad
0052 #include <G4ThreeVector.hh>  // for G4ThreeVector
0053 #include <G4VUserDetectorConstruction.hh>  // for G4VUserDetectorConstruction
0054 #include <G4ios.hh>  // for G4endl, G4cout
0055 #include <algorithm>  // for max
0056 #include <numeric>  // for accumulate
0057 #include <ostream>  // for operator<<, basic_ostream
0058 #include <string>  // for allocator, char_traits
0059 class G4VPhysicalVolume;
0060 #ifdef USE_INFERENCE
0061 #  include "Par04MLFastSimModel.hh"
0062 #endif
0063
0064 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0065
0066 Par04DetectorConstruction::Par04DetectorConstruction() : G4VUserDetectorConstruction()
0067 {
0068   fDetectorMessenger = new Par04DetectorMessenger(this);
0069   G4NistManager* nistManager = G4NistManager::Instance();
0070   fAbsorberMaterial[0] = nistManager->FindOrBuildMaterial("G4_PbWO4");
0071 }
0072
0073 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0074
0075 Par04DetectorConstruction::~Par04DetectorConstruction() = default;
0076
0077 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0078
0079 G4VPhysicalVolume* Par04DetectorConstruction::Construct()
0080 {
0081   //--------- Material definition ---------
0082   G4NistManager* nistManager = G4NistManager::Instance();
0083   G4Material* air = nistManager->FindOrBuildMaterial("G4_AIR");
0084
0085   //--------- Derived dimensions ---------
0086   G4double full2Pi = 2. * CLHEP::pi * rad;
0087   G4double layerThickness =
0088     std::accumulate(fAbsorberThickness.begin(), fAbsorberThickness.end(), 0.);
0089   G4double detectorOuterRadius = fDetectorInnerRadius + fNbOfLayers * layerThickness;
0090   G4double worldSizeXY = detectorOuterRadius * 4.;
0091   G4double worldSizeZ = fDetectorLength * 2;
0092   // check number of materials: (1 = homogeneous calo, 2 = sampling calo)
0093   G4int nbOfMaterials = 0;
0094   for (const auto material : fAbsorberMaterial) {
0095     if (material != nullptr) nbOfMaterials++;
0096   }
0097
0098   //--------- World ---------
0099   auto fSolidWorld = new G4Box("World",  // name
0100                                worldSizeXY / 2.,  // half-width in X
0101                                worldSizeXY / 2.,  // half-width in Y
0102                                worldSizeZ / 2.);  // half-width in Z
0103   auto fLogicWorld = new G4LogicalVolume(fSolidWorld,  // solid
0104                                          air,  // material
0105                                          "World");  // name
0106   auto fPhysicWorld = new G4PVPlacement(0,  // no rotation
0107                                         G4ThreeVector(),  // at (0,0,0)
0108                                         fLogicWorld,  // logical volume
0109                                         "World",  // name
0110                                         0,  // mother volume
0111                                         false,  // not used
0112                                         99999,  // copy number
0113                                         false);  // check overlaps
0114   fLogicWorld->SetVisAttributes(G4VisAttributes::GetInvisible());
0115
0116   //--------- Detector envelope ---------
0117   auto fSolidDetector = new G4Tubs("Detector",  // name
0118                                    fDetectorInnerRadius,  // inner radius
0119                                    detectorOuterRadius,  // outer radius
0120                                    fDetectorLength / 2.,  // half-width in Z
0121                                    0,  // start angle
0122                                    full2Pi);  // delta angle
0123   auto fLogicDetector = new G4LogicalVolume(fSolidDetector,  // solid
0124                                             air,  // material
0125                                             "Detector");  // name
0126   new G4PVPlacement(0,  // no rotation
0127                     G4ThreeVector(0, 0, 0),  // detector centre at (0,0,0)
0128                     fLogicDetector,  // logical volume
0129                     "Detector",  // name
0130                     fLogicWorld,  // mother volume
0131                     false,  // not used
0132                     999,  // copy number
0133                     false);  // check overlaps
0134
0135   // Region for fast simulation
0136   auto detectorRegion = new G4Region("DetectorRegion");
0137   detectorRegion->AddRootLogicalVolume(fLogicDetector);
0138
0139   //--------- Detector layers: material ---------
0140   std::array<G4VisAttributes, 2> attribs;
0141   attribs[0].SetColour(G4Colour(0, 0, 1, 0.1));
0142   attribs[0].SetForceSolid(true);
0143   attribs[1].SetColour(G4Colour(1, 0, 0, 0.1));
0144   attribs[1].SetForceSolid(true);
0145   /// useful variable
0146   G4double innerRadius = fDetectorInnerRadius;
0147   for (G4int iLayer = 0; iLayer < fNbOfLayers; iLayer++) {
0148     for (G4int iMaterial = 0; iMaterial < nbOfMaterials; iMaterial++) {
0149       auto fSolidLayer = new G4Tubs("Layer",  // name
0150                                     innerRadius,  // inner radius
0151                                     innerRadius + fAbsorberThickness[iMaterial],  // outer radius
0152                                     fDetectorLength / 2.,  // half-width in Z
0153                                     0,  // start angle
0154                                     full2Pi);  // delta angle
0155       G4LogicalVolume* logical = new G4LogicalVolume(fSolidLayer,  // solid
0156                                                      fAbsorberMaterial[iMaterial],  // material
0157                                                      "Layer");  // name
0158       new G4PVPlacement(0,  // no rotation
0159                         G4ThreeVector(),  // place at centre of mother volume
0160                         logical,  // logical volume
0161                         "Layer",  // name
0162                         fLogicDetector,  // mother volume
0163                         false,  // not used
0164                         iLayer * nbOfMaterials + iMaterial,  // copy number
0165                         false);  // check overlaps
0166       logical->SetVisAttributes(attribs[iMaterial]);
0167       innerRadius += fAbsorberThickness[iMaterial];
0168       if (fAbsorberSensitivity[iMaterial]) {
0169         fLayerLogical.push_back(logical);
0170       }
0171     }
0172   }
0173
0174   Print();
0175   return fPhysicWorld;
0176 }
0177
0178 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0179
0180 void Par04DetectorConstruction::ConstructSDandField()
0181 {
0182   Par04SensitiveDetector* caloSD =
0183     new Par04SensitiveDetector("sensitiveDetector", fMeshNbOfCells, fMeshSizeOfCells);
0184   G4SDManager::GetSDMpointer()->AddNewDetector(caloSD);
0185   for (const auto logical : fLayerLogical) {
0186     SetSensitiveDetector(logical, caloSD);
0187   }
0188
0189   auto detectorRegion = G4RegionStore::GetInstance()->GetRegion("DetectorRegion");
0190   // Par04DefineMeshModel needs to be first model to call
0191   new Par04DefineMeshModel("defineMesh", detectorRegion);
0192 #ifdef USE_INFERENCE
0193   new Par04MLFastSimModel("inferenceModel", detectorRegion);
0194 #endif
0195 }
0196
0197 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0198
0199 void Par04DetectorConstruction::Print() const
0200 {
0201   G4cout << "\n------------------------------------------------------"
0202          << "\n--- Detector length:\t" << G4BestUnit(fDetectorLength, "Length")
0203          << "\n--- Detector inner radius:\t" << G4BestUnit(fDetectorInnerRadius, "Length")
0204          << "\n--- Number of layers:\t" << fNbOfLayers << G4endl << "\n--- 1st layer: \t"
0205          << G4BestUnit(fAbsorberThickness[0], "Length") << " of "
0206          << (fAbsorberSensitivity[0] ? "active " : "passive ") << fAbsorberMaterial[0]->GetName()
0207          << G4endl;
0208   if (fAbsorberMaterial[1] != nullptr)
0209     G4cout << "--- 2nd layer: \t" << G4BestUnit(fAbsorberThickness[1], "Length") << " of "
0210            << (fAbsorberSensitivity[1] ? "active " : "passive ") << fAbsorberMaterial[1]->GetName()
0211            << G4endl;
0212   G4cout << "-----------------------------------------------------" << G4endl;
0213 }
0214
0215 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0216
0217 void Par04DetectorConstruction::SetAbsorberMaterial(const std::size_t aLayer, const G4String& aName)
0218 {
0219   // search material by its name
0220   G4Material* material = G4NistManager::Instance()->FindOrBuildMaterial(aName);
0221   if (material)
0222     fAbsorberMaterial[aLayer] = material;
0223   else
0224     G4Exception("Par04DetectorConstruction::SetAbsorberMaterial()", "InvalidSetup", FatalException,
0225                 ("Unknown material name: " + aName).c_str());
0226   G4RunManager::GetRunManager()->PhysicsHasBeenModified();
0227 }
0228
0229 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0230
0231 void Par04DetectorConstruction::SetAbsorberThickness(const std::size_t aLayer,
0232                                                      const G4double aThickness)
0233 {
0234   if (aLayer < fAbsorberThickness.size())
0235     fAbsorberThickness[aLayer] = aThickness;
0236   else
0237     G4Exception("Par04DetectorConstruction::SetAbsorberThickness()", "InvalidSetup", FatalException,
0238                 ("Requested layer " + std::to_string(aLayer)
0239                  + " is larger than number of available layers ("
0240                  + std::to_string(fAbsorberThickness.size()) + ").")
0241                   .c_str());
0242   G4RunManager::GetRunManager()->PhysicsHasBeenModified();
0243 }
0244
0245 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0246
0247 void Par04DetectorConstruction::SetAbsorberSensitivity(const std::size_t aLayer,
0248                                                        const G4bool aSensitivity)
0249 {
0250   if (aLayer < fAbsorberSensitivity.size())
0251     fAbsorberSensitivity[aLayer] = aSensitivity;
0252   else
0253     G4Exception(
0254       "Par04DetectorConstruction::SetAbsorberSensitivity()", "InvalidSetup", FatalException,
0255       ("Requested layer " + std::to_string(aLayer) + " is larger than number of available layers ("
0256        + std::to_string(fAbsorberSensitivity.size()) + ").")
0257         .c_str());
0258   G4RunManager::GetRunManager()->PhysicsHasBeenModified();
0259 }
0260
0261 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0262
0263 void Par04DetectorConstruction::SetInnerRadius(G4double aRadius)
0264 {
0265   fDetectorInnerRadius = aRadius;
0266 }
0267
0268 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0269
0270 void Par04DetectorConstruction::SetLength(G4double aLength)
0271 {
0272   fDetectorLength = aLength;
0273 }