Geant4_10
DetectorConstruction.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
28 //
29 // $Id: DetectorConstruction.cc 67268 2013-02-13 11:38:40Z ihrivnac $
30 //
31 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
32 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
33 
34 #include "DetectorConstruction.hh"
35 #include "DetectorMessenger.hh"
36 
37 #include "G4NistManager.hh"
38 #include "G4Material.hh"
39 #include "G4Box.hh"
40 #include "G4LogicalVolume.hh"
41 #include "G4PVPlacement.hh"
42 #include "G4PVReplica.hh"
43 #include "G4UniformMagField.hh"
44 
45 #include "G4GeometryManager.hh"
46 #include "G4PhysicalVolumeStore.hh"
47 #include "G4LogicalVolumeStore.hh"
48 #include "G4SolidStore.hh"
49 
50 #include "G4UImanager.hh"
51 #include "G4UnitsTable.hh"
52 #include "G4PhysicalConstants.hh"
53 #include "G4SystemOfUnits.hh"
54 #include <iomanip>
55 
56 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
57 
60  fDefaultMaterial(0),fSolidWorld(0),fLogicWorld(0),fPhysiWorld(0),
61  fSolidCalor(0),fLogicCalor(0),fPhysiCalor(0),
62  fSolidLayer(0),fLogicLayer(0),fPhysiLayer(0),
63  fMagField(0),fDetectorMessenger(0)
64 {
65  // default parameter values of the calorimeter
66  fNbOfAbsor = 2;
67  fAbsorThickness[1] = 2.3*mm;
68  fAbsorThickness[2] = 5.7*mm;
69  fNbOfLayers = 50;
70  fCalorSizeYZ = 40.*cm;
71  ComputeCalorParameters();
72 
73  // materials
74  DefineMaterials();
75  SetWorldMaterial("Galactic");
76  SetAbsorMaterial(1,"Lead");
77  SetAbsorMaterial(2,"liquidArgon");
78 
79  // create commands for interactive definition of the calorimeter
80  fDetectorMessenger = new DetectorMessenger(this);
81 }
82 
83 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
84 
86 {
87  delete fDetectorMessenger;
88 }
89 
90 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
91 
93 {
94  return ConstructCalorimeter();
95 }
96 
97 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
98 
99 void DetectorConstruction::DefineMaterials()
100 {
101  // This function illustrates the possible ways to define materials using
102  // G4 database on G4Elements
104  manager->SetVerbose(0);
105  //
106  // define Elements
107  //
108  G4double z,a;
109 
110  G4Element* H = manager->FindOrBuildElement(1);
111  G4Element* C = manager->FindOrBuildElement(6);
112  G4Element* N = manager->FindOrBuildElement(7);
113  G4Element* O = manager->FindOrBuildElement(8);
114  G4Element* Si = manager->FindOrBuildElement(14);
115  G4Element* Ge = manager->FindOrBuildElement(32);
116  G4Element* Sb = manager->FindOrBuildElement(51);
117  G4Element* I = manager->FindOrBuildElement(53);
118  G4Element* Cs = manager->FindOrBuildElement(55);
119  G4Element* Pb = manager->FindOrBuildElement(82);
120  G4Element* Bi = manager->FindOrBuildElement(83);
121 
122  //
123  // define an Element from isotopes, by relative abundance
124  //
125  G4int iz, n; //iz=number of protons in an isotope;
126  // n=number of nucleons in an isotope;
127  G4int ncomponents;
128  G4double abundance;
129 
130  G4Isotope* U5 = new G4Isotope("U235", iz=92, n=235, a=235.01*g/mole);
131  G4Isotope* U8 = new G4Isotope("U238", iz=92, n=238, a=238.03*g/mole);
132 
133  G4Element* U = new G4Element("enriched Uranium", "U", ncomponents=2);
134  U->AddIsotope(U5, abundance= 90.*perCent);
135  U->AddIsotope(U8, abundance= 10.*perCent);
136 
137  //
138  // define simple materials
139  //
141 
142  new G4Material("liquidH2", z=1., a= 1.008*g/mole, density= 70.8*mg/cm3);
143  new G4Material("Aluminium", z=13., a= 26.98*g/mole, density= 2.700*g/cm3);
144  new G4Material("Titanium", z=22., a= 47.867*g/mole, density= 4.54*g/cm3);
145  new G4Material("Iron", z=26., a= 55.85*g/mole, density= 7.870*g/cm3);
146  new G4Material("Copper", z=29., a= 63.55*g/mole, density= 8.960*g/cm3);
147  new G4Material("Tungsten", z=74., a= 183.85*g/mole, density= 19.30*g/cm3);
148  new G4Material("Gold", z=79., a= 196.97*g/mole, density= 19.32*g/cm3);
149  new G4Material("Uranium", z=92., a= 238.03*g/mole, density= 18.95*g/cm3);
150 
151  //
152  // define a material from elements. case 1: chemical molecule
153  //
154  G4int natoms;
155 
156  G4Material* H2O =
157  new G4Material("Water", density= 1.000*g/cm3, ncomponents=2);
158  H2O->AddElement(H, natoms=2);
159  H2O->AddElement(O, natoms=1);
161  H2O->SetChemicalFormula("H_2O");
162 
163  G4Material* CH =
164  new G4Material("Polystyrene", density= 1.032*g/cm3, ncomponents=2);
165  CH->AddElement(C, natoms=1);
166  CH->AddElement(H, natoms=1);
167 
168  G4Material* Sci =
169  new G4Material("Scintillator", density= 1.032*g/cm3, ncomponents=2);
170  Sci->AddElement(C, natoms=9);
171  Sci->AddElement(H, natoms=10);
172 
173  Sci->GetIonisation()->SetBirksConstant(0.126*mm/MeV);
174 
175  G4Material* Lct =
176  new G4Material("Lucite", density= 1.185*g/cm3, ncomponents=3);
177  Lct->AddElement(C, 59.97*perCent);
178  Lct->AddElement(H, 8.07*perCent);
179  Lct->AddElement(O, 31.96*perCent);
180 
181  G4Material* Sili =
182  new G4Material("Silicon", density= 2.330*g/cm3, ncomponents=1);
183  Sili->AddElement(Si, natoms=1);
184 
185  G4Material* SiO2 =
186  new G4Material("quartz", density= 2.200*g/cm3, ncomponents=2);
187  SiO2->AddElement(Si, natoms=1);
188  SiO2->AddElement(O , natoms=2);
189 
190  G4Material* G10 =
191  new G4Material("NemaG10", density= 1.700*g/cm3, ncomponents=4);
192  G10->AddElement(Si, natoms=1);
193  G10->AddElement(O , natoms=2);
194  G10->AddElement(C , natoms=3);
195  G10->AddElement(H , natoms=3);
196 
197  G4Material* CsI =
198  new G4Material("CsI", density= 4.534*g/cm3, ncomponents=2);
199  CsI->AddElement(Cs, natoms=1);
200  CsI->AddElement(I , natoms=1);
202 
203  G4Material* BGO =
204  new G4Material("BGO", density= 7.10*g/cm3, ncomponents=3);
205  BGO->AddElement(O , natoms=12);
206  BGO->AddElement(Ge, natoms= 3);
207  BGO->AddElement(Bi, natoms= 4);
208 
209  //SiNx
210  density= 3.1 *g/cm3;
211  G4Material* SiNx= new G4Material("SiNx", density, ncomponents=3);
212  SiNx-> AddElement(Si, 300);
213  SiNx-> AddElement(N, 310);
214  SiNx-> AddElement(H, 6);
215 
216  //
217  // define gaseous materials using G4 NIST database
218  //
219  G4double fractionmass;
220 
221  G4Material* Air = manager->FindOrBuildMaterial("G4_AIR");
222  manager->ConstructNewGasMaterial("Air20","G4_AIR",293.*kelvin,1.*atmosphere);
223 
224  G4Material* lAr = manager->FindOrBuildMaterial("G4_lAr");
225  G4Material* lArEm3 = new G4Material("liquidArgon", density= 1.390*g/cm3, ncomponents=1);
226  lArEm3->AddMaterial(lAr, fractionmass=1.0);
227 
228  //
229  // define a material from elements and others materials (mixture of mixtures)
230  //
231 
232  G4Material* Lead = new G4Material("Lead", density= 11.35*g/cm3, ncomponents=1);
233  Lead->AddElement(Pb, fractionmass=1.0);
234 
235  G4Material* LeadSb = new G4Material("LeadSb", density= 11.35*g/cm3, ncomponents=2);
236  LeadSb->AddElement(Sb, fractionmass=4.*perCent);
237  LeadSb->AddElement(Pb, fractionmass=96.*perCent);
238 
239  G4Material* Aerog = new G4Material("Aerogel", density= 0.200*g/cm3, ncomponents=3);
240  Aerog->AddMaterial(SiO2, fractionmass=62.5*perCent);
241  Aerog->AddMaterial(H2O , fractionmass=37.4*perCent);
242  Aerog->AddElement (C , fractionmass= 0.1*perCent);
243 
244  //
245  // examples of gas in non STP conditions
246  //
247  G4double temperature, pressure;
248 
249  G4Material* CO2 =
250  new G4Material("CarbonicGas", density= 27.*mg/cm3, ncomponents=2,
251  kStateGas, temperature= 325.*kelvin, pressure= 50.*atmosphere);
252  CO2->AddElement(C, natoms=1);
253  CO2->AddElement(O, natoms=2);
254 
255  G4Material* steam =
256  new G4Material("WaterSteam", density= 1.0*mg/cm3, ncomponents=1,
257  kStateGas, temperature= 273*kelvin, pressure= 1*atmosphere);
258  steam->AddMaterial(H2O, fractionmass=1.);
259 
260  new G4Material("ArgonGas", z=18, a=39.948*g/mole, density= 1.782*mg/cm3,
261  kStateGas, 273.15*kelvin, 1*atmosphere);
262  //
263  // examples of vacuum
264  //
265 
266  density = universe_mean_density; //from PhysicalConstants.h
267  pressure = 3.e-18*pascal;
268  temperature = 2.73*kelvin;
269  new G4Material("Galactic", z=1., a=1.008*g/mole, density,
270  kStateGas,temperature,pressure);
271 
272  density = 1.e-5*g/cm3;
273  pressure = 2.e-2*bar;
274  temperature = STP_Temperature; //from PhysicalConstants.h
275  G4Material* beam =
276  new G4Material("Beam", density, ncomponents=1,
277  kStateGas,temperature,pressure);
278  beam->AddMaterial(Air, fractionmass=1.);
279 
280  // G4cout << *(G4Material::GetMaterialTable()) << G4endl;
281 }
282 
283 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
284 
285 void DetectorConstruction::ComputeCalorParameters()
286 {
287  // Compute derived parameters of the calorimeter
288  fLayerThickness = 0.;
289  for (G4int iAbs=1; iAbs<=fNbOfAbsor; iAbs++) {
290  fLayerThickness += fAbsorThickness[iAbs];
291  }
292  fCalorThickness = fNbOfLayers*fLayerThickness;
293  fWorldSizeX = 1.2*fCalorThickness;
294  fWorldSizeYZ = 1.2*fCalorSizeYZ;
295 }
296 
297 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
298 
299 G4VPhysicalVolume* DetectorConstruction::ConstructCalorimeter()
300 {
301  // complete the Calor parameters definition
302  ComputeCalorParameters();
303 
304  // Cleanup old geometry
309 
310  //
311  // World
312  //
313 
314  fSolidWorld = new G4Box("World", //its name
315  fWorldSizeX/2,fWorldSizeYZ/2,fWorldSizeYZ/2); //its size
316 
317  fLogicWorld = new G4LogicalVolume(fSolidWorld, //its solid
318  fDefaultMaterial, //its material
319  "World"); //its name
320 
321  fPhysiWorld = new G4PVPlacement(0, //no rotation
322  G4ThreeVector(), //at (0,0,0)
323  fLogicWorld, //its fLogical volume
324  "World", //its name
325  0, //its mother volume
326  false, //no boolean operation
327  0); //copy number
328  //
329  // Calorimeter
330  //
331 
332  fSolidCalor = new G4Box("Calorimeter", //its name
333  fCalorThickness/2,fCalorSizeYZ/2,fCalorSizeYZ/2);//size
334 
335  fLogicCalor = new G4LogicalVolume(fSolidCalor, //its solid
336  fDefaultMaterial, //its material
337  "Calorimeter"); //its name
338 
339  fPhysiCalor = new G4PVPlacement(0, //no rotation
340  G4ThreeVector(), //at (0,0,0)
341  fLogicCalor, //its fLogical volume
342  "Calorimeter", //its name
343  fLogicWorld, //its mother volume
344  false, //no boolean operation
345  0); //copy number
346 
347  //
348  // Layers
349  //
350 
351  fSolidLayer = new G4Box("Layer", //its name
352  fLayerThickness/2,fCalorSizeYZ/2,fCalorSizeYZ/2); //size
353 
354  fLogicLayer = new G4LogicalVolume(fSolidLayer, //its solid
355  fDefaultMaterial, //its material
356  "Layer"); //its name
357  if (fNbOfLayers > 1)
358  fPhysiLayer = new G4PVReplica("Layer", //its name
359  fLogicLayer, //its fLogical volume
360  fLogicCalor, //its mother
361  kXAxis, //axis of replication
362  fNbOfLayers, //number of replica
363  fLayerThickness); //witdth of replica
364  else
365  fPhysiLayer = new G4PVPlacement(0, //no rotation
366  G4ThreeVector(), //at (0,0,0)
367  fLogicLayer, //its fLogical volume
368  "Layer", //its name
369  fLogicCalor, //its mother volume
370  false, //no boolean operation
371  0); //copy number
372 
373  //
374  // Absorbers
375  //
376 
377  G4double xfront = -0.5*fLayerThickness;
378  for (G4int k=1; k<=fNbOfAbsor; k++) {
379  fSolidAbsor[k] = new G4Box("Absorber", //its name
380  fAbsorThickness[k]/2,fCalorSizeYZ/2,fCalorSizeYZ/2);
381 
382  fLogicAbsor[k] = new G4LogicalVolume(fSolidAbsor[k], //its solid
383  fAbsorMaterial[k], //its material
384  fAbsorMaterial[k]->GetName());
385 
386  G4double xcenter = xfront+0.5*fAbsorThickness[k];
387  xfront += fAbsorThickness[k];
388  fPhysiAbsor[k] = new G4PVPlacement(0, //no rotation
389  G4ThreeVector(xcenter,0.,0.), //its position
390  fLogicAbsor[k], //its logical volume
391  fAbsorMaterial[k]->GetName(), //its name
392  fLogicLayer, //its mother
393  false, //no boulean operat
394  k); //copy number
395 
396  }
397 
398 
400 
401  //always return the fPhysical World
402  //
403  return fPhysiWorld;
404 }
405 
406 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
407 
409 {
410  G4cout << "\n-------------------------------------------------------------"
411  << "\n ---> The calorimeter is " << fNbOfLayers << " layers of:";
412  for (G4int i=1; i<=fNbOfAbsor; i++)
413  {
414  G4cout << "\n \t" << std::setw(12) << fAbsorMaterial[i]->GetName() <<": "
415  << std::setw(6) << G4BestUnit(fAbsorThickness[i],"Length");
416  }
417  G4cout << "\n-------------------------------------------------------------\n";
418 
419  G4cout << "\n" << fDefaultMaterial << G4endl;
420  for (G4int j=1; j<=fNbOfAbsor; j++)
421  G4cout << "\n" << fAbsorMaterial[j] << G4endl;
422 
423  G4cout << "\n-------------------------------------------------------------\n";
424 }
425 
426 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
427 
429 {
430  // search the material by its name
431  G4Material* pttoMaterial =
433  if (pttoMaterial) fDefaultMaterial = pttoMaterial;
434 }
435 
436 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
437 
439 {
440  // set the number of Layers
441  //
442  if (ival < 1)
443  { G4cout << "\n --->warning from SetfNbOfLayers: "
444  << ival << " must be at least 1. Command refused" << G4endl;
445  return;
446  }
447  fNbOfLayers = ival;
448 }
449 
450 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
451 
453 {
454  // set the number of Absorbers
455  //
456  if (ival < 1 || ival > (MaxAbsor-1))
457  { G4cout << "\n ---> warning from SetfNbOfAbsor: "
458  << ival << " must be at least 1 and and most " << MaxAbsor-1
459  << ". Command refused" << G4endl;
460  return;
461  }
462  fNbOfAbsor = ival;
463 }
464 
465 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
466 
468 {
469  // search the material by its name
470  //
471  if (ival > fNbOfAbsor || ival <= 0)
472  { G4cout << "\n --->warning from SetAbsorMaterial: absor number "
473  << ival << " out of range. Command refused" << G4endl;
474  return;
475  }
476 
477  G4Material* pttoMaterial =
479  if (pttoMaterial) fAbsorMaterial[ival] = pttoMaterial;
480 }
481 
482 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
483 
485 {
486  // change Absorber thickness
487  //
488  if (ival > fNbOfAbsor || ival <= 0)
489  { G4cout << "\n --->warning from SetAbsorThickness: absor number "
490  << ival << " out of range. Command refused" << G4endl;
491  return;
492  }
493  if (val <= DBL_MIN)
494  { G4cout << "\n --->warning from SetAbsorThickness: thickness "
495  << val << " out of range. Command refused" << G4endl;
496  return;
497  }
498  fAbsorThickness[ival] = val;
499 }
500 
501 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
502 
504 {
505  // change the transverse size
506  //
507  if (val <= DBL_MIN)
508  { G4cout << "\n --->warning from SetfCalorSizeYZ: thickness "
509  << val << " out of range. Command refused" << G4endl;
510  return;
511  }
512  fCalorSizeYZ = val;
513 }
514 
515 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
516 
517 #include "G4FieldManager.hh"
519 
521 {
522  //apply a global uniform magnetic field along Z axis
523  //
524  G4FieldManager* fieldMgr
526 
527  if(fMagField) delete fMagField; //delete the existing magn field
528 
529  if(fieldValue!=0.) // create a new one if non nul
530  { fMagField = new G4UniformMagField(G4ThreeVector(0.,0.,fieldValue));
531  fieldMgr->SetDetectorField(fMagField);
532  fieldMgr->CreateChordFinder(fMagField);
533  } else {
534  fMagField = 0;
535  fieldMgr->SetDetectorField(fMagField);
536  }
537 }
538 
539 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
540 
541 #include "G4RunManager.hh"
542 
544 {
545  G4RunManager::GetRunManager()->DefineWorldVolume(ConstructCalorimeter());
546 }
547 
548 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
G4IonisParamMat * GetIonisation() const
Definition: G4Material.hh:224
G4Material * FindOrBuildMaterial(const G4String &name, G4bool isotopes=true, G4bool warning=false)
G4Material * Air
Definition: TRTMaterials.hh:57
void SetChemicalFormula(const G4String &chF)
Definition: G4Material.hh:171
tuple a
Definition: test.py:11
CLHEP::Hep3Vector G4ThreeVector
G4bool SetDetectorField(G4Field *detectorField)
void AddMaterial(G4Material *material, G4double fraction)
Definition: G4Material.cc:450
int universe_mean_density
Definition: hepunit.py:307
Definition: G4Box.hh:63
void SetMeanExcitationEnergy(G4double value)
const G4String & GetName() const
Definition: G4Material.hh:176
G4VPhysicalVolume * Construct()
void SetBirksConstant(G4double value)
int atmosphere
Definition: hepunit.py:151
void SetWorldMaterial(const G4String &)
static void Clean()
Definition: G4SolidStore.cc:79
const G4int MaxAbsor
#define G4BestUnit(a, b)
#define G4_USE_G4BESTUNIT_FOR_VERBOSE 1
float STP_Temperature
Definition: hepunit.py:302
int G4int
Definition: G4Types.hh:78
static G4NistManager * Instance()
virtual void DefineWorldVolume(G4VPhysicalVolume *worldVol, G4bool topologyIsChanged=true)
static G4PhysicalVolumeStore * GetInstance()
string material
Definition: eplot.py:19
G4double density
Definition: TRTMaterials.hh:39
function g(Y1, Y2, PT2)
Definition: hijing1.383.f:5205
Char_t n[5]
G4GLOB_DLL std::ostream G4cout
void SetVerbose(G4int)
void SetAbsorMaterial(G4int, const G4String &)
G4double iz
Definition: TRTMaterials.hh:39
void AddIsotope(G4Isotope *isotope, G4double RelativeAbundance)
Definition: G4Element.cc:151
#define pascal
G4Material * Si
Definition: TRTMaterials.hh:78
static G4LogicalVolumeStore * GetInstance()
static G4SolidStore * GetInstance()
static G4GeometryManager * GetInstance()
static G4TransportationManager * GetTransportationManager()
G4FieldManager * GetFieldManager() const
static G4RunManager * GetRunManager()
Definition: G4RunManager.cc:74
#define DBL_MIN
Definition: templates.hh:75
tuple z
Definition: test.py:28
float perCent
Definition: hepunit.py:239
#define G4endl
Definition: G4ios.hh:61
void SetAbsorThickness(G4int, G4double)
**D E S C R I P T I O N
Definition: HEPEvtcom.cc:77
void OpenGeometry(G4VPhysicalVolume *vol=0)
void AddElement(G4Element *element, G4int nAtoms)
Definition: G4Material.cc:345
G4Material * ConstructNewGasMaterial(const G4String &name, const G4String &nameNist, G4double temp, G4double pres, G4bool isotopes=true)
double G4double
Definition: G4Types.hh:76
void CreateChordFinder(G4MagneticField *detectorMagField)
G4Element * FindOrBuildElement(G4int Z, G4bool isotopes=true)
Air AddElement(elN,.7)
G4Material * CO2
Definition: TRTMaterials.hh:81