electronScattering2/src/ElectronBenchmarkDetector.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 //
0027 /// \file medical/electronScattering2/src/ElectronBenchmarkDetector.cc
0028 /// \brief Implementation of the ElectronBenchmarkDetector class
0029
0030 #include "ElectronBenchmarkDetector.hh"
0031
0032 #include "ElectronBenchmarkDetectorMessenger.hh"
0033
0034 #include "G4Colour.hh"
0035 #include "G4GeometryManager.hh"
0036 #include "G4LogicalVolume.hh"
0037 #include "G4LogicalVolumeStore.hh"
0038 #include "G4Material.hh"
0039 #include "G4MultiFunctionalDetector.hh"
0040 #include "G4NistManager.hh"
0041 #include "G4PSCellFlux.hh"
0042 #include "G4PSPopulation.hh"
0043 #include "G4PVPlacement.hh"
0044 #include "G4PVReplica.hh"
0045 #include "G4PhysicalVolumeStore.hh"
0046 #include "G4RunManager.hh"
0047 #include "G4SDManager.hh"
0048 #include "G4SDParticleFilter.hh"
0049 #include "G4SolidStore.hh"
0050 #include "G4SystemOfUnits.hh"
0051 #include "G4Tubs.hh"
0052 #include "G4UImanager.hh"
0053 #include "G4VPrimitiveScorer.hh"
0054 #include "G4VisAttributes.hh"
0055
0056 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0057
0058 ElectronBenchmarkDetector::ElectronBenchmarkDetector() : G4VUserDetectorConstruction()
0059 {
0060   // Exit Window
0061   fPosWindow0 = 0.000000 * cm;
0062   fPosWindow1 = 0.004120 * cm;
0063
0064   // Scattering Foil
0065   fPosPrimFoil = 2.650000 * cm;
0066   fHalfThicknessPrimFoil = 0.0 * cm;
0067
0068   // Monitor Chamber
0069   fPosMon0 = 5.000000 * cm;
0070   fPosMon1 = 5.011270 * cm;
0071
0072   // Helium Bag
0073   fPosBag0 = 6.497500 * cm;
0074   fPosHelium0 = 6.500000 * cm;
0075   fPosHelium1 = 116.500000 * cm;
0076   fPosBag1 = 116.502500 * cm;
0077   fThicknessRing = 1.4 * cm;
0078
0079   // Scoring Plane
0080   fPosScorer = 118.200000 * cm;
0081   fThicknessScorer = 0.001 * cm;
0082   fWidthScorerRing = 0.1 * cm;
0083
0084   // Radii
0085   fRadOverall = 23.3 * cm;
0086   fRadRingInner = 20.0 * cm;
0087
0088   // Extra space remaining in world volume around apparatus
0089   fPosDelta = 1. * cm;
0090   fRadDelta = 0.1 * cm;
0091
0092   fMessenger = new ElectronBenchmarkDetectorMessenger(this);
0093   DefineMaterials();
0094 }
0095
0096 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0097
0098 ElectronBenchmarkDetector::~ElectronBenchmarkDetector()
0099 {
0100   delete fMessenger;
0101
0102   delete fWorldVisAtt;
0103   delete fWindowVisAtt;
0104   delete fPrimFoilVisAtt;
0105   delete fMonVisAtt;
0106   delete fBagVisAtt;
0107   delete fHeliumVisAtt;
0108   delete fRingVisAtt;
0109   delete fScorerVisAtt;
0110 }
0111
0112 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0113
0114 G4VPhysicalVolume* ElectronBenchmarkDetector::Construct()
0115 {
0116   return CreateGeometry();
0117 }
0118
0119 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0120
0121 void ElectronBenchmarkDetector::DefineMaterials()
0122 {
0123   // Use NIST database for elements and materials whereever possible.
0124   G4NistManager* man = G4NistManager::Instance();
0125   man->SetVerbose(1);
0126
0127   // Take all elements and materials from NIST
0128   man->FindOrBuildMaterial("G4_He");
0129   man->FindOrBuildMaterial("G4_Be");
0130   man->FindOrBuildMaterial("G4_Al");
0131   man->FindOrBuildMaterial("G4_Ti");
0132   man->FindOrBuildMaterial("G4_Ta");
0133   man->FindOrBuildMaterial("G4_AIR");
0134   man->FindOrBuildMaterial("G4_MYLAR");
0135
0136   G4Element* C = man->FindOrBuildElement("C");
0137   G4Element* Cu = man->FindOrBuildElement("Cu");
0138   G4Element* Au = man->FindOrBuildElement("Au");
0139   G4Element* Ti = man->FindOrBuildElement("Ti");
0140   G4Element* Al = man->FindOrBuildElement("Al");
0141   G4Element* V = man->FindOrBuildElement("V");
0142
0143   // Define materials not in NIST.
0144   // While the NIST database does contain default materials for C, Cu and Au,
0145   // those defaults have different densities than the ones used in the
0146   // benchmark specification.
0147   G4double density;
0148   G4int ncomponents;
0149   G4double fractionmass;
0150
0151   G4Material* G4_C = new G4Material("G4_C", density = 2.18 * g / cm3, ncomponents = 1);
0152   G4_C->AddElement(C, fractionmass = 1.00);
0153
0154   G4Material* G4_Cu = new G4Material("G4_Cu", density = 8.92 * g / cm3, ncomponents = 1);
0155   G4_Cu->AddElement(Cu, fractionmass = 1.00);
0156
0157   G4Material* G4_Au = new G4Material("G4_Au", density = 19.30 * g / cm3, ncomponents = 1);
0158   G4_Au->AddElement(Au, fractionmass = 1.00);
0159
0160   G4Material* TiAlloy = new G4Material("TiAlloy", density = 4.42 * g / cm3, ncomponents = 3);
0161   TiAlloy->AddElement(Ti, fractionmass = 0.90);
0162   TiAlloy->AddElement(Al, fractionmass = 0.06);
0163   TiAlloy->AddElement(V, fractionmass = 0.04);
0164
0165   // Print materials table
0166   G4cout << *(G4Material::GetMaterialTable()) << G4endl;
0167 }
0168
0169 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0170
0171 G4VPhysicalVolume* ElectronBenchmarkDetector::CreateGeometry()
0172 {
0173   if (fPhysiWorld) return fPhysiWorld;
0174
0175   // Instantiate the world
0176   fPhysiWorld = CreateWorld();
0177   fLogWorld = fPhysiWorld->GetLogicalVolume();
0178
0179   // Instantiate the geometry
0180   CreateExitWindow(fLogWorld);
0181   CreatePrimaryFoil(fLogWorld);
0182   CreateMonitor(fLogWorld);
0183   CreateHeliumBag(fLogWorld);
0184
0185   // Create the scorers
0186   CreateScorer(fLogWorld);
0187
0188   return fPhysiWorld;
0189 }
0190
0191 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0192
0193 G4VPhysicalVolume* ElectronBenchmarkDetector::CreateWorld()
0194 {
0195   G4double halfLengthWorld = fPosScorer / 2. + fPosDelta;
0196   G4double radWorld = fRadOverall + fRadDelta;
0197   G4VSolid* worldSolid =
0198     new G4Tubs("WorldSolid", 0. * cm, radWorld, halfLengthWorld, 0. * deg, 360. * deg);
0199   G4LogicalVolume* worldLog =
0200     new G4LogicalVolume(worldSolid, G4Material::GetMaterial("G4_AIR"), "WorldLog");
0201
0202   fWorldVisAtt = new G4VisAttributes(G4Colour(1.0, 1.0, 1.0));
0203   worldLog->SetVisAttributes(fWorldVisAtt);
0204
0205   G4VPhysicalVolume* worldPhys =
0206     new G4PVPlacement(0, G4ThreeVector(0., 0., 0.), worldLog, "World", 0, false, 0);
0207
0208   return worldPhys;
0209 }
0210
0211 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0212
0213 void ElectronBenchmarkDetector::CreateExitWindow(G4LogicalVolume* worldLog)
0214 {
0215   G4double halfLengthWorld = fPosScorer / 2.;
0216   G4double halfThicknessWindow = fPosWindow1 / 2.;
0217   G4VSolid* windowSolid =
0218     new G4Tubs("windowSolid", 0. * cm, fRadOverall, halfThicknessWindow, 0. * deg, 360. * deg);
0219   G4LogicalVolume* windowLog =
0220     new G4LogicalVolume(windowSolid, G4Material::GetMaterial("TiAlloy"), "windowLog");
0221
0222   fWindowVisAtt = new G4VisAttributes(G4Colour(0.5, 1.0, 0.5));
0223   windowLog->SetVisAttributes(fWindowVisAtt);
0224
0225   new G4PVPlacement(0, G4ThreeVector(0., 0., halfThicknessWindow - halfLengthWorld), windowLog,
0226                     "ExitWindow", worldLog, false, 0);
0227 }
0228
0229 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0230
0231 void ElectronBenchmarkDetector::CreatePrimaryFoil(G4LogicalVolume* worldLog)
0232 {
0233   G4double halfLengthWorld = fPosScorer / 2.;
0234
0235   // For some energies, we have no Primary Foil.
0236   if (fHalfThicknessPrimFoil == 0.) return;
0237
0238   fSolidPrimFoil =
0239     new G4Tubs("PrimFoilSolid", 0. * cm, fRadOverall, fHalfThicknessPrimFoil, 0. * deg, 360. * deg);
0240   fLogPrimFoil = new G4LogicalVolume(fSolidPrimFoil, fMaterialPrimFoil, "PrimFoilLog");
0241
0242   fPrimFoilVisAtt = new G4VisAttributes(G4Colour(0.5, 1.0, 0.5));
0243   fLogPrimFoil->SetVisAttributes(fPrimFoilVisAtt);
0244
0245   new G4PVPlacement(0,
0246                     G4ThreeVector(0., 0., fPosPrimFoil + fHalfThicknessPrimFoil - halfLengthWorld),
0247                     fLogPrimFoil, "ScatteringFoil", worldLog, false, 0);
0248 }
0249
0250 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0251
0252 void ElectronBenchmarkDetector::CreateMonitor(G4LogicalVolume* worldLog)
0253 {
0254   G4double halfLengthWorld = fPosScorer / 2.;
0255   G4double halfThicknessMon = (fPosMon1 - fPosMon0) / 2.;
0256   G4VSolid* monSolid =
0257     new G4Tubs("monSolid", 0. * cm, fRadOverall, halfThicknessMon, 0. * deg, 360. * deg);
0258   G4LogicalVolume* monLog =
0259     new G4LogicalVolume(monSolid, G4Material::GetMaterial("G4_MYLAR"), "monLog");
0260
0261   fMonVisAtt = new G4VisAttributes(G4Colour(0.5, 1.0, 0.5));
0262   monLog->SetVisAttributes(fMonVisAtt);
0263
0264   new G4PVPlacement(0, G4ThreeVector(0., 0., fPosMon0 + halfThicknessMon - halfLengthWorld), monLog,
0265                     "MonitorChamber", worldLog, false, 0);
0266 }
0267
0268 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0269
0270 void ElectronBenchmarkDetector::CreateHeliumBag(G4LogicalVolume* worldLog)
0271 {
0272   G4double halfLengthWorld = fPosScorer / 2.;
0273
0274   // Construct cylinder of Mylar
0275   G4double halfThicknessBag = (fPosBag1 - fPosBag0) / 2.;
0276   G4VSolid* bagSolid =
0277     new G4Tubs("bagSolid", 0. * cm, fRadOverall, halfThicknessBag, 0. * deg, 360. * deg);
0278   G4LogicalVolume* bagLog =
0279     new G4LogicalVolume(bagSolid, G4Material::GetMaterial("G4_MYLAR"), "bagLog");
0280
0281   fBagVisAtt = new G4VisAttributes(G4Colour(0.5, 1.0, 0.5));
0282   bagLog->SetVisAttributes(fBagVisAtt);
0283
0284   new G4PVPlacement(0, G4ThreeVector(0., 0., fPosBag0 + halfThicknessBag - halfLengthWorld), bagLog,
0285                     "HeliumBag", worldLog, false, 0);
0286
0287   // Insert cylinder of Helium into the Cylinder of Mylar
0288   G4double halfThicknessHelium = (fPosHelium1 - fPosHelium0) / 2.;
0289   G4VSolid* heliumSolid =
0290     new G4Tubs("heliumSolid", 0. * cm, fRadOverall, halfThicknessHelium, 0. * deg, 360. * deg);
0291   G4LogicalVolume* heliumLog =
0292     new G4LogicalVolume(heliumSolid, G4Material::GetMaterial("G4_He"), "heliumLog");
0293
0294   fHeliumVisAtt = new G4VisAttributes(G4Colour(0.5, 1.0, 0.5));
0295   heliumLog->SetVisAttributes(fHeliumVisAtt);
0296
0297   new G4PVPlacement(0, G4ThreeVector(0., 0., 0.), heliumLog, "Helium", bagLog, false, 0);
0298
0299   // Insert two rings of Aluminum into the Cylinder of Helium
0300   G4double halfThicknessRing = fThicknessRing / 2.;
0301   G4VSolid* ringSolid =
0302     new G4Tubs("ringSolid", fRadRingInner, fRadOverall, halfThicknessRing, 0. * deg, 360. * deg);
0303   G4LogicalVolume* ring0Log =
0304     new G4LogicalVolume(ringSolid, G4Material::GetMaterial("G4_Al"), "ring0Log");
0305   G4LogicalVolume* ring1Log =
0306     new G4LogicalVolume(ringSolid, G4Material::GetMaterial("G4_Al"), "ring1Log");
0307
0308   fRingVisAtt = new G4VisAttributes(G4Colour(0.5, 1.0, 0.5));
0309   ring0Log->SetVisAttributes(fRingVisAtt);
0310   ring1Log->SetVisAttributes(fRingVisAtt);
0311
0312   new G4PVPlacement(0, G4ThreeVector(0., 0., -halfThicknessHelium + halfThicknessRing), ring0Log,
0313                     "Ring0", heliumLog, false, 0);
0314
0315   new G4PVPlacement(0, G4ThreeVector(0., 0., halfThicknessHelium - halfThicknessRing), ring1Log,
0316                     "Ring1", heliumLog, false, 0);
0317 }
0318
0319 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0320
0321 void ElectronBenchmarkDetector::CreateScorer(G4LogicalVolume* worldLog)
0322 {
0323   G4double halfLengthWorld = fPosScorer / 2.;
0324   G4double halfThicknessScorer = fThicknessScorer / 2.;
0325
0326   G4VSolid* scorerSolid =
0327     new G4Tubs("scorerSolid", 0. * cm, fRadOverall, halfThicknessScorer, 0. * deg, 360. * deg);
0328   G4LogicalVolume* scorerLog =
0329     new G4LogicalVolume(scorerSolid, G4Material::GetMaterial("G4_AIR"), "scorerLog");
0330
0331   fScorerVisAtt = new G4VisAttributes(G4Colour(0.5, 1.0, 0.5));
0332   scorerLog->SetVisAttributes(fScorerVisAtt);
0333   new G4PVPlacement(0, G4ThreeVector(0., 0., halfLengthWorld - halfThicknessScorer), scorerLog,
0334                     "Scorer", worldLog, false, 0);
0335
0336   G4VSolid* scorerRingSolid =
0337     new G4Tubs("scorerRingSolid", 0. * cm, fRadOverall, halfThicknessScorer, 0. * deg, 360. * deg);
0338   fScorerRingLog =
0339     new G4LogicalVolume(scorerRingSolid, G4Material::GetMaterial("G4_AIR"), "scorerRingLog");
0340   new G4PVReplica("ScorerRing", fScorerRingLog, scorerLog, kRho,
0341                   G4int(fRadOverall / fWidthScorerRing), fWidthScorerRing);
0342 }
0343
0344 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0345
0346 // Note that this method is called both at start of job and again after
0347 // any command causes a change to detector geometry
0348 void ElectronBenchmarkDetector::ConstructSDandField()
0349 {
0350   G4SDManager::GetSDMpointer()->SetVerboseLevel(1);
0351
0352   // G4Cache mechanism is necessary for multi-threaded operation
0353   // as it allows us to store separate detector pointer per thread
0354   G4MultiFunctionalDetector*& sensitiveDetector = fSensitiveDetectorCache.Get();
0355
0356   if (!sensitiveDetector) {
0357     sensitiveDetector = new G4MultiFunctionalDetector("MyDetector");
0358
0359     G4VPrimitiveScorer* primitive;
0360
0361     G4SDParticleFilter* electronFilter = new G4SDParticleFilter("electronFilter", "e-");
0362
0363     primitive = new G4PSCellFlux("cell flux");
0364     sensitiveDetector->RegisterPrimitive(primitive);
0365
0366     primitive = new G4PSCellFlux("e cell flux");
0367     primitive->SetFilter(electronFilter);
0368     sensitiveDetector->RegisterPrimitive(primitive);
0369
0370     primitive = new G4PSPopulation("population");
0371     sensitiveDetector->RegisterPrimitive(primitive);
0372
0373     primitive = new G4PSPopulation("e population");
0374     primitive->SetFilter(electronFilter);
0375     sensitiveDetector->RegisterPrimitive(primitive);
0376   }
0377   G4SDManager::GetSDMpointer()->AddNewDetector(sensitiveDetector);
0378   fScorerRingLog->SetSensitiveDetector(sensitiveDetector);
0379 }
0380
0381 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0382
0383 void ElectronBenchmarkDetector::SetPrimFoilMaterial(const G4String& matname)
0384 {
0385   G4Material* material = G4NistManager::Instance()->FindOrBuildMaterial(matname);
0386
0387   if (material && material != fMaterialPrimFoil) {
0388     fMaterialPrimFoil = material;
0389     if (fLogPrimFoil) {
0390       fLogPrimFoil->SetMaterial(fMaterialPrimFoil);
0391     }
0392     G4RunManager::GetRunManager()->PhysicsHasBeenModified();
0393   }
0394 }
0395
0396 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0397
0398 void ElectronBenchmarkDetector::SetPrimFoilThickness(G4double thicknessPrimFoil)
0399 {
0400   fHalfThicknessPrimFoil = thicknessPrimFoil / 2.;
0401 }
0402
0403 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......