Back to home page

EIC code displayed by LXR
 
 

    

File indexing completed on 2025-01-18 09:17:05 

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 //
0028 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0029 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0030 
0031 #include "DetectorConstruction.hh"
0032 
0033 #include "G4AutoDelete.hh"
0034 #include "G4Box.hh"
0035 #include "G4Colour.hh"
0036 #include "G4Ellipsoid.hh"
0037 #include "G4GeometryManager.hh"
0038 #include "G4GlobalMagFieldMessenger.hh"
0039 #include "G4LogicalVolume.hh"
0040 #include "G4Material.hh"
0041 #include "G4NistManager.hh"
0042 #include "G4PVPlacement.hh"
0043 #include "G4PhysicalConstants.hh"
0044 #include "G4PhysicalVolumeStore.hh"
0045 #include "G4RunManager.hh"
0046 #include "G4Sphere.hh"
0047 #include "G4SystemOfUnits.hh"
0048 #include "G4UniformMagField.hh"
0049 #include "G4UnitsTable.hh"
0050 #include "G4UserLimits.hh"
0051 #include "G4VisAttributes.hh"
0052 
0053 #include "DetectorMessenger.hh"
0054 
0055 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0056 
0057 DetectorConstruction::DetectorConstruction()
0058   : G4VUserDetectorConstruction(),
0059     fAbsorberMaterial(nullptr),
0060     fWorldMaterial(nullptr),
0061     fSolidWorld(nullptr),
0062     fLogicWorld(nullptr),
0063     fPhysiWorld(nullptr),
0064     fSolidAbsorber(nullptr),
0065     fLogicAbsorber(nullptr),
0066     fPhysiAbsorber(nullptr),
0067     material1(nullptr),
0068     material2(nullptr),
0069     material3(nullptr),
0070     material4(nullptr),
0071     material5(nullptr),
0072     material6(nullptr),
0073     material_GDP(nullptr),
0074     ellipse1(nullptr),
0075     logicEllipse1(nullptr),
0076     physiEllipse1(nullptr),
0077     ellipse2(nullptr),
0078     logicEllipse2(nullptr),
0079     physiEllipse2(nullptr),
0080     ellipse3(nullptr),
0081     logicEllipse3(nullptr),
0082     physiEllipse3(nullptr),
0083     ellipse4(nullptr),
0084     logicEllipse4(nullptr),
0085     physiEllipse4(nullptr),
0086     ellipse5(nullptr),
0087     logicEllipse5(nullptr),
0088     physiEllipse5(nullptr),
0089     ellipse6(nullptr),
0090     logicEllipse6(nullptr),
0091     physiEllipse6(nullptr),
0092     solid_GDP(nullptr),
0093     logic_GDP(nullptr),
0094     physi_GDP(nullptr),
0095     fDetectorMessenger(nullptr)
0096 {
0097   phantom_type = 1;
0098 
0099   DefineMaterials();
0100 
0101   if (phantom_type == 1) {
0102     // default parameter values of the box
0103     fAbsorberThickness = 40. * um;
0104     fAbsorberSizeYZ = 40. * um;
0105 
0106     //    SetAbsorberMaterial("G4_Cu");
0107     SetAbsorberMaterial("Body_10times");
0108     //    SetAbsorberMaterial("Body_real");
0109   }
0110   fXposAbs = 0. * cm;
0111   ComputeGeomParameters();
0112 
0113   // materials
0114   SetWorldMaterial("G4_Galactic");
0115 
0116   // create commands for interactive definition of the box
0117   fDetectorMessenger = new DetectorMessenger(this);
0118 }
0119 
0120 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0121 
0122 DetectorConstruction::~DetectorConstruction()
0123 {
0124   delete fDetectorMessenger;
0125 }
0126 
0127 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0128 
0129 void DetectorConstruction::DefineMaterials()
0130 {
0131   //
0132   // define Elements
0133   //
0134   auto H = new G4Element("Hydrogen", "H", 1, 1.01 * g / mole);
0135   auto C = new G4Element("Carbon", "C", 6, 12.01 * g / mole);
0136   auto N = new G4Element("Nitrogen", "N", 7, 14.01 * g / mole);
0137   auto O = new G4Element("Oxygen", "O", 8, 16.00 * g / mole);
0138   auto P = new G4Element("Phosphorus", "P", 15, 30.97 * g / mole);
0139   auto S = new G4Element("Sulphur", "S", 16, 32.07 * g / mole);
0140   auto Cl = new G4Element("Chlorine", "Cl", 17, 35.45 * g / mole);
0141   auto K = new G4Element("Potassium", "K", 19, 39.10 * g / mole);
0142   auto Ca = new G4Element("Calcium", "Ca", 20, 40.08 * g / mole);
0143   auto Ti = new G4Element("Titanium", "Ti", 22, 47.87 * g / mole);
0144   auto Ge = new G4Element("Germanium", "Ge", 32., 72.59 * g / mole);
0145 
0146   //
0147   // define biological materials
0148   //
0149   auto BioMaterial = new G4Material("C10H17O3N2", 1.218 * g / cm3, 4);
0150   BioMaterial->AddElement(C, 10);
0151   BioMaterial->AddElement(H, 17);
0152   BioMaterial->AddElement(O, 3);
0153   BioMaterial->AddElement(N, 2);
0154 
0155   // FIRST ELLIPSOID    : external surface of the worm ("skin")
0156   // Original (real) composition in dry mass (the worm is dehydrated by
0157   // lyophilization) Each element content is defined by its mass fraction (must
0158   // be normalised to 1)
0159   auto Skin_real = new G4Material("Skin_real", 0.49729 * g / cm3, 5);
0160   Skin_real->AddElement(P, 0.47 * perCent);
0161   Skin_real->AddElement(S, 0.21 * perCent);
0162   Skin_real->AddElement(Cl, 0.05 * perCent);
0163   Skin_real->AddElement(K, 0.48 * perCent);
0164   Skin_real->AddMaterial(BioMaterial, 98.79 * perCent);
0165 
0166   // Skin1: mass fractions of mineral elements 1O times more than original
0167   // (real) mass fractions
0168   auto Skin_10times = new G4Material("Skin_10times", 0.49729 * g / cm3, 5);
0169   Skin_10times->AddElement(P, 4.71 * perCent);
0170   Skin_10times->AddElement(S, 2.10 * perCent);
0171   Skin_10times->AddElement(Cl, 0.46 * perCent);
0172   Skin_10times->AddElement(K, 4.77 * perCent);
0173   Skin_10times->AddMaterial(BioMaterial, 87.96 * perCent);
0174 
0175   // SECOND ELLIPSOID   : "body" of the worm, limited by the internal surface of
0176   // the "skin" Original (real) composition in dry mass (the worm is dehydrated
0177   // by lyophilization) Each element content is defined by its mass fraction
0178   // (must be normalised to 1)
0179   auto Body_real = new G4Material("Body_real", 0.40853 * g / cm3, 2);
0180   Body_real->AddElement(P, 0.72 * perCent);
0181   Body_real->AddMaterial(BioMaterial, 99.28 * perCent);
0182 
0183   // Body_10times: mass fractions of mineral elements 1O times more than
0184   // original (real) mass fractions
0185   auto Body_10times = new G4Material("Body_10times", 0.40853 * g / cm3, 2);
0186   Body_10times->AddElement(P, 7.23 * perCent);
0187   Body_10times->AddMaterial(BioMaterial, 92.77 * perCent);
0188 
0189   // THIRD ELLIPSOID : Nucleus of cell 1 (contained inside ellipsoid 2)
0190   // Original (real) composition in dry mass (the worm is dehydrated by
0191   // lyophilization) Each element content is defined by its mass fraction (must
0192   // be normalised to 1)
0193   auto Cell1_real = new G4Material("Cell1_real", 0.66262 * g / cm3, 3);
0194   Cell1_real->AddElement(P, 0.57 * perCent);
0195   Cell1_real->AddElement(Ca, 1.55 * perCent);
0196   Cell1_real->AddMaterial(BioMaterial, 97.88 * perCent);
0197 
0198   // Cell1_10times: mass fractions of mineral elements 1O times more than
0199   // original (real) mass fractions
0200   auto Cell1_10times = new G4Material("Cell1_10times", 0.66262 * g / cm3, 3);
0201   Cell1_10times->AddElement(P, 5.66 * perCent);
0202   Cell1_10times->AddElement(Ca, 15.49 * perCent);
0203   Cell1_10times->AddMaterial(BioMaterial, 78.85 * perCent);
0204 
0205   // FOURTH ELLIPSOID : Nucleus of cell 2 (contained inside ellipsoid 2)
0206   // Original (real) composition in dry mass (the worm is dehydrated by
0207   // lyophilization) Each element content is defined by its mass fraction (must
0208   // be normalised to 1)
0209   auto Cell2_real = new G4Material("Cell2_real", 0.60227 * g / cm3, 3);
0210   Cell2_real->AddElement(P, 0.61 * perCent);
0211   Cell2_real->AddElement(Ca, 1.44 * perCent);
0212   Cell2_real->AddMaterial(BioMaterial, 97.95 * perCent);
0213 
0214   // Cell2_10times: mass fractions of mineral elements 1O times more than
0215   // original (real) mass fractions
0216   auto Cell2_10times = new G4Material("Cell2_10times", 0.60227 * g / cm3, 3);
0217   Cell2_10times->AddElement(P, 6.10 * perCent);
0218   Cell2_10times->AddElement(Ca, 14.45 * perCent);
0219   Cell2_10times->AddMaterial(BioMaterial, 79.45 * perCent);
0220 
0221   // FIFTH ELLIPSOID: intestine of the worm
0222   // Original (real) composition in dry mass (the worm is dehydrated by
0223   // lyophilization) Each element content is defined by its mass fraction (must
0224   // be normalised to 1)
0225   auto Intestine_real = new G4Material("Intestine_real", 0.54058 * g / cm3, 4);
0226   Intestine_real->AddElement(S, 0.12 * perCent);
0227   Intestine_real->AddElement(Cl, 0.02 * perCent);
0228   Intestine_real->AddElement(K, 0.20 * perCent);
0229   Intestine_real->AddMaterial(BioMaterial, 99.66 * perCent);
0230 
0231   // Intestine_10times: mass fractions of mineral elements 1O times more than
0232   // original (real) mass fractions
0233   auto Intestine_10times = new G4Material("Intestine_10times", 0.54058 * g / cm3, 4);
0234   Intestine_10times->AddElement(S, 1.17 * perCent);
0235   Intestine_10times->AddElement(Cl, 0.23 * perCent);
0236   Intestine_10times->AddElement(K, 1.99 * perCent);
0237   Intestine_10times->AddMaterial(BioMaterial, 96.61 * perCent);
0238 
0239   // SIXTH ELLIPSOID: Nucleus of cell Ti (contained inside ellipsoid 5)
0240   // Original (real) composition in dry mass (the worm is dehydrated by
0241   // lyophilization) Each element content is defined by its mass fraction (must
0242   // be normalised to 1)
0243   auto CellTi_real = new G4Material("CellTi_real", 0.75138 * g / cm3, 5);
0244   CellTi_real->AddElement(S, 0.11 * perCent);
0245   CellTi_real->AddElement(Cl, 0.02 * perCent);
0246   CellTi_real->AddElement(K, 0.18 * perCent);
0247   CellTi_real->AddElement(Ti, 0.05 * perCent);
0248   CellTi_real->AddMaterial(BioMaterial, 99.64 * perCent);
0249 
0250   // CellTi_200times: mass fractions of mineral elements 1O times more than
0251   // original (real) mass fractions except for Ti which is multiplied by 200
0252   // times
0253   auto CellTi_200times = new G4Material("CellTi_200times", 0.75138 * g / cm3, 5);
0254   CellTi_200times->AddElement(S, 1.05 * perCent);
0255   CellTi_200times->AddElement(Cl, 0.22 * perCent);
0256   CellTi_200times->AddElement(K, 1.79 * perCent);
0257   CellTi_200times->AddElement(Ti, 10.89 * perCent);  // Ti density 200 times
0258   CellTi_200times->AddMaterial(BioMaterial, 86.05 * perCent);
0259 
0260   // CellTi_1000times: mass fractions of mineral elements 10 times more than
0261   // original (real) mass fractions except for Ti which is multiplied by 1000
0262   // times
0263   auto CellTi_1000times = new G4Material("CellTi_1000times", 0.75138 * g / cm3, 5);
0264   CellTi_1000times->AddElement(S, 1.05 * perCent);
0265   CellTi_1000times->AddElement(Cl, 0.22 * perCent);
0266   CellTi_1000times->AddElement(K, 1.79 * perCent);
0267   CellTi_1000times->AddElement(Ti, 54.47 * perCent);  // Ti density 1000 times
0268   CellTi_1000times->AddMaterial(BioMaterial, 42.47 * perCent);
0269 
0270   //
0271   // define simple materials  Cl-doped polystyrene capsules (PS)
0272   //
0273   // Polystyrene
0274   auto Cl_Polystyrene = new G4Material("Cl-doped Polystyrene", 1.18425 * g / cm3, 3);
0275   Cl_Polystyrene->AddElement(C, 97);
0276   Cl_Polystyrene->AddElement(H, 97);
0277   Cl_Polystyrene->AddElement(Cl, 6);
0278 
0279   // Ge-doped glow discharge polymer (GDP)
0280   auto GDP = new G4Material("GDP", 1.08 * g / cm3, 3);
0281   GDP->AddElement(C, 76430);
0282   GDP->AddElement(H, 107002);
0283   GDP->AddElement(Ge, 744);
0284 
0285   //
0286   // define a material from elements.
0287   //
0288   auto H2O = new G4Material("Water", 1.000 * g / cm3, 2);
0289   H2O->AddElement(H, 2);
0290   H2O->AddElement(O, 1);
0291   H2O->GetIonisation()->SetMeanExcitationEnergy(78 * eV);
0292 
0293   auto Air = new G4Material("Air", 1.290 * mg / cm3, 2);
0294   Air->AddElement(N, 0.7);
0295   Air->AddElement(O, 0.3);
0296 
0297   //
0298   // example of vacuum
0299   //
0300   G4double density = universe_mean_density;  // from PhysicalConstants.h
0301   G4double pressure = 3.e-18 * pascal;
0302   G4double temperature = 2.73 * kelvin;
0303   new G4Material("Galactic", 1, 1.01 * g / mole, density, kStateGas, temperature, pressure);
0304 }
0305 
0306 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0307 
0308 void DetectorConstruction::ComputeGeomParameters()
0309 {
0310   if (phantom_type == 1) {
0311     // Compute derived parameters of the box
0312     fXstartAbs = fXposAbs - 0.5 * fAbsorberThickness;
0313     fXendAbs = fXposAbs + 0.5 * fAbsorberThickness;
0314 
0315     G4double xmax = std::max(std::abs(fXstartAbs), std::abs(fXendAbs));
0316     fWorldSizeX = 4 * xmax;
0317     fWorldSizeYZ = 2 * fAbsorberSizeYZ;
0318 
0319     if (nullptr != fPhysiWorld) {
0320       ChangeGeometry();
0321     }
0322   }
0323   else if (phantom_type == 2) {
0324     fWorldSizeX = 1 * mm;
0325     fWorldSizeYZ = 1 * mm;
0326   }
0327   else if (phantom_type == 3) {
0328     fWorldSizeX = 1.5 * mm;
0329     fWorldSizeYZ = 1.5 * mm;
0330   }
0331 }
0332 
0333 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0334 
0335 G4VPhysicalVolume* DetectorConstruction::Construct()
0336 {
0337   if (nullptr != fPhysiWorld) {
0338     return fPhysiWorld;
0339   }
0340   // World
0341   //
0342   fSolidWorld = new G4Box("World",  // its name
0343                           fWorldSizeX / 2, fWorldSizeYZ / 2, fWorldSizeYZ / 2);  // its size
0344 
0345   fLogicWorld = new G4LogicalVolume(fSolidWorld,  // its solid
0346                                     fWorldMaterial,  // its material
0347                                     "World");  // its name
0348 
0349   fPhysiWorld = new G4PVPlacement(nullptr,  // no rotation
0350                                   G4ThreeVector(0., 0., 0.),  // at (0,0,0)
0351                                   fLogicWorld,  // its logical volume
0352                                   "World",  // its name
0353                                   nullptr,  // its mother  volume
0354                                   false,  // no boolean operation
0355                                   0);  // copy number
0356   if (phantom_type == 1) {
0357     Construct_Phantom1();
0358   }
0359   else if (phantom_type == 2) {
0360     Construct_Phantom2();
0361   }
0362   else if (phantom_type == 3) {
0363     Construct_Phantom3();
0364   }
0365 
0366   //    fLogicWorld->SetVisAttributes (G4VisAttributes::GetInvisible());
0367 
0368   // always return the physical World
0369   //
0370   return fPhysiWorld;
0371 }
0372 
0373 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0374 
0375 void DetectorConstruction::PrintGeomParameters()
0376 {
0377   if (phantom_type == 1) {
0378     G4cout << "\n" << fWorldMaterial << G4endl;
0379     G4cout << "\n" << fAbsorberMaterial << G4endl;
0380 
0381     G4cout << "\n The  WORLD   is made of " << G4BestUnit(fWorldSizeX, "Length") << " of "
0382            << fWorldMaterial->GetName();
0383     G4cout << ". The transverse size (YZ) of the world is " << G4BestUnit(fWorldSizeYZ, "Length")
0384            << G4endl;
0385     G4cout << " The ABSORBER is made of " << G4BestUnit(fAbsorberThickness, "Length") << " of "
0386            << fAbsorberMaterial->GetName();
0387     G4cout << ". The transverse size (YZ) is " << G4BestUnit(fAbsorberSizeYZ, "Length") << G4endl;
0388     G4cout << " X position of the middle of the absorber " << G4BestUnit(fXposAbs, "Length");
0389     G4cout << G4endl;
0390   }
0391   else if (phantom_type == 2) {
0392     G4cout << "The upper part of C.elegans is made of 6 ellipsoids" << G4endl;
0393   }
0394   else if (phantom_type == 3) {
0395     G4cout << "A microsphere inertial confinement fusion target is contructed, made of: " << G4endl;
0396     G4cout << "\n" << material_GDP << G4endl;
0397   }
0398 }
0399 
0400 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0401 
0402 void DetectorConstruction::SetAbsorberMaterial(const G4String& materialChoice)
0403 {
0404   // search the material by its name
0405   G4Material* pttoMaterial = G4NistManager::Instance()->FindOrBuildMaterial(materialChoice);
0406 
0407   if (pttoMaterial && fAbsorberMaterial != pttoMaterial) {
0408     fAbsorberMaterial = pttoMaterial;
0409     if (fLogicAbsorber) {
0410       fLogicAbsorber->SetMaterial(fAbsorberMaterial);
0411     }
0412     G4RunManager::GetRunManager()->PhysicsHasBeenModified();
0413   }
0414 }
0415 
0416 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0417 
0418 void DetectorConstruction::SetWorldMaterial(const G4String& materialChoice)
0419 {
0420   // search the material by its name
0421   G4Material* pttoMaterial = G4NistManager::Instance()->FindOrBuildMaterial(materialChoice);
0422 
0423   if (pttoMaterial && fWorldMaterial != pttoMaterial) {
0424     fWorldMaterial = pttoMaterial;
0425     if (fLogicWorld) {
0426       fLogicWorld->SetMaterial(fWorldMaterial);
0427     }
0428     G4RunManager::GetRunManager()->PhysicsHasBeenModified();
0429   }
0430 }
0431 
0432 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0433 
0434 void DetectorConstruction::SetAbsorberThickness(G4double val)
0435 {
0436   fAbsorberThickness = val;
0437   ComputeGeomParameters();
0438 }
0439 
0440 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0441 
0442 void DetectorConstruction::SetAbsorberSizeYZ(G4double val)
0443 {
0444   fAbsorberSizeYZ = val;
0445   ComputeGeomParameters();
0446 }
0447 
0448 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0449 
0450 void DetectorConstruction::SetWorldSizeX(G4double val)
0451 {
0452   fWorldSizeX = val;
0453   ComputeGeomParameters();
0454 }
0455 
0456 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0457 
0458 void DetectorConstruction::SetWorldSizeYZ(G4double val)
0459 {
0460   fWorldSizeYZ = val;
0461   ComputeGeomParameters();
0462 }
0463 
0464 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0465 
0466 void DetectorConstruction::SetAbsorberXpos(G4double val)
0467 {
0468   fXposAbs = val;
0469 }
0470 
0471 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
0472 
0473 void DetectorConstruction::ConstructSDandField()
0474 {
0475   if (fFieldMessenger.Get() == nullptr) {
0476     // Create global magnetic field messenger.
0477     // Uniform magnetic field is then created automatically if
0478     // the field value is not zero.
0479     G4ThreeVector fieldValue = G4ThreeVector();
0480     auto msg = new G4GlobalMagFieldMessenger(fieldValue);
0481     // msg->SetVerboseLevel(1);
0482     G4AutoDelete::Register(msg);
0483     fFieldMessenger.Put(msg);
0484   }
0485 }
0486 
0487 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0488 
0489 void DetectorConstruction::ChangeGeometry()
0490 {
0491   if (phantom_type == 1) {
0492     fSolidWorld->SetXHalfLength(fWorldSizeX * 0.5);
0493     fSolidWorld->SetYHalfLength(fWorldSizeYZ * 0.5);
0494     fSolidWorld->SetZHalfLength(fWorldSizeYZ * 0.5);
0495 
0496     fSolidAbsorber->SetXHalfLength(fAbsorberThickness * 0.5);
0497     fSolidAbsorber->SetYHalfLength(fAbsorberSizeYZ * 0.5);
0498     fSolidAbsorber->SetZHalfLength(fAbsorberSizeYZ * 0.5);
0499   }
0500 }
0501 
0502 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
0503 
0504 void DetectorConstruction::SetPhantomType(G4int value)
0505 {
0506   phantom_type = value;
0507   ComputeGeomParameters();
0508 }
0509 
0510 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
0511 
0512 void DetectorConstruction::Construct_Phantom1()
0513 {
0514   fSolidAbsorber =
0515     new G4Box("Absorber", fAbsorberThickness / 2, fAbsorberSizeYZ / 2, fAbsorberSizeYZ / 2);
0516 
0517   fLogicAbsorber = new G4LogicalVolume(fSolidAbsorber,  // its solid
0518                                        fAbsorberMaterial,  // its material
0519                                        "Absorber");  // its name
0520 
0521   G4RotationMatrix rm;
0522   G4double angle = 0 * deg;
0523   rm.rotateZ(angle);
0524 
0525   fPhysiAbsorber = new G4PVPlacement(G4Transform3D(rm, G4ThreeVector(fXposAbs * um, 0., 0.)),
0526                                      fLogicAbsorber,  // its logical volume
0527                                      "Absorber",  // its name
0528                                      fLogicWorld,  // its mother
0529                                      false,  // no boolean operation
0530                                      0);  // copy number
0531 
0532   auto logVisAtt = new G4VisAttributes(G4Colour(0.670588, 0.333333, 0.862745));
0533   logVisAtt->SetForceSolid(true);
0534   fLogicAbsorber->SetVisAttributes(logVisAtt);
0535 
0536   PrintGeomParameters();
0537 }
0538 
0539 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
0540 
0541 void DetectorConstruction::Construct_Phantom2()
0542 {
0543   G4NistManager* nist = G4NistManager::Instance();
0544   // PHANTOM 3: C. ELEGANS WORM **********************************
0545   // made of 6 ellipsoids
0546 
0547   // FIRST ELLIPSOID    : external surface of the worm ("skin")
0548 
0549   material1 = nist->FindOrBuildMaterial("Skin_10times");
0550   G4ThreeVector center1 =
0551     G4ThreeVector(fXposAbs, 0,
0552                   0);  // position of the center of the first ellipsoid is (fXposAbs,0,0)
0553   // if we change fXposAbs, all the phantom will be translated along X axis
0554   // as the position of the center of all other ellipsoids is defined according
0555   // to fXposAbs
0556   G4double halfx1 = 20.61 * um;  // Semiaxis in X
0557   G4double halfy1 = 21.42 * um;  // Semiaxis in Y
0558   G4double halfz1 = 187.82 * um;  // Semiaxis in Z
0559   G4double pzTopCut1 = 187.82 * um;  // here the top of the ellipsoid is not cut away
0560   G4double pzBottomCut1 = 0 * um;  // here the bottom of the ellipsoid is cut
0561                                    // away below z = 0 (original configuration)
0562   // G4double pzBottomCut1 = 18.07*um;  // test with position of the slice of
0563   // interest for vis.mac visualisation of the slice z = 18.07 µm
0564 
0565   ellipse1 = new G4Ellipsoid("Ellipse1", halfx1, halfy1, halfz1, pzBottomCut1,
0566                              pzTopCut1);  // cut ellipsoid (original configuration)
0567 
0568   logicEllipse1 = new G4LogicalVolume(ellipse1,  // its solid
0569                                       material1,  // its material
0570                                       "Ellipse1");  // its name
0571 
0572   physiEllipse1 = new G4PVPlacement(nullptr,  // no rotation
0573                                     center1,  // its position
0574                                     logicEllipse1,  // its logical volume
0575                                     "Ellipse1",  // its name
0576                                     fLogicWorld,  // its mother
0577                                     false,  // no boolean operation
0578                                     0);  // copy number
0579 
0580   // SECOND ELLIPSOID   : "body" of the worm, limited by the internal surface of
0581   // the "skin"
0582   material2 = nist->FindOrBuildMaterial("Body_10times");
0583   G4ThreeVector center2 = G4ThreeVector(-0.39 * um, 0, 0);  // position according to ellipsoid 1
0584   G4double halfx2 = 18.61 * um;  // Semiaxis in X
0585   G4double halfy2 = 19.01 * um;  //  Semiaxis in Y
0586   G4double halfz2 = 186.64 * um;  // Semiaxis in Z
0587 
0588   G4double pzTopCut2 = 186.64 * um;  // here the top of the ellipsoid is not cut away
0589   G4double pzBottomCut2 = 0 * um;  // here the bottom of the ellipsoid is cut
0590                                    // away below z = 0 (original configuration)
0591   // G4double pzBottomCut2 = 18.07*um; // test with position of the slice of
0592   // interest for vis.mac visualisation of the slice z = 18.07 µm
0593 
0594   // Ellipse2 = new G4Ellipsoid("Ellipse2",halfx2,halfy2,halfz2); // test with
0595   // entire ellipsoid
0596   ellipse2 = new G4Ellipsoid("Ellipse2", halfx2, halfy2, halfz2, pzBottomCut2,
0597                              pzTopCut2);  // cut ellipsoid
0598 
0599   logicEllipse2 = new G4LogicalVolume(ellipse2,  // its solid
0600                                       material2,  // its material
0601                                       "Ellipse2");  // its name
0602 
0603   physiEllipse2 = new G4PVPlacement(nullptr,  // no rotation
0604                                     center2,  // its position
0605                                     logicEllipse2,  // its logical volume
0606                                     "Ellipse2",  // its name
0607                                     logicEllipse1,  // its mother
0608                                     false,  // no boolean operation
0609                                     0);  // copy number
0610 
0611   // THIRD ELLIPSOID : Nucleus of cell 1 (contained inside ellipsoid 2)
0612   material3 = nist->FindOrBuildMaterial("Cell1_10times");
0613   G4ThreeVector center3 =
0614     G4ThreeVector(2.36 * um, -7.09 * um, 18.07 * um);  // position according to ellipsoid 2
0615   G4double halfx3 = 1.95 * um;  // Semiaxis in X
0616   G4double halfy3 = 3.23 * um;  // Semiaxis in Y
0617   G4double halfz3 = 4.32 * um;  // Semiaxis in Z
0618 
0619   // if we want to visualize slice at z = 18.07 µm using vis.mac,
0620   // we should cut all the ellipsoids (cut away all parts below 18.07 µm)
0621   // Visualization test: the lower part of ellipse 3 is cut off
0622   // G4double pzTopCut3 =4.32 *um;  // Visualization test
0623   // G4double pzBottomCut3 = 0*um;  // Visualization test: cut at zero before
0624   // moving ellipsoid 3
0625   // // to position defined just above => it will be cut at 18.07 µm after
0626   // moving
0627 
0628   ellipse3 = new G4Ellipsoid("Ellipse3", halfx3, halfy3, halfz3);  // entire ellipsoid
0629   // ellipse3 = new
0630   // G4Ellipsoid("Ellipse3",halfx3,halfy3,halfz3,pzBottomCut3,pzTopCut3); //
0631   // Visualization test
0632 
0633   logicEllipse3 = new G4LogicalVolume(ellipse3,  // its solid
0634                                       material3,  // its material
0635                                       "Ellipse3");  // its name
0636 
0637   physiEllipse3 = new G4PVPlacement(nullptr,  // no rotation
0638                                     center3,  // its position
0639                                     logicEllipse3,  // its logical volume
0640                                     "Ellipse3",  // its name
0641                                     logicEllipse2,  // its mother
0642                                     false,  // no boolean operation
0643                                     0);  // copy number
0644 
0645   // FOURTH ELLIPSOID : Nucleus of cell 2 (contained inside ellipsoid 2)
0646   material4 = nist->FindOrBuildMaterial("Cell2_10times");
0647   G4ThreeVector center4 =
0648     G4ThreeVector(8.66 * um, -3.15 * um, 18.07 * um);  // position according to ellipsoid 2
0649   G4double halfx4 = 2.08 * um;  // Semiaxis in X
0650   G4double halfy4 = 2.46 * um;  // Semiaxis in Y
0651   G4double halfz4 = 4.32 * um;  // Semiaxis in Z
0652 
0653   // means the lower part of ellipse 4 was cut off
0654   // G4double pzTopCut4 = 4.32*um;
0655   // G4double pzBottomCut4 = 0*um;// Visualization test: cut at zero before
0656   // moving ellipsoid 4
0657   // // to position defined just above => it will be cut at 18.07 µm after
0658   // moving
0659 
0660   ellipse4 = new G4Ellipsoid("Ellipse4", halfx4, halfy4, halfz4);  // entire ellipsoid
0661   // Ellipse4 = new
0662   // G4Ellipsoid("Ellipse4",halfx4,halfy4,halfz4,pzBottomCut4,pzTopCut4); //
0663   // Visualization test
0664 
0665   logicEllipse4 = new G4LogicalVolume(ellipse4,  // its solid
0666                                       material4,  // its material
0667                                       "Ellipse4");  // its name
0668 
0669   physiEllipse4 = new G4PVPlacement(nullptr,  // no rotation
0670                                     center4,  // its position
0671                                     logicEllipse4,  // its logical volume
0672                                     "Ellipse4",  // its name
0673                                     logicEllipse2,  // its mother
0674                                     false,  // no boolean operation
0675                                     0);  // copy number
0676 
0677   // FIFTH ELLIPSOID: intestine of the worm
0678   // Original (real) composition in dry mass (the worm is dehydrated by
0679   // lyophilization) Each element content is defined by its mass fraction (must
0680   // be normalised to 1)
0681   material5 = nist->FindOrBuildMaterial("Intestine_10times");
0682   G4ThreeVector center5 =
0683     G4ThreeVector(1.64 * um, 0.61 * um, 0 * um);  // position according to ellipsoid 2
0684   G4RotationMatrix rm5;
0685   G4double angle5 = -58.986 * deg;  // ellipse 5 is rotated around z axis
0686                                     // (rotation defined according to ellipse 2)
0687   // G4double angle5=0*deg;
0688   rm5.rotateZ(angle5);
0689 
0690   G4double halfx5 = 3.67 * um;  // Semiaxis in X
0691   G4double halfy5 = 16.48 * um;  // Semiaxis in Y
0692   G4double halfz5 = 28.68 * um;  // Semiaxis in Z
0693 
0694   G4double pzTopCut5 = 28.68 * um;  // here the top of the ellipsoid is not cut away
0695   G4double pzBottomCut5 = 0 * um;  // here the bottom of the ellipsoid is cut
0696                                    // away below z = 0 // original configuration
0697   // G4double pzBottomCut5 =18.07*um; // position of the slice of interest    //
0698   // Visualization test
0699 
0700   // ellipse5 = new G4Ellipsoid("Ellipse5",halfx5,halfy5,halfz5); // entire
0701   // ellipsoid // Visualization test
0702   ellipse5 = new G4Ellipsoid("Ellipse5", halfx5, halfy5, halfz5, pzBottomCut5,
0703                              pzTopCut5);  // original configuration
0704 
0705   logicEllipse5 = new G4LogicalVolume(ellipse5,  // its solid
0706                                       material5,  // its material
0707                                       "Ellipse5");  // its name
0708 
0709   physiEllipse5 = new G4PVPlacement(G4Transform3D(rm5, center5),  // rotation
0710                                     logicEllipse5,  // its logical volume
0711                                     "Ellipse5",  // its name
0712                                     logicEllipse2,  // its mother
0713                                     false,  // no boolean operation
0714                                     0);  // copy number
0715 
0716   // SIXTH ELLIPSOID: Nucleus of cell Ti (contained inside ellipsoid 5)
0717   // Original (real) composition in dry mass (the worm is dehydrated by
0718   // lyophilization) Each element content is defined by its mass fraction (must
0719   // be normalised to 1)
0720 
0721   material6 = nist->FindOrBuildMaterial("CellTi_1000times");
0722   G4ThreeVector center6 =
0723     G4ThreeVector(0.005 * um, 5.83 * um, 18.07 * um);  // position according to ellipsoid 5
0724   G4double halfx6 = 1.62 * um;  // Semiaxis in X
0725   G4double halfy6 = 1.95 * um;  // Semiaxis in Y
0726   G4double halfz6 = 1.62 * um;  // Semiaxis in Z
0727 
0728   // means the lower part of ellipse 6 was cut off
0729   // G4double pzTopCut6 = 1.62*um;
0730   // G4double pzBottomCut6 = 0*um;  // Visualization test: cut at zero before
0731   // moving ellipsoid 6
0732   // // to position defined just above => it will be cut at 18.07 µm after
0733   // moving. As ellipsoid 6 is contained in ellipsoid 5, it has been rotated at
0734   // the same time as ellipsoid 5. We want to counterbalance this rotation and
0735   // put ellipsoid 6 in the right direction (not rotated at all). So we must
0736   // rotate ellipsoid 6 of the opposite angle as ellipsoid 5
0737   G4RotationMatrix rm6;
0738   G4double angle6 = 58.986 * deg;  // ellipsoid 6 is rotated around z axis
0739                                    // (rotation defined according to ellipse 5)
0740   // G4double angle6= 0*deg;
0741   rm6.rotateZ(angle6);
0742   ellipse6 = new G4Ellipsoid("Ellipse6", halfx6, halfy6, halfz6);  // entire ellipsoid
0743   // ellipse6 = new
0744   // G4Ellipsoid("Ellipse6",halfx6,halfy6,halfz6,pzBottomCut6,pzTopCut6);    //
0745   // Visualization test
0746   logicEllipse6 = new G4LogicalVolume(ellipse6,  // its solid
0747                                       material6,  // its material
0748                                       "Ellipse6");  // its name
0749 
0750   physiEllipse6 =
0751     new G4PVPlacement(G4Transform3D(rm6, center6),  // So now ellipsoid 6 is vertical again
0752                                                     // (appears not rotated at all)
0753                       logicEllipse6,  // its logical volume
0754                       "Ellipse6",  // its name
0755                       logicEllipse5,  // its mother
0756                       false,  // no boolean operation
0757                       0);  // copy number
0758 
0759   /* logicEllipse1->SetUserLimits(new G4UserLimits(0.5*um,DBL_MAX,DBL_MAX,0));
0760 
0761   logicEllipse2->SetUserLimits(new G4UserLimits(0.5*um,DBL_MAX,DBL_MAX,0)); */
0762 
0763   // As the material composing the C. elegans worm is very transparent to
0764   // photons, the default steps in Geant4 (distance between pre-step and
0765   // post-step point) are too long. To perform a more precise simulation,
0766   // maximal values are defined for step length in each volume, according to its
0767   // size and configuration.
0768 
0769   logicEllipse1->SetUserLimits(new G4UserLimits(0.2 * um));
0770   logicEllipse2->SetUserLimits(new G4UserLimits(0.3 * um));
0771   logicEllipse3->SetUserLimits(new G4UserLimits(0.4 * um));
0772   logicEllipse4->SetUserLimits(new G4UserLimits(0.4 * um));
0773   logicEllipse5->SetUserLimits(new G4UserLimits(0.7 * um));
0774   logicEllipse6->SetUserLimits(new G4UserLimits(0.3 * um));
0775   // fLogicWorld->SetUserLimits(new G4UserLimits(1*um)); // world does not have
0776   // to be simulated as precisely
0777 
0778   //    fLogicWorld ->SetVisAttributes(new
0779   //    G4VisAttributes(G4Colour(1.0,1.0,1.0)));
0780   logicEllipse1->SetVisAttributes(new G4VisAttributes(G4Colour(1.0, 0.8, 1.0, 0.5)));
0781   logicEllipse2->SetVisAttributes(new G4VisAttributes(G4Colour(1.0, 0.5, 1.0, 0.5)));
0782   //    logicEllipse3 ->SetVisAttributes(new
0783   //    G4VisAttributes(G4Colour(0.0,1.0,1.0))); logicEllipse4
0784   //    ->SetVisAttributes(new G4VisAttributes(G4Colour(0.5,0.3,0.5)));
0785   logicEllipse5->SetVisAttributes(new G4VisAttributes(G4Colour(0.5, 0.5, 1.0, 0.5)));
0786   //    logicEllipse6 ->SetVisAttributes(new
0787   //    G4VisAttributes(G4Colour(0.5,0.8,0.2)));
0788 
0789   auto logVisAtt3 = new G4VisAttributes(G4Colour(0.0, 1.0, 1.0));
0790   logVisAtt3->SetForceSolid(true);
0791   logicEllipse3->SetVisAttributes(logVisAtt3);
0792 
0793   auto logVisAtt4 = new G4VisAttributes(G4Colour(0.5, 0.3, 0.5));
0794   logVisAtt4->SetForceSolid(true);
0795   logicEllipse4->SetVisAttributes(logVisAtt4);
0796 
0797   auto logVisAtt6 = new G4VisAttributes(G4Colour(0.5, 0.8, 0.2));
0798   logVisAtt6->SetForceSolid(true);
0799   logicEllipse6->SetVisAttributes(logVisAtt6);
0800 
0801   PrintGeomParameters();
0802 }
0803 
0804 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
0805 
0806 void DetectorConstruction::Construct_Phantom3()
0807 {
0808   G4NistManager* nist = G4NistManager::Instance();
0809   material_GDP = nist->FindOrBuildMaterial("GDP");
0810 
0811   G4double pRmin = 171 * um;  // thickness = 25 um
0812   G4double pRmax = 196 * um;
0813 
0814   //    G4double pRmin = 146.8*um;  // thickness = 10 um
0815   //    G4double pRmax = 156.8*um;
0816 
0817   //    G4double pRmin = 293.6*um;  // thickness = 20 um
0818   //    G4double pRmax = 313.6*um;
0819 
0820   solid_GDP = new G4Sphere("GDP",  // name
0821                            pRmin, pRmax, 0., 2 * pi, 0., pi);
0822 
0823   logic_GDP = new G4LogicalVolume(solid_GDP,  // its solid
0824                                   material_GDP,  // its material
0825                                   "GDP");  // its name
0826 
0827   physi_GDP = new G4PVPlacement(nullptr,  // no rotation
0828                                 G4ThreeVector(fXposAbs, 0., 0.),  // its position
0829                                 logic_GDP,  // its logical volume
0830                                 "GDP",  // its name
0831                                 fLogicWorld,  // its mother
0832                                 false,  // no boolean operation
0833                                 0);  // copy number
0834 
0835   logic_GDP->SetUserLimits(new G4UserLimits(1 * um));
0836 
0837   auto logVisAtt = new G4VisAttributes(G4Colour(0.670588, 0.333333, 0.862745, 0.5));
0838   //    new G4VisAttributes(G4Colour(0.670588, 0.333333, 0.862745));
0839   logVisAtt->SetForceSolid(true);
0840   //  logVisAtt->SetVisibility(true);
0841 
0842   logic_GDP->SetVisAttributes(logVisAtt);
0843 
0844   PrintGeomParameters();
0845 }
0846 
0847 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!