Gabriele Cosmo
394 lines
16 KiB
C++
394 lines
16 KiB
C++
//
|
|
// ********************************************************************
|
|
// * License and Disclaimer *
|
|
// * *
|
|
// * The Geant4 software is copyright of the Copyright Holders of *
|
|
// * the Geant4 Collaboration. It is provided under the terms and *
|
|
// * conditions of the Geant4 Software License, included in the file *
|
|
// * LICENSE and available at http://cern.ch/geant4/license . These *
|
|
// * include a list of copyright holders. *
|
|
// * *
|
|
// * Neither the authors of this software system, nor their employing *
|
|
// * institutes,nor the agencies providing financial support for this *
|
|
// * work make any representation or warranty, express or implied, *
|
|
// * regarding this software system or assume any liability for its *
|
|
// * use. Please see the license in the file LICENSE and URL above *
|
|
// * for the full disclaimer and the limitation of liability. *
|
|
// * *
|
|
// * This code implementation is the result of the scientific and *
|
|
// * technical work of the GEANT4 collaboration. *
|
|
// * By using, copying, modifying or distributing the software (or *
|
|
// * any work based on the software) you agree to acknowledge its *
|
|
// * use in resulting scientific publications, and indicate your *
|
|
// * acceptance of all terms of the Geant4 Software license. *
|
|
// ********************************************************************
|
|
//
|
|
/// \file DetectorConstruction.cc
|
|
/// \brief Implementation of the B5::DetectorConstruction class
|
|
|
|
#include "DetectorConstruction.hh"
|
|
|
|
#include "CellParameterisation.hh"
|
|
#include "Constants.hh"
|
|
#include "DriftChamberSD.hh"
|
|
#include "EmCalorimeterSD.hh"
|
|
#include "HadCalorimeterSD.hh"
|
|
#include "HodoscopeSD.hh"
|
|
#include "MagneticField.hh"
|
|
|
|
#include "G4Box.hh"
|
|
#include "G4Colour.hh"
|
|
#include "G4FieldManager.hh"
|
|
#include "G4GenericMessenger.hh"
|
|
#include "G4LogicalVolume.hh"
|
|
#include "G4Material.hh"
|
|
#include "G4NistManager.hh"
|
|
#include "G4PVParameterised.hh"
|
|
#include "G4PVPlacement.hh"
|
|
#include "G4PVReplica.hh"
|
|
#include "G4RotationMatrix.hh"
|
|
#include "G4RunManager.hh"
|
|
#include "G4SDManager.hh"
|
|
#include "G4SystemOfUnits.hh"
|
|
#include "G4ThreeVector.hh"
|
|
#include "G4Tubs.hh"
|
|
#include "G4UserLimits.hh"
|
|
#include "G4VisAttributes.hh"
|
|
|
|
namespace B5
|
|
{
|
|
|
|
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
|
|
|
|
G4ThreadLocal MagneticField* DetectorConstruction::fMagneticField = nullptr;
|
|
G4ThreadLocal G4FieldManager* DetectorConstruction::fFieldMgr = nullptr;
|
|
|
|
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
|
|
|
|
DetectorConstruction::DetectorConstruction()
|
|
{
|
|
fArmRotation = new G4RotationMatrix();
|
|
fArmRotation->rotateY(fArmAngle);
|
|
|
|
// define commands for this class
|
|
DefineCommands();
|
|
}
|
|
|
|
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
|
|
|
|
DetectorConstruction::~DetectorConstruction()
|
|
{
|
|
delete fArmRotation;
|
|
delete fMessenger;
|
|
}
|
|
|
|
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
|
|
|
|
G4VPhysicalVolume* DetectorConstruction::Construct()
|
|
{
|
|
// Construct materials
|
|
ConstructMaterials();
|
|
auto air = G4Material::GetMaterial("G4_AIR");
|
|
// auto argonGas = G4Material::GetMaterial("_Ar");
|
|
auto argonGas = G4Material::GetMaterial("G4_Ar");
|
|
auto scintillator = G4Material::GetMaterial("G4_PLASTIC_SC_VINYLTOLUENE");
|
|
auto csI = G4Material::GetMaterial("G4_CESIUM_IODIDE");
|
|
auto lead = G4Material::GetMaterial("G4_Pb");
|
|
|
|
// Option to switch on/off checking of volumes overlaps
|
|
//
|
|
G4bool checkOverlaps = true;
|
|
|
|
// geometries --------------------------------------------------------------
|
|
// experimental hall (world volume)
|
|
auto worldSolid = new G4Box("worldBox", 10. * m, 3. * m, 10. * m);
|
|
auto worldLogical = new G4LogicalVolume(worldSolid, air, "worldLogical");
|
|
auto worldPhysical = new G4PVPlacement(nullptr, G4ThreeVector(), worldLogical, "worldPhysical",
|
|
nullptr, false, 0, checkOverlaps);
|
|
|
|
// Tube with Local Magnetic field
|
|
|
|
auto magneticSolid = new G4Tubs("magneticTubs", 0., 1. * m, 1. * m, 0., 360. * deg);
|
|
|
|
fMagneticLogical = new G4LogicalVolume(magneticSolid, air, "magneticLogical");
|
|
|
|
// placement of Tube
|
|
|
|
auto fieldRot = new G4RotationMatrix();
|
|
fieldRot->rotateX(90. * deg);
|
|
new G4PVPlacement(fieldRot, G4ThreeVector(), fMagneticLogical, "magneticPhysical", worldLogical,
|
|
false, 0, checkOverlaps);
|
|
|
|
// set step limit in tube with magnetic field
|
|
auto userLimits = new G4UserLimits(1 * m);
|
|
fMagneticLogical->SetUserLimits(userLimits);
|
|
|
|
// first arm
|
|
auto firstArmSolid = new G4Box("firstArmBox", 1.5 * m, 1. * m, 3. * m);
|
|
auto firstArmLogical = new G4LogicalVolume(firstArmSolid, air, "firstArmLogical");
|
|
new G4PVPlacement(nullptr, G4ThreeVector(0., 0., -5. * m), firstArmLogical, "firstArmPhysical",
|
|
worldLogical, false, 0, checkOverlaps);
|
|
|
|
// second arm
|
|
auto secondArmSolid = new G4Box("secondArmBox", 2. * m, 2. * m, 3.5 * m);
|
|
auto secondArmLogical = new G4LogicalVolume(secondArmSolid, air, "secondArmLogical");
|
|
auto x = -5. * m * std::sin(fArmAngle);
|
|
auto z = 5. * m * std::cos(fArmAngle);
|
|
fSecondArmPhys = new G4PVPlacement(fArmRotation, G4ThreeVector(x, 0., z), secondArmLogical,
|
|
"fSecondArmPhys", worldLogical, false, 0, checkOverlaps);
|
|
|
|
// hodoscopes in first arm
|
|
auto hodoscope1Solid = new G4Box("hodoscope1Box", 5. * cm, 20. * cm, 0.5 * cm);
|
|
fHodoscope1Logical = new G4LogicalVolume(hodoscope1Solid, scintillator, "hodoscope1Logical");
|
|
|
|
for (auto i = 0; i < kNofHodoscopes1; i++) {
|
|
G4double x1 = (i - kNofHodoscopes1 / 2) * 10. * cm;
|
|
new G4PVPlacement(nullptr, G4ThreeVector(x1, 0., -1.5 * m), fHodoscope1Logical,
|
|
"hodoscope1Physical", firstArmLogical, false, i, checkOverlaps);
|
|
}
|
|
|
|
// drift chambers in first arm
|
|
auto chamber1Solid = new G4Box("chamber1Box", 1. * m, 30. * cm, 1. * cm);
|
|
auto chamber1Logical = new G4LogicalVolume(chamber1Solid, argonGas, "chamber1Logical");
|
|
|
|
for (auto i = 0; i < kNofChambers; i++) {
|
|
G4double z1 = (i - kNofChambers / 2) * 0.5 * m;
|
|
new G4PVPlacement(nullptr, G4ThreeVector(0., 0., z1), chamber1Logical, "chamber1Physical",
|
|
firstArmLogical, false, i, checkOverlaps);
|
|
}
|
|
|
|
// "virtual" wire plane
|
|
auto wirePlane1Solid = new G4Box("wirePlane1Box", 1. * m, 30. * cm, 0.1 * mm);
|
|
fWirePlane1Logical = new G4LogicalVolume(wirePlane1Solid, argonGas, "wirePlane1Logical");
|
|
new G4PVPlacement(nullptr, G4ThreeVector(0., 0., 0.), fWirePlane1Logical, "wirePlane1Physical",
|
|
chamber1Logical, false, 0, checkOverlaps);
|
|
|
|
// hodoscopes in second arm
|
|
auto hodoscope2Solid = new G4Box("hodoscope2Box", 5. * cm, 20. * cm, 0.5 * cm);
|
|
fHodoscope2Logical = new G4LogicalVolume(hodoscope2Solid, scintillator, "hodoscope2Logical");
|
|
|
|
for (auto i = 0; i < kNofHodoscopes2; i++) {
|
|
G4double x2 = (i - kNofHodoscopes2 / 2) * 10. * cm;
|
|
new G4PVPlacement(nullptr, G4ThreeVector(x2, 0., 0.), fHodoscope2Logical, "hodoscope2Physical",
|
|
secondArmLogical, false, i, checkOverlaps);
|
|
}
|
|
|
|
// drift chambers in second arm
|
|
auto chamber2Solid = new G4Box("chamber2Box", 1.5 * m, 30. * cm, 1. * cm);
|
|
auto chamber2Logical = new G4LogicalVolume(chamber2Solid, argonGas, "chamber2Logical");
|
|
|
|
for (auto i = 0; i < kNofChambers; i++) {
|
|
G4double z2 = (i - kNofChambers / 2) * 0.5 * m - 1.5 * m;
|
|
new G4PVPlacement(nullptr, G4ThreeVector(0., 0., z2), chamber2Logical, "chamber2Physical",
|
|
secondArmLogical, false, i, checkOverlaps);
|
|
}
|
|
|
|
// "virtual" wire plane
|
|
auto wirePlane2Solid = new G4Box("wirePlane2Box", 1.5 * m, 30. * cm, 0.1 * mm);
|
|
fWirePlane2Logical = new G4LogicalVolume(wirePlane2Solid, argonGas, "wirePlane2Logical");
|
|
new G4PVPlacement(nullptr, G4ThreeVector(0., 0., 0.), fWirePlane2Logical, "wirePlane2Physical",
|
|
chamber2Logical, false, 0, checkOverlaps);
|
|
|
|
// CsI calorimeter
|
|
auto emCalorimeterSolid = new G4Box("EMcalorimeterBox", 1.5 * m, 30. * cm, 15. * cm);
|
|
auto emCalorimeterLogical = new G4LogicalVolume(emCalorimeterSolid, csI, "EMcalorimeterLogical");
|
|
new G4PVPlacement(nullptr, G4ThreeVector(0., 0., 2. * m), emCalorimeterLogical,
|
|
"EMcalorimeterPhysical", secondArmLogical, false, 0, checkOverlaps);
|
|
|
|
// EMcalorimeter cells
|
|
auto cellSolid = new G4Box("cellBox", 7.5 * cm, 7.5 * cm, 15. * cm);
|
|
fCellLogical = new G4LogicalVolume(cellSolid, csI, "cellLogical");
|
|
G4VPVParameterisation* cellParam = new CellParameterisation();
|
|
new G4PVParameterised("cellPhysical", fCellLogical, emCalorimeterLogical, kXAxis, kNofEmCells,
|
|
cellParam);
|
|
|
|
// hadron calorimeter
|
|
auto hadCalorimeterSolid = new G4Box("HadCalorimeterBox", 1.5 * m, 30. * cm, 50. * cm);
|
|
auto hadCalorimeterLogical =
|
|
new G4LogicalVolume(hadCalorimeterSolid, lead, "HadCalorimeterLogical");
|
|
new G4PVPlacement(nullptr, G4ThreeVector(0., 0., 3. * m), hadCalorimeterLogical,
|
|
"HadCalorimeterPhysical", secondArmLogical, false, 0, checkOverlaps);
|
|
|
|
// hadron calorimeter column
|
|
auto HadCalColumnSolid = new G4Box("HadCalColumnBox", 15. * cm, 30. * cm, 50. * cm);
|
|
auto HadCalColumnLogical = new G4LogicalVolume(HadCalColumnSolid, lead, "HadCalColumnLogical");
|
|
new G4PVReplica("HadCalColumnPhysical", HadCalColumnLogical, hadCalorimeterLogical, kXAxis,
|
|
kNofHadColumns, 30. * cm);
|
|
|
|
// hadron calorimeter cell
|
|
auto HadCalCellSolid = new G4Box("HadCalCellBox", 15. * cm, 15. * cm, 50. * cm);
|
|
auto HadCalCellLogical = new G4LogicalVolume(HadCalCellSolid, lead, "HadCalCellLogical");
|
|
new G4PVReplica("HadCalCellPhysical", HadCalCellLogical, HadCalColumnLogical, kYAxis, kNofHadRows,
|
|
30. * cm);
|
|
|
|
// hadron calorimeter layers
|
|
auto HadCalLayerSolid = new G4Box("HadCalLayerBox", 15. * cm, 15. * cm, 2.5 * cm);
|
|
auto HadCalLayerLogical = new G4LogicalVolume(HadCalLayerSolid, lead, "HadCalLayerLogical");
|
|
new G4PVReplica("HadCalLayerPhysical", HadCalLayerLogical, HadCalCellLogical, kZAxis,
|
|
kNofHadCells, 5. * cm);
|
|
|
|
// scintillator plates
|
|
auto HadCalScintiSolid = new G4Box("HadCalScintiBox", 15. * cm, 15. * cm, 0.5 * cm);
|
|
fHadCalScintiLogical =
|
|
new G4LogicalVolume(HadCalScintiSolid, scintillator, "HadCalScintiLogical");
|
|
new G4PVPlacement(nullptr, G4ThreeVector(0., 0., 2. * cm), fHadCalScintiLogical,
|
|
"HadCalScintiPhysical", HadCalLayerLogical, false, 0, checkOverlaps);
|
|
|
|
// visualization attributes ------------------------------------------------
|
|
|
|
G4VisAttributes invisible(G4VisAttributes::GetInvisible());
|
|
G4VisAttributes invisibleBlue(false, G4Colour::Blue());
|
|
G4VisAttributes invisibleGreen(false, G4Colour::Green());
|
|
G4VisAttributes invisibleYellow(false, G4Colour::Yellow());
|
|
G4VisAttributes blue(G4Colour::Blue());
|
|
G4VisAttributes cgray(G4Colour::Gray());
|
|
G4VisAttributes green(G4Colour::Green());
|
|
G4VisAttributes red(G4Colour::Red());
|
|
G4VisAttributes yellow(G4Colour::Yellow());
|
|
|
|
worldLogical->SetVisAttributes(invisible);
|
|
firstArmLogical->SetVisAttributes(invisible);
|
|
secondArmLogical->SetVisAttributes(invisible);
|
|
|
|
fMagneticLogical->SetVisAttributes(cgray);
|
|
fHodoscope1Logical->SetVisAttributes(red);
|
|
fHodoscope2Logical->SetVisAttributes(red);
|
|
|
|
chamber1Logical->SetVisAttributes(green);
|
|
chamber2Logical->SetVisAttributes(green);
|
|
fWirePlane1Logical->SetVisAttributes(invisibleGreen);
|
|
fWirePlane2Logical->SetVisAttributes(invisibleGreen);
|
|
|
|
emCalorimeterLogical->SetVisAttributes(invisibleYellow);
|
|
fCellLogical->SetVisAttributes(yellow);
|
|
|
|
hadCalorimeterLogical->SetVisAttributes(blue);
|
|
HadCalColumnLogical->SetVisAttributes(invisibleBlue);
|
|
HadCalCellLogical->SetVisAttributes(invisibleBlue);
|
|
HadCalLayerLogical->SetVisAttributes(invisibleBlue);
|
|
fHadCalScintiLogical->SetVisAttributes(invisibleBlue);
|
|
|
|
// return the world physical volume ----------------------------------------
|
|
|
|
return worldPhysical;
|
|
}
|
|
|
|
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
|
|
|
|
void DetectorConstruction::ConstructSDandField()
|
|
{
|
|
// sensitive detectors -----------------------------------------------------
|
|
auto sdManager = G4SDManager::GetSDMpointer();
|
|
G4String SDname;
|
|
|
|
auto hodoscope1 = new HodoscopeSD(SDname = "/hodoscope1");
|
|
sdManager->AddNewDetector(hodoscope1);
|
|
fHodoscope1Logical->SetSensitiveDetector(hodoscope1);
|
|
|
|
auto hodoscope2 = new HodoscopeSD(SDname = "/hodoscope2");
|
|
sdManager->AddNewDetector(hodoscope2);
|
|
fHodoscope2Logical->SetSensitiveDetector(hodoscope2);
|
|
|
|
auto chamber1 = new DriftChamberSD(SDname = "/chamber1");
|
|
sdManager->AddNewDetector(chamber1);
|
|
fWirePlane1Logical->SetSensitiveDetector(chamber1);
|
|
|
|
auto chamber2 = new DriftChamberSD(SDname = "/chamber2");
|
|
sdManager->AddNewDetector(chamber2);
|
|
fWirePlane2Logical->SetSensitiveDetector(chamber2);
|
|
|
|
auto emCalorimeter = new EmCalorimeterSD(SDname = "/EMcalorimeter");
|
|
sdManager->AddNewDetector(emCalorimeter);
|
|
fCellLogical->SetSensitiveDetector(emCalorimeter);
|
|
|
|
auto hadCalorimeter = new HadCalorimeterSD(SDname = "/HadCalorimeter");
|
|
sdManager->AddNewDetector(hadCalorimeter);
|
|
fHadCalScintiLogical->SetSensitiveDetector(hadCalorimeter);
|
|
|
|
// magnetic field ----------------------------------------------------------
|
|
fMagneticField = new MagneticField();
|
|
fFieldMgr = new G4FieldManager();
|
|
fFieldMgr->SetDetectorField(fMagneticField);
|
|
fFieldMgr->CreateChordFinder(fMagneticField);
|
|
G4bool forceToAllDaughters = true;
|
|
fMagneticLogical->SetFieldManager(fFieldMgr, forceToAllDaughters);
|
|
}
|
|
|
|
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
|
|
|
|
void DetectorConstruction::ConstructMaterials()
|
|
{
|
|
auto nistManager = G4NistManager::Instance();
|
|
|
|
// Air
|
|
nistManager->FindOrBuildMaterial("G4_AIR");
|
|
|
|
// Argon gas
|
|
nistManager->FindOrBuildMaterial("G4_Ar");
|
|
// With a density different from the one defined in NIST
|
|
// G4double density = 1.782e-03*g/cm3;
|
|
// nistManager->BuildMaterialWithNewDensity("_Ar","G4_Ar",density);
|
|
// !! cases segmentation fault
|
|
|
|
// Scintillator
|
|
// (PolyVinylToluene, C_9H_10)
|
|
nistManager->FindOrBuildMaterial("G4_PLASTIC_SC_VINYLTOLUENE");
|
|
|
|
// CsI
|
|
nistManager->FindOrBuildMaterial("G4_CESIUM_IODIDE");
|
|
|
|
// Lead
|
|
nistManager->FindOrBuildMaterial("G4_Pb");
|
|
|
|
// Vacuum "Galactic"
|
|
// nistManager->FindOrBuildMaterial("G4_Galactic");
|
|
|
|
// Vacuum "Air with low density"
|
|
// auto air = G4Material::GetMaterial("G4_AIR");
|
|
// G4double density = 1.0e-5*air->GetDensity();
|
|
// nistManager
|
|
// ->BuildMaterialWithNewDensity("Air_lowDensity", "G4_AIR", density);
|
|
|
|
G4cout << G4endl << "The materials defined are : " << G4endl << G4endl;
|
|
G4cout << *(G4Material::GetMaterialTable()) << G4endl;
|
|
}
|
|
|
|
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
|
|
|
|
void DetectorConstruction::SetArmAngle(G4double val)
|
|
{
|
|
if (!fSecondArmPhys) {
|
|
G4cerr << "Detector has not yet been constructed." << G4endl;
|
|
return;
|
|
}
|
|
|
|
fArmAngle = val;
|
|
*fArmRotation = G4RotationMatrix(); // make it unit vector
|
|
fArmRotation->rotateY(fArmAngle);
|
|
auto x = -5. * m * std::sin(fArmAngle);
|
|
auto z = 5. * m * std::cos(fArmAngle);
|
|
fSecondArmPhys->SetTranslation(G4ThreeVector(x, 0., z));
|
|
|
|
// tell G4RunManager that we change the geometry
|
|
G4RunManager::GetRunManager()->GeometryHasBeenModified();
|
|
}
|
|
|
|
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
|
|
|
|
void DetectorConstruction::DefineCommands()
|
|
{
|
|
// Define /B5/detector command directory using generic messenger class
|
|
fMessenger = new G4GenericMessenger(this, "/B5/detector/", "Detector control");
|
|
|
|
// armAngle command
|
|
auto& armAngleCmd = fMessenger->DeclareMethodWithUnit(
|
|
"armAngle", "deg", &DetectorConstruction::SetArmAngle, "Set rotation angle of the second arm.");
|
|
armAngleCmd.SetParameterName("angle", true);
|
|
armAngleCmd.SetRange("angle>=0. && angle<180.");
|
|
armAngleCmd.SetDefaultValue("30.");
|
|
}
|
|
|
|
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
|
|
|
|
} // namespace B5
|