ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
DetectorConstruction.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file DetectorConstruction.cc
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 //
30 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
31 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
32 
33 #include "DetectorConstruction.hh"
34 #include "DetectorMessenger.hh"
35 
36 #include "G4Material.hh"
37 #include "G4Box.hh"
38 #include "G4LogicalVolume.hh"
39 #include "G4PVPlacement.hh"
40 #include "G4UniformMagField.hh"
41 
42 #include "G4GeometryManager.hh"
43 #include "G4PhysicalVolumeStore.hh"
44 #include "G4LogicalVolumeStore.hh"
45 #include "G4SolidStore.hh"
46 
47 #include "G4UnitsTable.hh"
48 #include "G4NistManager.hh"
49 #include "G4RunManager.hh"
50 
51 #include "G4PhysicalConstants.hh"
52 #include "G4SystemOfUnits.hh"
53 
54 #include "G4GlobalMagFieldMessenger.hh"
55 #include "G4AutoDelete.hh"
56 
57 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
58 
59 DetectorConstruction::DetectorConstruction()
60  : G4VUserDetectorConstruction(),
61  fAbsorberMaterial(nullptr),fWorldMaterial(nullptr),fDefaultWorld(true),
62  fSolidWorld(nullptr),fLogicWorld(nullptr),fPhysiWorld(nullptr),
63  fSolidAbsorber(nullptr),fLogicAbsorber(nullptr),fPhysiAbsorber(nullptr),
64  fDetectorMessenger(nullptr)
65 {
66  // default parameter values of the calorimeter
67  fAbsorberThickness = 1.*cm;
68  fAbsorberSizeYZ = 2.*cm;
69  fXposAbs = 0.*cm;
70  ComputeGeomParameters();
71 
72  // materials
73  DefineMaterials();
74  SetWorldMaterial ("G4_Galactic");
75  SetAbsorberMaterial("G4_Si");
76 
77  // create commands for interactive definition of the calorimeter
78  fDetectorMessenger = new DetectorMessenger(this);
79 }
80 
81 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
82 
83 DetectorConstruction::~DetectorConstruction()
84 {
85  delete fDetectorMessenger;
86 }
87 
88 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
89 
90 void DetectorConstruction::DefineMaterials()
91 {
92  //This function illustrates the possible ways to define materials
93 
94  G4String symbol; //a=mass of a mole;
95  G4double a, z, density; //z=mean number of protons;
96 
97  G4int ncomponents, natoms;
98  G4double fractionmass;
99  G4double temperature, pressure;
100 
101  //
102  // define Elements
103  //
104 
105  G4Element* H = new G4Element("Hydrogen",symbol="H", z= 1, a= 1.01*g/mole);
106  G4Element* C = new G4Element("Carbon", symbol="C", z= 6, a= 12.01*g/mole);
107  G4Element* N = new G4Element("Nitrogen",symbol="N", z= 7, a= 14.01*g/mole);
108  G4Element* O = new G4Element("Oxygen", symbol="O", z= 8, a= 16.00*g/mole);
109  G4Element* Na = new G4Element("Sodium", symbol="Na", z=11, a= 22.99*g/mole);
110  G4Element* Ar = new G4Element("Argon", symbol="Ar", z=18, a= 39.95*g/mole);
111  G4Element* I = new G4Element("Iodine", symbol="I" , z=53, a= 126.90*g/mole);
112  G4Element* Xe = new G4Element("Xenon", symbol="Xe", z=54, a= 131.29*g/mole);
113 
114  //
115  // define simple materials
116  //
117 
118  new G4Material("H2Liq" , z= 1, a= 1.01*g/mole, density= 70.8*mg/cm3);
119  new G4Material("Beryllium", z= 4, a= 9.01*g/mole, density= 1.848*g/cm3);
120  new G4Material("Aluminium", z=13, a=26.98*g/mole, density= 2.700*g/cm3);
121  new G4Material("Silicon" , z=14, a=28.09*g/mole, density= 2.330*g/cm3);
122 
123  G4Material* lAr =
124  new G4Material("liquidArgon", density= 1.390*g/cm3, ncomponents=1);
125  lAr->AddElement(Ar, natoms=1);
126 
127  new G4Material("Iron", z=26, a= 55.85*g/mole, density= 7.870*g/cm3);
128  new G4Material("Copper", z=29, a= 63.55*g/mole, density= 8.960*g/cm3);
129  new G4Material("Germanium",z=32, a= 72.61*g/mole, density= 5.323*g/cm3);
130  new G4Material("Silver", z=47, a=107.87*g/mole, density= 10.50*g/cm3);
131  new G4Material("Tungsten", z=74, a=183.85*g/mole, density= 19.30*g/cm3);
132  new G4Material("Gold", z=79, a=196.97*g/mole, density= 19.32*g/cm3);
133  new G4Material("Lead", z=82, a=207.19*g/mole, density= 11.35*g/cm3);
134 
135  //
136  // define a material from elements. case 1: chemical molecule
137  //
138 
139  G4Material* H2O = new G4Material("Water",density= 1.000*g/cm3,ncomponents=2);
140  H2O->AddElement(H, natoms=2);
141  H2O->AddElement(O, natoms=1);
142  H2O->GetIonisation()->SetMeanExcitationEnergy(78*eV);
143 
144  G4Material* CH = new G4Material("Plastic",density= 1.04*g/cm3,ncomponents=2);
145  CH->AddElement(C, natoms=1);
146  CH->AddElement(H, natoms=1);
147 
148  G4Material* NaI = new G4Material("NaI", density= 3.67*g/cm3, ncomponents=2);
149  NaI->AddElement(Na, natoms=1);
150  NaI->AddElement(I , natoms=1);
151  NaI->GetIonisation()->SetMeanExcitationEnergy(452*eV);
152 
153  //
154  // define a material from elements. case 2: mixture by fractional mass
155  //
156 
157  G4Material* Air = new G4Material("Air", density= 1.290*mg/cm3, ncomponents=2);
158  Air->AddElement(N, fractionmass=0.7);
159  Air->AddElement(O, fractionmass=0.3);
160 
161  G4Material* Air20 =
162  new G4Material("Air20", density= 1.205*mg/cm3, ncomponents=2,
163  kStateGas, 293.*kelvin, 1.*atmosphere);
164  Air20->AddElement(N, fractionmass=0.7);
165  Air20->AddElement(O, fractionmass=0.3);
166 
167  //Graphite
168  //
169  G4Material* Graphite =
170  new G4Material("Graphite", density= 1.7*g/cm3, ncomponents=1);
171  Graphite->AddElement(C, fractionmass=1.);
172 
173  //Havar
174  //
175  G4Element* Cr = new G4Element("Chrome", "Cr", z=24, a= 51.996*g/mole);
176  G4Element* Fe = new G4Element("Iron" , "Fe", z=26, a= 55.845*g/mole);
177  G4Element* Co = new G4Element("Cobalt", "Co", z=27, a= 58.933*g/mole);
178  G4Element* Ni = new G4Element("Nickel", "Ni", z=28, a= 58.693*g/mole);
179  G4Element* W = new G4Element("Tungsten","W", z=74, a= 183.850*g/mole);
180 
181  G4Material* Havar =
182  new G4Material("Havar", density= 8.3*g/cm3, ncomponents=5);
183  Havar->AddElement(Cr, fractionmass=0.1785);
184  Havar->AddElement(Fe, fractionmass=0.1822);
185  Havar->AddElement(Co, fractionmass=0.4452);
186  Havar->AddElement(Ni, fractionmass=0.1310);
187  Havar->AddElement(W , fractionmass=0.0631);
188 
189  //
190  // examples of gas
191  //
192  new G4Material("ArgonGas", z=18, a=39.948*g/mole, density= 1.782*mg/cm3,
193  kStateGas, 273.15*kelvin, 1*atmosphere);
194 
195  new G4Material("XenonGas", z=54, a=131.29*g/mole, density= 5.458*mg/cm3,
196  kStateGas, 293.15*kelvin, 1*atmosphere);
197 
198  G4Material* CO2 =
199  new G4Material("CarbonicGas", density= 1.977*mg/cm3, ncomponents=2);
200  CO2->AddElement(C, natoms=1);
201  CO2->AddElement(O, natoms=2);
202 
203  G4Material* ArCO2 =
204  new G4Material("ArgonCO2", density= 1.8223*mg/cm3, ncomponents=2);
205  ArCO2->AddElement (Ar, fractionmass=0.7844);
206  ArCO2->AddMaterial(CO2, fractionmass=0.2156);
207 
208  //another way to define mixture of gas per volume
209  G4Material* NewArCO2 =
210  new G4Material("NewArgonCO2", density= 1.8223*mg/cm3, ncomponents=3);
211  NewArCO2->AddElement (Ar, natoms=8);
212  NewArCO2->AddElement (C, natoms=2);
213  NewArCO2->AddElement (O, natoms=4);
214 
215  G4Material* ArCH4 =
216  new G4Material("ArgonCH4", density= 1.709*mg/cm3, ncomponents=3);
217  ArCH4->AddElement (Ar, natoms=93);
218  ArCH4->AddElement (C, natoms=7);
219  ArCH4->AddElement (H, natoms=28);
220 
221  G4Material* XeCH =
222  new G4Material("XenonMethanePropane", density= 4.9196*mg/cm3, ncomponents=3,
223  kStateGas, 293.15*kelvin, 1*atmosphere);
224  XeCH->AddElement (Xe, natoms=875);
225  XeCH->AddElement (C, natoms=225);
226  XeCH->AddElement (H, natoms=700);
227 
228  G4Material* steam =
229  new G4Material("WaterSteam", density= 1.0*mg/cm3, ncomponents=1);
230  steam->AddMaterial(H2O, fractionmass=1.);
231  steam->GetIonisation()->SetMeanExcitationEnergy(71.6*eV);
232 
233  G4Material* rock1 = new G4Material("StandardRock",
234  2.65*CLHEP::g/CLHEP::cm3, 1, kStateSolid);
235  rock1->AddElement(Na, 1);
236 
237  //
238  // example of vacuum
239  //
240  density = universe_mean_density; //from PhysicalConstants.h
241  pressure = 3.e-18*pascal;
242  temperature = 2.73*kelvin;
243  new G4Material("Galactic", z=1, a=1.01*g/mole,density,
244  kStateGas,temperature,pressure);
245 }
246 
247 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
248 
249 void DetectorConstruction::ComputeGeomParameters()
250 {
251  // Compute derived parameters of the calorimeter
252  fXstartAbs = fXposAbs-0.5*fAbsorberThickness;
253  fXendAbs = fXposAbs+0.5*fAbsorberThickness;
254 
255  G4double xmax = std::max(std::abs(fXstartAbs), std::abs(fXendAbs));
256 
257  // change world size by the flag or if the absorber is large
258  if (fDefaultWorld || 2*xmax >= fWorldSizeX ||
259  fAbsorberSizeYZ >= fWorldSizeYZ)
260  {
261  fWorldSizeX = 3*xmax;
262  fWorldSizeYZ= 1.2*fAbsorberSizeYZ;
263  }
264 }
265 
266 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
267 
268 G4VPhysicalVolume* DetectorConstruction::Construct()
269 {
270  if(fPhysiWorld) { return fPhysiWorld; }
271  // World
272  //
273  fSolidWorld = new G4Box("World", //its name
274  fWorldSizeX/2,fWorldSizeYZ/2,fWorldSizeYZ/2); //its size
275 
276  fLogicWorld = new G4LogicalVolume(fSolidWorld, //its solid
277  fWorldMaterial, //its material
278  "World"); //its name
279 
280  fPhysiWorld = new G4PVPlacement(0, //no rotation
281  G4ThreeVector(0.,0.,0.), //at (0,0,0)
282  fLogicWorld, //its logical volume
283  "World", //its name
284  0, //its mother volume
285  false, //no boolean operation
286  0); //copy number
287 
288  // Absorber
289  //
290  fSolidAbsorber = new G4Box("Absorber",
291  fAbsorberThickness/2,fAbsorberSizeYZ/2,fAbsorberSizeYZ/2);
292 
293  fLogicAbsorber = new G4LogicalVolume(fSolidAbsorber, //its solid
294  fAbsorberMaterial, //its material
295  "Absorber"); //its name
296 
297  fPhysiAbsorber = new G4PVPlacement(0, //no rotation
298  G4ThreeVector(fXposAbs,0.,0.), //its position
299  fLogicAbsorber, //its logical volume
300  "Absorber", //its name
301  fLogicWorld, //its mother
302  false, //no boulean operat
303  0); //copy number
304 
305  PrintGeomParameters();
306 
307  //always return the physical World
308  //
309  return fPhysiWorld;
310 }
311 
312 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
313 
314 void DetectorConstruction::PrintGeomParameters()
315 {
316  G4cout << "\n" << fWorldMaterial << G4endl;
317  G4cout << "\n" << fAbsorberMaterial << G4endl;
318 
319  G4cout << "\n The WORLD is made of " << G4BestUnit(fWorldSizeX,"Length")
320  << " of " << fWorldMaterial->GetName();
321  G4cout << ". The transverse size (YZ) of the world is "
322  << G4BestUnit(fWorldSizeYZ,"Length") << G4endl;
323  G4cout << " The ABSORBER is made of "
324  <<G4BestUnit(fAbsorberThickness,"Length")
325  << " of " << fAbsorberMaterial->GetName();
326  G4cout << ". The transverse size (YZ) is "
327  << G4BestUnit(fAbsorberSizeYZ,"Length") << G4endl;
328  G4cout << " X position of the middle of the absorber "
329  << G4BestUnit(fXposAbs,"Length");
330  G4cout << G4endl;
331 }
332 
333 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
334 
335 void DetectorConstruction::SetAbsorberMaterial(const G4String& materialChoice)
336 {
337  // search the material by its name
338  G4Material* pttoMaterial =
339  G4NistManager::Instance()->FindOrBuildMaterial(materialChoice);
340 
341  if (pttoMaterial && fAbsorberMaterial != pttoMaterial) {
342  fAbsorberMaterial = pttoMaterial;
343  if(fLogicAbsorber) { fLogicAbsorber->SetMaterial(fAbsorberMaterial); }
344  G4RunManager::GetRunManager()->PhysicsHasBeenModified();
345  }
346 }
347 
348 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
349 
350 void DetectorConstruction::SetWorldMaterial(const G4String& materialChoice)
351 {
352  // search the material by its name
353  G4Material* pttoMaterial =
354  G4NistManager::Instance()->FindOrBuildMaterial(materialChoice);
355 
356  if (pttoMaterial && fWorldMaterial != pttoMaterial) {
357  fWorldMaterial = pttoMaterial;
358  if(fLogicWorld) { fLogicWorld->SetMaterial(fWorldMaterial); }
359  G4RunManager::GetRunManager()->PhysicsHasBeenModified();
360  }
361 }
362 
363 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
364 
365 void DetectorConstruction::SetAbsorberThickness(G4double val)
366 {
367  fAbsorberThickness = val;
368  ComputeGeomParameters();
369  if(fPhysiWorld) { ChangeGeometry(); }
370 }
371 
372 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
373 
374 void DetectorConstruction::SetAbsorberSizeYZ(G4double val)
375 {
376  fAbsorberSizeYZ = val;
377  ComputeGeomParameters();
378  if(fPhysiWorld) { ChangeGeometry(); }
379 }
380 
381 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
382 
383 void DetectorConstruction::SetWorldSizeX(G4double val)
384 {
385  fWorldSizeX = val;
386  fDefaultWorld = false;
387  ComputeGeomParameters();
388  if(fPhysiWorld) { ChangeGeometry(); }
389 }
390 
391 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
392 
393 void DetectorConstruction::SetWorldSizeYZ(G4double val)
394 {
395  fWorldSizeYZ = val;
396  fDefaultWorld = false;
397  ComputeGeomParameters();
398  if(fPhysiWorld) { ChangeGeometry(); }
399 }
400 
401 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
402 
403 void DetectorConstruction::SetAbsorberXpos(G4double val)
404 {
405  if(!fPhysiWorld) { fXposAbs = val; }
406 }
407 
408 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
409 
410 void DetectorConstruction::ConstructSDandField()
411 {
412  if ( fFieldMessenger.Get() == 0 ) {
413  // Create global magnetic field messenger.
414  // Uniform magnetic field is then created automatically if
415  // the field value is not zero.
416  G4ThreeVector fieldValue = G4ThreeVector();
417  G4GlobalMagFieldMessenger* msg =
418  new G4GlobalMagFieldMessenger(fieldValue);
419  //msg->SetVerboseLevel(1);
420  G4AutoDelete::Register(msg);
421  fFieldMessenger.Put( msg );
422  }
423 }
424 
425 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
426 
427 void DetectorConstruction::ChangeGeometry()
428 {
429  fSolidWorld->SetXHalfLength(fWorldSizeX*0.5);
430  fSolidWorld->SetYHalfLength(fWorldSizeYZ*0.5);
431  fSolidWorld->SetZHalfLength(fWorldSizeYZ*0.5);
432 
433  fSolidAbsorber->SetXHalfLength(fAbsorberThickness*0.5);
434  fSolidAbsorber->SetYHalfLength(fAbsorberSizeYZ*0.5);
435  fSolidAbsorber->SetZHalfLength(fAbsorberSizeYZ*0.5);
436 }
437 
438 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
439