source: trunk/examples/extended/analysis/A01/src/A01DetectorConstruction.cc @ 1341

Last change on this file since 1341 was 1337, checked in by garnier, 14 years ago

tag geant4.9.4 beta 1 + modifs locales
File size: 17.6 KB
Line 
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//
26// $Id: A01DetectorConstruction.cc,v 1.10 2009/11/21 00:22:55 perl Exp $
27// --------------------------------------------------------------
28//
29
30#include "A01DetectorConstruction.hh"
31
32#include "G4FieldManager.hh"
33#include "G4TransportationManager.hh"
34#include "G4Mag_UsualEqRhs.hh"
35
36#include "G4Material.hh"
37#include "G4Element.hh"
38#include "G4MaterialTable.hh"
39#include "G4NistManager.hh"
40
41#include "G4VSolid.hh"
42#include "G4Box.hh"
43#include "G4Tubs.hh"
44#include "G4LogicalVolume.hh"
45#include "G4VPhysicalVolume.hh"
46#include "G4PVPlacement.hh"
47#include "G4PVParameterised.hh"
48#include "G4UserLimits.hh"
49
50#include "G4SDManager.hh"
51#include "G4VSensitiveDetector.hh"
52#include "G4RunManager.hh"
53
54#include "G4VisAttributes.hh"
55#include "G4Colour.hh"
56
57#include "G4ios.hh"
58
59#include "A01DetectorConstMessenger.hh"
60#include "A01MagneticField.hh"
61#include "A01CellParameterisation.hh"
62#include "A01Hodoscope.hh"
63#include "A01DriftChamber.hh"
64#include "A01EmCalorimeter.hh"
65
66#include "G4PVReplica.hh"
67#include "A01HadCalorimeter.hh"
68
69A01DetectorConstruction::A01DetectorConstruction()
70 : air(0), argonGas(0), scintillator(0), CsI(0), lead(0),
71   worldVisAtt(0), magneticVisAtt(0),
72   armVisAtt(0), hodoscopeVisAtt(0), chamberVisAtt(0),
73   wirePlaneVisAtt(0), EMcalorimeterVisAtt(0), cellVisAtt(0),
74   HadCalorimeterVisAtt(0), HadCalorimeterCellVisAtt(0),
75   armAngle(30.*deg), secondArmPhys(0)
76
77{
78  messenger = new A01DetectorConstMessenger(this);
79  magneticField = new A01MagneticField();
80  fieldMgr = new G4FieldManager();
81  armRotation = new G4RotationMatrix();
82  armRotation->rotateY(armAngle);
83}
84
85A01DetectorConstruction::~A01DetectorConstruction()
86{
87  delete armRotation;
88  delete magneticField;
89  delete fieldMgr;
90  delete messenger;
91
92  DestroyMaterials();
93
94  delete worldVisAtt;
95  delete magneticVisAtt;
96  delete armVisAtt;
97  delete hodoscopeVisAtt;
98  delete chamberVisAtt;
99  delete wirePlaneVisAtt;
100  delete EMcalorimeterVisAtt;
101  delete cellVisAtt;
102  delete HadCalorimeterVisAtt;
103  delete HadCalorimeterCellVisAtt;
104}
105
106G4VPhysicalVolume* A01DetectorConstruction::Construct()
107{
108  // All managed (deleted) by SDManager
109  G4VSensitiveDetector* hodoscope1;
110  G4VSensitiveDetector* hodoscope2;
111  G4VSensitiveDetector* chamber1;
112  G4VSensitiveDetector* chamber2;
113  G4VSensitiveDetector* EMcalorimeter;
114//@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
115  G4VSensitiveDetector* HadCalorimeter;
116//@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
117  ConstructMaterials();
118
119  // Local Magnetic Field
120  static G4bool fieldIsInitialized = false;
121
122  if(!fieldIsInitialized)
123  {
124    fieldMgr->SetDetectorField(magneticField);
125    fieldMgr->CreateChordFinder(magneticField);
126    fieldIsInitialized = true;
127  }
128
129  // geometries --------------------------------------------------------------
130  // experimental hall (world volume)
131  G4VSolid* worldSolid = new G4Box("worldBox",10.*m,3.*m,10.*m);
132  G4LogicalVolume* worldLogical
133    = new G4LogicalVolume(worldSolid,air,"worldLogical",0,0,0);
134  G4VPhysicalVolume* worldPhysical
135    = new G4PVPlacement(0,G4ThreeVector(),worldLogical,"worldPhysical",0,0,0);
136
137  // Tube with Local Magnetic field
138   
139  G4VSolid* magneticSolid = new G4Tubs("magneticTubs",0.,1.*m,1.*m,0.,360.*deg);
140  G4NistManager* man = G4NistManager::Instance();
141  G4Material* G4_Galactic = man->FindOrBuildMaterial("G4_Galactic");
142 
143   G4LogicalVolume* magneticLogical
144    = new G4LogicalVolume(magneticSolid,G4_Galactic,"magneticLogical",fieldMgr,0,0);
145  //                                                                  ********
146 
147  // placement of Tube
148 
149  G4RotationMatrix* fieldRot = new G4RotationMatrix();
150  fieldRot->rotateX(90.*deg);
151  new G4PVPlacement(fieldRot,G4ThreeVector(),magneticLogical,
152                    "magneticPhysical",worldLogical,0,0);
153
154  // set "user limits" for drawing smooth curve
155  G4UserLimits* userLimits = new G4UserLimits(5.0*cm);
156  magneticLogical->SetUserLimits(userLimits);
157
158  // first arm
159  G4VSolid* firstArmSolid = new G4Box("firstArmBox",1.5*m,1.*m,3.*m);
160  G4LogicalVolume* firstArmLogical
161    = new G4LogicalVolume(firstArmSolid,air,"firstArmLogical",0,0,0);
162  new G4PVPlacement(0,G4ThreeVector(0.,0.,-5.*m),firstArmLogical,
163                    "firstArmPhysical",worldLogical,0,0);
164
165  // second arm
166  G4VSolid* secondArmSolid = new G4Box("secondArmBox",2.*m,2.*m,3.5*m);
167  G4LogicalVolume* secondArmLogical
168    = new G4LogicalVolume(secondArmSolid,air,"secondArmLogical",0,0,0);
169  G4double x = -5.*m * std::sin(armAngle);
170  G4double z = 5.*m * std::cos(armAngle);
171  secondArmPhys
172    = new G4PVPlacement(armRotation,G4ThreeVector(x,0.,z),secondArmLogical,
173                        "secondArmPhys",worldLogical,0,0);
174
175  // hodoscopes in first arm
176  G4VSolid* hodoscope1Solid = new G4Box("hodoscope1Box",5.*cm,20.*cm,0.5*cm);
177  G4LogicalVolume* hodoscope1Logical
178    = new G4LogicalVolume(hodoscope1Solid,scintillator,"hodoscope1Logical",0,0,0);
179  for(int i1=0;i1<15;i1++)
180  {
181    G4double x1 = (i1-7)*10.*cm;
182    new G4PVPlacement(0,G4ThreeVector(x1,0.,-1.5*m),hodoscope1Logical,
183                      "hodoscope1Physical",firstArmLogical,0,i1);
184  }
185
186  // drift chambers in first arm
187  G4VSolid* chamber1Solid = new G4Box("chamber1Box",1.*m,30.*cm,1.*cm);
188  G4LogicalVolume* chamber1Logical
189    = new G4LogicalVolume(chamber1Solid,argonGas,"chamber1Logical",0,0,0);
190  for(int j1=0;j1<5;j1++)
191  {
192    G4double z1 = (j1-2)*0.5*m;
193    new G4PVPlacement(0,G4ThreeVector(0.,0.,z1),chamber1Logical,
194                      "chamber1Physical",firstArmLogical,0,j1);
195  }
196
197  // "virtual" wire plane
198  G4VSolid* wirePlane1Solid = new G4Box("wirePlane1Box",1.*m,30.*cm,0.1*mm);
199  G4LogicalVolume* wirePlane1Logical
200    = new G4LogicalVolume(wirePlane1Solid,argonGas,"wirePlane1Logical",0,0,0);
201  new G4PVPlacement(0,G4ThreeVector(0.,0.,0.),wirePlane1Logical,
202                    "wirePlane1Physical",chamber1Logical,0,0);
203
204  // hodoscopes in second arm
205  G4VSolid* hodoscope2Solid = new G4Box("hodoscope2Box",5.*cm,20.*cm,0.5*cm);
206  G4LogicalVolume* hodoscope2Logical
207    = new G4LogicalVolume(hodoscope2Solid,scintillator,"hodoscope2Logical",0,0,0);
208  for(int i2=0;i2<25;i2++)
209  {
210    G4double x2 = (i2-12)*10.*cm;
211    new G4PVPlacement(0,G4ThreeVector(x2,0.,0.),hodoscope2Logical,
212                      "hodoscope2Physical",secondArmLogical,0,i2);
213  }
214
215  // drift chambers in second arm
216  G4VSolid* chamber2Solid = new G4Box("chamber2Box",1.5*m,30.*cm,1.*cm);
217  G4LogicalVolume* chamber2Logical
218    = new G4LogicalVolume(chamber2Solid,argonGas,"chamber2Logical",0,0,0);
219  for(int j2=0;j2<5;j2++)
220  {
221    G4double z2 = (j2-2)*0.5*m - 1.5*m;
222    new G4PVPlacement(0,G4ThreeVector(0.,0.,z2),chamber2Logical,
223                      "chamber2Physical",secondArmLogical,0,j2);
224  }
225
226  // "virtual" wire plane
227  G4VSolid* wirePlane2Solid = new G4Box("wirePlane2Box",1.5*m,30.*cm,0.1*mm);
228  G4LogicalVolume* wirePlane2Logical
229    = new G4LogicalVolume(wirePlane2Solid,argonGas,"wirePlane2Logical",0,0,0);
230  new G4PVPlacement(0,G4ThreeVector(0.,0.,0.),wirePlane2Logical,
231                    "wirePlane2Physical",chamber2Logical,0,0);
232
233  // CsI calorimeter
234  G4VSolid* EMcalorimeterSolid = new G4Box("EMcalorimeterBox",1.5*m,30.*cm,15.*cm);
235  G4LogicalVolume* EMcalorimeterLogical
236    = new G4LogicalVolume(EMcalorimeterSolid,CsI,"EMcalorimeterLogical",0,0,0);
237  new G4PVPlacement(0,G4ThreeVector(0.,0.,2.*m),EMcalorimeterLogical,
238                    "EMcalorimeterPhysical",secondArmLogical,0,0);
239
240  // EMcalorimeter cells
241  G4VSolid* cellSolid = new G4Box("cellBox",7.5*cm,7.5*cm,15.*cm);
242  G4LogicalVolume* cellLogical
243    = new G4LogicalVolume(cellSolid,CsI,"cellLogical",0,0,0);
244  G4VPVParameterisation* cellParam = new A01CellParameterisation();
245  new G4PVParameterised("cellPhysical",cellLogical,EMcalorimeterLogical,
246                         kXAxis,80,cellParam);
247
248//@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
249  // hadron calorimeter
250  G4VSolid* HadCalorimeterSolid
251    = new G4Box("HadCalorimeterBox",1.5*m,30.*cm,50.*cm);
252  G4LogicalVolume* HadCalorimeterLogical
253    = new G4LogicalVolume(HadCalorimeterSolid,lead,"HadCalorimeterLogical",0,0,0);
254  new G4PVPlacement(0,G4ThreeVector(0.,0.,3.*m),HadCalorimeterLogical,
255                    "HadCalorimeterPhysical",secondArmLogical,0,0);
256
257  // hadron calorimeter column
258  G4VSolid* HadCalColumnSolid
259    = new G4Box("HadCalColumnBox",15.*cm,30.*cm,50.*cm);
260  G4LogicalVolume* HadCalColumnLogical
261    = new G4LogicalVolume(HadCalColumnSolid,lead,"HadCalColumnLogical",0,0,0);
262  new G4PVReplica("HadCalColumnPhysical",HadCalColumnLogical,
263                   HadCalorimeterLogical,kXAxis,10,30.*cm);
264
265  // hadron calorimeter cell
266  G4VSolid* HadCalCellSolid
267    = new G4Box("HadCalCellBox",15.*cm,15.*cm,50.*cm);
268  G4LogicalVolume* HadCalCellLogical
269    = new G4LogicalVolume(HadCalCellSolid,lead,"HadCalCellLogical",0,0,0);
270  new G4PVReplica("HadCalCellPhysical",HadCalCellLogical,
271                   HadCalColumnLogical,kYAxis,2,30.*cm);
272
273  // hadron calorimeter layers
274  G4VSolid* HadCalLayerSolid
275    = new G4Box("HadCalLayerBox",15.*cm,15.*cm,2.5*cm);
276  G4LogicalVolume* HadCalLayerLogical
277    = new G4LogicalVolume(HadCalLayerSolid,lead,"HadCalLayerLogical",0,0,0);
278  new G4PVReplica("HadCalLayerPhysical",HadCalLayerLogical,
279                  HadCalCellLogical,kZAxis,20,5.*cm);
280
281  // scintillator plates
282  G4VSolid* HadCalScintiSolid
283    = new G4Box("HadCalScintiBox",15.*cm,15.*cm,0.5*cm);
284  G4LogicalVolume* HadCalScintiLogical
285    = new G4LogicalVolume(HadCalScintiSolid,scintillator,"HadCalScintiLogical",0,0,0);
286  new G4PVPlacement(0,G4ThreeVector(0.,0.,2.*cm),HadCalScintiLogical,
287                    "HadCalScintiPhysical",HadCalLayerLogical,0,0);
288//@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
289
290  // sensitive detectors -----------------------------------------------------
291  G4SDManager* SDman = G4SDManager::GetSDMpointer();
292  G4String SDname;
293
294  hodoscope1 = new A01Hodoscope(SDname="/hodoscope1");
295  SDman->AddNewDetector(hodoscope1);
296  hodoscope1Logical->SetSensitiveDetector(hodoscope1);
297  hodoscope2 = new A01Hodoscope(SDname="/hodoscope2");
298  SDman->AddNewDetector(hodoscope2);
299  hodoscope2Logical->SetSensitiveDetector(hodoscope2);
300
301  chamber1 = new A01DriftChamber(SDname="/chamber1");
302  SDman->AddNewDetector(chamber1);
303  wirePlane1Logical->SetSensitiveDetector(chamber1);
304  chamber2 = new A01DriftChamber(SDname="/chamber2");
305  SDman->AddNewDetector(chamber2);
306  wirePlane2Logical->SetSensitiveDetector(chamber2);
307
308  EMcalorimeter = new A01EmCalorimeter(SDname="/EMcalorimeter");
309  SDman->AddNewDetector(EMcalorimeter);
310  cellLogical->SetSensitiveDetector(EMcalorimeter);
311
312  HadCalorimeter = new A01HadCalorimeter(SDname="/HadCalorimeter");
313  SDman->AddNewDetector(HadCalorimeter);
314  HadCalScintiLogical->SetSensitiveDetector(HadCalorimeter);
315
316  // visualization attributes ------------------------------------------------
317
318  worldVisAtt = new G4VisAttributes(G4Colour(1.0,1.0,1.0));
319  worldVisAtt->SetVisibility(false);
320  worldLogical->SetVisAttributes(worldVisAtt);
321
322  magneticVisAtt = new G4VisAttributes(G4Colour(0.9,0.9,0.9));   // LightGray
323  magneticLogical->SetVisAttributes(magneticVisAtt);
324
325  armVisAtt = new G4VisAttributes(G4Colour(1.0,1.0,1.0));
326  armVisAtt->SetVisibility(false);
327  firstArmLogical->SetVisAttributes(armVisAtt);
328  secondArmLogical->SetVisAttributes(armVisAtt);
329
330  hodoscopeVisAtt = new G4VisAttributes(G4Colour(0.8888,0.0,0.0));
331  hodoscope1Logical->SetVisAttributes(hodoscopeVisAtt);
332  hodoscope2Logical->SetVisAttributes(hodoscopeVisAtt);
333
334  chamberVisAtt = new G4VisAttributes(G4Colour(0.0,1.0,0.0));
335  chamber1Logical->SetVisAttributes(chamberVisAtt);
336  chamber2Logical->SetVisAttributes(chamberVisAtt);
337
338  wirePlaneVisAtt = new G4VisAttributes(G4Colour(0.0,0.8888,0.0));
339  wirePlaneVisAtt->SetVisibility(false);
340  wirePlane1Logical->SetVisAttributes(wirePlaneVisAtt);
341  wirePlane2Logical->SetVisAttributes(wirePlaneVisAtt);
342
343  EMcalorimeterVisAtt = new G4VisAttributes(G4Colour(0.8888,0.8888,0.0));
344  EMcalorimeterVisAtt->SetVisibility(false);
345  EMcalorimeterLogical->SetVisAttributes(EMcalorimeterVisAtt);
346
347  cellVisAtt = new G4VisAttributes(G4Colour(0.9,0.9,0.0));
348  cellLogical->SetVisAttributes(cellVisAtt);
349
350  HadCalorimeterVisAtt = new G4VisAttributes(G4Colour(0.0, 0.0, 0.9));
351  HadCalorimeterLogical->SetVisAttributes(HadCalorimeterVisAtt);
352  HadCalorimeterCellVisAtt = new G4VisAttributes(G4Colour(0.0, 0.0, 0.9));
353  HadCalorimeterCellVisAtt->SetVisibility(false);
354  HadCalColumnLogical->SetVisAttributes(HadCalorimeterCellVisAtt);
355  HadCalCellLogical->SetVisAttributes(HadCalorimeterCellVisAtt);
356  HadCalLayerLogical->SetVisAttributes(HadCalorimeterCellVisAtt);
357  HadCalScintiLogical->SetVisAttributes(HadCalorimeterCellVisAtt);
358
359  // return the world physical volume ----------------------------------------
360
361  G4cout << G4endl << "The geometrical tree defined are : " << G4endl << G4endl;
362  DumpGeometricalTree(worldPhysical);
363
364  return worldPhysical;
365}
366
367void A01DetectorConstruction::ConstructMaterials()
368{
369  G4double a;
370  G4double z;
371  G4double density;
372  G4double weightRatio;
373  G4String name;
374  G4String symbol;
375  G4int nElem;
376
377  // Argon gas
378  a = 39.95*g/mole;
379  density = 1.782e-03*g/cm3;
380  argonGas = new G4Material(name="ArgonGas", z=18., a, density);
381
382  // elements for mixtures and compounds
383  a = 1.01*g/mole;
384  G4Element* elH = new G4Element(name="Hydrogen", symbol="H", z=1., a);
385  a = 12.01*g/mole;
386  G4Element* elC = new G4Element(name="Carbon", symbol="C", z=6., a);
387  a = 14.01*g/mole;
388  G4Element* elN = new G4Element(name="Nitrogen", symbol="N", z=7., a);
389  a = 16.00*g/mole;
390  G4Element* elO = new G4Element(name="Oxigen", symbol="O", z=8., a);
391  a = 126.9*g/mole;
392  G4Element* elI = new G4Element(name="Iodine", symbol="I", z=53., a);
393  a = 132.9*g/mole;
394  G4Element* elCs= new G4Element(name="Cesium", symbol="Cs", z=55., a);
395
396  // Air
397  density = 1.29*mg/cm3;
398  air = new G4Material(name="Air", density, nElem=2);
399  air->AddElement(elN, weightRatio=.7);
400  air->AddElement(elO, weightRatio=.3);
401
402  // Scintillator
403  density = 1.032*g/cm3;
404  scintillator = new G4Material(name="Scintillator", density, nElem=2);
405  scintillator->AddElement(elC, 9);
406  scintillator->AddElement(elH, 10);
407
408  // CsI
409  density = 4.51*g/cm3;
410  CsI = new G4Material(name="CsI", density, nElem=2);
411  CsI->AddElement(elI, weightRatio=.5);
412  CsI->AddElement(elCs,weightRatio=.5);
413
414  // Lead
415  a = 207.19*g/mole;
416  density = 11.35*g/cm3;
417  lead = new G4Material(name="Lead", z=82., a, density);
418
419  G4cout << G4endl << "The materials defined are : " << G4endl << G4endl;
420  G4cout << *(G4Material::GetMaterialTable()) << G4endl;
421}
422
423void A01DetectorConstruction::DestroyMaterials()
424{
425  // Destroy all allocated elements and materials
426  size_t i;
427  G4MaterialTable* matTable = (G4MaterialTable*)G4Material::GetMaterialTable();
428  for(i=0;i<matTable->size();i++)
429  { delete (*(matTable))[i]; }
430  matTable->clear();
431  G4ElementTable* elemTable = (G4ElementTable*)G4Element::GetElementTable();
432  for(i=0;i<elemTable->size();i++)
433  { delete (*(elemTable))[i]; }
434  elemTable->clear();
435}
436
437void A01DetectorConstruction::DumpGeometricalTree(G4VPhysicalVolume* aVolume,G4int depth)
438{
439  for(int isp=0;isp<depth;isp++)
440  { G4cout << "  "; }
441  G4cout << aVolume->GetName() << "[" << aVolume->GetCopyNo() << "] "
442         << aVolume->GetLogicalVolume()->GetName() << " "
443         << aVolume->GetLogicalVolume()->GetNoDaughters() << " "
444         << aVolume->GetLogicalVolume()->GetMaterial()->GetName();
445  if(aVolume->GetLogicalVolume()->GetSensitiveDetector())
446  {
447    G4cout << " " << aVolume->GetLogicalVolume()->GetSensitiveDetector()
448                            ->GetFullPathName();
449  }
450  G4cout << G4endl;
451  for(int i=0;i<aVolume->GetLogicalVolume()->GetNoDaughters();i++)
452  { DumpGeometricalTree(aVolume->GetLogicalVolume()->GetDaughter(i),depth+1); }
453}
454
455void A01DetectorConstruction::SetArmAngle(G4double val)
456{
457  if(!secondArmPhys)
458  {
459    G4cerr << "Detector has not yet been constructed." << G4endl;
460    return;
461  }
462
463  armAngle = val;
464  *armRotation = G4RotationMatrix();  // make it unit vector
465  armRotation->rotateY(armAngle);
466  G4double x = -5.*m * std::sin(armAngle);
467  G4double z = 5.*m * std::cos(armAngle);
468  secondArmPhys->SetTranslation(G4ThreeVector(x,0.,z));
469
470  // tell G4RunManager that we change the geometry
471  G4RunManager::GetRunManager()->GeometryHasBeenModified();
472}
Note: See TracBrowser for help on using the repository browser.